A majority of m6A residues are in the last exons, allowing the potential for 3' UTR regulation

Abstract

We adapted UV CLIP (cross-linking immunoprecipitation) to accurately locate tens of thousands of m(6)A residues in mammalian mRNA with single-nucleotide resolution. More than 70% of these residues are present in the 3'-most (last) exons, with a very sharp rise (sixfold) within 150-400 nucleotides of the start of the last exon. Two-thirds of last exon m(6)A and >40% of all m(6)A in mRNA are present in 3' untranslated regions (UTRs); contrary to earlier suggestions, there is no preference for location of m(6)A sites around stop codons. Moreover, m(6)A is significantly higher in noncoding last exons than in next-to-last exons harboring stop codons. We found that m(6)A density peaks early in the 3' UTR and that, among transcripts with alternative polyA (APA) usage in both the brain and the liver, brain transcripts preferentially use distal polyA sites, as reported, and also show higher proximal m(6)A density in the last exons. Furthermore, when we reduced m6A methylation by knocking down components of the methylase complex and then examined 661 transcripts with proximal m6A peaks in last exons, we identified a set of 111 transcripts with altered (approximately two-thirds increased proximal) APA use. Taken together, these observations suggest a role of m(6)A modification in regulating proximal alternative polyA choice.

Keywords: alternative polyadenylation; last exon; m6A-CLIP/IP; microRNA.

Figure 1.

Single-nucleotide resolution achieved by m⁶A-CLIP/immunoprecipitation (IP). (A) m⁶A-CLIP/IP illustration. (B) The m⁶A PS (defined in the text) maps m⁶A sites at its 0 position. “Enrichment of RRACU or RAC” is the fold enrichment of RRACU (red) or RAC (blue) motif density at that position compared with the background motif density (horizontal gray dot line). The background motif density was calculated as the motif density of RRACU or RAC at regions 400 nt upstream of or downstream from PSs where the density value was flattened. The vertical black dotted line represents the PS 0 position. The motif logo with its P-value represents the de novo motif. (C) Substitution of cross-linking-induced mutation sites (CIMSs) maps m⁶A sites at its −1 position. (D) m⁶A-induced truncation sites (MITSs) map m⁶A sites at its +1 position. (E) Summary of m⁶A sites identified by different approaches. “All RACs” is the total number of RACs within m⁶A peak regions. (F) Precisely mapped m⁶A sites are conserved in vertebrates regardless of the mapping approaches. Box plot bars are the first quartile to third quartile; the black line in the middle is the median value. (**) P < 10⁻⁶; (*) P < 0.01, Wilcoxon rank sum test. (G) m⁶A-CLIP/IP mapped a precise m⁶A site in TUG1 RNA by two independent means: at −2 of a deletion CIMS and also at +1 of a MITS. “m⁶A-IP enrich” is m⁶A-IP reads normalized to its input. “m⁶A-CLIP/IP” is by m⁶A-CLIP/IP (the A marked in red). SCARLET (site-specific cleavage and radioactive labeling followed by ligation-assisted extraction and thin-layer chromatography) is by Liu et al. (2013) at the underlined motif. (H) m⁶A-CLIP/IP mapped a precise m⁶A site in BSG mRNA that is located at the +1 position of MITSs. (I) m⁶A-CLIP/IP mapped precise m⁶A sites in TPT1 mRNA. From left to right, these sites were located at −1 of substitution CIMSs, +1 of MITSs, and 0 of PSs and were co-mapped by −1 of substitution CIMSs and +1 of MITSs. See also Supplemental Figures 1 and 2.

Figure 2.

m⁶A is enriched when entering last exons but not at stop codons. (A) Entering the last exon, m⁶A density increased sharply (bottom panel) in contrast to its lagging increase when approaching the stop codon (top panel). Data are from mouse brains. “m⁶A peak density” was calculated as the number of m⁶A peak regions in a 10-nt interval divided by the total number of mRNAs that contained this position. (B) m⁶A is enriched in the last exons but not around the stop codons. The top panel is the m⁶A peak distribution in the region around the stop codons; the two bottom panels are the distribution of the m⁶A peak regions and stop codons around the last exon start. mRNAs were grouped according to their stop codon locations to the last exon start. (C) m⁶A is enriched in the last exons but not around the stop codon when the stop codon is not in the last exon. (Top and bottom panels) mRNA with noncoding last exons anchored at the last exon start or the stop codon. (D) Most exonic m⁶As locate in last exon (mouse brains). The three pie graphs show the relative proportions of m⁶A peaks and the RAC and RRACU motifs in the last exon and other exons, with 100% representing all m⁶A peaks and RAC and RRACU sequences in mRNA. See also Supplemental Figure 3 for human and previously published data.

Figure 3.

Short last exons have a lower m⁶A density (mouse brains). The m⁶A density in the last exons were compared between longer and shorter last exons (long: >850 nt, 5916 genes; short: ≤850 nt, 3624 genes). All of the genes considered are multiexon-coding genes with RPKM (reads per kilobase per million mapped reads) ≥ 1 in the mouse brain. The results showing more proximal m⁶A in transcripts using distal polyA sites are normalized to transcript abundance: “m⁶A peak density” was calculated as the number of m⁶A peak regions in a 10-nt interval divided by the total number of mRNAs that contained this position (the same definition is used throughout the text). The error bar indicates the standard error of the mean at each position. See also Supplemental Figure 4 for human data and RAC motif controls.

Figure 4.

Higher m⁶A levels in the last exons correlates with more distal polyA site usage. (A) m⁶A is more abundant near proximal sites that are used less frequently (black line; left panel); when distal sites are used more, the m⁶A level is higher until the polyA site is reached (black line; right panel). The error bar indicates the standard error of the mean at each position. The shaded area highlights the major area of difference between two groups. (**) P < 10⁻³⁵, Fisher's exact test. “Used more“ means ≥60% usage of all polyA sites in the last exon, and “used less” means ≤40% usage. Other cutoffs, including 50% or 70%, produce the same findings. The expression of mRNAs in regions between proximal and distal polyA sites is adequate for m⁶A detection (RPKM ≥ 1). Mouse brain data are presented here; see also Supplemental Figure 5A for human data. (B) Distal APA sites in the last exons are used more often in brains for coexpressed mRNAs. APA site pairs with statistically significant differences in usage are highlighted in orange (higher usage of distal sites in the brain) or blue (higher usage of distal sites in the liver). FDR = 20%, Fisher's exact test. The same conclusion is true for FDR = 5% (Supplemental Fig. 5B). (C) m⁶A peak regions in the last exons have a higher enrichment value in the brain. “m⁶A peak enrich. value” is the m⁶A-IP reads normalized to the input reads for each m⁶A peak region. m⁶A peak regions with statistically significant differences in “m⁶A peak enrich. value” are highlighted in orange (higher in the brain) or blue (higher in the liver). FDR = 5%, Fisher's exact test. (D) Positional plot of m⁶A peaks around APA sites in the last exons. “m⁶A peak density (brain–liver)” for each position was calculated as the number of m⁶A peak regions that are significantly higher in the brain (orange points in C) in a 10-nt interval divided by the total number of mRNAs that contained this position. (Orange lines) More distal usage; (gray lines) no change in the brain. The expression of mRNAs in regions between the proximal and distal polyA sites is adequate for m⁶A detection in both tissues (RPKM ≥ 1). The error bar indicates the standard error of the mean. The shaded area highlights the major area of difference between two groups. (**) P < 10⁻²⁰, Fisher's exact test.

Figure 5.

The effects of global limitation of m⁶A on alternative polyadenylation. (A) Triple knockdown of METTL3, METTL14, and WTAP reduced m⁶A of polyA⁺ RNAs to 26% in human A549 cells. The three different controls are biological replicates of control siRNA. The four different knockdown samples are as follows: KD1a and KD1b are biological replicates of one siRNA set, KD2a and KD2b are biological replicates of the other siRNA set, and the two siRNA sets generate essentially the same result. (**) P < 10⁻⁷, Student's t-test. (B) m⁶A may contribute to the inhibition of proximal polyA sites in the last exons. (C) The global limitation of m⁶A changed the choice of APA sites in the last exons. APA site pairs with statistically significant differences in usage are highlighted in orange (higher usage of proximal sites in knockdown) or blue (higher usage of proximal sites in the control). FDR = 20%, Fisher's exact test (the same conclusion is true for FDR = 5%, as shown in Supplemental Fig. S6). The mRNAs using two APA sites in the last exon and with m⁶A present around the proximal polyA sites are summarized. See also Supplemental Figures S6 and S7.

Figure 6.

Overlap of m⁶A and Ago-binding sites in 3′ UTRs of mouse brains. (% m⁶A sites in overlap) The percentage of m⁶A sites (±200-nt flank) overlapped with the RNA-binding protein-binding site in mouse brains (mapped by published HITS-CLIP data) (see the Materials and Methods for reference details). The control was the randomly picked RAC motif (±200-nt flank) in the same 3′ UTR as m⁶A. Ten independent controls were performed to estimate the statistics. The error bar indicates the standard error of the mean. (***) P < 10⁻⁵⁰; (**) P < 10⁻¹⁰, one sample t-test.