Abstract
N6-methyladenosine (m6A) is the most prevalent internal modification of mammalian messenger RNAs (mRNAs) and long non-coding RNAs. The biological functions of this reversible RNA modification can be interpreted by cytoplasmic and nuclear ‘m6A reader’ proteins to fine-tune gene expression, such as mRNA degradation and translation initiation. Here we profiled transcriptome-wide m6A sites in adult mouse cerebral cortex, underscoring that m6A is a widespread epitranscriptomic modification in brain. Interestingly, the mRNA targets of fragile X mental retardation protein (FMRP), a selective RNA-binding protein, are enriched for m6A marks. Loss of functional FMRP leads to Fragile X syndrome (FXS), the most common inherited form of intellectual disability. Transcriptome-wide gene expression profiling identified 2035 genes differentially expressed in the absence of FMRP in cortex, and 92.5% of 174 downregulated FMRP targets are marked by m6A. Biochemical analyses indicate that FMRP binds to the m6A sites of its mRNA targets and interacts with m6A reader YTHDF2 in an RNA-independent manner. FMRP maintains the stability of its mRNA targets while YTHDF2 promotes the degradation of these mRNAs. These data together suggest that FMRP regulates the stability of its m6A-marked mRNA targets through YTHDF2, which could potentially contribute to the molecular pathogenesis of FXS.
Introduction
Epigenetic regulation has been shown to play pivotal roles in the development and functions of central nervous system. N6-methyladenosine (m6A) is the most prevalent internal modification of mammalian messenger RNAs (mRNAs) and long non-coding RNAs (lncRNAs) (1). A multi-component methyltransferase complex containing several subunits (METTL3, METTL14, WTAP, KIAA1429/VIRMA, RBM15 and RBM15B) installs m6A (2–6). Two demethylases of the AlkB dioxygenase family, FTO and ALKBH5, remove m6A (7,8). The biological functions of this dynamic RNA modification can be interpreted by many cytoplasmic and nuclear proteins (known as ‘m6A readers’) to fine-tune multiple steps in RNA life cycle, such as nuclear pre-mRNA splicing (YTHDC1 and HNRNPA2/B1), nuclear export of mRNA (YTHDC1), microRNA generation (HNRNPA2/B1), translation initiation (YTHDF1, YTHDF3, YTHDC2 and eIF3) and mRNA degradation (YTHDF2 and YTHDC2) (9–19). In addition, m6A can change the accessibility of its surrounding sequence by altering the local structures in mRNAs and lncRNAs, thus modulate the interaction between RNA and RNA binding proteins (a mechanism known as ‘m6A-switch’) (20,21). m6A is particularly abundant in mammalian brain and the methylome is conserved in mammals (22,23). We and others have revealed prominent functions of m6A methylase and demethylase in the proliferation and differentiation of adult neural stem cells, cortical neurogenesis, hippocampal learning and memory, midbrain dopaminergic signaling and circadian clock control (24–28).
Fragile X syndrome (FXS, OMIM#300624) is the most common inherited form of intellectual disability and a leading genetic cause of autism spectrum disorders (ASD) (29–31). FXS is caused by loss of functional fragile X mental retardation protein (FMRP) (32). FMRP, encoded by FMR1 gene, is a selective RNA-binding protein associated with translating polyribosomes (polysome) (33–35). It represses and regulates translation of a set of synaptic plasticity-related transcripts. Coinciding with m6A, FMRP is also particularly abundant in brain with high expression in neurons. Intriguingly, the consensus sequence for m6A methylation, RGACH (R = A/G; H = A or C/U) largely overlaps with FMRP-bound YGGA (Y = C/U) and GAC motifs and is found in 57% of FMRP consensus binding sites (22,23,36–38). Furthermore, FMRP has recently been shown to bind an m6A-containing RNA probe (GGACU)4 (39).
In this study, we determined the transcriptome-wide m6A sites in adult mouse cerebral cortex using anti-m6A RNA immunoprecipitation (RIP) followed by high-throughput sequencing (m6A-seq). We found that the mRNA targets of FMRP are significantly enriched for m6A marks, and the loss of FMRP alters both transcriptome and epitranscriptome in cortex. Surprisingly, among FMRP targets, at least 92% of Fmr1 KO-downregulated mRNAs were marked by m6A [with an odds ratio (OR) near 7], suggesting that FMRP could have interplay with m6A-mediated regulation of mRNA stability. Our biochemical analyses indicated that FMRP could bind to m6A sites in its targets and interact with YTHDF2, an ‘m6A reader’, in an RNA-independent manner. Functionally, we showed that FMRP maintains the stability of its mRNA targets while YTHDF2 promotes the degradation of these transcripts. Together, our works demonstrate that FMRP regulates the stability of its mRNA targets marked by m6A, most likely through YTHDF2. This mechanism represents an additional layer of gene regulation mediated by FMRP, which could potentially contribute to the molecular pathogenesis of FXS.
Results
m6A RNA methylation is a widespread mRNA modification in mouse cerebral cortex
To determine the landscape of m6A in mammalian brain, we profiled the transcriptome-wide m6A sites in the cerebral cortex of 6 week old male mice using m6A-seq (Supplementary Material, Fig. S1A) (28,40). Poly(A)(+) RNA isolated from cortex tissue of six mice were pooled for each assay. m6A-seq identified about 22000 genome-based peaks (false discovery rate, FDR, < 0.001) shared between two biological replicates as high-confidence m6A peaks (or ‘HC peaks’) and 7587 mRNAs and ncRNAs as m6A-marked transcripts (Supplementary Material, Fig. S1B and C and Table S1), suggesting that m6A is a widespread regulatory mechanism for transcripts in adult cortex. In total 93% of HC peaks occur within genic regions (Fig. 1A), which are abundant in exons (75%), and highly enriched for 3’ untranslated region (UTR) (38-fold), coding sequences (CDS) (36-fold) and 5’ UTR (25-fold). In addition, a deduced consensus motif, RGACH (R = A/G; H = A or C/U), was greatly enriched in the HC peaks (Fig. 1B). Two variants of this consensus motif, WGGACH (W = A or U) and RGACU, were both centrally enriched within peaks (Supplementary Material, Fig. S1D). Intriguingly, FMRP target transcripts are significantly enriched for m6A marks. About three quarters of 822 stably expressed FMRP targets (36) (Fragments Per Kilobase of transcript per Million mapped reads, or FPKM, ≥ 0.2) contained at least one HC peak (Fig. 1C, and Supplementary Material, Tables S1 and S2), insinuating that FMRP may directly bind to the m6A sites in its targets as recently suggested by two research groups (39,41), and/or m6A pathway and FMRP may potentially regulate the same set of genes. Overall, m6A marks moderately expressed transcripts (23) and FMRP target transcripts marked by m6A were expressed at lower levels than FMRP target transcripts with no m6A peak (Fig. 1D), suggesting that the significant overlap between m6A-marked transcripts identified by m6A-seq and FMRP targets identified by high-throughput sequencing of RNAs isolated by crosslinking immunoprecipitation (HITS-CLIP) is not predominantly based on transcript abundance. Moreover, gene ontology (GO) analyses indicate that m6A-marked FMRP targets are significantly enriched for GO terms related to nervous system development, synaptic transmission, intracellular signal transduction, learning and behaviors (Fig. 1E).
The mRNA targets of FMRP are significantly enriched for m6A marks in cerebral cortex. (A) Pie chart illustrating the region distribution of HC peaks identified by m6A-seq in 6 week old WT male mouse cortex; promoter-TSS, from −1 kb to +100 bp of transcription starting site (TSS); TTS, from −100 bp to + 1 kb of TTS; 5’ UTR, from +101 bp of TSS to −1 bp of start codon; CDS from start codon to stop codon; 3’ UTR, from +1 bp of stop codon to −101 bp of TTS; ncRNA, non-coding RNA; pseudo, pseudogene. (B) Top two enriched motifs identified by DREME tools with HC peaks from WT cortex. The un-erased E-value of each motif is shown. Out of top 988 scored peaks, 504 and 331 sites were found under WT motifs 1 and 2, respectively. (C) A significant overlap of genes containing HC peaks in WT cortex and genes encoding FMRP targets. Total number of stably expressed genes (FPKM ≥ 0.2) in WT cortex was 15 705. Chi-square test. (D) Cumulative distribution function plot of the expression levels of stably expressed cortical transcripts in different sectors shown in (C). Transcripts were grouped and color-coded based on their sectors. The numbers of transcripts in each sector are shown in parentheses. Two-tailed Mann–Whitney U-test. (E) GO analysis of stably expressed FMRP mRNA targets containing m6A sites in WT cortex. Top 20 enriched biological processes are shown.
Loss of FMRP alters m6A landscape in mouse cerebral cortex
Given the potential link between FMRP and m6A pathway, we then investigated potential role(s) of FMRP in cortical m6A epitranscriptome using Fmr1 knockout (KO) mice. We performed m6A-seq on two biological replicates of pooled poly(A)(+) RNA isolated from the cortex tissue of 6 week old males. m6A-seq in KO cortex identified about 21000 HC peaks as m6A sites, and 7590 mRNAs and ncRNAs as m6A-marked transcripts (Supplementary Material, Fig. S2A and B, and Table S1), close to those numbers in wild-type (WT) cortex. In addition, top two enriched motifs, WGGAC and AVAGACU (V = A/G or C), in the HC peaks from KO cortex were both centrally enriched within peaks and overlapped greatly with the WT motifs (Supplementary Material, Fig. S2C and D).
Next, we compared the m6A peaks identified in WT and Fmr1 KO cortex. About 14000 HC peaks were present in all four m6A-seq datasets from either WT or KO cortex (defined as common peaks) (Fig. 2A and Supplementary Material, Table S1). However, 2858 peaks were found only in WT cortex but not in any replicate of KO cortex (defined as WT-specific peaks), while 2769 peaks were found only in KO but not in any replicate of WT (defined as KO-specific peaks) (Fig. 2A−C, Supplementary Material, Fig. S2E−G and Table S1). Like most RIP-based experiments, m6A signal detection relies on transcripts expressed at sufficient level. We compared the expression levels of 2016 transcripts marked by WT-specific peaks in WT and KO cortex and found that only < 14% of these transcripts were significantly downregulated in KO cortex (Supplementary Material, Fig. S2H, and Tables S1 and S2). Similarly, among 1143 transcripts containing KO-specific peaks, < 14% were significantly upregulated by Fmr1 KO in cortex (Supplementary Material, Fig. S2H and Tables S1 and S2). These findings indicate that the majority of WT- and KO-specific peaks identified by our m6A-seq experiments were not due to altered transcript abundance. In addition, we noticed that the expression levels of transcripts marked by WT-specific peaks and by common peaks were overall very close to each other in both WT and KO cortex, but the transcripts marked by KO-specific peaks in general were expressed at lower level than the transcripts marked by common peaks in both WT and KO cortex (Fig. 2D). Finally, there were 2016 and 1143 genes whose transcripts were marked by WT- and KO-specific m6A peaks, respectively (Fig. 2E and F). These two groups of genes overlapped significantly with stably expressed FMRP targets (Supplementary Material, Fig. S2I). Taken together, these results suggested that the loss of FMRP in cortex could alter m6A epitranscriptome.
![Loss of FMRP leads to altered m6A landscape in mouse cerebral cortex. (A) The number of HC peaks found in both WT and KO cortex (common), or peaks found in both replicates of one genotype but in neither replicate of another (specific). HC peaks of one genotype that were found in only one replicate of another genotype were not shown. (B and C) Examples of WT-specific (B) or KO-specific (C) m6A peak in 6-week-old male mouse cortex. IGV tracks displaying overlays of normalized and combined read coverage from two biological replicates of m6A-seq. Read clusters from WT IP samples, FMR1 KO IP samples and corresponding non-IP (input) samples are shown in blue, red and grey, respectively. The legends below KO track show genomic features. Line, intron; thin box, UTR of exon; thick box, coding region of exon; and arrow, direction of reading frame. The lowermost track marks the location of m6A motif [GAC] or its reverse complementary [GTC]. The data ranges were indicated on the right side of each panel. The blue or red bar below WT track indicates WT-specific or KO-specific m6A sites, respectively. Last: last exon or intron in the transcript. See Supplementary Material, Fig. S2F and G for more examples. (D) Violin and bar plots depicting the expression levels of genes whose transcripts contained common HC peaks or HC genotype-specific m6A peaks in WT or KO cortex. Lower/middle/upper position in bar plots indicate 25/50/75% quantile, respectively. *, P-value < 0.05; ****, P-value < 0.0001; n.s., not significant, two-tailed Mann-Whitney U-test. (E and F) GO analyses of genes whose transcripts contained at least one WT-specific (E) or KO-specific (F) m6A peak in 6 week old cortex.](https://cdn.statically.io/img/oup.silverchair-cdn.com/oup/backfile/Content_public/Journal/hmg/27/22/10.1093_hmg_ddy292/2/m_ddy292f2.jpeg?Expires=1724365936&Signature=ZIrehK4YMedL~3opRcpc2T44ozsYu-fzN1xkCKZVNVLL5T4Iz4oXjYUkXxdo-Q~shbOkA9~vNPeHO0rZUqkMQGG7J2prQ7ROFg0Y6rDfz-QL9PKIifs-bhig5x~6sdQrSWXIp8XILvlO9PbCpouiGmQZc12LQa8WeAn5V9ZeKgHl7mKLbHa6jvc6Cd7Cfu25WT5Vken9x~KNr-hLIcJ~j8Q0aCxDftodQgk0VrJ5aQe2TuWC5W-mxVUjd94LoH4fQDZQRWuj3XjM~JBx~vQHDeBpO72Dj-7tGfU9aFak1pXidq6WpckHVVo8cQ-Q3wZ7QY17V4yN02OKaOzDE3RK5g__&Key-Pair-Id=APKAIE5G5CRDK6RD3PGA)
Loss of FMRP leads to altered m6A landscape in mouse cerebral cortex. (A) The number of HC peaks found in both WT and KO cortex (common), or peaks found in both replicates of one genotype but in neither replicate of another (specific). HC peaks of one genotype that were found in only one replicate of another genotype were not shown. (B and C) Examples of WT-specific (B) or KO-specific (C) m6A peak in 6-week-old male mouse cortex. IGV tracks displaying overlays of normalized and combined read coverage from two biological replicates of m6A-seq. Read clusters from WT IP samples, FMR1 KO IP samples and corresponding non-IP (input) samples are shown in blue, red and grey, respectively. The legends below KO track show genomic features. Line, intron; thin box, UTR of exon; thick box, coding region of exon; and arrow, direction of reading frame. The lowermost track marks the location of m6A motif [GAC] or its reverse complementary [GTC]. The data ranges were indicated on the right side of each panel. The blue or red bar below WT track indicates WT-specific or KO-specific m6A sites, respectively. Last: last exon or intron in the transcript. See Supplementary Material, Fig. S2F and G for more examples. (D) Violin and bar plots depicting the expression levels of genes whose transcripts contained common HC peaks or HC genotype-specific m6A peaks in WT or KO cortex. Lower/middle/upper position in bar plots indicate 25/50/75% quantile, respectively. *, P-value < 0.05; ****, P-value < 0.0001; n.s., not significant, two-tailed Mann-Whitney U-test. (E and F) GO analyses of genes whose transcripts contained at least one WT-specific (E) or KO-specific (F) m6A peak in 6 week old cortex.
Loss of FMRP leads to downregulation of long m6A-marked FMRP mRNA targets. (A) Volcano plot showing transcriptome-wide log2 fold change in gene expression and corresponding P-value. Fmr1 KO-upregulated and -downregulated genes are highlighted in red and blue, respectively. (B and C) 3-D Venn diagrams showing the overlaps between stably expressed m6A genes, FMRP targets and Fmr1 KO downregulated (B) or upregulated (C) genes. Total number of stably expressed genes in WT and Fmr1 KO cortex combined was 16 003. OR, odds ratio. Fisher’s exact test (contingency within II) and Chi-square test (I ∩ III and II ∩ III). (D and E) Cumulative distribution function plots of m6A-marked FMRP mRNA targets based on either the length of CDS (D) or the level of expression (E). Transcripts were grouped and coded based on indicated bins. If a gene has multiple transcript variants, the longest one in RefSeq was used (D). Sig down, transcripts significantly downregulated in cortex by Fmr1 KO. WT, in WT cortex; KO, in Fmr1 KO cortex. nt, nucleotide. Two-tailed Mann–Whitney U-test. (F) Bar plots showing that among 590 m6A-marked FMRP targets, the 161 Fmr1 KO-downregulated ones contained more m6A peaks per transcript than the rest 429 in both WT and Fmr1 KO cortex. (G) Metagene plot showing the distribution of m6A peaks and FMRP-CLIP tags across the length of m6A-marked FMRP target mRNAs. 5’ UTR, CDS and 3’ UTR of Ensembl mRNAs were individually scaled, and the center of exome-based m6A peaks in WT cortex or the center of FMRP-CLIP tags (GSM1098057) was used to calculate the moving average of density. The median widths of exome-based m6A peaks are comparable between the two indicated groups (201 versus 200 nt). If a gene has multiple transcript variants, the longest one in Ensembl was used.
Loss of FMRP leads to the downregulation of long m6A-marked FMRP targets
We next investigated how the loss of FMRP could impact on its targets marked by m6A. Since FMRP is a polyribosome-associated RNA-binding protein (33,36), we applied MeTPeak (42) and MeTDiff (43) to identify exome-based m6A peaks in mature transcripts (taking exon–exon junctions and reads’ variances among replicates into consideration) and to detect differential methylation (DM) among these m6A sites (Fmr1 KO versus WT). In total the transcripts of 5963 genes were found to harbor at least one exome-based m6A peak in both WT and KO cortex (Supplementary Material, Table S3). The MeRIP tag enrichment of exome-based m6A peaks within FMRP target transcripts overall showed little bias towards WT or KO cortex (Supplementary Material, Fig. S3A), and only 3 out of 590 m6A-marked FMRP targets contained differentially methylated m6A peaks (KO versus WT, FDR < 0.05; Supplementary Material, Table S3). But on average the expression level of m6A-marked FMRP targets was greatly reduced by Fmr1 KO (Supplementary Material, Fig. S3B). We postulated that the loss of FMRP could alter transcriptome, so pairwise differential gene expression analysis was performed (KO versus WT). Indeed, a considerable number of differentially expressed genes were identified (Fig. 3A and Supplementary Material, Fig. S3C and Table S2). GO analyses showed that 952 Fmr1 KO-upregulated genes were mainly enriched for GO terms related to translation and energy (ATP) production, while 1083 Fmr1 KO-downregulated genes were enriched for biological processes like transcription, intracellular signaling transduction, epigenetic modification and ion transport (Supplementary Material, Fig. S3D and E). Surprisingly, only those downregulated genes, but not the upregulated ones, showed significant overlaps with both FMRP target genes (II ∩ III) and with genes whose transcripts were marked by m6A (I ∩ III) (Fig. 3B and C). More importantly, more than 92% of 174 Fmr1 KO-downregulated FMRP target mRNAs were marked by m6A with an OR of 6.31, meaning that among stably expressed FMRP targets in KO cortex, m6A-marked ones were significantly downregulated 6.31 times more often than the ones with no m6A (Fig. 3B). Oppositely, if they were marked by m6A, FMRP targets were less likely to be upregulated in KO cortex than in WT cortex (Fig. 3C). To investigate whether these m6A-marked FMRP targets were also downregulated at protein level in Fmr1 KO cortex, we sampled two genes, Ntrk3 (neurotrophic receptor tyrosine kinase 3) and Kcnq2 (potassium voltage-gated channel subfamily Q member 2). Ntrk3 and Kcnq2 were downregulated by 61 and 47% at mRNA level, respectively, by Fmr1 KO (Supplementary Material, Table S2). Western blotting of cortical samples of WT and KO male littermates validated that the proteins encoded by Ntrk3 and Kcnq2 were also downregulated in Fmr1 KO cortex (Supplementary Material, Fig. S3F). So at least for some FMRP targets in brain, the small overall increase in protein synthesis in the absence of functional FMRP (36) was not able to mask the effect of decreased mRNA level.
Next, we scrutinized the 161 genes whose transcripts were m6A-marked, targeted by FMRP and downregulated by Fmr1 KO in cortex. GO analysis showed that these genes were enriched for specific GO terms (Supplementary Material, Fig. S3G) and pathways (Supplementary Material, Fig. S3H). Moreover, among the 590 FMRP targets containing exome-based m6A peaks, we compared the 161 downregulated in Fmr1 KO cortex with the 429 not downregulated to determine the differences in various aspects of mRNAs and m6A modification. Interestingly, the 161 mRNAs downregulated by Fmr1 KO overall are much longer and contain more exons than the 429 mRNAs which were not downregulated (Supplementary Material, Fig. S4A and B). The difference in mRNA length is largely due to the longer CDS in Fmr1 KO-downregulated FMRP targets (Fig. 3D and Supplementary Material, Fig. S4C and D). Furthermore, even in WT cortex, these downregulated mRNAs on average were expressed at lower level than those not downregulated (Fig. 3E). Also, in either WT or KO cortex, the 161 downregulated mRNAs in general harbored more exome-based m6A peaks than the 429 mRNAs not downregulated (Fig. 3F), but tag enrichment (IP/non-IP ratio within m6A peaks) between these two groups showed no significant difference (Supplementary Material, Fig. S4E). Lastly, we compared the distribution patterns of exome-based m6A peaks across the length of these two groups of mRNAs. Metagene plot revealed that m6A peaks from 161 Fmr1 KO-downregulated FMRP targets were more enriched in the CDS than m6A peaks from 429 FMRP targets that were not downregulated by Fmr1 KO. But the distribution patterns of FMRP-CLIP tags (36) in these two groups of mRNAs were similar (predominantly in the CDS) (Fig. 3G). In summary, we discovered that the loss of FMRP in cortex led to significant downregulation of a subset of m6A-marked FMRP target mRNAs sharing characteristics like having long CDS marked by m6A peaks and being moderately expressed. These findings suggested a mechanistic link between the m6A modification on this subset of FMRP targets and the reduced levels of these mRNAs when FMRP was missing.
FMRP binds to m6A sites and is associated with an m6A reader, YTHDF2
The role for m6A in regulating mRNA stability has been well established. YTHDF2, a cytoplasmic m6A reader protein of the YT521-B homology (YTH) domain family, expedites the decay of m6A-marked transcripts by recruiting the CCR4–NOT deadenylase complex (44) and by targeting its cognate mRNAs for degradation in processing bodies (18). Notably, the more m6A/YTHDF2-binding sites one transcript contains, the larger destabilizing effect YTHDF2 will bring to that transcript (18). Since the consensus binding motif of YTHDF2 (18) is almost identical to the m6A motifs found in cortex (Fig. 1B and Supplementary Material, Fig. S2C), and these m6A motifs greatly overlap with FMRP consensus binding motifs (38), we hypothesized that FMRP, as an RNA-binding protein, might bind to the m6A sites in its target mRNAs directly. The first evidence comes from an integrated analysis of our m6A-seq results and previously published matched HITS-CLIP data for FMRP and ELAVL1/HuR/HuA, two selective neuronal RNA binding proteins, in mouse brain (36) (Fig. 4A). FMRP binding is clearly enriched at exome-based m6A sites, consistent with a previous study (39). But ELAVL1 binding is not. About 54% of m6A sites on FMRP target transcripts overlap with FMRP-CLIP tags, while only 28% of m6A sites on ELAVL1 target transcripts overlap with ELAVL1-CLIP tags. Subsequently, the direct binding between FMRP and m6A-marked RNA sequences was tested with an in vitro binding assay using purified c-MYC-tagged FMRP and synthetic RNA probes with methylated or unmethylated m6A motifs. The 40 nt long 5’-biotinylated RNA probes were designed based on an m6A peak within the last exon (E14) of Ntrk3 transcript variant 2 (NM_182809), which was significantly downregulated in Fmr1 KO cortex (Supplementary Material, Table S2). The positions of m6A modification in the probes were designed by considering the summit of the exome-based m6A peak in cortex (Supplementary Material, Table S3) and the location of m6A consensus motifs within the peak, and by referring to an m6A-CLIP (cross-linking immunoprecipitation) site in whole brain (GSM1828593) (45). Anti-FMRP blotting after streptavidin pull down showed that FMRP could directly bind to the m6A-marked RNA probe. However, FMRP also bounds to the same RNA sequence with no methylation on the two indicated adenosines, and its binding affinities for A and m6A were indistinguishable (Fig. 4B and Supplementary Material, Fig. S5A). Furthermore, to investigate if FMRP binding depends on sequence context instead, we mutated the sequence of the biotinylated RNA probe to remove all known FMRP-binding motifs. Although the m6A modifications were unchanged, the binding affinity of FMRP for the 40 nt long mutant probe has been reduced (Fig. 4B).
FMRP binds to m6A sites and is associated with YTHDF2. (A) Density plot depicting the distances from exome-based m6A peaks to their nearest FMRP- or ELAVL1-CLIP tags in WT mouse brain. In total, 68 178 and 98 832 unique FMRP- and ELAVL1-CLIP tags from the same P13 mouse brain lysate (GSM1098057 and GSM1098059), respectively, were used. In total, 1592 m6A peaks in FMRP target transcripts (blue) were analyzed. The distances from 10 305 m6A peaks in transcripts not targeted by FMRP to their nearest FMRP-CLIP tags (green) and the distances from 340 m6A peaks in ELAVL1 target transcripts to their nearest ELAVL1-CLIP tags (red) are shown as controls for specificity. Median size of FMRP-CLIP tags is 46 bp. If an m6A peak overlaps with any fraction of an FMRP-CLIP tag, the distance between them is considered as 0 bp. (B) Western blot showing recombinant FMRP pulled down by 40-nucleotide biotinylated RNA probes. The design of synthetic RNA probes was based on an m6A-marked region within exon 14 of Ntrk3 mRNA (transcript variant 2, NM_182809). The sequences of the probes are provided at the bottom, with two m6A motifs underlined and the other two FMRP binding sequences highlighted in cyan. Red, m6A; blue, mutated bases. (C) Western blotting showing endogenous FMRP co-immunoprecipitated with FLAG-tagged YTHDF2 by anti-FLAG antibodies in N2A cells. The migrations of protein mass standards are indicated on the left. Input, non-IP sample. (D) Anti-FLAG western blotting showing that FLAG-tagged YTHDF2 were co-immunoprecipitated with EGFP-FLAG-tagged FMRP by ant-GFP antibodies in 293T cells. Blots in (B−D) are representative images of two (B) or three (C and D) biological replicates and have been cropped to improve the clarity and conciseness of the presentation.
We further examined whether FMRP could directly interact with YTHDF2. In a co-immunoprecipitation (co-IP) assay, FLAG-tagged YTHDF2 was transfected into mouse Neuro 2a (N2A) cell line. Anti-FLAG IP followed by anti-FMRP blotting demonstrated that endogenous FMRP interacted with YTHDF2 (Fig. 4C). The interaction between these two RNA-binding proteins was not RNA-dependent since the RNase treatment did not abolish their association (Fig. 4C). Furthermore, EGFP-FLAG-tagged FMRP and FLAG-tagged YTHDF2 were co-transfected into human HEK293T cells for reciprocal IP. Anti-GFP IP followed by anti-FLAG blotting confirmed that FMRP interacts with YTHDF2 in an RNA-independent manner (Fig. 4D). Taken together, these results demonstrate that FMRP could bind to m6A-marked RNA and interact with the m6A reader YTHDF2.
FMRP modulates the stability of its mRNA targets marked by m6A through YTHDF2
Given that FMRP could bind to m6A sites and associate with YTHDF2, we hypothesized that the loss of FMRP could accelerate the degradation of those m6A-marked FMRP targets through YTHDF2. To test this hypothesis, we first investigated if FMRP targets in brain are also bound by FMRP in N2A cells. RNA immunoprecipitation (RIP) was performed for FMRP in N2A cell lysate (Supplementary Material, Fig. S5B), and gene-specific reverse transcription quantitative polymerase chain reaction (RT-qPCR) was used to measure the relative level of transcripts co-immunoprecipitated with FMRP. Six m6A-marked mRNAs (Ntrk3; serine/arginine repetitive matrix 2, Srrm2; Kcnq2; tau tubulin kinase 2, Ttbk2; DOT1-like histone lysine methyltransferase, Dot1l; and tetratricopeptide repeat, ankyrin repeat and coiled-coil containing 2, Tanc2), which are FMRP targets in brain and were among top Fmr1 KO-downregulated mRNAs in cortex, were found to be enriched by RIP and associated with FMRP to a similar level as two validated FMRP targets (microtubule associated protein 1B, Map1b; and calcium/calmodulin-dependent protein kinase II α, Camk2a) (29), indicating that these six mRNAs are also targets of FMRP in N2A cells (Supplementary Material, Fig. S5C and D). We then applied RNA interference to N2A cells. Knockdown (KD) Fmr1 using small interfering RNA (siRNA) downregulated the level of the six m6A-marked FMRP target mRNAs (Ntrk3, Srrm2, Kcnq2, Ttbk2, Dot1l and Tanc2) and two FMRP target mRNAs without m6A peak (GRAM domain containing 1B, Gramd1b; and ankyrin repeat domain 52, Ankrd52) in N2A cells, recapitulating the finding in 6 week old cerebral cortex (Fig. 5A and Supplementary Material, Fig. S5E).
FMRP regulates the stability of its m6A-marked mRNA targets. (A) Gene-specific RT-qPCR results showing relative mRNA levels in N2A cells with Fmr1 KD (si-Fmr1). Rpl15 mRNA served as an internal control. Individual data points and mean ± SD of three biological replicates (each with three technique replicates) are shown. *, P-value < 0.05; **, P-value < 0.01, two-tailed one-sample t-test. (B and E) Representative plots for calculating the decay rates of transcripts in N2A cells under si-control or si-Mettl14 (B) and under si-control, si-Fmr1, si-Ythdf2 or si-Fmr1 plus si-Ythdf2 (double KD) (E). Transcription was inhibited by Actinomycin D treatment. Dot1l mRNA abundance was normalized to a spike-in control (luciferase mRNA). Mean ± SD of three technique replicates in one of the three biological replicates are shown. Decay rate equals to the slope of the linear regression line. See Supplementary Material, Fig. S6 for more examples of FMRP targets marked by m6A and negative controls. (C, D and F) Transcript decay rates in N2A cells under si-control or si-Mettl14 (C), under si-control or si-Ythdf1 (D) and under si-control, si-Fmr1, si-Ythdf2 or si-Fmr1 plus si-Ythdf2 (F). Actinomycin D was used to block global transcription. Individual data points and mean ± SEM of two to three biological replicates (each with three technique replicates) are shown. *, P-value < 0.05; **, P-value < 0.01, two-tailed t-test.
Next, we measured the decay rates of mRNAs in N2A cells transfected with siRNA(s). Global transcription in N2A cells was blocked by actinomycin D treatment. The mRNA level in relation to a spike-in mRNA was measured at different time points after transcription inhibition, and the first-order decay rate was determined by linear regression. METTL14 plays a critical role for substrate recognition in the multi-component m6A ‘writer’ complex (3,4). KD or conditional KO of METTL14/Mettl14 greatly reduces m6A levels (2,27,46). As expected, Mettl14 KD in N2A cells significantly reduced the decay rates of the six m6A-marked FMRP targets compared with the scrambled control, but the decay rates of Rpl15 (ribosomal protein L15) and Hprt (hypoxanthine phosphoribosyltransferase) mRNAs, which lack m6A, were unchanged, suggesting that the degradation of those m6A-marked FMRP targets was regulated by m6A (Fig. 5B and C and Supplementary Material, Fig. S5F and G and Fig. S6A−F). Furthermore, on average about 80% KD of Fmr1 significantly increased the decay rates of the six Fmr1 KO-downregulated FMRP targets marked by m6A, in comparison to the scrambled control (Fig. 5E and F and Supplementary Material, Fig. S6G−J), indicating that expedited degradation played an evident role in reducing the steady-state level of these mRNAs during Fmr1 downregulation. However, Fmr1 KD did not accelerate the decay of the two Fmr1 KO-downregulated FMRP targets with no m6A site (Gramd1b and Ankrd52) (Fig. 5F and Supplementary Material, Fig. S6K), suggesting that the stabilizing effect of FMRP for its targets was mediated by m6A. It is worth noting that we have examined three additional Fmr1 KO-downregulated FMRP target mRNAs without m6A peak (mindbomb E3 ubiquitin protein ligase 1, Mib1; chromodomain helicase DNA binding protein 3, Chd3; and FCH domain only 1, Fcho1), but their decay rates could not be determined because none of them showed appreciable amount of degradation at 6 h after actinomycin D treatment.
Lastly, we investigated if the degradation of these m6A-marked FMRP targets was mediated by YTHDF2. On average about three-quarter KD of Ythdf2 reduced the decay rates of the six m6A-marked FMRP targets by at least 50% in comparison to the scrambled control (Fig. 5E and F and Supplementary Material, Fig. S5J, S5K and S6G−J). On the contrary, downregulation of YTHDF1, another cytoplasmic m6A reader protein and a homolog of YTHDF2, did not alter the decay rates of those mRNAs (Fig. 5D and Supplementary Material, Fig. S5H and S5I). Furthermore, Ythdf2 KD abolished the expedited degradation of these m6A-marked FMRP targets induced by Fmr1 KD (Fig. 5E and F and Supplementary Material, Fig. S6G−J). As a negative control, the stability of Hprt mRNA, which was neither an FMRP target (36) nor an YTHDF2 target (18), was not affected by si-Fmr1 and/or si-Ythdf2 (Supplementary Material, Fig. S6L). Taken together, these results suggested that the degradation of the Fmr1 KO-downregulated and m6A-marked FMRP targets was mainly through YTHDF2, and FMRP could stabilize this subset of its targets.
Discussion
FXS is the most common inherited form of intellectual disability and a leading genetic cause of ASD. FXS is caused by the loss of functional FMRP (29–32). FMRP is a selective RNA-binding protein associated with polyribosome (33–35). It has been well established that FMRP could regulate the translation of its mRNA targets. Here we show that in brain a significant portion of FMRP mRNA targets are marked by m6A, the most prevalent internal RNA modification of mammalian mRNAs. Biochemically FMRP could bind to the m6A sites of its mRNA targets and interact with YTHDF2, an m6A reader, in an RNA-independent manner. Mechanistically we demonstrated that FMRP could maintain the stability of its mRNA targets while YTHDF2 promotes the degradation of these mRNAs. These data together provide the first link between FMRP and mRNA stability regulated through m6A, revealing an additional layer of gene regulation mediated by FMRP.
The observation that the consensus sequence for m6A methylation (RGACH) largely overlaps with FMRP-bound YGGA and GAC motifs raised the possibility that FMRP could directly bind to the m6A sites in its targets. It has been shown that the YTH domain of m6A reader proteins, YTHDF1 and YTHDF2, exhibits 16- to 20-fold higher binding affinity for an m6A-marked RNA probe than for an unmethylated one with the same sequence (18). FMRP is a context-dependent RNA-binding protein (47) and three categories of RNA binding motifs are present in FMRP: two Agenet motifs, three KH (heterogeneous nuclear RNP K homology) motifs and an RGG (Arg-Gly-Gly) box (48). However, FMRP does not carry YTH domain. Using RNA-binding assay, we found that FMRP did not display preferentially binding to an m6A-marked RNA probe (based on an in vivo m6A site) over an unmethylated one with the same sequence, albeit a recent study reported that FMRP displayed about 2-fold preference for a 20 nt m6A-marked RNA probe with sequence (GGACU)4 over an unmethylated one (39). Nevertheless, it is still possible that some ‘m6A-switches’ (20,21) in vivo can alter local RNA structures (like stem-loops or G-quartets) to modulate the accessibility of FMRP binding sites. In addition, there is overt difference in the distribution of FMRP binding sites and m6A sites across the length of FMRP target mRNAs. It has been shown that FMRP predominantly binds the CDS with no specific position relative to the start and stop codons (Fig. 3G) (36,38), while m6A sites are abundant not only in the CDS but also in the 3’ UTR and are highly enriched near the stop codon (22). Moreover, for individual FMRP target mRNAs, an even distribution of FMRP binding along the CDS has been observed (36), but on average there are only ∼3−5 m6A per mRNA which most likely will restrict m6A readers to very specific binding sites (49). Taken together, this evidence supports that FMRP can bind to m6A-marked mRNA, but unlike m6A reader proteins, its binding most likely depends on sequence context rather than methylation status.
We found that the loss of FMRP could alter m6A epitranscriptome in cortex, leading to both gain and loss of m6A peaks in transcripts enriched for transcriptional and synaptic functions. Surprisingly, a significant portion of those downregulated mRNAs, but not the upregulated ones, are FMRP targets and almost all of them are marked by m6A.
We also noticed that the m6A peaks from Fmr1 KO-downregulated FMRP targets were more likely to locate in the CDS than the m6A peaks from FMRP targets that were not downregulated by Fmr1 KO. Given the fact that FMRP binding sites located predominantly in the CDS, this finding implies higher level of interplay between FMRP and m6A reader in regulating those Fmr1 KO-downregulated FMRP targets. Using mouse neuronal cells, we demonstrated that FMRP maintains the stability of its long mRNA targets marked by m6A, while m6A marks and YTHDF2 promotes their degradation, suggesting a competition between FMRP and YTHDF2 for binding to the m6A sites on these mRNAs. In agreement with the literature (13), YTHDF1 did not seem to regulate the stability of those mRNAs. In addition, we showed that FMRP interacts with YTHDF2 in an RNA-independent manner. Thus, it is possible that FMRP not only competes with m6A readers for binding to overlapped sites in certain neural mRNAs, but also modulates m6A-mediated regulation of mRNA homeostasis directly via interacting with m6A reader proteins (Supplementary Material, Fig. S7).
FMRP is known to repress the translation of a set of synaptic plasticity-related mRNAs by causing translating ribosomes to stall and accumulate on the mRNA (36). This repression can be relieved by dephosphorylation and ubiquitination of FMRP, presumably in response to synaptic activation, resulting in rapid protein production with precise timing (50,51). Recently, a competition between FMRP and YTHDF1 for binding to m6A sites has been speculated (39). If confirmed, this competition will suggest an antagonism between FMRP and m6A pathway in the control of translational efficiency of mRNAs. In addition, a role for FMRP in the control of mRNA stability has been proposed. For example, FMRP stabilizes one of its targets, Dlg4 (PSD-95) mRNA, in hippocampus but not in cortex (52). FMRP has also been shown to destabilize Nxf1 mRNA in N2A cells (53). However, neither of these two mRNAs contain exome-based m6A peak in cortex (Supplementary Material, Table s3). Herein, we showed that FMRP stabilizes m6A-marked neural mRNAs containing long CDS most likely by impeding m6A/YTHDF2-mediated degradation, and the proteins encoded by these mRNAs are significantly involved in synapses, signaling pathways and long-term potentiation. This novel mechanism, together with m6A-regulated local translation (54) and potential interplay among FMRP and other m6A readers in YTH domain family, could further ensure timely local protein synthesis from existing FMRP target mRNAs upon the stimulation.
Our m6A profiling in adult mouse cerebral cortex revealed more than 47% of stably expressed cortical genes whose transcripts contain high-confidence m6A sites, corroborating that m6A is a widespread epitranscriptomic modification in cortex. FXS not only represents the most common inherited form of intellectual disability but also is the most common monogenic cause of ASD (31). FXS patients often display autistic behaviors, with about 70% of them meeting diagnostic criteria for ASD (55). A significant overlap of FMRP targets and 117 autism-related genes has been reported (36). To date, among 870 mouse homologs of autism candidate genes in the SFARI database, 172 (20%) of them are regulated by FMRP (a significant overlap, P-value = 1.9E-97, Chi-square test) with 84% (144 genes) marked by m6A in this study (Supplementary Material, Tables S1 and S2). Notably, within 144 ASD-related and m6A-marked FMRP targets, one-third (48 genes) were significantly downregulated by Fmr1 KO, in contrast to only 2% (3 genes) upregulated (Supplementary Material, Tables S2). This finding supports that FMRP-mediated regulation on mRNA stability may be relevant to the autistic features in FXS. Since a stringent set of 842 FMRP target transcripts and stringent m6A peak calling algorithms were used in our study, the true number of mRNAs whose decay is modulated by both FMRP and m6A pathway may be underestimated.
In summary, we show that m6A is a widespread epitranscriptomic modification in mammalian cortex, and FMRP modulates the stability of its mRNA targets in concert with m6A pathway. Also, we revealed dysregulated transcriptome and epitranscriptome in the absence of FMRP, providing a resource for studying molecular pathogenesis of FXS. The work presented here represents the first link between FMRP and m6A-regulated neural mRNA stability, which could contribute to the pathogenesis of FXS.
Materials and Methods
Animals
Control WT mice (C57BL/6J, Jackson Laboratory, Bar Harbor, ME, stock # 000664) and Fmr1 KO mice (generated on the same background, C57BL/6J Fmr1tm1Cgr, available from Jackson Laboratory, stock # 003025) were housed, maintained and euthanized according to the Emory University Institutional Animal Care and Use Committee guidelines. Fmr1 KO males and their WT male littermates were obtained by crossing WT males with Fmr1tm1Cgr heterozygous females. Genotyping was determined by standard PCR using a mutant reverse primer (5′-TTCGACCACCAAGCGAAACATC-3′) or a WT reverse primer (5′-TCCAGCTTGATCTTATGGAAAG-3′), with a common forward primer (5′-GTGATAGAATATGCAGCATGTG-3′). Hemizygous mutant males and WT males generate ∼750 bp and ∼120 bp PCR products, respectively.
m6A-seq
Six week old male mice were euthanized using CO2. Total cellular RNA was isolated from fresh mouse cortex using TRIzol Reagent (Thermo Fisher, Waltham, MA) and Phase Lock Gel, Heavy (QuantaBio, Beverly, MA) per manufacturers’ protocols. For each biological replicate, total RNA samples isolated from the cortex of six 6 week old male mice were pooled. mRNA isolation and clean up, m6A-RIP, RNA-seq library preparation and next-generation sequencing were performed as previously described (28).
Mammalian cell culture and siRNA transfection
N2A mouse neuroblastoma cells (ATCC, Manassas, VA, CCL-131) were grown in Minimum Essential Medium (Corning, Corning, NY, 10-009-CV), supplemented with 10% (vol/vol) fetal bovine serum (FBS) (Atlanta Biologicals, Flowery Branch, GA) and 100 units/mL penicillin plus 100 μg/mL streptomycin (Thermo Fisher). 293T human embryonic kidney epithelial cells expressing SV40 T antigen (ATCC, CRL-3216) were grown in Dulbecco’s Modified Eagle Medium (DMEM) (Thermo Fisher, 11995065), supplemented with 10% FBS, 2 mm L-glutamine (Thermo Fisher) and 100 units/mL penicillin plus 100 μg/mL streptomycin. The genotypes of transformed cell lines were verified by ATCC, and only cells with low passage number (< 20) were used. All cells were maintained in a humidified incubator at 37°C in an atmosphere of 5% CO2. Dharmacon ON-TARGETplus SMARTpool mouse Fmr1 siRNA (L-045448-00), mouse Ythdf1 siRNA (L-050937-01) and mouse Ythdf2 siRNA (L-058271-01) were used to downregulate Fmr1, Ythdf1 and Ythdf2 in N2A cells, respectively. Qiagen negative control siRNA (1027310) was used as the scrambled control. For RNA interference, 25 nM siRNA oligos were transfected into N2A cells (for example, 0.5 million cells plated about 16 h ago in 2 mL of growth medium) using Lipofectamine RNAiMAX (Thermo Fisher) per manufacturer’s protocol. One day later, transfected cells were split (1:3) with fresh growth medium.
Western blotting
Six week old male mouse cortex was homogenized and lysed in RIPA buffer (25 mm Tris•HCl pH 7.6, 150 mm NaCl, 1% NP-40, 1% sodium deoxycholate and 0.1% SDS) supplemented with cOmplete Protease Inhibitor Cocktail (Roche, Basel, Switzerland). The protein concentration of cleared lysate was determined by Pierce BCA Protein Assay Kit (Thermo Fisher) and a Synergy H4 Multi-Detection Microplate Reader (BioTek, Winooski, VT) per manufacturer’s instructions. For RNAi experiments, cell lysate was prepared at 72 h after transfection. Lysate samples were denatured in 2× Laemmli Sample Buffer (Bio-Rad, Hercules, CA) at 95°C, separated by SDS–polyacrylamide gel electrophoresis (PAGE) on Mini-PROTEAN TGX Precast Mini Gels (Bio-Rad), transferred onto 0.2 μm PVDF membrane by Trans-Blot Turbo Blotting System and Mini Transfer Pack (Bio-Rad) and probed with primary antibodies against FMRP (1:2000, rabbit pAb, abcam, Cambridge, United Kingdom, ab17722), TrkC (1:1000, rabbit mAb (C44H5), Cell Signaling Technology, Danvers, MA, 3376), Kv7.2 (1:500, rabbit mAb (D9L5S), Cell Signaling Technology, 14752), METTL14 (1:2000, rabbit pAb, Sigma-Aldrich, St. Louis, MO, HPA038002), YTHDF1 (1:1000, rabbit pAb, proteintech, Chicago, IL, 17479-1-AP), YTHDF2 (1:1000, rabbit pAb, abcam, ab170118), FLAG tag (1:2000, mouse mAb (M2), Sigma-Aldrich, F3165 or rabbit pAb, Sigma-Aldrich, SAB1306078) and GAPDH (1:10000, mouse mAb (6C5), Santa Cruz Biotechnology, Dallas, TX, sc-32233) in 5% of nonfat dry milk in 1× TBST (0.1% Tween 20) at 4°C overnight. Primary antibodies were labeled by HRP-linked secondary antibodies (1:4000, Cell Signaling Technology). In western blotting analyses of anti-FMRP RIP, the primary anti-FMRP rabbit pAb was labeled by TrueBlot Anti-Rabbit IgG HRP (1:1000, Rockland Immunochemicals, Limerick, PA) to avoid the detection of denatured anti-FMRP antibody used for IP. The enhanced chemiluminescent signals were detected using HyGLO Quick Spray Chemiluminescent HRP Antibody Detection Reagent (Thomas Scientific, Swedesboro, NJ) and were captured by either ChemiDoc Touch Imaging System (Bio-Rad) or Hyblot CL Autoradiography Film (Denville Scientific). Densitometric quantification of western blotting images was performed using ImageJ/Fiji (56). No pixel used for quantification was overexposed.
Quantitative RT-PCR
Total cellular RNA was purified from N2A cells (72 h after transfection) using TRIzol Reagent and Direct-zol RNA MiniPrep Kit (Zymo Research, Irvine, CA) per manufacturers’ protocols. RNA concentration was determined by a NanoDrop 2000 Spectrophotometer (Thermo Fisher), and 1 μg of total RNA was used to synthesize first-strand cDNA. SuperScript III First-Strand Synthesis System for RT-PCR (Thermo Fisher) and random hexamer priming were used. SYBR Green qPCR was performed as previously described (57). TaqMan gene expression assay was used for Neat1: 1 μL of sequence-specific TaqMan gene expression primers (20×; Thermo Fisher, Mm01720914_g1), 1 μL of reverse transcription product, 10 μL of TaqMan Universal PCR Master Mix II No UNG (2×; Thermo Fisher) and 8 μL of nuclease-free water. All qPCR reactions were performed in triplicates, and the results were analyzed using the ΔΔCT method. The internal control for siRNA-transfected N2A cells was Rpl15. All SYBR Green qPCR primers were designed to span exon–exon junctions (see Supplementary Material, Table S4 for their sequences).
RNA decay
N2A cells transfected with siRNA(s) were re-plated into 12-well plates at 24 h after transfection and were cultured to ∼80% confluence at 72 h after transfection. 5 μm of actinomycin D (Sigma-Aldrich, A9415) was added to inhibit transcription. At 15 min, 2 or 4 h after actinomycin D treatment, cells were washed once with PBS and then lysed in 0.35 mL of Buffer RLT (Qiagen, Hilden, Germany) supplemented with 1% β-mercaptoethanol. Cell lysate was passed through a QIAshredder spin column (Qiagen), and 10 μL of homogenized cell lysate was taken for the quantitation of DNA content using Quant-iT PicoGreen dsDNA Kit (Thermo Fisher) and a Synergy H4 Multi-Detection Microplate Reader (BioTek) per manufacturer’s instructions. Before RNA purification, 0.5 ng of luciferase control RNA (Promega, Madison, WI) per 1 μg of DNA content was added to each sample. Total cellular RNA was extracted using RNeasy Micro Kit (Qiagen), and RT-qPCR was performed. For each transcript, ΔCT was calculated by subtracting the average of CT numbers of the spike-in luciferase control RNA from the CT number of the transcript in each sample. And then ΔΔCT for each sample was calculated by subtracting the average of ΔCT for the total RNA isolated at 15 min after actinomycin D treatment. The log2 relative levels of the transcript at the three time points (15 min, 2 and 4 h) were plotted against time, and a linear regression was used to determine the slope, which is the decay rate (k).
co-IP
Human full-length FMR1 tagged with EGFP and FLAG at the N terminus in the FUGW plasmid (58) and pcDNA 3.0 with N-FLAG-tagged human YTHDF2 (18) were prepared with ZymoPURE Plasmid Maxiprep Kit (Zymo Research). For co-IP experiments, N2A cells (in a 10 cm culture dish) were transfected with 20 μg of N-Flag-YTHDF2 plasmid, and 293T cells (in a 10 cm culture dish) were co-transfected with 20 μg of EGFP-FLAG-FMR1 plasmid and 20 μg of N-FLAG-YTHDF2 plasmid, using Lipofectamine 3000 (Thermo Fisher) per manufacturer’s protocol. At 48 h after transfection, cells were lysed using 1 mL of co-IP Lysis Buffer (50 mm Tris•HCl pH 7.5, 300 mm NaCl, 30 mm EDTA and 0.5% Triton X-100) supplemented with protease inhibitor cocktail at 4°C with rotation for 30 min. Cell lysates were cleared by centrifugation at 20000× g at 4°C for 10 min and were incubated with or without RNase A (Sigma-Aldrich, R6148) and with either 10 μg of anti-GFP rabbit polyclonal antibody (abcam, ab290)/anti-FLAG mouse monoclonal antibody (Sigma-Aldrich, F3165) or 10 μg of normal rabbit/mouse IgG (as a control; Millipore, Burlington, MA) at 4°C with rotation overnight. The immunocomplex was captured by 0.9 mg of Dynabeads Protein G (Thermo Fisher) and was washed three times with co-IP Lysis Buffer. Proteins bound to the beads were eluted by 2× Laemmli Sample Buffer for SDS-PAGE followed by western blotting with anti-FLAG mouse monoclonal antibody (Sigma-Aldrich, F3165).
Pull-down assay
5’-Biotinylated RNA probes were custom synthesized by Integrated DNA Technologies (Coralville, IA) and were purified by RNase-free HPLC. 30 μL of Dynabeads M-280 Streptavidin (Thermo Fisher) was washed three times with 1 mL of immobilization buffer (10 mm Tris•HCl pH 7.5, 1 M NaCl, 1 mm EDTA and 0.01% Tween-20). 10 μm of RNA probes were denatured at 68°C for 6 min and then immediately chilled on ice. 50 pmol of denatured RNA probes were immobilized to 30 μL of washed Dynabeads M-280 Streptavidin in 0.2 mL of immobilization buffer at room temperature for 15 min with rotation. After washing beads twice with 1 mL of binding/washing buffer, 0.5 mg of bovine serum albumin (Sigma-Aldrich) in 0.2 mL of binding/washing buffer (20 mm HEPES pH 7.9, 150 mm NaCl, 2 mm MgCl2 and 0.2% NP-40) supplemented with 1 U/μL of SUPERase•In RNase Inhibitor (Thermo Fisher) was added to the beads. After incubation at 4°C for 1 h with rotation, beads were washed twice with 1 mL of binding/washing buffer and then 0.2 mL of binding/washing buffer supplemented with 1 U/μL of SUPERase•In RNase Inhibitor and indicated amounts of C-terminal MYC/FLAG-tagged FMRP recombinant protein (OriGene, Rockville, MD) were added to the beads. In no probe control (mock), 0.2 mL of binding/washing buffer with RNase inhibitor and 6 pmol of recombinant FMRP were added. After binding at 4°C for 1 h with rotation, beads were washed six times with 1 mL of binding/washing buffer (5 min for each wash) and the bound FMRP was eluted by 2× Laemmli Sample Buffer for SDS-PAGE followed by western blotting.
RIP
RIP procedure was adapted from a previous report (59). Two days after plating 1.7 million cells per 10 cm dish, N2A cells in two 10 cm dishes were washed twice with ice-cold PBS, collected by cell scraper into ice-cold PBS and pelleted by centrifuge at 500 g, 4°C for 5 min. The fresh cell pellet was re-suspended with 2 mL of RIP lysis buffer (300 mm NaCl, 50 mm Tris•HCl pH 7.4, 30 mm EDTA and 0.5% Triton X-100) supplemented with protease inhibitor cocktail and 200 U/ml RNaseOUT Recombinant Ribonuclease Inhibitor (Thermo Fisher) and incubated on ice for 10 min. The cell lysate was passed through a small needle (e.g. 26G) several times and centrifuged at 20000 g, 4°C for 15 min. 80 μL of supernatant was saved as input (10%) for RNA extraction, and another 80 μL of supernatant was set aside for protein analysis. 50 μl of Dynabeads Protein G (Thermo Fisher) was washed twice with 0.7 mL of RIP lysis buffer and re-suspended with 150 μL of RIP lysis buffer. 10 μg of rabbit anti-FMRP polyclonal antibody (abcam, ab17722) or 10 μg of normal rabbit IgG (as a control; Cell Signaling Technology) was added to the re-suspended magnetic beads, and the mixture was rotated continuously at room temperature for 30 min. After incubation, beads-antibody complex was washed four times with 0.7 mL of RIP lysis buffer, mixed with 0.8 mL of cleared cell lysate and then rotated continuously at 4°C for 3 h. The tubes were centrifuged briefly and placed on a magnetic separator. 80 μL of supernatant was set aside as flow-through (10%) for protein analysis; the beads were collected, washed five times with 0.7 mL of RIP lysis buffer at 4°C and re-suspended with 0.5 mL of RIP lysis buffer. 50 μL of mixture was saved as IP (10%) for protein analysis. The rest of beads (90%) were transferred into a new tube, collected and re-suspended with 150 μL of proteinase K digestion solution [125 μL of NT2 buffer (150 mm NaCl, 50 mm Tris•HCl pH 7.4, 1 mm MgCl2 and 0.05% NP40), 15 μL of 10% SDS and 10 μL of proteinase K (20 mg/ml; Thermo Scientific)]. In parallel, 25 μL of nuclease-free water, 15 μL of 10% SDS and 10 μL of proteinase K were added to the input sample thawed on ice. Proteinase K digestion was carried out at 55°C for 30 min with constant shaking in a thermomixer (1200 rpm). After the digestion, the tubes were centrifuged briefly and placed on a magnetic separator. The supernatant was transferred into a 2 ml tube of Phase Lock Gel Heavy. 250 μL of NT2 buffer was added to the tube, followed by 400 μL of phenol:chloroform:isoamyl alcohol (25:24:1, pH 6.5; Sigma-Aldrich). Tubes were vortexed for 10 s and centrifuged at 15000 g, 4°C for 10 min. 400 μL of the aqueous phase was transferred into a new tube. RNA was recovered from the aqueous phase by overnight ethanol precipitation with Pellet Paint Co-Precipitant (Millipore) per manufacturer’s instructions. Purified input and immunoprecipitated RNA were used in RT-qPCR. In addition, input, flow-through and IP samples were analyzed by SDS-PAGE and western blotting.
Bioinformatics analyses
High-throughput sequencing data in this study are summarized in (Supplementary Material, Table S5). Genome-based m6A peaks were identified and annotated as previously described (28), except updated versions of Bowtie (v2.2.4) (60) and Homer (v4.7.2) (61) software were used. When the distance between two peak centers is within 100 bp (about the average size of RNA fragments used in m6A-IP), the two peaks are considered as overlapping (Supplementary Material, Fig. S1C). Integrative Genomics Viewer (IGV) (62,63) was used for creating tiled data files (.tdf) from BAM files (.bam). Genome-wide coverage data were scaled by 1000000/total read count and visualized in IGV. DREME (v4.12.0) (64) was used to discover enriched motifs with top 1000 scored HC peaks. Shuffled input sequences were used as control sequences, and only given strand was searched. CentriMo (v4.12.0) (65) was used to search centrally enriched motifs in all HC peaks. Only given strand was searched, and weighted moving average with a 30 nt window was used to smooth plots. Exome-based m6A peaks (fold enrichment ≥ 2) and DM (FDR ≤ 0.05) in m6A-seq data of WT and Fmr1 KO cortex were identified and annotated by R using MeTPeak (v1.1) (42) and MeTDiff (v1.1.0) (43) packages from RNA-seq reads aligned to mouse genome (build mm9) by Tophat v2.0.13 (66). DM sites are defined as the ones with FDR ≤ 0.05 and fold enrichment of the peaks on the hypermethylated side being consistently larger than 2. The FPKM values for gene expression were calculated by Cufflinks (v2.2.1) (67) using UCSC annotation mm9. Stably expressed genes were defined as FPKM ≥ 0.2. Differentially expressed genes were identified by Cuffdiff (v2.2.1) (68) using UCSC mm9 and defined as q-value ≤ 0.05. R was used for statistical processing and visualizing data. MetaPlotR (69), a Perl/R pipeline was used to create the metagene plot. BEDTools closestBED (70) was used to calculate the distance between an exome-based m6A peak and the nearest FMRP- or ELAVL1-CLIP tag. GSM1098057 in the dataset GSE45148 was chosen for FMRP-CLIP, because it contains significantly more tags obtained by Illumina sequencing (36). GSM1098059 is from the matched ELAVL1-CLIP experiment (36). The list of ELAVL1 targets in mouse brain was obtained from the CLIPdb (71) using Piranha peak caller (v1.2.0) (72). GO analysis was performed by the Database for Annotation, Visualization and Integrated Discovery v6.8 (73). GO mapping of biological process was directly annotated by the source database. No parent terms were shown. Significantly enriched GO terms were defined as FDR (Benjamini–Hochberg) ≤ 0.05. Pathway analysis was performed by the WEB-based Gene SeT AnaLysis Toolkit (WebGestalt) (74) update 2017 with KEGG data source (release 1 October 2016). Significantly enriched pathways were defined as adjusted P-value ≤ 0.05. Other bioinformatics analyses were performed in the R computational environment.
Statistical analyses
Investigators were not blinded to in vitro experiments. The differences between two groups were determined with at least three replicates, and equal variances were assumed. Details of statistical analyses are included in figure legends. An alpha level of 0.05 was used for significance. Two-tailed t-test, Mann–Whitney U-test and Wilcoxon matched-pairs signed rank test, and linear regression were performed in Prism 7 (GraphPad Software). Fisher’s exact test, Chi-square test and Pearson’s r test were performed in R. P-values are indicated in figures and figure legends. When calculating statistical significance of overlapping genes sets, all genes stably expressed in our RNA-seq dataset were used as the total set of genes, and the numbers are indicated in figure legend.
Data and code availability
The m6A-seq data reported in this paper have been submitted to Gene Expression Omnibus (http://www.ncbi.nlm.nih.gov/geo/) with accession number GSE107434. Other data including m6A peak calling and annotation, FMRP targets, autism spectrum disorder-associated genes, differential gene expression analyses and overlap analyses are shown in (Supplementary Material, Tables S1–3). All other relevant data are available from the corresponding author. Command used to run Homer findPeaks is shown in (Supplementary Material, Table S1). Commands used to run Bowtie, IGV tools, Tophat, MeTPeak, MeTDiff, Cufflinks, Cuffdiff or MetaPlotR are included in the manuals of corresponding software, which are available on the internet. R code used to perform bioinformatics analyses is available upon request.
Acknowledgements
We thank Nirav Patel at Nonhuman Primate Genomics Core of Emory University/Yerkes National Primate Research Center for the assistance with high-throughput sequencing. pcDNA 3.0 with N-FLAG-tagged human YTHDF2 was a gift from Dr. Chuan He at the University of Chicago.
Conflict of Interest statement. None declared.
Funding
National Institutes of Health (NS051630, NS091859, NS079625 and MH102690 to P.J., HG008935 and NS097206 to P.J. and H.S., NS047344 to H.S.; GM122083 to H.W.); the Simons Foundation Autism Research Initiative (239320 to P.J. and 308988 to H.S.).
References
