1. Developmental Biology
Download icon

m6A RNA methylation impacts fate choices during skin morphogenesis

  1. Linghe Xi
  2. Thomas Carroll
  3. Irina Matos
  4. Ji-Dung Luo
  5. Lisa Polak
  6. H Amalia Pasolli
  7. Samie R Jaffrey
  8. Elaine Fuchs  Is a corresponding author
  1. Howard Hughes Medical Institute, Robin Chemers Neustein Laboratory of Mammalian Cell Biology and Development, The Rockefeller University, United States
  2. Bioinformatics Resource Center, The Rockefeller University, United States
  3. Electron Microscopy Resource Center, The Rockefeller University, United States
  4. Department of Pharmacology, Weill Cornell Medicine, Cornell University, United States
Research Article
  • Cited 5
  • Views 2,734
  • Annotations
Cite this article as: eLife 2020;9:e56980 doi: 10.7554/eLife.56980

Abstract

N6-methyladenosine is the most prominent RNA modification in mammals. Here, we study mouse skin embryogenesis to tackle m6A’s functions and physiological importance. We first landscape the m6A modifications on skin epithelial progenitor mRNAs. Contrasting with in vivo ribosomal profiling, we unearth a correlation between m6A modification in coding sequences and enhanced translation, particularly of key morphogenetic signaling pathways. Tapping physiological relevance, we show that m6A loss profoundly alters these cues and perturbs cellular fate choices and tissue architecture in all skin lineages. By single-cell transcriptomics and bioinformatics, both signaling and canonical translation pathways show significant downregulation after m6A loss. Interestingly, however, many highly m6A-modified mRNAs are markedly upregulated upon m6A loss, and they encode RNA-methylation, RNA-processing and RNA-metabolism factors. Together, our findings suggest that m6A functions to enhance translation of key morphogenetic regulators, while also destabilizing sentinel mRNAs that are primed to activate rescue pathways when m6A levels drop.

Introduction

Sophisticated gene expression regulatory machineries are necessary to achieve proper fate choices during mammalian embryogenesis. While transcriptional regulation has been studied in great depth, much less is known about how post-transcriptional regulation comes into play in tissue formation. Being an important element of post-transcriptional regulation of gene expression, chemical modifications on RNAs have been demonstrated to affect a wide range of RNA bio-activities. Among the over 100 types of RNA modifications identified thus far, m6A is the most abundant modification on mRNAs. Previous studies have identified m6A methyltransferases (writers), m6A demethylases (erasers) and also factors that recognize m6A-modified RNAs (readers) (Alarcón et al., 2015; Bokar et al., 1997; Liu et al., 2014; Meyer et al., 2015; Meyer and Jaffrey, 2017; Patil et al., 2016; Zheng et al., 2013Figure 1A).

Figure 1 with 1 supplement see all
miCLIP and ribosomal profiling analyses of the mouse skin epithelial progenitors.

(A) Schematic depicting the major factors involved in regulating the cellular dynamics of m6A modification. (B) Schematic depicting embryonic development of mammalian epithelial skin progenitors. WNThi/lo implies cells that show strong WNT or low WNT signaling as judged by Axin2-LacZ transgene expression (Matos et al., 2020). HF morphogenesis occurs in temporal waves, with mature HFs emerging shortly after birth. (C) Schematic depicting the process for enzymatically isolating epithelial progenitors from mouse skin at age P0. For miCLIP, cells were subjected to FACS purification as described in the Materials and methods. (D) Consensus sequence motif enriched around miCLIP-identified m6A sites in P0 skin progenitor mRNAs. (E) Metagene plots depicting the distribution of miCLIP-identified m6A sites along mRNAs. Data from three independent replicates are shown. (F) Schematic depicts comparison of the miCLIP data, which measures m6A modification to the ribosome profiling data, which landscapes bound ribosomes on neonatal skin progenitor mRNAs. The empirical cumulative distribution function (ECDF) plots compare the relative mRNA translation efficiency of the top 20% and bottom 20% of m6A-modified mRNAs. The data reveal that transcripts with higher levels of m6A modification (assessed by the sum of normalized-to-input uTPM value of m6A along the full-length transcript) tend to have higher levels of translation efficiency. The correlation between translation efficiency and the sum of normalized-to-input uTPM value of m6A at different regions of the mRNAs (5’ UTR, coding sequence, 3’ UTR) shows that the coding sequence m6A gives the best correlation to translation efficiency. The p values were calculated through Wilcoxon rank sum test. (G) GSEA of the overlap between mRNAs that are the top 20% heavily m6A-modified in coding sequence (assessed by the sum of normalized-to-input uTPM value of m6A along coding sequence) and the top 20% most efficiently translated mRNAs (assessed by ribosome profiling). Shown are the top eight enriched KEGG signaling pathways, each of which has a p value <0.05 and >5 enriched mRNAs. (H) Pathways known to play essential roles in regulating skin lineage specification.

A number of studies have focused on dissecting m6A’s functions in cell culture systems, where transcriptome-wide m6A profiles and expression patterns of the writers, erasers and readers were found to vary among cell types (An et al., 2020; Delaunay and Frye, 2019; Roundtree et al., 2017). In recent years, models have been made to conditionally disrupt specific m6A writers, erasers or readers in mouse tissues (Cheng et al., 2019; Geula et al., 2015; Hsu et al., 2017; Ivanova et al., 2017; Lee et al., 2019; Li et al., 2017; Li et al., 2018; Lin et al., 2017; Shi et al., 2018; Wu et al., 2018; Xu et al., 2017; Yoon et al., 2017; Zhang et al., 2018). While uniformly underscoring the physiological importance of m6A, these studies have also unveiled a profound complexity in m6A’s biological functions (Du et al., 2016; Kennedy et al., 2016; Meyer et al., 2015; Wang et al., 2014a; Wang et al., 2015; Zaccara and Jaffrey, 2020). Thus in some cell types and contexts, m6A modification promotes differentiation (Geula et al., 2015; Lee et al., 2019), while in others, it blocks differentiation (Li et al., 2017; Vu et al., 2017). Sometimes, as in the case of Myc, such opposing effects can be observed in different cell types sharing a common m6A target (Lee et al., 2019; Vu et al., 2017; Yao et al., 2018; Zou et al., 2019). In other cases, m6A targets can differ based upon cellular contexts and expression. Finally, m6A can exert context-dependent effects on distinct steps of RNA metabolism, compounding the unpredictability of whether a modification will impart to a given target increased or decreased expression.

The skin epithelium is an excellent model to begin to unravel the mechanisms underlying m6A’s influences on various morphogenetic processes. Embryonic skin epithelium begins as a single layer of multipotent epithelial progenitors, which will develop into three strikingly distinct tissues: epidermis, hair follicles (HFs) and sebaceous glands (Figure 1B). Launching tissue diversification is WNT signaling, which begins heterogeneously within the progenitor population.

WNTs play a critical role in specifying the HF fate (Andl et al., 2002; Gat et al., 1998; van Genderen et al., 1994; Huelsken et al., 2001Xu et al., 2015). Within the basal epidermal plane, WNThi embryonic progenitors cluster into hair placodes (Ahtiainen et al., 2016), which emerge in three waves beginning at embryonic day 14.5 (E14.5) and develop into fully mature HFs by birth (Duverger and Morasso, 2009; Figure 1B). The first divisions in WNThi placode cells are asymmetric, generating WNTlo daughters that will give rise to hair follicle stem cells (HFSCs) and outer root sheath, and WNThi daughters that will generate the inner root sheath and hair shaft of the follicles (Ouspenskaia et al., 2016). Near birth, oil-rich sebaceous glands spawn from the upper portion of each HF (Figure 1B). Sebaceous gland development is favored over HFs in mice whose skin expresses ΔNLEF1, a dominant-negative disrupter of WNT signaling (Merrill et al., 2001; Niemann et al., 2003), indicating that WNT levels play an important role in balancing these fate choices.

Lower levels of WNT signaling also favor epidermal over HF fates, as in the developing plane of embryonic progenitors, those with lower levels of WNT signaling become fated to become epidermal progenitors, which fuel the production of upward columns of terminally differentiated cells that form the skin’s barrier that excludes pathogens and retains body fluids. This differentiation program is regulated by NOTCH signaling and MYC activation (Blanpain et al., 2006; Frye et al., 2003; Gandarillas and Watt, 1997; Moriyama et al., 2008; Rangarajan et al., 2001; Watt et al., 2008a; Watt et al., 2008b).

While the levels of external signals and their downstream transcriptional effectors function critically in these epithelial fate choices within the skin, post-transcriptional regulation, such as translational control, has received more emphasis on balancing proper homeostasis once the tissue fate has been selected. Similar to the hematopoietic system (Buszczak et al., 2014; Signer et al., 2014), skin progenitors maintain reduced levels of protein synthesis relative to their differentiating progeny (Blanco et al., 2016; Liakath-Ali et al., 2018; Sendoel et al., 2017). In differentiating skin cells, 5-methylcytosine (m5C) is enhanced, which protects tRNAs from cleavage and facilitates protein translation (Blanco et al., 2011). Under conditions of stress in the epidermis, EIF2α is phosphorylated, dampening protein translation, while a less efficient, more promiscuous initiator, EIF2A, translates key proteins that allow progenitor survival in a harsher environment (Sendoel et al., 2017). Interestingly, however, when translational control is compromised in the skin, epidermal, HF and sebaceous gland stem cell compartments respond differently, with epidermal cells implementing a feedback mechanism to increase global translation (Liakath-Ali et al., 2018). While context again surfaces as being key in triggering these varied responses, their purpose appears to be rooted in restoring homeostasis to the tissue.

How signaling, transcriptional regulation and cellular fate choices are linked to translational control and tissue homeostasis has remained elusive. Here, we show that m6A modifications of mRNAs may be involved in this connection. We began by mapping m6A at single-nucleotide resolution in skin epithelia. Exploiting our prior in vivo epidermal ribosomal profiling data (Sendoel et al., 2017), we then interrogated how m6A correlates with overall translational efficiency. Probing the physiological relevance of our findings, we employed genetics to conditionally ablate the deposition of m6A on RNAs. After analyzing the striking and differential consequences of m6A loss to the fates of the skin epithelial progenitors in vivo, we turned to single-cell RNA sequencing and functional studies to gain insights into the major cellular activities that are targeted by this regulatory machinery, and the significance of m6A to these three diverse programs of epithelial morphogenesis in the skin.

Results

Skin transcripts most highly modified by m6A are involved in HF morphogenesis

To assess the location and extent of m6A modification on the mRNAs expressed in mouse epidermal cells under physiological conditions, we performed m6A individual-nucleotide-resolution cross-linking and immunoprecipitation (miCLIP) (Grozhik et al., 2017; Linder et al., 2015) on poly(A)+ RNAs directly isolated from basal skin progenitors in vivo (Figure 1C). Adapting a procedure previously used for in vivo ribosomal profiling of these cells (Sendoel et al., 2017), we applied dispase treatment on newborn (postnatal day 0, P0) skin, which effectively removes the dermis, elongated hair pegs and HFs (Rhee et al., 2006; Figure 1—figure supplement 1A). This enabled us to apply fluorescence-activated cell sorting (FACS) for integrin α6 (CD49f) to enrich for progenitors of epidermis and hair placodes/germs (Figure 1C; Figure 1—figure supplement 1B). Sequencing libraries were then constructed from the m6A-immunoprecipitated mRNAs as well as the corresponding input RNA samples (Figure 1—figure supplement 1C).

Based upon the crosslinking-induced mutations detected in the miCLIP libraries, 11,420 transcripts displayed a total of 157,641 sites that were modified by m6A and typified by a ‘DRACH’ motif (Fu et al., 2014; Linder et al., 2015Figure 1D; Supplementary file 1). The data were highly reproducible across three independent replicates, and as previously noted, m6A density was highest around mRNA stop codons (Figure 1E). To quantify the extent of m6A modification at each site, we calculated the unique tag counts per million (uTPM) as described in Patil et al., 2016 around each site in the miCLIP data. We then normalized to the uTPM around the same site in the corresponding input data. This normalized-to-input uTPM value was also highly consistent across biological replicates (Figure 1—figure supplement 1D).

Exploiting our prior in vivo ribosome profiling data of mouse skin epidermal progenitors (Sendoel et al., 2017), we next examined in an unbiased fashion whether the degree of m6A modification of an mRNA correlated with its translational efficiency (Figure 1F, schematic). For every m6A-modified mRNA, we calculated the sum of normalized-to-input uTPM at each m6A site identified on full-length, 5’ UTR, coding sequence, and 3’ UTR. For each, we then compared the translational efficiencies of the top 20% most highly and bottom 20% least m6A-modified mRNAs from each way of calculation. Intriguingly, mRNAs with high levels of m6A modification tended to be more highly translated than those with low modification (Figure 1F). This correlation was most striking when we focused on m6A modifications in the coding sequence.

We then ranked the m6A-modified mRNAs according to the sum of normalized-to-input uTPM at each m6A site within the coding sequence (Supplementary file 2). Gene set enrichment analyses (GSEA) upon that value featured ‘basal cell carcinoma’, dominated by WNT and SHH signaling factors, and ‘NOTCH signaling’ among the top three categories (Figure 1—figure supplement 1E,F, Supplementary file 2). Interestingly, these top three categories tended to be more efficiently translated than other pathways (Figure 1—figure supplement 1G). Overall, these findings pointed to a potential importance of m6A modifications specifically within the coding sequence in promoting translation of an mRNA.

We then performed GSEA to interrogate the KEGG pathways enriched with the top 20% most highly modified and highly translated mRNAs. We focused on pathways with p values <0.05 and >5 enriched mRNAs/category. Although we do not discount the potential importance of any individual mRNA, it was easier to predict this potential for mRNAs with high coding sequence m6A levels, high translation and residing within larger enrichment categories.

The top categories sharing these criteria were signaling pathways including ECM-receptor interactions, basal cell carcinoma (dominated by mRNAs encoding WNT and SHH signaling factors) and NOTCH signaling (Figure 1G; Supplementary file 2). These pathways play essential roles in HF morphogenesis (Figure 1H). Notably, WNT signaling, which prompts basal progenitors to organize into hair placodes, precedes SHH signaling (Woo et al., 2012), and WNT signaling is required to maintain SHH pathway-driven basal cell carcinomas in mice and in humans (Sánchez-Danés et al., 2018).

Given the large number of mRNAs that were modified by m6A, we did not expect nor did we see a single category with an extraordinarily high enrichment score. However, the fact that the WNT and SHH pathways surfaced in multiple top enrichment categories was promising and suggested that m6A might function in promoting the HF fate.

Conditional ablation of Mettl3 in epidermal progenitors results in a marked defect in HF morphogenesis

The METTL3-METTL14 writer complex is responsible for adding the vast majority of m6A onto RNAs (Figure 1A). Upon ablation of Mettl3 or Mettl14, global m6A is dramatically diminished in a wide variety of cell types (Batista et al., 2014; Cheng et al., 2019; Cui et al., 2017; Geula et al., 2015; Lin et al., 2017; Liu et al., 2014; Vu et al., 2017; Wang et al., 2014b; Weng et al., 2018; Yao et al., 2018; Yoon et al., 2017). Thus, to investigate the impact of m6A deficiency on skin morphogenesis, we generated Mettl3 conditional knockout (cKO) mice harboring the following alleles: Mettl3fl/fl, Rosa26-YFPfl/+, Krt14-Cre (Figure 2—figure supplement 1). Active exclusively in the epidermal progenitors beginning at E13.5 (Vasioukhin et al., 1999), Krt14-Cre had excised the stop codon in the Rosa26 locus by E14.5 (Figure 2A). Whole-mount immunofluorescence imaging revealed that by E16.5, expression of METTL3 was largely abolished within the targeted YFP+ cells in the cKO skin epithelium (Figure 2B). And consistent with its known location in other cell types (Kwon et al., 2019; Liu et al., 2014; Schöller et al., 2018; Zhang et al., 2018), METTL3 was predominantly nuclear in control (Ctrl) epithelial progenitors.

Figure 2 with 2 supplements see all
Krt14-Cre driven conditional Mettl3 knockout mice display severe defects in HF morphogenesis.

(A) Representative pictures of Krt14-Cre-/-, Rosa26-YFPfl/+ and Krt14-Cre+/-, Rosa26-YFPfl/+ littermate embryos demonstrating the onset of uniform YFP expression in the E14.5 skin epithelium (scale bars: 2 mm). (B) Confocal images of E16.5 whole-mount back skin immunolabeled for P-cadherin (PCAD), METTL3 and YFP (scale bars: 20 µm). Note that nuclear METTL3 immunofluorescence is selectively depleted from the YFP+ cells in cKO skin. White dashed lines denote the dermal-epidermal border. (C) Left panel: representative pictures of thin layer chromatography (TLC) on Poly(A)+ RNA samples isolated from E16.5 skin epithelial cells. Right panel: quantification of m6A levels based on TLC (error bars: standard deviation, for Mettl3+/+ n = 2 biological replicates, for each of the other conditions n = 3 biological replicates, **p<0.01 by unpaired two-tailed Student’s t-test). (D) Representative pictures of E16.5, P0 control and cKO littermates (scale bars: 1 cm). (E) Hematoxylin and eosin (H and E) stained back skin sagittal sections at indicated time points (scale bars: 100 µm). (F) Confocal images of whole-mount back skin immunolabeled for PCAD at indicated time points (scale bars: 100 µm) and quantification of HFs at different developmental stages (for E16.5, n = 3 biological replicates ×10 images per replicate and for P0, n = 18 images from four biological replicates, **p<0.01 and ***p<0.001 by unpaired two-tailed Student’s t-test).

We confirmed the quantitative loss of m6A by thin-layer chromatography (TLC), which we performed on poly(A)+ RNA isolated from the FACS-purified YFP+ epidermal cells of Mettl3 wild-type, heterozygous and null skin epithelium (Figure 2C). Quantification of the m6A/A ratio based on the radioactive signals revealed a comparable m6A content on mRNAs from wild-type Mettl3+/+ and heterozygous Mettl3+/- cells, indicating that one allele of Mettl3 was enough to sustain a normal m6A level on mRNAs. Thus from this point forward, we used heterozygous Mettl3 fl/+, Rosa26-YFPfl/+, Krt14-Cre mice as our control.

Throughout embryogenesis, control and cKO mice were grossly comparable in appearance (Figure 2D). From postnatal day 2 (P2) onward, however, phenotypic abnormalities of the cKO animals began to emerge. Most striking was the progressive decline in body sizes and weights compared to their control littermates (Figure 2—figure supplement 2A,B). Most cKO animals did not survive beyond postnatal day 6 (P6).

To investigate potential causes of the lethality, we assessed the epidermal barrier function by measuring the trans-epidermal water loss (TEWL) rate of the control and cKO skin at P0 and P6. No significant TEWL differences were observed, ruling against dehydration as a cause of death (Figure 2—figure supplement 2C). Examination of oral tissues indicated that both teeth and tongue formed, but tongue filiform papillae, necessary for suckling and food intake at this stage, were notably diminished (Figure 2—figure supplement 2D,E). Back skin HF morphogenesis was also severely altered in the cKO animals. This began to surface soon after Mettl3 ablation where the numbers of de novo placodes began to decline, and progressed thereafter as judged by the paucity of hair germs and pegs at birth (Figure 2E,F). Analyses of the few pups still alive by P6 revealed that with the exception of the sparse, large guard hairs, which form prior to Krt14-Cre activity and do not use WNT signaling in their specification, HF downgrowth was largely impaired and HFs failed to mature beyond the peg stage (Figure 2—figure supplement 2F).

Based upon these data, the defects did not appear to reflect developmental delays, but rather altered morphogenesis. In this regard, the severity of defects in the tongue filiform papillae and skin HFs was particularly intriguing as both require WNT and SHH signaling for their morphogenesis (DasGupta and Fuchs, 1999; Iwatsuki et al., 2007; Järvinen et al., 2006). At the top of these morphogenetic cascades is LEF1, the transcription factor which binds to WNT effector β-catenin and translocates to the nucleus, where it is required to mediate WNT signaling and launch cellular fates (Adam et al., 2018; van Genderen et al., 1994; Ouspenskaia et al., 2016; Zhou et al., 1995).

Immunofluorescence revealed that at E17.5, nuclear LEF1 was significantly diminished in the basal (placode) progenitors within the epidermal plane as well as WNThi cells in the developing HFs of cKO skin (Figure 3A). At P0, LEF1 levels remained low in the WNThi cells. It was also diminished in the underlying dermal condensates (precursors to the dermal papilla, DP), whose WNT signaling and nuclear LEF1 is known to be dependent upon epithelial WNT signaling (Mok et al., 2019; Figure 3B). Notably, immunohistochemistry revealed that in contrast to control HFs, which showed appreciable anti-β-catenin nuclear staining at the follicle:DP interface, in Mettl3 null HFs, β-catenin was mostly at the intercellular borders, reflective of its WNT-independent role at adherens junctions (Figure 3C). Additionally, these HF progenitors were perturbed in engulfing the DP, a feature seen when WNT signaling and/or SHH signaling are inhibited (Heitman et al., 2020; Matos et al., 2020).

Figure 3 with 1 supplement see all
Loss of m6A results in diminished WNT signaling and signs of perturbed HF fate.

(A) Left panel: representative images from E17.5 sagittal sections immunolabeled for LEF1 and PCAD (scale bars: 50 µm). White solid lines denote skin surface and dashed lines denote dermal-epidermal border. Right panel: quantification of LEF1 immunofluorescence at indicated compartments (for each condition n = 3 biological replicates ×7 images per replicate, *p<0.05 and ***p<0.001 by unpaired two-tailed Student’s t-test). (B) Left panel: representative images from P0 sagittal sections immunolabeled for LEF1 and PCAD (scale bars: 50 µm). Solid and dashed lines as in (A). Right panel: quantifications of LEF1 immunofluorescence in indicated compartments (for each condition n = 3 biological replicates ×7 images per replicate, *p<0.05 and ***p<0.001 by unpaired two-tailed Student’s t-test). (C) Representative images from P0 sagittal sections with immunohistochemistry staining of β-catenin, counter-stained with hematoxylin (scale bars: 50 µm). (D) Representative images from P0 sagittal sections stained with X-gal and nuclear fast red to examine the expression of the Gli1-LacZ transgene, a proxy for SHH signaling (scale bars: 100 µm). (E) Representative images from P6 sagittal sections stained with Oil Red O and counterstained with hematoxylin to visualize signs of HF to sebocyte fate switching within the epithelium (scale bars: 100 µm). The staining in the control dermis reflects adipocyte-derived lipids, missing in the cKO pups, which lose weight after birth. (F) Representative pictures of control (Ctrl) and cKO back skins engrafted onto Nude (Nu/Nu) mice and analyzed 15 days later. (G) Representative Hematoxylin and eosin (H and E) stained sagittal sections of 15-dayengrafted back skins (scale bars: 100 µm). (H) Representative sagittal sections as in (G) and counterstained with Oil Red O and hematoxylin (scale bars: 100 µm).

Figure 3—source data 1

Quantification of LEF1 immunofluorescence signals at E17.5 in (A).

https://cdn.elifesciences.org/articles/56980/elife-56980-fig3-data1-v2.xlsx
Figure 3—source data 2

Quantification of LEF1 immunofluorescence signals at P0 in (B).

https://cdn.elifesciences.org/articles/56980/elife-56980-fig3-data2-v2.xlsx

Since SHH signaling is dependent upon WNT signaling, we next tested for the activation of this pathway. To do so, we mated our mice to Gli1-lacZ mice, and performed X-gal histochemical staining on skin sections. While SHH signaling was not blocked in the absence of METTL3, it was aberrant, with an overall reduction, particularly in the epithelium (Figure 3D). These results were intriguing in light of the knowledge that laminin 5-1-1 is required for HF downgrowth and SHH signaling (Fleger-Weckmann et al., 2016; Gao et al., 2008), and that Lama5 and Lamac1 were among the most highly modified and efficiently translated mRNAs in skin (Figure 1G).

As an appendage of HFs, sebaceous glands use a different hedgehog pathway and require reduced WNT signaling for their formation (Merrill et al., 2001; Niemann et al., 2003). Consistent with diminished WNT signaling in the Mettl3 cKO skin, an increase in Oil Red O staining was seen in the aberrant HFs (Figure 3E). The paucity of lipids in the P6 dermis was expected due to the failure to suckle and weight loss likely caused by the tongue defects. These phenotypes bore a striking resemblance to the mice expressing a dominant negative form of LEF1 in their skin epithelium.

To interrogate the long-term consequences of epidermal m6A depletion, we engrafted the back skin of P0 control and cKO littermates onto Nude mice. By the 15th day post-engraftment, the epithelium of both control and cKO skins were YFP+, indicating that they had survived the engraftment procedure (Figure 3—figure supplement 1). However, whereas control skin had generated a thick hair coat, the cKO skin showed a striking paucity of hair growth (Figure 3F,G). Oil Red O staining further corroborated our P6 analyses and suggested that cKO progenitors follow the sebocyte versus HF path of morphogenesis (Figure 3H).

Many of these features were confirmed at the ultrastructural level, including the arrested growth of HFs and the presence of Mettl3 null HFs that failed to engulf the DP (Figure 4A). Although we did not see gross signs of junctional defects, we did see morphological signs of apoptosis in P0 HFs, which we substantiated by activated (cleaved) Caspase-3 and TUNEL staining (Figure 4—figure supplement 1A), suggesting that the defects went beyond mere signaling perturbations. We also noted Mettl3 null HF cells that were positive for TUNEL, but not activated Caspase-3. While this is a feature of sebocytes (Fischer et al., 2017), it could also be a sign of general cell degradation. Moreover, the organization of HF cells also appeared to be perturbed. In the hair placodes that did initiate, spindle orientations appeared to be largely unaltered and remained perpendicular to the underlying basement membrane (Figure 4—figure supplement 1B). However, as judged by phalloidin staining to visualize F-actin, apical cytoskeletal polarization of progenitors facing the DP appeared to be perturbed (Figure 4B).

Figure 4 with 1 supplement see all
Perturbations in WNT-driven dermal papilla engulfment and in actin-mediated cellular polarity within Mettl3 cKO HFs.

(A) Ultrastructure of HF in control and Mettl3 cKO mice in P0 back skin. Matrix (Mx) cells engulf the dermal papilla (DP, colored in pink) in the control HF but often fail to do so in the cKO. Dashed line indicates the boundary between matrix and dermal papilla. The boundary between dermal papilla cell #1 and matrix cell #2 is magnified in the panel below. BM, basement membrane. Scale bars: 10 µm (upper panel), 600 nm (lower panel). Ap, apoptotic bodies. (B) Representative images from P0 sagittal sections labeled for F-actin by phalloidin (scale bars: 25 µm). White solid lines denote skin surface and dashed lines denote dermal-epidermal border. Note apical polarization of F-actin in control HFs, often missing in Mettl3 cKO follicles.

Loss of m6A results in perturbations within the single-cell transcriptomics of the HF lineage

To further interrogate how loss of METTL3 affects the skin epithelium, we first addressed whether we might be able to obtain sufficient Mettl3 null progenitors to perform in vivo ribosomal profiling and proteomics/immunoblot analyses. We first tried dispase treatment to enrich for in vivo skin epithelium, but in contrast to control skin, cKO hair placodes/germs were left behind in the dermal fraction as early as E17.5 (Figure 5—figure supplement 1A,B). When we turned to in vitro culture, we learned that despite the ability of cKO skin progenitors to survive in vivo until birth, Mettl3 null skin epithelial progenitors failed to form colonies (Figure 5—figure supplement 1C). Although the underlying basis for the puzzling differential sensitivity to dispase and inability to survive in vitro are likely complex and beyond the scope of the present study, it compromised our ability to perform ribosomal profiling or proteomics/immunoblot analyses on Mettl3 cKO progenitors.

Although these limitations precluded high throughput analyses to interrogate the direct consequences of m6A loss to translation and protein production, we were still able to ascertain the consequences to global gene expression through single-cell RNA sequencing (scRNA-seq), an analysis requiring only a few thousand cells. We focused on E17, a time when METTL3 and transcriptome-wide m6A were lost, but before marked differences in the skin epithelium or its classical cell identity markers had surfaced (Figure 5—figure supplement 1D). We FACS-purified YFP+ cells from trypsinized skins of control and Mettl3 cKO embryos, and then binned intact YFP+ cells according to their exclusion of 4’, 6-diamidino-2-phenylindole (DAPI) staining and their relative levels of integrins (α6, β1) (Figure 5—figure supplement 1E). We then processed these isolated cells through the 10X Genomics Next GEM system for scRNA-seq, from which we obtained high-quality sequencing data from 4443 control cells and 3455 Mettl3 null cells.

Unbiased clustering of the cells based on their gene expression profiles indicated seven major clusters as revealed by t-distributed stochastic neighbor embedding (t-SNE) plots whose identities were ascertained based upon established markers (Figure 5A, Figure 5—figure supplement 1F). Three clusters were related to basal epidermal progenitors (Epi basal) and were typified by high Krt5 and Krt14 expression. We classified one cluster as the differentiating suprabasal epidermal keratinocytes (Epi suprabasal), based upon their high Krt1 and Krt10 expression. As expected, some of these Krt1+Krt10+ cells were still cycling, a feature of the asymmetric cell divisions that occur at this stage of embryonic development and which place daughter cells from proliferative parents into the suprabasal layers (Williams et al., 2011). The cluster with elevated expression of Shh, Lef1 and Ctnnb1 (encoding WNT-activated β-catenin) was identified as the WNThi cells of specified HFs, while the two clusters with high expression of Krt17 and Sox9 and little or no Lef1 were typical of WNTlo cells.

Figure 5 with 1 supplement see all
Single-cell transcriptomics of Mettl3 cKO compared to control skin epithelial lineages.

(A) YFP+ progenitors were FACS-enriched from E17 control and Mettl3 cKO whole back skins and subject to scRNA-seq as described in the Materials and methods. Shown are the data from unbiased clustering of the single cells projected by t-SNE. (Clustering was performed on all cells from control and cKO samples pooled together.) Seven major clusters and three minor clusters were identified, which are present in both samples. (B) Heat map illustrating how the seven populations were identified based upon their expression of established marker genes. Each column depicts data from one single cell belonging to the cluster that is indicated at the top. Each row illustrates cluster-wide expression data for the specific mRNA listed at the right. Color coding of each cluster is according to the color scheme in the t-SNE plots in (A). (C) Pseudotime plots showing the estimated lineage trajectories of each of the clusters as derived from the scRNA-seq data. The p values were calculated through Pearson's chi-square test.

A heat map of the expression of these markers and their clustering into these seven groups is presented in Figure 5B. The molecular distinctions among the three basal epidermal progenitor populations and among the two WNTlo populations appeared to at least partially reside within cell cycle stage differences, as judged by expression of cell cycle mRNAs such as Cdk1 and Mki67 (Figure 5B, Figure 5—figure supplement 1F). However, since markers such as SCA1 (Ly6a), that molecularly distinguish adult progenitors of the interfollicular epidermis from those of the upper outer root sheath were not yet expressed at this time, other differences in these signatures could represent early signs of their specification. Regarding the WNTlo populations, the Lhx2-expressing cells showed features more characteristic of early HF stem cells (Rhee et al., 2006), while the Krt79-expressing cells were likely to be early sebaceous gland progenitors (Veniaminova et al., 2019).

To assess lineage differentiation trajectories, we performed pseudotime analyses on both control and cKO transcriptomes (Figure 5C). As expected, for control skin, most of the cells fell along two main branches: a hair germ branch composed largely of WNThi (Krt17+Lhx2+Sox9neg) and WNTlo (Krt17+Sox9+Tgfbr2+) cells; and an interfollicular epidermal branch of Krt5hiKrt14hi basal progenitors and their suprabasal Krt1+Krt10+ progeny. In agreement with the phenotypic abnormalities, the population of WNThi cells along the HF branch was markedly diminished in the cKO, and some WNThi cells were abnormally positioned along an extra lineage branch. These findings were consistent with the reduced density of hair placodes seen at this time (Figure 2F).

Signaling pathways and canonical translation initiation factors are among the mRNAs that are downregulated upon METTL3 loss

To gain further insights into m6A’s function, we next evaluated how steady-state levels of each mRNA within a skin cell population changed upon METTL3 depletion (Supplementary file 3), and then correlated this change to m6A modification levels. In addition to the sum of normalized-to-input uTPM parameter that we used earlier, we also tested other parameters to assess m6A levels (Supplementary file 4). The degree of correlation with normalized-to-input uTPM was higher than merely counting m6A site numbers. Correlations were also higher when we normalized to RNA length (i.e. per nucleotide, per nt) and focused on the density of m6A within the coding sequence rather than the 5’ or 3’ UTRs or full-length mRNA (Figure 6—figure supplement 1). Therefore, in the following analysis with the scRNA-seq data, the per nt sum of normalized-to-input uTPM in coding sequence (coding sequence SN-uTPM per nt) parameter was used to evaluate the importance of m6A levels with regards to mRNA level changes arising from METTL3 loss.

Given our data thus far, we began by focusing on those RNAs that were downregulated upon METTL3 loss. Since HFs are specified within the basal epidermal plane and are WNThi, we primarily focused on the most significantly downregulated mRNAs in these categories. Although the dynamic range of fold changes for this RNA cohort was low, 2555 genes were downregulated with a cut-off Z score of <−1.96 (all dots on the left of the dotted line in the plots of Figure 6A; Supplementary file 3). Intriguingly, many of these mRNAs fell into signaling pathways that had surfaced when we discovered a correlation between high m6A with efficient translation (Figure 6B; compare with Figure 1G). Moreover, although the mRNAs within these categories differed, some of the most significantly downregulated mRNAs upon Mettl3 ablation were ones that were amongst the most highly m6A-modified in wild-type cells (blue dots in Figure 6A; blue-highlighted mRNAs in Figure 6B). Of additional note, canonical translation initiation complex factors were also included in the cohort of significantly downregulated mRNAs, suggesting an additional layer of effects on global translation through m6A.

Figure 6 with 2 supplements see all
Investigation of RNAs whose levels diminish upon METTL3 loss.

(A) scRNA-seq data from Figure 5A were binned according to four major classifications: Epi basal, Epi suprabasal, HF WNThi and HF WNTlo. Scatter plots of mRNAs in these cells were then analyzed according to their expression changes (cKO/Ctrl) assessed by scRNA-seq Z score and to their coding sequence m6A density in wild-type (WT) skin epithelium assessed by the miCLIP SN-uTPM per nt value. Dots on the left of the dashed line in each plot indicate RNAs which from scRNA-seq have a Z score (cKO/Ctrl)<−1.96 (Supplementary file 3). Among those, blue dots denote mRNAs with m6A coding sequence SN-uTPM per nt among the top 20% (Supplementary file 4). (B) Major pathways and their associated RNAs that were downregulated upon METTL3 loss. Shown are the data for the basal epidermal progenitors, which at E17 contained both epidermal and hair placode cells. mRNAs highlighted in blue correspond to blue dots in (A), and were among the most significantly downregulated upon METTL3 loss but heavily m6A-modified in wild-type. Note that many of these pathways also corresponded to those whose heavily m6A-modified mRNAs were also efficiently translated. (C) Violin plots illustrating the relative expression levels of Lef1 mRNA in the Mettl3 cKO versus control basal epidermal progenitors and WNThi progenitors. Z score assessment of expressional difference between cKO and control [Z (cKO/Ctrl)] and false discovery rate (FDR) is calculated by MAST. The down-regulation was verified with qPCR on total RNA samples extracted from YFP+ skin epithelial cells FACS isolated from E16.5 embryos with Tbp mRNA as internal control (error bars: standard deviation, for each condition n = 3 biological replicates, **p<0.01 by unpaired two-tailed Student’s t-test). (D) Confocal images of E16.5 whole-mount back skin immunolabeled for HES1 and YFP (scale bars: 20 µm). HES1 expression was quantified in the stratified layers of skin epithelium (middle line corresponds to the mean; for each condition, the data are from two biological replicates; ****p<0.0001 by unpaired two-tailed Student’s t-test). (E) Pulse-chase assay examining the rate of epidermal cell flux from basal to suprabasal layers. Control and cKO animals were pulsed at E18.5 with EdU and the signal was then chased until P1. Before tissue collection, the P1 pups were treated with a short (1 hr) BrdU pulse. P1 back skin sagittal sections were subjected to immunofluorescence to examine the EdU and BrdU labeling in the basal versus suprabasal layers of the epidermis (scale bars: 25 µm). White solid lines denote skin surface and dashed lines denote dermal-epidermal border. Cell flux rates were quantified based on the ratio of EdU+ cells in the suprabasal layer to BrdU+ basal cells (for each condition n = 4 biological replicates ×10 images per replicate, *p<0.05 by unpaired two-tailed Student’s t-test). (F) Radial histograms depicting the division orientation of epidermal basal cells during anaphase/telophase at E17.5 and P0, assessed by IF staining of Survivin, integrin β4 (CD104) and PCAD as described in Williams et al., 2011. For each condition, three biological replicates were analyzed and n indicates the total number of anaphase/telophase cells examined from the embryos. (G) Ultrastructure of epidermis in control and Mettl3 cKO P0 back skin. Ba, basal layer, colored in green; Sp, spinous layer, colored in greenish yellow; Gr, granular layer, colored in yellow; SC, stratum corneum, colored in orange. Note the increased numbers of cells in the spinous layer of cKO, and the presence of nuclei in many cells of the granular layer. The boundary between dermis (Der) and the basal layer is shown in the middle panel. KF, keratin filaments; HD, hemidesmosomes. The border between cell #1 (basal) and cell #2 (suprabasal) is shown in the lower panel. Intercellular membranes are sealed in the control. Note small gaps (arrow) are present at the intercellular border, more frequently in cKO than in control. Scale bars: 10 µm (upper panel), 600 nm (middle and lower panel).

Given that a number of phenotypic features of METTL3 loss are recapitulated in WNT pathway mutants, we naturally gravitated toward the WNT signaling mRNAs that were significantly downregulated. Although Lef1 was not within the top 20% of m6A-modified mRNAs, it was within the top cohort of significantly downregulated mRNAs upon Mettl3 targeting and was confirmed by quantitative PCR (qPCR) (Figure 6C). Lef1 mRNAs were reduced in both basal epidermal and WNThi progenitors, further suggesting that some basal epidermal progenitors were unable to progress to form HFs without m6A modification. Many other factors required for WNT-mediated cell fate specification were also significantly downregulated (Figure 6B). Although changes were relatively modest at this early stage, levels of WNT signaling are known to profoundly impact fate outcomes and proper tissue morphogenesis (Buechling and Boutros, 2011; Tortelote et al., 2017).

It was also notable that the NOTCH signaling pathway also resurfaced in our loss of function analyses of basal epidermal progenitors. Both Jag1 and Maml2 were among the most highly m6A-modified mRNAs and among the most significantly downregulated (Figure 6B). NOTCH signaling functions suprabasally at the transition between the basal and spinous layer, where it is typified by expression of classical NOTCH target, HES1 (Blanpain et al., 2006). Indeed, as judged by whole-mount immunofluorescence, HES1 was diminished (Figure 6D).

While HES1 downregulation correlated with the changes in gene expression revealed by scRNA-seq, it was also possible that it was reflective of a change in the flux of basal cells into the spinous layers. We were particularly intrigued by this latter possibility since components of hemidesmosomes, for example BPAG1 (Col17A1) and integrin β4 (Itgb4), which are responsible for the bulk of basal cell adhesion to the basement membrane (Dowling et al., 1996), were among the most significantly downregulated mRNAs upon METTL3 loss (Figure 6B). We therefore performed a pulse-chase experiment, by first administering a short pulse of deoxythymidine analog 5-ethynyl-2′-deoxyuridine (EdU) on E18.5 pups, and then treating with bromodeoxyuridine (BrdU) just prior to analyses (Figure 6E). Although labeling was somewhat higher in the basal cells of cKO compared to control skin, BrdU incorporation was limited to the basal progenitors in both. Most notably, the percentage of EdU-suprabasal:BrdU-basal cells was greater in the cKO skin, indicating that Mettl3 null basal progenitors fluxed at a higher rate than normal into the suprabasal layers.

Additional signs of basal:suprabasal perturbations in the epidermis were revealed by immunofluorescence imaging for the basal progenitor cadherin, P-cadherin (PCAD), relative to E-cadherin (ECAD), a pan-epidermal marker. While fluorescence intensities showed that both cadherins were elevated at intercellular borders, the PCAD:ECAD ratio was significantly lower than normal (Figure 6—figure supplement 2A). PCAD is often associated with cells that undergo migration and invasion, while ECAD is typically associated with more static intercellular contacts (Kümper and Ridley, 2010). Thus, the reduction in PCAD:ECAD ratio, along with the diminished actin regulators, were further consistent with enhanced departure of Mettl3 null cells from the basal layer, where cells have traction.

The spindle orientations of Mettl3 null cell divisions were largely unaltered and remained within the plane of the basal layer (Figure 6F). However, at P0, a 45 min pulse with EdU revealed signs of elevated proliferation in the cKO basal cells (Figure 6—figure supplement 2B). This was surprising as pups did not grow in size appreciably post-birth (Figure 2—figure supplement 2A). As we show below, this could also not be accounted for by an appreciable rise in apoptosis. Rather, we attribute this change to reflect replacements due to the increased departure of basal cells into the suprabasal layers. This was accompanied by a marked increase in cellularity, particularly within the spinous layers (Figure 6G, Figure 6—figure supplement 2B). Suprabasal cells were also smaller than normal (Figure 6—figure supplement 2C). The increase of cell number and decrease of cell size in the suprabasal layers was also confirmed by FACS quantifications, in which we found the numbers of K10+ cells were greater, with somewhat lower values on the forward scatter projection (Figure 6—figure supplement 2D).

Although the four stages of terminal differentiation in Mettl3 cKO epidermis were still recognizable by morphology and by molecular markers (Figure 6G, Figure 6—figure supplement 2E), perturbations reverberated into the later stages of differentiation. In normal skin, when spinous cells transition to the granular layer, they initially remain transcriptionally active, but then a destructive phase ensues as they lose nuclei and other organelles and flatten to become dead ‘squames’ (Quiroz et al., 2020). This is normally accompanied by DNA destruction, detected by TUNEL staining (Figure 6—figure supplement 2F). In contrast, Mettl3 cKO epidermis displayed an appreciable number of nuclei in their granular layers (Figure 6G) and showed little or no signs of TUNEL positive cells (Figure 6—figure supplement 2F).

Additionally, although desmosomes and hemidesmosomes were still present, small gaps were often noted ultrastructurally at intercellular borders in the Mettl3 cKO epidermis (arrow in Figure 6G, middle frame). We also observed occasional signs of spinous cell cytolysis, which was also reflected in the increased incidence of early suprabasal cells positive for TUNEL but negative for activated Caspase-3 (Figure 6—figure supplement 2F). Given the many perturbations seen in the Mettl3 null epidermis, it was remarkable that the skin barrier of neonatal mice was still functional as judged by its ability to retain fluids (Figure 2—figure supplement 2C).

Signs of compensatory mechanisms and RNA metabolism revealed in the mRNAs that are upregulated upon METTL3 loss

Finally, we turned to the mRNAs that were significantly upregulated upon METTL3 loss. We began by examining the expression status of Mettl3 cKO versus control E17 skin mRNAs whose coding sequence SN-uTPM per nt values from our miCLIP data had been among the top 20% of heavily m6A-modified transcripts. Interestingly, a considerable cohort surfaced whose mRNA levels were significantly upregulated upon m6A loss, and which had a higher dynamic range of fold-changes than the downregulated mRNAs (Figure 7A, red dots in the right quadrant; Supplementary file 3, 4 with a cut-off Z score >1.96). In fact, mRNAs that showed high m6A modifications within their coding sequences were generally more up- than downregulated upon METTL3 loss (Figure 6—figure supplement 1 and Figure 7—figure supplement 1A). These findings were consistent with prior reports that m6A can promote RNA degradation (Du et al., 2016; Wang et al., 2014a; Zaccara and Jaffrey, 2020).

Figure 7 with 1 supplement see all
Investigation of RNAs upregulated upon METTL3 loss.

(A) Scatter plots of mRNAs in the indicated groups of cells (as in Figure 5A) based on expression changes (cKO/Ctrl) assessed by scRNA-seq Z score and correlated with the coding sequence m6A density in wild-type (WT) skin epithelium assessed by the miCLIP SN-uTPM per nt value. Dots on the right of the dashed line in each plot indicate RNAs which from scRNA-seq have a Z score (cKO/Ctrl)>1.96 (Supplementary file 3). Among those, red dots denote mRNAs whose m6A coding sequence SN-uTPM per nt was among the top 20% (Supplementary file 4). (B) GSEA of transcripts with Z score (cKO/Ctrl)>1.96, FDR <0.05 in scRNA-seq and m6A coding sequence SN-uTPM per nt among the top 20% indicating the GO terms enriched. (C) Examination of the expression of selected upregulated mRNAs on the t-SNE plots. mRNA names are color-coded according to the GO terms they belong to in (B). Those mRNAs whose names are in black are known to be involved in translation regulation, and scored as elevated significantly upon METTL3 loss. *denotes mRNAs encoding factors related to m6A dynamics. (D) qPCR validation of increased levels of mRNAs in cKO/Ctrl. qPCR is on total RNA samples extracted from YFP+ skin epithelial cells FACS isolated from E16.5 embryos with Tbp mRNA as internal control (error bars: standard deviation, for each condition n = 3 biological replicates, *p<0.05 **p<0.01 ***p<0.001 by unpaired two-tailed Student’s t-test). (E) Proposed model summarizing the effects of m6A loss on mRNA translation and degradation in the skin epithelia, and the consequences to their integrity and fate choices. The ovals represent all mRNAs that were either downregulated or upregulated significantly upon METTL3 loss. Pink circle: by first examining mRNAs that were heavily modified by m6A and also efficiently translated (Figure 1), and then independently identifying mRNAs that were among the most significantly downregulated upon METTL3 loss, we discovered considerable overlap in their pathways, suggestive of a translation block. The finding that a number of mRNAs involved in canonical translation were also significantly downregulated added to this notion. Factors involved in WNT signaling, NOTCH signaling and adhesion were featured prominently, in agreement with the morphogenetic defects observed. Green circle: these mRNAs were significantly upregulated upon METTL3 loss and were heavily m6A-modified in wild-type skin, but they did not correlate with translation efficiency. Rather, they encompassed mRNAs indicative of translational rescue pathways as well as feedback mechanisms.

In all four populations, the major enriched GO term categories of highly modified mRNAs in wild-type that were significantly upregulated upon m6A loss revealed changes in various aspects of RNA metabolism, including RNA binding, RNA processing, RNA modification, ribosome and translation, as well as ribonucleoprotein complex (Figure 7B; Supplementary file 5). Examples of the mRNAs within specific categories (color coded by category) are shown in the t-SNE plots in Figure 7C and validations by qPCR are shown in Figure 7D. The categories of mRNAs enriched upon METTL3 loss and highly m6A-modified in wild-type conditions raised the possibility that the upregulated mRNAs might represent a global response to m6A loss.

Interestingly, CAP-independent translation was among the upregulated GO terms, consistent with the reduced mRNA levels that we had detected for canonical initiation. RNA methylation was also among the GO terms. Additionally, of note, mRNAs encoding all three of the major m6A reader proteins—Ythdf1, Ythdf2 and Ythdf3 (Wang et al., 2015; Wang et al., 2014a)were not only heavily m6A-modified in their coding sequences (Supplementary file 4), but also significantly upregulated upon m6A loss (Z-scores ranging from 3.8 to 16) (Supplementary file 3). These findings were suggestive of a possible feedback mechanism to sense m6A paucity and compensate. The upregulation of Pelota (Pelo) was also intriguing (Figure 7C). Pelota is a component of the evolutionarily conserved quality control machinery that rescues stalled ribosomes, and it is also essential for skin epidermal integrity (Liakath-Ali et al., 2018).

Of additional note, MYC has long been linked to ribosomal biogenesis and translational control (Piazzi et al., 2019; van Riggelen et al., 2010), and MYC overexpression in embryonic epidermal progenitors is also known to favor epidermal and sebocyte differentiation at the expense of the HF lineage (Berta et al., 2010). In this regard, it is also noteworthy that both Myc and Bmyc mRNAs were not only highly modified by m6A within their coding sequences, but also by scRNA-seq, were significantly elevated within the Mettl3 null epidermal progenitors of E17 skin (Figure 7—figure supplement 1B,C). Immunofluorescence microscopy was consistent with these differences (Figure 7—figure supplement 1D). Although modest, these increases followed the right trend to be a contributing factor to the Mettl3 cKO phenotype.

Discussion

Despite the increasing knowledge of the molecular mechanisms that underlie the roles of m6A modifications on mRNA degradation and translation, it has remained elusive how this abundant modification affects tissue biology. Our current study uncovered several new insights into m6A’s actions and made headway in tying together and advancing our knowledge of certain recurrent themes that have emerged from prior in vivo and in vitro studies. As importantly, our findings place m6A as a key integrator of stem cell fate choices and translational control (Figure 7E).

In gaining new insights, we were aided by having previously performed in vivo ribosomal profiling on embryonic murine skin epithelium (Sendoel et al., 2017). By performing high-resolution miCLIP on this tissue, we were uniquely poised to tease out, on a global in vivo scale, a prominent correlation between the degree of m6A modifications that an mRNA has within its coding sequence and its translation efficiency. This correlation was particularly pronounced in mRNAs encoding signaling pathways known to affect skin progenitor fates. Moreover, by performing scRNA-seq on control and Mettl3 cKO skin progenitors, we were further able to appreciate that many of the most significantly downregulated mRNAs upon m6A loss were associated with these same pathways.

The most striking phenotypic defects surfacing upon m6A loss were in HF morphogenesis, where both specification and down growth were markedly curtailed. These processes require WNT signaling, and mRNAs encoding a number of WNT regulators were significantly downregulated. Moreover, WNT-effectors nuclear LEF1 and β-catenin, were both diminished at the sites where WNT signaling is known to be essential. Further signs that WNT signaling was deleteriously affected included the promotion of sebaceous gland fate characteristics at the expense of HF fate ones (Merrill et al., 2001; Niemann et al., 2003), the perturbations in SHH signaling and also the difficulty that Mettl3 null HFs seemed to face in engulfing the DP.

Recent studies on other cell types and cancers have also reported connections between WNT signaling and m6A levels ( Bai et al., 2019Cui et al., 2020; Han et al., 2020; Liu et al., 2019; Miao et al., 2019; Tang et al., 2020; Zhang et al., 2019). In some of these cases, m6A loss resulted in an elevation of WNT signaling, while in others, it was diminished, and a variety of mechanisms have been proposed to explain these seemingly context-dependent effects. Our findings show that not only are WNT pathway mRNAs among the most heavily m6A-modified and efficiently translated, but they are also among the most significantly downregulated early after m6A loss. When coupled with the physiological consequences to skin development, our findings place the impact of m6A alterations at the helm of this pathway.

For the skin, the link between WNT signaling and m6A was particularly intriguing because it offered a hitherto elusive connection between mRNA regulation and cellular fate choices. Additionally, given that overall levels of WNT signaling can profoundly impact fate specification, our studies reveal how even subtle changes in this pathway can have a profound impact on fate outcomes. That said, despite the importance of WNT signaling in skin progenitor fate specification and the consistency of the METTL3 loss of function phenotype with a decrease in this pathway, m6A modifications in the skin epithelium went well beyond this pathway. Indeed, our comprehensive analyses of the Mettl3 cKO phenotype in the skin provided graphic appreciation that m6A RNA modifications function broadly and diversely in regulating cellular responses.

In addition to WNT signaling, the pathways affected when highly m6A-modified and efficiently translated mRNAs experienced depletion of m6A encompassed actin regulators, cell polarity, ECM-receptor interaction and NOTCH signaling. A number of cellular changes and phenotypic correlations were consistent with these differences. Thus for example, the perturbations we found in apical cytoskeletal polarity seem likely to contribute to the decreased ability of HF cells to engulf the DP. Additionally, the decrease in HES1 in early spinous cells was reflective of diminished NOTCH signaling and consistent with mRNA downregulation of NOTCH ligand JAG2 and NOTCH co-activator Mastermind like 2. And although we did not observe overt morphological differences in hemidesmosomes, the heightened departure of progenitors from the basal layer and their inability to form colonies in vitro were in agreement with possible defects in Mettl3 null basal progenitor:ECM adhesion.

The high number of mRNAs upregulated significantly upon METTL3 loss were equally informative but in a strikingly different way. Although many of these transcripts were among the most highly m6A-modified in the wild-type state, their GO terms were not among the top m6A-modified mRNAs that were efficiently translated. On the contrary, the marked increase in their mRNA levels upon m6A ablation was suggestive that this mark normally functions to destabilize this cohort.

In contemplating the physiological relevance of the upregulated mRNAs in Mettl3 cKO embryos, we were drawn to several curious possibilities. First, it was both intriguing and paradoxical to find that that Myc and Bmyc were upregulated upon m6A loss. In cancer cells, MYC has been implicated as a key positive regulator of CAP-dependent mRNA translation (Cargnello and Topisirovic, 2019), a process that was predicted to be downregulated upon Mettl3 ablation in embryonc skin. That said, MYC overexpression has also been shown to promote differentiation in the epidermis and sebaceous glands at the expense of HF morphogenesis (Arnold and Watt, 2001; Cottle et al., 2013; Watt et al., 2008b), phenocopying some of the key features we observed upon m6A loss.

Notably, upregulation of RNA processing factors has been observed in the skin progenitors of both MYC-overexpressing skin (Arnold and Watt, 2001; Frye et al., 2003), and Mettl3 cKO skin. As global rates in RNA metabolism can profoundly affect stem cells' activities as well as how they respond to stressful environments (Blanco et al., 2016; Sampath et al., 2008; Sendoel et al., 2017; Signer et al., 2014; Starck et al., 2016; Zismanov et al., 2016), it seems likely that the block in the degradation of this selected cohort of RNAs might impact the translation of the signaling pathways that orchestrate lineage specification.

A second major curiosity was that by miCLIP, Ythdf1, Ythdf2 and Ythdf3 fell within the top 20% of mRNAs whose coding sequence in basal skin progenitors was highly modified, and in Mettl3 null basal progenitors, these mRNAs were among the most significantly increased in levels. YTHDF1-3 function as m6A readers. While YTHDF1 and 3 have been proposed to increase the stability and translation of mRNAs (Shi et al., 2017; Wang et al., 2015), recent studies suggest that all three may function in localizing m6A-modified mRNAs to the decay machinery (Wang et al., 2014a; Zaccara and Jaffrey, 2020). Irrespective of their precise reader roles, their collective m6A modification in control progenitors and upregulation upon m6A loss is suggestive of a key rescue juncture for the cell.

Indeed many of the other upregulated pathways were also suggestive of potential feedback mechanisms, activated upon m6A loss. Among them were CAP-independent translation genes, RNA methylation genes, RNA processing pathways and ribosome-associated mRNA quality control genes. When coupled with the downregulation in mRNAs that are normally m6A-modified and highly translated, and with the reduction in canonical translation transcripts, these upregulated mRNA perturbations pointed not only to further signs of an altered protein translation state, but also to a potential built-in mechanism to boost alternative translational routes and RNA metabolism lost upon Mettl3 ablation. Thus by having a cohort of key regulatory mRNAs that are heavily modified by m6A, progenitors are able to sense when the modification is lost and respond by rapidly enhancing their mRNA pool and calling into action rapid responses and rescue pathways to control the damage involved.

A model based on our findings and summarizing this speculation is presented in Figure 7E. Future studies will be needed to probe the significance of m6A modifications of individual mRNAs and to unearth the mechanisms underlying why some m6A-modified mRNAs were downregulated upon m6A loss while others were upregulated. It will also take further investigation to sift through the myriad of mRNAs and pathways that were affected both negatively and positively by m6A loss, and which exposed an underlying complexity in m6A's roles in translational regulation. The problem becomes even more daunting when considering that we have not yet examined the consequences of METTL3 loss to miRNAs in the skin. In this regard, in a breast cancer cell line, METTL3 was found to mark pri-miRNAs for recognition and processing (Alarcón et al., 2015), and in skin, METTL3 depletion shares certain similarities with global loss of miRNAs (Yi et al., 2008; Yi and Fuchs, 2011).

While many areas remain to be explored, our results show that the two arms of the pathways for down- and up-regulated mRNAs are intricately interwoven. In this regard, our many findings reported here shed new light on the previously identified roles of m6A in translational regulation. Finally, and intriguingly, the circuitry funneled to progenitor fate choices, each of which uses a distinct repertoire of signaling pathways and has distinct needs for protein synthesis. Our studies on m6A through the lens of skin morphogenesis now place m6A on the center stage of this arena.

Materials and methods

Mouse strains

Request a detailed protocol

The Mettl3fl mouse strain was a gift from the Brüning Lab at Max Planck Institute for Metabolism Research and Policlinic for Endocrinology, Diabetes and Preventive Medicine. It was generated by inserting LoxP sites at both sides of the fourth exon of the Mettl3 gene through homologous recombination with the constructed template from the Knockout Mouse Project Repository (Cheng et al., 2019). Sequences of primers for genotyping of the Mettl3fl allele are listed in Supplementary file 6. The Mettl3fl strain was then crossed with Krt14-Cre, Rosa26-YFPfl and Gli1-LacZ strains to generate breeders. The Gli1-LacZ strain (Swiss Webster background) is from the Joyner Lab at Memorial Sloan Kettering Cancer Center. Nude(Nu/Nu) mice for grafts are from Charles Rivers Laboratories Strain 088. The wild-type mice for miCLIP and the Mettl3fl, Rosa26-YFPfl, Krt14-Cre mouse strains are under the C57BL/6 background. All animals were maintained in an American Association for the Accreditation of Laboratory Animal Care Internationally approved Comparative Bio-Science Center at the Rockefeller University and procedures were performed using Institutional Animal Care and Use Committee-approved protocols that ad-here to the standards of the National Institutes of Health.

Back skin engraftment

Request a detailed protocol

Back skin dissected from P0 pups was incubated on a moist surface with phosphate-buffered saline (PBS) at 4°C overnight and then grafted to Nu/Nu female mice (6–8 weeks old) with the control and cKO skin on each side of the back. After 15 days, bandages were removed, and the grafted skin pieces were harvested for analysis.

EdU and BrdU labeling

Request a detailed protocol

For EdU pulse labeling at P0, EdU was applied on the pups at 25 µg per gram body weight through intraperitoneal injection at 45 min before sacrifice. For EdU chase-BrdU pulse, EdU was applied on pregnant female mice bearing E18.5 embryos at 25 µg per gram body weight through intraperitoneal injection, and BrdU was applied on the born P1 pups at 50 µg per gram body weight through intraperitoneal injection at 1 hr before sacrifice.

Histology analysis, X-gal staining, Oil Red O staining and immunofluorescence on sagittal sections

Request a detailed protocol

For sagittal sections, back skin or whole embryos were embedded in OCT compound and cut into 10–14 µm sections on a Leica cryostat. The sections were mounted on SuperFrost Plus Adhesion slides (VWR), air-dried and fixed in PBS with 4% paraformaldehyde (PFA) at room temperature for 15 min.

Hematoxylin and eosin staining was performed with an adapted protocol from the University of California San Diego (http://mousepheno.ucsd.edu/hematoxylin.shtml).

X-gal staining was performed for intracellular β-galactosidase activity assessment. Slides were first stained with X-gal staining buffer (100 mM sodium phosphate buffer pH = 7.3, 1.3 mM MgCl2, 3 mM K3Fe(CN)6, 3 mM K4Fe(CN)6, 1 mg/ml X-gal) at 37°C for 3 hr and then counter stained with nuclear fast red (Sigma-Aldrich, N3020) at room temperature for 5 min.

Oil Red O staining was performed in 2 µg/ml Oil Red O-isopropanol solution at room temperature for 1 hr followed with nuclei staining by hematoxylin.

For immunofluorescence, the slides were incubated in blocking buffer (1% fish gelatin, 5% normal donkey serum, 1% bovine serum albumin, 0.3% Triton X-100 in 1x PBS) at room temperature for 1 hr. Primary antibody staining was performed in the blocking buffer at 4°C overnight and secondary staining was performed in the blocking buffer at room temperature for 1 hr. Primary antibodies used: rat anti-BrdU (1:200, Abcam, ab6326); rat anti-CD104 (1:500, BD Pharmingen, 553745); armenian hamster anti-CD29 (1:500, BioLegend, 102201); rat anti-CD49f (1:1000, BioLegend, 313602); rabbit anti-cleaved Caspase-3 (1:200, R and D Systems, AF835); rabbit anti-FLG (1:2000, Fuchs lab); rabbit anti-iNV (1:2000, BioLegend, 924401); guinea pig anti-K5 (1:500, Fuchs lab); rabbit anti-K10 (1:1000, Covance, PRB-159P-100); guinea pig anti-LEF1 (1:5000, Fuchs lab); rabbit anti-LHX2 (1:5000, Fuchs lab); rabbit anti-LOR (1:4000, Fuchs lab); rabbit anti-METTL3 (1:500, Abcam, ab195352); rabbit anti-MYC (1:100, Abcam, ab32072); goat anti-PCAD (1:300, R and D Systems, AF761); rabbit anti-SOX9 (1:1000, Millipore, AB5535); rabbit anti-Survivin (1:500, Cell Signaling, 2808); chicken anti-GFP/YFP (1:1000, Abcam, ab13970). Secondary antibodies conjugated to Alexa Fluor 488, Rhodamine Red-X and Alexa Fluor 647 are from Life Technologies. In Situ Cell Death Detection Kit, TMR red (Sigma-Aldrich, Roche-12156792910) was used for TUNEL staining. AF647 Click-iT EdU cell proliferation kit for imaging (Thermo Fisher Scientific, C10340) was used for EdU labeling. AF647 Phalloidin (Thermo Fisher Scientific, A22287) was used for F-actin staining. Slides were mounted in ProLong Gold Antifade Mountant with DAPI (Thermo Fisher Scientific, P36941) and imaged under a Zeiss Axio Observer.Z1 epifluorescence microscope equipped with a Hamamatsu ORCA-ER camera (Hamamatsu Photonics), and with an ApoTome.2 (Zeiss) slider that reduces the light scatter in the fluorescent samples controlled by Zen software (Zeiss). Images were processed through Fiji (ImageJ).

For quantitative analyses, different groups of cells (epi basal, suprabasal, HF WNThi and HF WNTlo) were identified according to PCAD staining. Average LEF1 and MYC immunofluorescence intensity was measured among cells identified in the same group from one image. Background was measured and subtracted for each channel. The fluorescence intensity was then normalized to the average values from the corresponding control samples to get relative signal values. Cell division angles were quantified as described in Williams et al., 2011. Generally, the division angle was between the spindle orientation determined by Survivin and DAPI staining and the basal membrane determined by integrin β4 staining. Other quantifications (HF numbers, EdU, BrdU, cell sizes) are based on counting/measuring specific features in the images. TUNEL and cleaved Caspase 3-positive cells were quantified per HF (placode, germ and peg) and per cell number in the stratified layers, where positive cells were observed (basal, suprabasal layer one and late granular).

To reduce any bias in data collection, all data from each group were not analyzed until all images were collected. No statistical method was used to predetermined sample size, randomization and experiment blinding was not used. Each experiment was repeated with at least two replicates and data presented is from three or more embryos, same age. Significance of p value was set at <0.05. Statistical details for each experiment, including the statistical test used, the sample size for each experiment and p value can be found in the corresponding figure legend.

Immunohistochemistry

Request a detailed protocol

Back skin dissected from neonates was fixed in PBS with 4% PFA at 4°C overnight. Then the tissue was washed in PBS with 0%, 35% and 70% EtOH gradually and sent to Histowiz for further processing. Generally, immunohistochemistry was performed on a Bond Rx autostainer (Leica Biosystems) with enzyme treatment (1:1000) using standard protocols. Antibodies used were rat monoclonal F4/80 primary antibody (1:200, eBioscience, 14–4801) and rabbit anti-rat secondary (1:100, Vector). Bond Polymer Refine Detection (Leica Biosystems) was used according to manufacturer’s protocol. After staining, sections were dehydrated and film coverslipped using a TissueTek-Prisma and Coverslipper (Sakura). Whole slide scanning (40x) was performed on an Aperio AT2 (Leica Biosystems).

Whole-mount immunofluorescence and confocal microscopy

Request a detailed protocol

For whole-mount immunofluorescence, embryos were fixed in PBS with 4% PFA for 1 hr, followed by extensive washing in PBS. Samples were then permeabilized for 3 hr in 0.3% Triton X-100 in PBS and treated with Gelatin Block (2.5% fish gelatin, 5% normal donkey serum, 3% BSA, 0.3% Triton, 1x PBS). The following primary antibodies were used: rabbit anti-ECAD (1:500, Cell Signaling, 3195); rabbit anti-HES1 (1:200; Fuchs lab); rabbit anti-METTL3 (1:100, Abcam, ab195352); goat anti-PCAD (1:600, R and D Systems, AF761); chicken anti-GFP/YFP (1:1200, Abcam, ab13970). Primary antibodies were incubated at 4°C for 36 hr. After washing with 0.1% Triton X-100 in PBS, samples were incubated overnight at 4°C with secondary antibodies conjugated to Alexa Fluor 488, Rhodamine Red-X and Alexa Fluor 647 (1:1000, Life Technologies). For E16.5 whole skin and P0 epidermis, samples were washed, counterstained with DAPI and mounted in SlowFade Diamond Antifade Mountant (Invitrogen). For P0 whole skin, samples were washed, counterstained with DAPI and processed though tissue clearing with ethyl cinnamate as described in Gur-Cohen et al., 2019. Confocal images of whole-mounts were acquired using a spinning disk confocal system (Andor Technology Ltd) equipped with an Andor Zyla 4.2 and a Yokogawa CSU-W1 (Yokogawa Electric, Tokyo) unit based on a Nikon TE2000-E inverted microscope. Four laser lines (405, 488, 561 and 625 nm) were used for near simultaneous excitation of DAPI, Alexa448, Rhodamine Red-X and Alexa647 fluorophores. The system was driven by Andor IQ3 software. 40x oil objective was used to acquire Z stacks of 0.5–1 µm steps.

HES1+ cells were quantified using Fiji (ImageJ). Briefly, 40x Z stacks of spinning disk confocal images of inter follicular epidermis 2500 µm2 regions were converted into composite images in which DAPI was in blue channel, YFP in the green channel and HES1 in the red channel. From each region, HES1+ and YFP+ cells were quantified per optical sections of basal, suprabasal 1 (S-1) and suprabasal 2 (S-2) layers. The numbers of HES1+ cells fwere recorded and the proportions calculated relative to the total YFP+ cells for each represented layer. A minimum of 25 regions were analyzed per embryo.

PCAD and ECAD immunofluorescence quantifications were performed using Fiji (ImageJ). Briefly, using spinning disk Z stacks of whole-mount 40x confocal images, we summed the intensity across the basal layer of 2500 µm2 regions. The integrated density of PCAD and ECAD immuolabeling through the perimeter of the cell was measured and recorded for 10 cells in each region to a minimum of 18 regions per embryo. Background was then measured and subtracted for each channel. The ratio PCAD:ECAD was calculated per cell.

To reduce any bias in data collection, all data from each group were not analyzed until all images were collected. No statistical method was used to predetermined sample size, randomization and experiment blinding was not used. Each experiment was repeated with at least two replicates and data presented is from two or more embryos, same age. Distributions were tested for normality using D’Agostino and Pearson test. To test significance, unpaired or paired two-tailed Student’s t-tests were used for normal distribution and nonparametric Mann-Whitney test when the distribution did not follow a normal distribution. Significance of p value was set at <0.05. Statistical details for each experiment, including the statistical test used, the sample size for each experiment and p values can be found in the corresponding figure legend. All graphs and statistics were produced using GraphPad Prism 8.2 for MAC, GraphPad Software, San Diego, California, www.graphpad.com.

Electron microscopy

Request a detailed protocol

Skin samples were fixed in 2% glutaraldehyde, 4% PFA, and 2 mM CaCl2 in 0.05 M sodium cacodylate buffer, pH 7.2, at room temperature for over 1 hr. Then the samples were post-fixed in 1% osmium tetroxide, and processed for Epon embedding. Ultrathin sections (60 to 65 nm) were counterstained with uranyl acetate and lead citrate. Images were taken with a transmission electron microscope (Tecnai G2-12; FEI) equipped with a digital camera (AMT BioSprint29).

Flow cytometry and cytospin analysis of skin epithelial cells

Request a detailed protocol

Preparation of embryonic and neonatal back skin for isolation or examination of epithelial cells by flow cytometry and staining protocols were performed as previously described (Asare et al., 2017; Sendoel et al., 2017). For E16.5 and E17 samples, back skin dissected from the embryos were treated with 0.25% collagenase and then 0.25% trypsin to get cell suspension. For P0 skin, back skin dissected from the neonates were treated with dispase first to separate the epidermis from dermis. Then the epidermis was digested with 0.25% trypsin to get cell suspension. Antibodies used: anti-CD29 PECy7 (1:1000, Thermo Fisher Scientific, 25–0291-82); anti-CD31 APC (1:1000, Thermo Fisher Scientific, 17–0311-82); anti CD45-APC (1:1000, Thermo Fisher Scientific, 17–0451-83); anti-CD49f PE (1:1000, Thermo Fisher Scientific, 12–0495-81); anti-CD49f PECy7 (1:1000, Thermo Fisher Scientific, 25–0495-82); anti CD117-APC (1:1000, Thermo Fisher Scientific, 17–1172-81); anti CD140a-APC (1:500, Thermo Fisher Scientific, 17–1401-81); rabbit anti-K10 (1:500, Biolegend, PRB-159P); mouse anti-rabbit IgG PE-Cy7 (1:300, Santa Cruz, sc-516721). CD31, CD45, CD117 and CD140a are lineage-negative markers. For FACS, dead cells were excluded by DAPI staining. The samples were proceeded to FACS on BD FACSAriaII sorter commended by the Diva software (BD Biosciences). For flow cytometry-based examination of K10+ suprabasal cells, dead cells were excluded by LIVE/DEAD Fixable Aqua Dead Cell Stain Kit (Thermo Fisher Scientific, L34965). Antibody staining was performed after cells were fixed with Fixation/Permeabilization Solution Kit (BD Biosciences, 554714). Samples were analyzed on BD LSRII flow cytometer commended by the BD FACSDiva software (BD Biosciences) and data were processed with FlowJo. For cytospin analysis, cells were spin onto SuperFrost Plus Adhesion slides (VWR) with a Cytospin4 unit (Thermo/Shandon), and immunofluorescence staining was performed afterwards.

Colony formation assay

Request a detailed protocol

YFP+ epidermal cells were FACS isolated from E16.5 embryos and cell viability is examined by Trypan blue staining. 75,000 living cells were plated in each well of 6-well plates covered with mitomycin C-treated dermal fibroblasts. Cells were cultured in E-media supplemented with 15% (vol/vol) fetal bovine serum and 300 µM Ca2+ in an incubator with 7% CO2 in the air for 10 days. The cells were then fixed with 4% PFA in PBS and stained with chicken anti-GFP/YFP (1:1000, Abcam, ab13970) primary antibody and secondary antibody conjugated to Alexa Fluor 488. Image of each well was scanned with a Cytation 5 Cell Imaging Multi-Mode Reader (BioTek Instruments).

Thin-layer chromatography

Request a detailed protocol

Measurement of m6A/A ratio in Poly(A)+ RNAs was performed as described in Kruse et al., 2011. YFP+ epidermal cells were FACS isolated from control and cKO E16.5 embryos and collected in TRIzol-LS (Sigma-Aldrich, T3934) for total RNA extraction. Two rounds of Poly(A)+ RNA enrichment were performed with the Dynabeads mRNA Purification Kit (Thermo Fisher Scientific, 61006). 100 ng of the purified Poly(A)+ RNA was used for each biological replicate (pool of embryos with the same genotype) to perform the assay. The signal on the membrane was detected with a Typhoon Trio PhosphorImager (GE Healthcare) and quantified with Fiji (ImageJ).

TEWL measurements

Request a detailed protocol

TEWL rate was assessed by Tewameter TM300 (Courage + Khazaka electronic GmbH as described in Quiroz et al., 2020). Basically neonates were sacrificed and their back skin was dissected and immediately spread over a clean, smooth surface. Over four TEWL measurements were collected on each piece of fully acclimatized skin. The values reported by the instrument were then normalized to the average values from the corresponding control samples to get relative water loss rate.

Quantitative PCR

Request a detailed protocol

YFP+ epidermal cells were FACS isolated from E16.5 embryos and collected in TRIzol-LS (Sigma-Aldrich, T3934). Total RNA were extracted with Direct-zol RNA Miniprep kit (Zymo Research, R2051) and treated with the RQ1 RNase-free DNase (Promega, M6101) to remove DNA contamination. cDNA was generated with the High-Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific, 4368814) and qPCR was performed with SYBR green PCR Master Mix (Thermo Fisher Scientific, 4367660) on a QuantStudio 6 Flex Real-Time PCR System (Thermo Fisher Scientific), with Tbp mRNA as internal loading control. The numbers were further normalized to the average values from the corresponding control samples to get relative mRNA levels. Sequences of the primers are listed in Supplementary file 6.

m6A individual-nucleotide-resolution cross-linking and immunoprecipitation

Request a detailed protocol

miCLIP was performed as described in Grozhik et al., 2017. Integrin α6-high epidermal cells were FACS isolated from wild-type P0 pups and collected in TRIzol-LS (Sigma-Aldrich, T3934) for total RNA extraction. Poly(A)+ RNA was then extracted with the Dynabeads mRNA Purification Kit (Thermo-Fisher, 61006) and treated with the RQ1 RNase-free DNase (Promega, M6101) to remove DNA contamination. 7–8 µg Poly(A)+ RNA was used for each biological replicate (pool of 3 litters of pups) to perform miCLIP. After fragmentation, 1/10 of the sample was saved as input to perform parallel library construction without CLIP. Libraries were sequenced on the Illumina Hi-seq platform to generate paired-ended 50 bp reads. Sequencing data was processed as described in Grozhik et al., 2017. Sequence composition nearby the m6A sites identified by CIMS analysis was analyzed with seqLogo/R package. Normalized-to-input uTPM over identified m6A sites (Lawrence et al., 2013) was counted using the GenomicRanges Bioconductor/R package. Basically a 21-nucleotide window centered at each identified m6A site was selected. The ratio between the total uTPM within the window counted from the miCLIP dataset and that from the corresponding input dataset was recorded as the normalized-to-input uTPM. The m6A sites were assigned to genes using the ChIPseeker package (Yu et al., 2015). All exon annotation as well as the nested coding sequence and UTR annotation of exons was extracted from our gene models in GTF format using the GenomicFeatures Bioconductor/R package (Lawrence et al., 2013). SN-uTPM for individual transcripts was then calculated as the sum of signal for all normalized-to-input uTPMs within the exons of a specific transcript. Visualization of m6A site distribution was performed as described in Olarerin-George and Jaffrey, 2017. Gene set enrichment analysis of m6A levels was performed using the fgsea Bioconductor package (Sergushichev, 2016) and the C2 (Curated pathways) and mSigDB C5 (GO)gene sets (Subramanian et al., 2005).

Correlation analysis between the miCLIP and ribosome profiling data

Request a detailed protocol

The m6A abundance levels for transcripts were calculated as the sum of the normalized m6A counts within a gene's identified m6A sites over the normalized paired RNA-seq counts within these sites. Genes were binned into quintiles based on these normalized m6A abundance levels. Translational efficiencies were retrieved from the lab’s prior in vivo ribosomal profiling of neonatal skin progenitors (Sendoel et al., 2017). The translational efficiency of the top and bottom 20% gene quintiles were compared visually by empirical cumulative distribution function (ECDF) plots and the significance of observed differences between these quintiles was assessed by Wilcoxon rank sum tests.

To capture putative functional enrichment for pathways in the genes which are both highly m6A-modified and have a high translational efficiency in wild-type cells, genes were binned again into quintiles by translational efficiency. Genes within the top 20% of translational efficiency and top 20% normalized m6A modification level were tested for functional enrichment of KEGG pathways using the GOseq R/Bioconductor package.

For the comparison of translational efficiency in functional gene sets, gene set members were compared to all genes outside the functional set by both visual inspection in ECDF plots and statistical testing by Wilcoxon rank sum tests.

Single-cell RNA sequencing

Request a detailed protocol

YFP+ epidermal cells were FACS isolated from E17 embryos and cell viability was examined by Trypan blue staining. For both control and cKO samples, the ratio of living cells was >90%. 8000 of single cells from each sample was processed with the Chromium Single Cell 3’ Library and Gel Bead Kit (V2) (10x Genomics, PN-120267) to prepare scRNA-seq libraries, which were then sequenced on Illumina NextSeq 500 sequenceras 26 × 57 × 8.

The sequencing data were analyzed and aggregated using Cellranger (version 2.1.1) count and Aggr function with default setting, respectively. The aggregated datasets were processed by Seurat (version 2.3.4) (Stuart et al., 2018). The cells expressing less than 1800 genes and the genes expressed in less than 10 cells were removed from further consideration. Counts of genes in each cell were normalized and log10 transformed. The cell cycle phase was estimated by Seurat with default setting. And the mitochondrial transcripts in each cell were also calculated. Then, the data was rescaled to remove the effects of cell cycle and mitochondrial transcripts. The top 2000 variable genes were used to principle component analysis. The first 10 principle components were selected for clustering analysis. The clustering results were projected in t-SNE plots. The cell types were assigned according to the well-known marker gene expressions, including Krt14, Krt15, Bmpr1b, Wnt3a, Fzd10, Krt1, Krt10, Sox9, Trps1, Shh, Lhx2, Dkk4, Fgfr1, Pthlh, Nfatc1, Tgfb2, Wif1, Krt17 and Sox2.

Pseudo-time estimation of lineage differentiation trajectories was performed with monocle (version 2.10.1) (Qiu et al., 2017a; Qiu et al., 2017b; Trapnell et al., 2014). In the beginning, the expression levels of Dkk4, Shh and Lhx2 were applied to cell type differentiation. Then, the cell clustering, pseudo-time and trajectory estimation was based on the user manual of monocle (http://cole-trapnell-lab.github.io/monocle-release/docs/#constructing-single-cell-trajectories).

Differential gene expression analysis and single-cell GSEA were processed by MAST (version 1.8.2, Finak et al., 2015). The comparison between control and cKO in each cluster was processed by default setting with 50 times bootstrap. The genes with p<0.05 and pathways with false discovery rates <0.25 were selected.

Visualization of m6A levels against scRNA-seq Z scores was performed using the ggplot2 R libraries. Functional enrichment analyses of genes with high m6A (top 20%) and upregulated in Mettl3 cKO conditions [Z score (cKO/Ctrl)>1.96] were performed using the weighted hypogeometric test in the GOseq Bioconductor package using the C5 and C2 mSigDB gene sets (Subramanian et al., 2005).

Statistical analyses

Request a detailed protocol

All statistical analyses are provided for each of the individual methods sections. Additionally, statistical and graphical data analyses were performed using Microsoft Excel and Prism 8 (GraphPad) software. For measurements, ≥2 biological replicates and two or more technical replicates were used, where applicable. To determine the significance between two groups, comparisons were made using an unpaired two-tailed Student’s t-test or analysis of variance, as appropriate. Multiple testing correction was done using the Benjamini–Hochberg method. In the box and whisker plots, the middle line represents the median, the upper and lower hinges correspond to the first and third quartiles, and the upper and lower whiskers display the full range of variation (minimum to maximum). Most experiments were repeated on ≥3 pairs of sample and control sets.

Appendix 1

Appendix 1—key resources table
Reagent type
(species) or resource
DesignationSource or referenceIdentifiersAdditional information
Genetic reagent (M. musculus)C57BL/6JJackson labStock #: 000664
Genetic reagent (M. musculus)Mettl3floxPMID:31412241
J. Brüning (Max Planck Institute for Metabolism Research and Policlinic for Endocrinology, Diabetes and Preventive Medicine)
Genetic reagent (M. musculus)Krt14-CrePMID:10411913
E. Fuchs (Rockefeller University)
MGI:1926500
Genetic reagent (M. musculus)Rosa26-YFPfloxJackson labStock #: 006148
Genetic reagent (M. musculus)Gli-LacZPMID:12361967
A. Joyner (MSKCC)
MGI:2449767
Genetic reagent (M. musculus)Nude (Nu/Nu)Charles River LaboratoriesStrain code: 088
AntibodyRat anti-BrdUAbcamCat. #: ab6326
RRID:AB_305426
IMF (1:200)
AntibodyRat anti-CD104BD PharmingenCat. #: 553745
RRID: AB_395027
IMF (1:500)
AntibodyArmenian hamster anti-CD29BioLegendCat. #: 102201
RRID:AB_312878
IMF (1:500)
AntibodyRat anti-CD49fBioLegendCat. #: 313602
RRID:AB_345296
IMF (1:1000)
AntibodyRabbit anti-Cleaved Caspase-3R&D SystemsCat. #: AF835
RRID:AB_2243952
IMF (1:200)
AntibodyRabbit anti-ECADCell SignalingCat. #: 3195
RRID:AB_2291471
IMF (1:500)
AntibodyRabbit anti-FLGE. Fuchs (Rockefeller University)IMF (1:2000)
AntibodyRabbit anti-HES1E. Fuchs (Rockefeller University)IMF (1:200)
AntibodyRabbit anti-iNVBioLegendCat. #: 924401
RRID:AB_2565452
IMF (1:2000)
AntibodyGuinea pig anti-K5E. Fuchs (Rockefeller University)IMF (1:500)
AntibodyRabbit anti-K10CovanceCat. #: PRB-159P-100
RRID:AB_291580
IMF (1:1000)
Flow (1:500)
AntibodyGuinea pig anti-LEF1E. Fuchs (Rockefeller University)IMF (1:5000)
AntibodyRabbit anti-LHX2E. Fuchs (Rockefeller University)IMF (1:5000)
AntibodyRabbit anti-LORE. Fuchs (Rockefeller University)IMF (1:4000)
AntibodyRabbit anti-METTL3AbcamCat. #: ab195352
RRID:AB_2721254
IMF (1:100-1:500)
AntibodyRabbit anti-MYCAbcamCat. #: ab32072
RRID:AB_731658
IMF (1:100)
AntibodyGoat anti-PCADR&D SystemsCat. #: AF761
RRID:AB_355581
IMF (1:300-1:600)
AntibodyRabbit anti-SOX9MilliporeCat. #: AB5535
RRID:AB_2239761
IMF (1:1000)
AntibodyRabbit anti-SurvivinCell SignallingCat. #: 2808
RRID:AB_2063948
IMF (1:500)
AntibodyChicken anti-GFP/YFPAbcamCat. #: ab13970
RRID:AB_300798
IMF (1:1000-1:1200)
Antibodyanti-CD29 PECy7Thermo Fisher ScientificCat. #: 25-0291-82
RRID: AB_1234962
FACS (1:1000)
Antibodyanti-CD31 APCThermo Fisher ScientificCat. #: 17-0311-82
RRID:AB_657735
FACS (1:1000)
Antibodyanti CD45-APCThermo Fisher ScientificCat. #: 17-0451-83
RRID:AB_469393
FACS (1:1000)
Antibodyanti-CD49f PEThermo Fisher ScientificCat. #: 12-0495-81
RRID:AB_891478
FACS (1:1000)
Antibodyanti-CD49f PECy7Thermo Fisher ScientificCat. #: 25-0495-82
RRID:AB_10804881
FACS (1:1000)
Antibodyanti CD117-APCThermo Fisher ScientificCat. #: 17-1172-81
RRID:AB_469432
FACS (1:1000)
Antibodyanti CD140a-APCThermo Fisher ScientificCat. #: 17-1401-81
RRID:AB_529482
FACS (1:500)
AntibodyMouse anti-rabbit IgG PE-Cy7Santa CruzCat. #: sc-516721
Flow (1:300)
Chemical compound, drugAlexa Fluor 647 PhalloidinThermo Fisher ScientificCat. #: A22287
RRID:AB_2620155
IMF (1:50)
Commercial assay or kitIn situ cell death detection kit, TMR redSigma-AldrichCat. #: 12156792910
Commercial assay or kitClick-iT EdU cell proliferation kit for imaging, Alexa Fluor 647 dyeThermo Fisher ScientificCat. #: C10340
Commercial assay or kitLIVE/DEADfixable Aqua dead cell stain kit, for 405 nm excitationThermo Fisher ScientificCat. #: L34957
Commercial assay or kitFixation/ permeabilization solution kitBD BiosciencesCat. #: 554714
Commercial assay or kitDynabeads mRNA purification kitThermo Fisher ScientificCat. #: 61006
Commercial assay or kitDirect-zol RNA miniprep kitZymo ResearchCat. #: R2051
Commercial assay or kitRQ1 RNase-free DNasePromegaCat. #: M6101
Commercial assay or kitHigh-capacity cDNA reverse transcription kitThermo Fisher ScientificCat. #: 4368814
Commercial assay or kitPower SYBR Green PCR Master MixThermo Fisher ScientificCat. #: 4367660
Commercial assay or kitChromium single cell 3' reagent kit10x GenomicsCat. #: PN-120267
Chemical compound, drugTRI Reagent
LS
Sigma-AldrichCat. #: T3934
Software, algorithmZenZeisshttps://www.zeiss.com/microscopy/int/products/microscope-software/zen.html 
Software, algorithmFiji ImageJImageJhttps://imagej.net/Fiji
Software, algorithmAndor IQ3Oxford Instrumentshttps://andor.oxinst.com/products/iq-live-cell-imaging-software/
Software, algorithmGraphPad Prism 8.2GraphPadhttps://www.graphpad.com/
Software, algorithmBD FACSDivaBD Bioscienceshttps://www.bdbiosciences.com/en-us/instruments/research-instruments/research-software/flow-cytometry-acquisition/facsdiva-software
Software, algorithmFlowJoFlowJo, LLChttps://www.flowjo.com/solutions/flowjoVersion 10
Software, algorithmRCRANhttps://cran.r-project.orgR Version 3.6.3
Bioconductor Version 3.10
Software, algorithmFlexbarPMID:24832523https://github.com/seqan/flexbarVersion 2.5
Software, algorithmpyCRACS. Granneman (SynthSys)http://sandergranneman.bio.ed.ac.uk/pycrac-softwareVersion 1.1.3
Software, algorithmNovoalignNovocrafthttp://www.novocraft.com/products/novoalign/Version 3.04.06
Software, algorithmCIMSC. Zhang (Columbia University)https://zhanglab.c2b2.columbia.edu/index.php/CTK_Documentation
Software, algorithmBedtoolsPMID:20110278
A. Quinlan (University of Utah)
https://github.com/arq5x/bedtools2.gitVersion 2.25.0
Software, algorithmSamtoolsPMID:19505943
H. Li (Harvard University)
http://samtools.sourceforge.netVersion 1.2
Software, algorithmCell Ranger10x Genomicshttps://support.10xgenomics.com/single-cell-gene-expression/software/overview/welcomeVersion 2.1.1
Software, algorithmSeuratPMID:29608179
R. Satija (New York University)
https://satijalab.org/seurat/
Version 2.3.4
Software, algorithmMonoclePMID:28114287
C. Trapnell (University of Washington)
http://cole-trapnell-lab.github.io/monocle-release/Version 2.10.1
Software, algorithmMASTPMID:26653891
A. McDavid (University of Rochester Medical Center)
https://www.bioconductor.org/packages/release/bioc/html/MAST.htmlVersion 1.8.2

Data availability

The miCLIP and scRNA-seq data that support the findings of this study have been deposited to the Gene Expression Omnibus (GEO) repository with the accession codes GSE147415, GSE147489, and GSE14749.

The following data sets were generated
    1. Xi L
    2. Fuchs E
    (2020) NCBI Gene Expression Omnibus
    ID GSE147415. Single-cell RNA-seq of embryonic day 17 (E17) mouse skin epithelial cells with or without Mettl3 knockout.
    1. Xi L
    2. Fuchs E
    (2020) NCBI Gene Expression Omnibus
    ID GSE147489. miCLIP-seq of postnatal day 0 (P0) normal mouse skin epithelial cells.
The following previously published data sets were used
    1. Sendoel A
    2. Fuchs E
    (2017) NCBI Gene Expression Omnibus
    ID GSE83332. Epidermis-specific ribosome profiling to describe the translational landscape of SOX2.

References

    1. Bokar JA
    2. Shambaugh ME
    3. Polayes D
    4. Matera AG
    5. Rottman FM
    (1997)
    Purification and cDNA cloning of the AdoMet-binding subunit of the human mRNA (N6-adenosine)-methyltransferase
    RNA 3:1233–1247.
  1. Book
    1. Buechling T
    2. Boutros M
    (2011) Chapter two - Wnt Signaling: Signaling at and Above the Receptor Level
    In: Birchmeier C, editors. Current Topics in Developmental Biology Growth Factors in Development. Academic Press. pp. 21–53.
    https://doi.org/10.1016/B978-0-12-385975-4.00008-5
    1. DasGupta R
    2. Fuchs E
    (1999)
    Multiple roles for activated LEF/TCF transcription complexes during hair follicle development and differentiation
    Development 126:4557–4568.
  2. Book
    1. Grozhik AV
    2. Linder B
    3. Olarerin-George AO
    4. Jaffrey SR
    (2017) Mapping m6A at Individual-Nucleotide Resolution Using Crosslinking and Immunoprecipitation (miCLIP)
    In: Lusser A, editors. RNA Methylation: Methods and Protocols, Methods in Molecular Biology. New York: Springer. pp. 55–78.
    https://doi.org/10.1007/978-1-4939-6807-7_5

Decision letter

  1. Valerie Horsley
    Reviewing Editor; Yale University, United States
  2. Marianne E Bronner
    Senior Editor; California Institute of Technology, United States
  3. Rui Yi
    Reviewer; University of Colorado Boulder, United States

In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.

Acceptance summary:

Overall, this manuscript provides informative data sets of m6A profiles and clear-cut evidence showing the consequences of m6A loss in skin epidermal progenitors. This study highlights the imperative role of m6A modification in modulating fate decision of skin epidermal progenitors.

Decision letter after peer review:

Thank you for submitting your article "m6A impacts fate choices during skin morphogenesis" for consideration by eLife. Your article has been reviewed by Marianne Bronner as the Senior Editor, a Reviewing Editor, and two reviewers. The following individuals involved in review of your submission have agreed to reveal their identity: Rui Yi (Reviewer #1).

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.

As the editors have judged that your manuscript is of interest, but as described below that additional experiments are required before it is published, we would like to draw your attention to changes in our revision policy that we have made in response to COVID-19 (https://elifesciences.org/articles/57162). First, because many researchers have temporarily lost access to the labs, we will give authors as much time as they need to submit revised manuscripts. We are also offering, if you choose, to post the manuscript to bioRxiv (if it is not already there) along with this decision letter and a formal designation that the manuscript is "in revision at eLife". Please let us know if you would like to pursue this option. (If your work is more suitable for medRxiv, you will need to post the preprint yourself, as the mechanisms for us to do so are still in development.)

Summary:

In this manuscript, Xi et al., studied m6A RNA modification and its role in embryonic skin development. They first profiled m6A modification in epithelial progenitors and investigated the correlation between m6A modification and translation efficiency as well as mRNA stability. By conditionally deleting Mettl3, a critical writer for m6A, in epithelial cells during skin development, they examined the physiological impact of Mettl3/m6A modification and observed compromised hair follicle morphogenesis.

Essential revisions:

1) The correlation between m6A and translation efficiency (TE) seems weak, and the way they show this in Figure 1E is confusing. In the Materials and methods section, I couldn't find the definition of "Relative TE" and the cut-off for "high", "med" and "low" m6A modification. As it's shown in Figure 1E, ~50% of "high" m6A mRNA have relative TE greater than 0 and ~50% less than 0. This doesn't support "high" m6A enhances TE. Even in the CDS panel, "med" and "low" m6A transcripts show only very mild reduction of overall TE. They should fully define the parameters used in this analysis; use non-modified transcripts as control; calculate the statistical significance of each comparison e.g. "high" vs "med", "med" vs "low" and all against non-modified transcripts. To truly support the claim that m6A enhances TE, they should perform ribosome profiling in WT and mettl3 ko epithelial cells and calculate the differential TE. And for a few examples, they should quantify protein levels by Western blotting and/or staining to support the increased TE indeed leads to a higher concentration of protein or the reduced protein levels in cKO.

2) They relied heavily on gene ontology analysis to identify pathways that are preferentially affected by m6A. However, this analysis is confusing. They stated they "ranked the mRNAs carrying the modification according to the sum of normalized-to-input uTPM at each m6A site within the coding sequence, and divided by CDS length" (subsection “Skin transcripts most highly modified by m6A are involved in hair follicle 129 morphogenesis”). However, I'm not sure if they used the ranking system at all when they searched for enriched KEGG pathways. It seems that they simply identified pathways with most transcripts with m6A but not necessarily highly modified. They need to show (1) which pathways the most highly m6A modified mRNAs are enriched in? (2) whether the top enriched pathways contain many highly m6A modified mRNAs? In addition, they should apply their TE analysis to these enriched pathways. Can they find reduced TE in these highly modified mRNAs that belong to the top KEGG pathways such as BCC and hedgehog signaling?

3) The phenotypical analysis should be strengthened by more careful analysis of cell proliferation such as colony formation in addition to EdU and Ki67 staining and cell adhesion analysis. The demonstrated phenotypes also seem similar to Dicer1 KO, when all miRNAs are depleted, including neonatal lethality due to the lack of weight gain, degenerated hair follicles and diminished filiform papillae formation.

4) They mentioned that WNT and SHH pathways surfaced in multiple pathways in their miCLIP data and suggested that m6A might function in hair follicle fate through these pathways. However, the authors did not show how these pathways behaved in WT vs cKO mice. Such analysis will strengthen their claims. Aberrant LEF and SHH signaling was associated with failure to form tight dermal condensates. With this data, they concluded that WNT signaling has been altered by m6A modification. However, there were not mechanistic studies to prove that WNT has been altered. The authors can do qPCR analysis of transcripts that are under WNT regulation to further show altered WNT signaling (e.g. axin) or IHC imaging of β-catenin etc. to future strengthen their conclusion.

5) They performed single-cell RNA-seq to measure the changed transcriptome and cellular state. However, their cell clustering analysis should be performed by using control and cKO samples together, rather than clustering separately and relying on marker genes to identify each population in control and cKO (Figure 5—figure supplement 1D, F). Are all corresponding cell populations in control and cKO clustered together? If not, they should explore the underlying reasons. The pseudotime analysis in Figure 4B should be strengthened by statistical analysis rather than eyeballing the population/lineage differences between Ctrl and cKO. They also should keep in mind they only have 1 pair of biological samples, which may reduce the robustness of the analysis. For differential gene analysis, they should independently confirm each data (Figure 4C, Figure 5D etc.) by qPCR or bulk RNA-seq from more sample pairs.

6) Based on ribosome profiling analysis, they suggest that m6A enhances translation. Based on RNA-seq, they suggest that m6A reduces mRNA stability. Globally, these are two opposite effects on mRNA's potential to make protein. It's important to distinguish (1) whether these two opposite effects of m6A are on different mRNAs? (2) if the same mRNAs are under the control of both effects, what is the net effect on protein synthesis? Are they cancelled out or one effect is more dominant than the other? How does this explain/correlate with the phenotypes?

7) The observation of prematurely detached basal cells is interesting. Are they the direct consequence of reduced cell-cell adhesion and/or cell-bm adhesion?

8) The reviewers were not sure about the point of OPP experiments. OPP label may be too low resolution to show any changes in Mettl3 cKO.

9) Figure 3E-F: it is clear that in the grafted skin Mettl3 cKO HFs underwent HF degeneration in accompany with sebocyte differentiation. However, it is not clear if this fate switch only occurred in Mettl3 cKO HFs upon wounding (grafting in this case). Oil red O staining with P6 ungrafted skin sections will partially address this question. If this fate switch only happens upon wounding, authors should discuss the effects of m6A loss in the differential circumstance as well as speculate the effect of wounding in m6A-mediated regulation.

10) Subsection “Conditional ablation of Mettl3 in epidermal progenitors results in a marked defect in HF morphogenesis”: "the engrafted cKO epidermis was hyperthickened, a feature which had also been evident in ungrafted cKO epidermis." To this reviewer, the hyperthicken epidermis could NOT be appreciated in ungrafted cKO epidermis based on the H&E staining shown in Figure 2G. The quantification in the thickness of ungrafted and grafted epidermis will be essential. In line with point #2, this hyperthicken epidermis in grafted cKO skin might be caused by wounding if hyperthickness is not detected in ungrafted Mettl3 cKO epidermis.

11) Figure 5E: Authors stated that MYC expression was elevated within Mettl3 cKO epidermal progenitors (subsection “A cohort of highly m6A-decorated mRNAs whose levels rise upon removal of their modification”). It is not clear to this reviewer how the quantification of MYC expression was done. Which cell markers did authors use to define "epidermal progenitors" (e.g. integrin α6 or K14…)? Did authors also include P-cad+ WNY-hi, WNT-lo HF cells, and epidermal suprabasal cells? Please clarify this by adding information in the figure legend as well as Materials and methods.

12) Figure 3—figure supplement 1: It is interesting to see high levels of PCAD expression in the basal layer of grafted Mettl3 cKO epidermis. Did authors also find this striking pattern in ungrafted Mettl3 cKO epidermis (e.g. P6 ungrafted skin)? What could be the explanation of this significant upregulation of P-Cadherin in the Mettl3 cKO epidermal basal cells, e.g. alteration in cell junction? Confusion in fate specification? It will be helpful if authors can check on ungrafted skin and speculate the possibilities.

https://doi.org/10.7554/eLife.56980.sa1

Author response

Essential revisions:

1) The correlation between m6A and translation efficiency (TE) seems weak, and the way they show this in Figure 1E is confusing. In the Materials and methods section, I couldn't find the definition of "Relative TE" and the cut-off for "high", "med" and "low" m6A modification. As it's shown in Figure 1E, ~50% of "high" m6A mRNA have relative TE greater than 0 and ~50% less than 0. This doesn't support "high" m6A enhances TE. Even in the CDS panel, "med" and "low" m6A transcripts show only very mild reduction of overall TE. They should fully define the parameters used in this analysis; use non-modified transcripts as control; calculate the statistical significance of each comparison e.g. "high" vs "med", "med" vs "low" and all against non-modified transcripts.

We thank the reviewers for these comments and agree completely! Specifically, we have now added a schematic in new Figure 1F to show how the ribosomal and miCLIP analyses were performed on a population of neonatal skin progenitors. In response to the excellent reviewer comment, we also reanalyzed the data, this time taking the miCLIP-measured top 20% and bottom 20% of m6A-modified mRNAs across (1) full-length mRNA, (2) 5’UTR mRNA, (3) coding sequence mRNA and (4) 3’UTR mRNA, and then analyzing them for their relative translational efficiency as measured by ribosomal profiling. All the words are now spelled out to make it easy for the reader to see.

We have also performed the appropriate statistical tests to evaluate shifts in translation efficiency for our newly defined m6A abundance bins. Visualization of the empirical cumulative distribution function (ECDF) of the translation efficiency of these top 20% and bottom 20% quintiles was used to illustrate the correlation of m6A modification level with translational efficiency. To evaluate the statistical significance of these m6A dependent shifts in translation efficiency we have performed Mann-Whitney tests between our successive quintiles. This new analysis confirms and extends our earlier observation that transcripts with high m6A modification levels are more efficiently translated. All these new analyses are currently demonstrated in Figure 1 and Figure 1—figure supplement 1 and described in detail in the Materials and methods section.

To truly support the claim that m6A enhances TE, they should perform ribosome profiling in WT and mettl3 ko epithelial cells and calculate the differential TE.

We discussed this point in our initial response to the reviewers, and our response was accepted and approved by the editors/reviewers back in April. We actually had already tried, but technical hurdles beyond our control prevented us from performing ribosomal profiling on Mettl3 KO epithelial cells either in vitro or in vivo. in vitro, the cells do not survive loss of Mettl3. Even when we culture the cells and then KO Mettl3 in culture, the cells do not tolerate LOF. in vivo, the HF defects result in the placodes and germs fractionating with the dispase-treated dermis and not the epidermis at E17.5 and as development proceeds, the reduced HFs makes the fractionation unequal between KO and control. We now mention these issues in the text and provide the problems (new Figure 5—figure supplement 1A-C) so that other readers will know.

And for a few examples, they should quantify protein levels by Western blotting and/or staining to support the increased TE indeed leads to a higher concentration of protein or the reduced protein levels in cKO.

The inability to obtain sufficient numbers of matched cKO and Control cells (see above) precluded us from not only ribosomal profiling but also immunoblotting. However, as suggested, we have been successful at immunostaining and now add LEF1, β-catenin and Gli1-lacZ for our WNT and SHH analyses, and we provide closeups to show that the HF-DP borders are affected in ways that are expected when WNT signaling is inhibited (Matos et al., 2020). The data are now shown in new Figure 3A-D. By focusing on the top 20% of modified genes that are efficiently translated, the new comparisons also unearthed NOTCH signaling genes as ones highly modified and efficiently translated. We add data on NOTCH but later in the manuscript (new Figure 6D), when we focus on the epidermis. Equally importantly, by examining how the 20% most heavily modified mRNAs fare in scRNA seq, we uncover hitherto unappreciated strong ties between the data in Figure 1 and the single-cell Mettl3 LOF data (new Figure 6).

2) They relied heavily on gene ontology analysis to identify pathways that are preferentially affected by m6A. However, this analysis is confusing. They stated they "ranked the mRNAs carrying the modification according to the sum of normalized-to-input uTPM at each m6A site within the coding sequence, and divided by CDS length" (subsection “Skin transcripts most highly modified by m6A are involved in hair follicle morphogenesis”). However, I'm not sure if they used the ranking system at all when they searched for enriched KEGG pathways. It seems that they simply identified pathways with most transcripts with m6A but not necessarily highly modified. They need to show (1) which pathways the most highly m6A modified mRNAs are enriched in? (2) whether the top enriched pathways contain many highly m6A modified mRNAs? In addition, they should apply their TE analysis to these enriched pathways. Can they find reduced TE in these highly modified mRNAs that belong to the top KEGG pathways such as BCC and hedgehog signaling?

We apologize for the confusion. We actually used the ranking system when we searched for enriched KEGG pathways, but we agree the way it was written and presented was confusing. We have now clarified this both in the figures and in the text. We expanded our description of the normalization used to emphasize that we have used this ranking of m6A modification over input for this analysis. We have added an analysis to demonstrate that KEGG pathways enriched for the highly modified transcripts tend to be more efficiently translated (new Figure 1—figure supplement 1G). We have also added a new analysis for functional enrichment for genes which are highly m6A modified (top 20%) and have a high translational efficiency (top 20%), as was suggested to in point (1) above. These KEGG pathways in Figure 1G now show the top 8 pathways that contain >5 genes within a category and that surface when we analyze the top 20% of highly modified mRNAs that are also among the top 20% of highly translated mRNAs. As you can see, they are enriched for signaling pathways important for skin lineages, similar to before, but will now be much more convincing and interesting to readers. The statistical parameters suggest the extent of the enrichment. We also add new Figure 1H for readers who don’t necessarily know how these pathways pertain to skin development.

We also greatly thank the reviewers for asking how these pathways for highly methylated, highly translated mRNAs fare upon METTL3 loss. We’ve now separated the scRNA data into significantly down as well as significantly upregulated mRNAs. As shown in new Figure 6A,B, the top 20% of genes that are downregulated in our scRNA seq data include signaling pathway genes (including WNT, SHH and NOTCH), actin regulators and cellular adhesion genes and among these are mRNAs that are heavily modified by m6A. Of particular intrigue, we discovered that many canonical translational initiators are among the 20% most significantly downregulated mRNAs upon Mettl3 LOF, offering a possible explanation as to why efficiently translated, heavily modified mRNAs tend to be downregulated (Figure 6B).

In our initial version, we did not include any information on the downregulated, heavily modified genes. These new analyses now beautifully tie together the data in Figure 1 with the scRNA seq data and the phenotypic data.

3) The phenotypical analysis should be strengthened by more careful analysis of cell proliferation such as colony formation in addition to EdU and Ki67 staining and cell adhesion analysis. The demonstrated phenotypes also seem similar to Dicer1 KO, when all miRNAs are depleted, including neonatal lethality due to the lack of weight gain, degenerated hair follicles and diminished filiform papillae formation.

The reviewer raises some excellent points. The cKO cells don’t survive in vitro cell culture, as alluded to above, but we now add the colony formation data (new Figure 5—figure supplement 1C). We now show the epidermal flux experiment as part of Figure 6 and with additional data in new Figure 6—figure supplement 2. As NOTCH signaling pathway surfaced among the top 20% of mRNAs modified and top 20% of mRNAs efficiently translated and among the top 20% most significantly downregulated upon METTL3 loss, we also confirm that HES1, the main NOTCH target, is downregulated in the differentiating epidermal cells (new Figure 6D). This is further interesting, as it is consistent with precocious outward flux as well as the increase in cellularity. We corroborate increased cellularity at the ultrastructural level (new Figure 6G). Finally, we show that apoptosis is not increased in the basal layer and thus, cannot account for why proliferation basally is up in the METTL3 LOF pups (new Figure 6—figure supplement 2F). Rather, the data are together, consistent with the increased flux. We did investigate cell adhesion. We see a decrease in P-cad:E-cad ratio, which correlates with a reduction in cell migration and increased departure from the basal layer (new Figure 6—figure supplement 2A). However, even though we saw downregulated Col17a1, Itgb4 at the mRNA level, we saw hemidesmosomes (Figure 6G). We leave open in the text now the possibility that a decrease in attachment could contribute to enhanced flux.

It is intriguing that to some extent, the phenotype resembles that of mice, whose skin lacks miRNAs. In 2015, Alcaron et al., (2015) used a breast cancer cell line to show that METTL3 methylates pri-miRNAs, marking them for recognition and processing by DGCR8, and that METTL3 depletion results in the global reduction of mature miRNAs and concomitant accumulation of unprocessed pri-miRNAs. While beyond the scope of our current study, our in vivo Mettl3 LOF coupled with prior Dgrcr8 and Dicer LOF in skin (Yi et al. 2008; 2011) point to an interesting parallel. That said, the skin phenotypes are different in the most fundamental way, and that is that when miRNAs are lost, HFs invaginate, something we don’t see in our mice. Since we don’t analyze miRNAs in our paper, it makes discussion a bit problematic. We have now added some text and references to this effect, in our last few paragraphs of the Discussion section.

4) They mentioned that WNT and SHH pathways surfaced in multiple pathways in their miCLIP data and suggested that m6A might function in hair follicle fate through these pathways. However, the authors did not show how these pathways behaved in WT vs cKO mice. Such analysis will strengthen their claims. Aberrant LEF and SHH signaling was associated with failure to form tight dermal condensates. With this data, they concluded that WNT signaling has been altered by m6A modification. However, there were not mechanistic studies to prove that WNT has been altered. The authors can do qPCR analysis of transcripts that are under WNT regulation to further show altered WNT signaling (e.g. axin) or IHC imaging of β-catenin etc. to future strengthen their conclusion.

This was an excellent point. To this end, we first add data for β-catenin IHC (new Figure 3C), higher magnification of the interface between DP and HF progenitors to show defects in both WNT and SHH signaling in New Figure 3C,D, qPCR results for Lef1 in new Figure 6C as requested and most importantly, a previously not included analysis of the downregulated genes from our single cell RNA seq analysis of control and Mettl3 cKO skin progenitors (new Figure 6A,B). This new analysis shows that among the most significantly downregulated genes are those involved in WNT and SHH signaling. Moreover, a number of these were both highly modified and highly translated in wild-type. Finally, we show that among the mRNAs downregulated are many canonical translation genes, providing further possibilities as to why genes highly modified and efficiently translated are likely down upon Mettl3 LOF. These data beautifully strengthen our conclusions and tie together the manuscript. We really thank the reviewers for raising this point!

5) They performed single-cell RNA-seq to measure the changed transcriptome and cellular state. However, their cell clustering analysis should be performed by using control and cKO samples together, rather than clustering separately and relying on marker genes to identify each population in control and cKO (Figure 5—figure supplement 1D, F). Are all corresponding cell populations in control and cKO clustered together? If not, they should explore the underlying reasons. The pseudotime analysis in Figure 4B should be strengthened by statistical analysis rather than eyeballing the population/lineage differences between Ctrl and cKO. They also should keep in mind they only have 1 pair of biological samples, which may reduce the robustness of the analysis. For differential gene analysis, they should independently confirm each data (Figure 4C, Figure 5D etc.) by qPCR or bulk RNA-seq from more sample pairs.

For the clustering of cell populations, we did analyze all cells in control and cKO together, so we agree on this point completely with the reviewers. We have now stated this more clearly in the text, as we agree the way we presented the story was confusing. We also now include more details in the Materials and methods section and also show the t-SNE plots with cells from the two conditions projected together in new Figure 5A. We have also performed statistical analysis for the pseudotime. Results are in new Figure 5C. qPCR and IMF verification for the top downregulated genes are substantiated in Figure 3, Figure 4 and Figure 6. Conversely, t-SNE plots, qPCR and IMF verification for upregulated genes of mRNAs highly modified in WT are now added in new Figure 6C, D and new Figure 7—figure supplement 1B-D.

6) Based on ribosome profiling analysis, they suggest that m6A enhances translation. Based on RNA-seq, they suggest that m6A reduces mRNA stability. Globally, these are two opposite effects on mRNA's potential to make protein. It's important to distinguish (1) whether these two opposite effects of m6A are on different mRNAs? (2) if the same mRNAs are under the control of both effects, what is the net effect on protein synthesis? Are they canceled out or one effect is more dominant than the other? How does this explain/correlate with the phenotypes?

Again, these points were incredibly helpful in our revising the manuscript so that others could clearly grasp the significance of our findings. As alluded to above, the comments made us realize that up until the scRNA seq data, we had focused on phenotypes corresponding to highly modified, efficiently translated mRNAs and which belonged to pathways that were suggested to be downregulated. However, in our initial version, we had not delved into the significantly downregulated mRNAs that surfaced upon scRNA seq of Mettl3 null vs Control transcripts. We’ve now analyzed the most significantly downregulated genes and find that excitingly, these mRNAs encode proteins that are involved in WNT and SHH signaling and actin/ECM/cell adhesion proteins! These data are now presented in a new Figure 6A, B and they add compelling evidence to show that heavily modified, efficiently translated mRNAs involve pathways that are downregulated upon loss of m6A. Additionally, we find that markedly downregulated genes include a number of mRNAs involved in canonical translation (Figure 6B, lower right), further providing a possible explanation for why efficiently translated, highly modified mRNAs might be selectively downregulated upon LOF.

Conversely, the upregulated mRNAs include those which are heavily modified in the wild-type and encode proteins involved in RNA methylation, CAP-independent translation and different aspects of RNA metabolism. So yes, their KEGG profiles and mRNAs are distinct! The data suggest that the upregulated mRNAs are likely involved in (1) rescue mechanisms to boost stalled ribosome machinery and activate alternative pathways for translation; and (2) built-in negative feedback mechanisms for key mRNAs, e.g. the m6A readers (Ythdf1-3) to sense when m6A is low and boost efficiency of reading. Overall, the data beautifully demonstrate that m6A functions in multifaceted ways, which cannot be categorized as simply regulating mRNA translation nor mRNA degradation. The data further show that the functions of down and up regulated mRNAs upon Mettl3 loss are markedly different—in one case affecting signaling and cell fate decisions and in the other case affecting pathways that become triggered when levels of m6A drop. We now highlight these data in the figures (new Figure 7C, D) and in the text, and add a new Figure 7E model to emphasize these points.

7) The observation of prematurely detached basal cells is interesting. Are they the direct consequence of reduced cell-cell adhesion and/or cell-bm adhesion?

Again, this is an excellent comment. We point out that ECM-receptor interactions are among the top 20% of mRNAs modified heavily and top 20% most efficiently translated. Upon loss of METTL3, Focal adhesion and intercellular adhesion are also categories among the most significantly downregulated in cKO and most heavily m6A modified in Control. This could explain the premature detachment. We provide new data that shows that cell polarity and cell organization are altered (new Figure 4) and that P-cad:E-cad is down in the cKO (Figure 6—figure supplement 2A), expected for a decrease in migration/actin dynamics normally seen in basal cells. Our new data point to NOTCH pathway genes as being highly m6A modified and efficiently translated and that NOTCH signaling is downregulated in the cKO epidermal progenitors as they exit the basal layer (new Figure 6D). Additionally, ECM-receptor interaction and actin regulators are also among the most significantly downregulated upon LOF. These data presented in new Figure 6B are supportive of the reviewers’ hypothesis. However, electron microscopy showed that although the epidermis is highly disordered and show increased cellularity and signs of basal features suprabasally, overall hemidesmosomes and intercellular junctions are largely intact, even at P0 (new Figure 6G). This is corroborated by TEWL, showing intact barrier function in newborn mice. Our data together point to the notion that precocious departure of basal progenitors is occurring and is reflected in basal features lingering suprabasally, but these defects are subtle. We attribute this to a remarkable ability of the epidermal progenitors to activate compensatory mechanisms that work for a while until eventually, the deleterious effects of m6A loss cannot be overcome. We thank the reviewers for encouraging us to explore this avenue in depth, and although our answers are not simple, they do provide valuable insights for future explorations.

8) The reviewers were not sure about the point of OPP experiments. OPP label may be too low resolution to show any changes in Mettl3 cKO.

We agree and felt that in the absence of ribosomal profiling, the data is best omitted, which we’ve done. However, our new data in Figure 6 show that canonical translational regulators are significantly downregulated upon METTL3 LOF as are signaling pathways whose mRNAs include among the most highly modified and most efficiently translated (new Figure 6B). Conversely, CAP-independent translational regulators are upregulated as are other Rna processing factors that suggest alternative pathways for translation (new Figure 7C-E). Finally, we also see signs of rescue mechanisms. First is Myc, known to function in enhanced translation (new Figure 7—figure supplement 1B-D). Second is YTHDF1, YTHDF2 and YTHDF3, all of which we show are highly methylated by m6A in WT progenitors and highly upregulated upon m6A loss. As these are key m6A readers—affecting translation (1 and 3) and RNA degradation (2) they too are suggestive of a means to boost m6A recognition in the face of m6A deficiency. These arguments are stronger and more interesting than the OPP data.

9) Figure 3E-F: it is clear that in the grafted skin Mettl3 cKO HFs underwent HF degeneration in accompany with sebocyte differentiation. However, it is not clear if this fate switch only occurred in Mettl3 cKO HFs upon wounding (grafting in this case). Oil red O staining with P6 ungrafted skin sections will partially address this question. If this fate switch only happens upon wounding, authors should discuss the effects of m6A loss in the differential circumstance as well as speculate the effect of wounding in m6A-mediated regulation.

Excellent point! We have performed the Oil red O staining on P6 ungrafted skin (new Figure 3E) and added discussion about the potential influences caused by wounding.

10) Subsection “Conditional ablation of Mettl3 in epidermal progenitors results in a marked defect in HF morphogenesis”: "the engrafted cKO epidermis was hyperthickened, a feature which had also been evident in ungrafted cKO epidermis." To this reviewer, the hyperthicken epidermis could NOT be appreciated in ungrafted cKO epidermis based on the H&E staining shown in Figure 2G. The quantification in the thickness of ungrafted and grafted epidermis will be essential. In line with point #2, this hyperthicken epidermis in grafted cKO skin might be caused by wounding if hyperthickness is not detected in ungrafted Mettl3 cKO epidermis.

We agree and feel that the point is best omitted from the text, even though we did quantify and confirmed the thickness. Our EM data showing increased cellularity is clear and supported by quantifications at the IMF level (new Figure 6G; Figure 6—figure supplement 2B).

11) Figure 5E: Authors stated that MYC expression was elevated within Mettl3 cKO epidermal progenitors (subsection “A cohort of highly m6A-decorated mRNAs whose levels rise upon removal of their modification”). It is not clear to this reviewer how the quantification of MYC expression was done. Which cell markers did authors use to define "epidermal progenitors" (e.g. integrin α6 or K14…)? Did authors also include P-cad+ WNY-hi, WNT-lo HF cells, and epidermal suprabasal cells? Please clarify this by adding information in the figure legend as well as Materials and methods section.

We now clarify this. We identified populations from our scRNA seq data according to known marker gene expressions (Figure 5B; now color coded as in the t-SNE plots in Figure 5A). For our substantiating qPCR, we used basal skin progenitors FACS-purified from whole E16.5 skins as shown in Figure 5—figure supplement 1E. We now state this in the text and methods. For MYC, we quantified the t-SNE data also by IMF and there we focus on basal progenitors at both E17.5 and P0. We moved the MYC data to the supplement as it is confusing to most readers not in the skin field. The Watt group showed that MYC functions in terminal differentiation of the epidermis while others have shown that Myc can act as a WNT target gene in e.g. colon. In the APC LOF skin mouse, Myc is actually downregulated even though WNT signaling is upregulated. Thus, although our data are consistent with the antagonistic effects of WNT signaling and MYC in the skin, it is not an essential point to our work here.

12) Figure 3—figure supplement 1: It is interesting to see high levels of PCAD expression in the basal layer of grafted Mettl3 cKO epidermis. Did authors also find this striking pattern in ungrafted Mettl3 cKO epidermis (e.g. P6 ungrafted skin)? What could be the explanation of this significant upregulation of P-Cadherin in the Mettl3 cKO epidermal basal cells, e.g. alteration in cell junction? Confusion in fate specification? It will be helpful if authors can check on ungrafted skin and speculate the possibilities.

First, we have selected a more representative image for the PCAD staining in the grafted skin (Figure 3—figure supplement 1). Second, we have examined the expression levels of PCAD and ECAD at E16.5 closely. Results are included in Figure 6—figure supplement 2A. Our data show clearly that P-cadherin:E-cadherin ratios shift in the Mettl3 cKO. We also provide data that there are some gaps, visible by EM in the adherens junctions of the cKO vs WT, but these are subtle (Figure 6G). We think that these differences, coupled with the marked downregulation of actin regulators, polarity genes and adhesion genes in the cKO (Figure 6B, Figure 4B) are all suggestive that there are subtle defects here. We show the data and discuss this possibility.

https://doi.org/10.7554/eLife.56980.sa2

Article and author information

Author details

  1. Linghe Xi

    Howard Hughes Medical Institute, Robin Chemers Neustein Laboratory of Mammalian Cell Biology and Development, The Rockefeller University, New York, United States
    Contribution
    Conceptualization, Data curation, Formal analysis, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing - original draft, Project administration, Writing - review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-1239-5316
  2. Thomas Carroll

    Bioinformatics Resource Center, The Rockefeller University, New York, United States
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Validation, Visualization, Methodology, Writing - review and editing
    Contributed equally with
    Irina Matos
    Competing interests
    No competing interests declared
  3. Irina Matos

    Howard Hughes Medical Institute, Robin Chemers Neustein Laboratory of Mammalian Cell Biology and Development, The Rockefeller University, New York, United States
    Contribution
    Conceptualization, Data curation, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing - review and editing
    Contributed equally with
    Thomas Carroll
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-6100-8020
  4. Ji-Dung Luo

    Bioinformatics Resource Center, The Rockefeller University, New York, United States
    Contribution
    Data curation, Software, Formal analysis, Visualization, Methodology
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-0150-1440
  5. Lisa Polak

    Howard Hughes Medical Institute, Robin Chemers Neustein Laboratory of Mammalian Cell Biology and Development, The Rockefeller University, New York, United States
    Contribution
    Investigation, Methodology, Acquisition of data and intellectual input
    Competing interests
    No competing interests declared
  6. H Amalia Pasolli

    Electron Microscopy Resource Center, The Rockefeller University, New York, United States
    Contribution
    Investigation, Visualization, Methodology, Data Analysis and Interpretation
    Competing interests
    No competing interests declared
  7. Samie R Jaffrey

    Department of Pharmacology, Weill Cornell Medicine, Cornell University, New York, United States
    Contribution
    Conceptualization, Resources, Supervision, Writing - review and editing
    Competing interests
    Scientific founder, advisor to, and owns equity in Gotham Therapeutics
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-3615-6958
  8. Elaine Fuchs

    Howard Hughes Medical Institute, Robin Chemers Neustein Laboratory of Mammalian Cell Biology and Development, The Rockefeller University, New York, United States
    Contribution
    Conceptualization, Supervision, Funding acquisition, Methodology, Writing - original draft, Project administration, Writing - review and editing
    For correspondence
    fuchslb@rockefeller.edu
    Competing interests
    Reviewing editor, eLife
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-0978-5137

Funding

National Institutes of Health (R01-AR27883)

  • Elaine Fuchs

National Institutes of Health (R01-CA186702)

  • Samie R Jaffrey

National Institutes of Health (R21-CA224391)

  • Samie R Jaffrey

Damon Runyon Cancer Research Foundation (DRG-2263-16)

  • Linghe Xi

HHMI

  • Elaine Fuchs

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

We thank Dr Jens Brüning and Dr Martin E. Hess at Max Planck Institute for Metabolism Research and Policlinic for Endocrinology, Diabetes and Preventive Medicine, Cologne, Germany for sharing the Mettl3fl mouse strain. We thank Rockefeller University’s Flow Cytometry Resource Center (Svetlana Mazel, director), the Genomics Resource Center (Connie Zhao, director) and the American Association for the Accreditation of Laboratory Animal Care-accredited comparative biology center (R Tolwani, director) for their services. We thank the Epigenomics Core Facility (Yushan Li, director) at Weill Cornell University for sample processing. We also thank members of the Fuchs lab, specifically John Levorse, Brian Hurwitz, Stephanie Ellis, Katherine S Stewart, Matthew Tierney, Felipe Garcia Quiroz, Vincent F Fiore, Shaopeng Yuan, Rachel Niec, Nicholas Gomez, Amma Asare, Aaron Mertz, Sanjee Baksh, Hanseul Yang, Melanie Laurin, Yejing Ge, Maria Nikolova, Ellen Wong, June Racelis, Megan Sribour and Lynette Hidalgo as well as those of the Jaffrey lab, specifically Anya V Grozhik, Brian Pickering, Sara Zacara, Anthony O Olarerin-George, Jan Mauer and Pierre Klein for thorough discussions and for assistance in setting up the various experiments performed in this study. Dr. Robert Roeder’s lab at Rockefeller University also offered us valuable technical support in setting up the miCLIP experiment. LX is the Dale F and Betty Ann Frey Fellow of the Damon Runyon Cancer Research Foundation, DRG-2263–16. EF is an investigator of the Howard Hughes Medical Institute. The work is also supported by grants from the National Institutes of Health (R01-AR27883 to EF, R01-CA186702 and R21-CA224391 to SRJ)

Ethics

Animal experimentation: All mouse strains were housed in an AAALAC-accredited facility and experiments were conducted according to the Rockefeller University's Institutional Animal Care and Use Committee (IACUC), and NIH guidelines for Animal Care and Use. All animal procedures used in this study are described in our #20012-H & #17091-H protocols, which have been previously reviewed and approved by the Rockefeller University IACUC.

Senior Editor

  1. Marianne E Bronner, California Institute of Technology, United States

Reviewing Editor

  1. Valerie Horsley, Yale University, United States

Reviewer

  1. Rui Yi, University of Colorado Boulder, United States

Publication history

  1. Received: March 17, 2020
  2. Accepted: August 25, 2020
  3. Accepted Manuscript published: August 26, 2020 (version 1)
  4. Version of Record published: October 5, 2020 (version 2)

Copyright

© 2020, Xi et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 2,734
    Page views
  • 460
    Downloads
  • 5
    Citations

Article citation count generated by polling the highest count across the following sources: PubMed Central, Crossref, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Developmental Biology
    2. Neuroscience
    Benoît Boulan et al.
    Research Article

    Neurodevelopmental axonal pathfinding plays a central role in correct brain wiring and subsequent cognitive abilities. Within the growth cone, various intracellular effectors transduce axonal guidance signals by remodeling the cytoskeleton. Semaphorin-3E (Sema3E) is a guidance cue implicated in development of the fornix, a neuronal tract connecting the hippocampus to the hypothalamus. Microtubule-Associated Protein 6 (MAP6) has been shown to be involved in the Sema3E growth-promoting signaling pathway. In this study, we identified the Collapsin Response Mediator Protein 4 (CRMP4) as a MAP6 partner and a crucial effector in Sema3E growth-promoting activity. CRMP4-KO mice displayed abnormal fornix development reminiscent of that observed in Sema3E-KO mice. CRMP4 was shown to interact with the Sema3E tripartite receptor complex within Detergent-Resistant Membrane (DRM) domains, and DRM domain integrity was required to transduce Sema3E signaling through the Akt/GSK3 pathway. Finally, we showed that the cytoskeleton-binding domain of CRMP4 is required for Sema3E's growth-promoting activity, suggesting that CRMP4 plays a role at the interface between Sema3E receptors, located in DRM domains, and the cytoskeleton network. As the fornix is affected in many psychiatric diseases, such as schizophrenia, our results provide new insights to better understand the neurodevelopmental components of these diseases.

    1. Developmental Biology
    2. Neuroscience
    Miguel Turrero García et al.
    Research Article

    The septum is a ventral forebrain structure known to regulate innate behaviors. During embryonic development, septal neurons are produced in multiple proliferative areas from neural progenitors following transcriptional programs that are still largely unknown. Here, we use a combination of single cell RNA sequencing, histology and genetic models to address how septal neuron diversity is established during neurogenesis. We find that the transcriptional profiles of septal progenitors change along neurogenesis, coinciding with the generation of distinct neuron types. We characterize the septal eminence, an anatomically distinct and transient proliferative zone composed of progenitors with distinctive molecular profiles, proliferative capacity and fate potential compared to the rostral septal progenitor zone. We show that Nkx2.1-expressing septal eminence progenitors give rise to neurons belonging to at least three morphological classes, born in temporal cohorts that are distributed across different septal nuclei in a sequential fountain-like pattern. Our study provides insight into the molecular programs that control the sequential production of different neuronal types in the septum, a structure with important roles in regulating mood and motivation.