Widespread Exon Junction Complex Footprints in the RNA Degradome Mark mRNA Degradation before Steady State Translation

5′P Ends Are Enriched in the Canonical EJC Region

If EJCs are able to block 5′ to 3′ decay of RNA and leave footprints in the RNA degradome as ribosomes do, it is expected that enrichment of 5′P ends in the canonical EJC region will be observed. To test this possibility, we analyzed the previously reported 5′P end data of five evolutionarily distant species. Consistent with this prediction, metagene analyses of both Arabidopsis PARE data (German et al., 2008) and rice (Oryza sativa) Degradome-Seq data (Li et al., 2010) revealed that the relative occurrence of 5′P ends was higher in a region 25 to 30 nucleotides upstream of the 3′ end of the exon (Figure 1A). As these two data sets only contain the 5′P ends of polyadenylated degradation fragments, we also examined Arabidopsis PARE data for nonpolyadenylated degradation fragments generated by Nagarajan et al. (2019). Interestingly, 5′P end enrichment was also detected in the same region of nonpolyadenylated degradation fragments but to a lesser degree compared with that in polyadenylated degradation fragments (Supplemental Figure 1A). Analysis of the Arabidopsis 5′P ends obtained using the GMUCT method (Willmann et al., 2014) gave an enrichment pattern in this region similar to that observed for the PARE data (German et al., 2008; Supplemental Figure 1B). Similarly, a 5′P peak spanning the same region was also detected when analyzing the GMUCT data for human HEK293 cells (Willmann et al., 2014; Figure 1A). Analysis of Degradome-Seq data for worm (Park et al., 2013), in which EJCs are not required for triggering NMD (Longman et al., 2007), also revealed an enrichment of 5′P ends in a similar region with a 1 nucleotide offset (Figure 1A). However, analysis of the PARE data for budding yeast (Harigaya and Parker, 2012), which lacks EJCs and has only a few intron-containing genes (Bannerman et al., 2018), revealed no enrichment of 5′P ends in this region (Figure 1A). The common patterns observed in the 5′P end data sets of four evolutionarily distant species indicate that this phenomenon is truly conserved across species possessing EJCs. Given that the canonical deposition site of human EJCs is centered at a position 24 nucleotides upstream of the exon-exon junction (Saulière et al., 2012; Singh et al., 2012) and that deposition of an EJC results in the protection of an 8-nucleotide fragment from complete RNase digestion (Le Hir et al., 2000; Figure 1A, top), the 5′P peaks spanning the −25 to −30 region correspond to the 5′ edges of EJCs deposited in the canonical region. Therefore, these 5′P peaks are very likely to be EJC-protected 5′ termini of RNA degradation intermediates naturally produced in cells.

• Download figure
• Open in new tab
• Download powerpoint

Figure 1.

The Canonical EJC Region Is a Conserved Hot Spot for 5′P Ends.

(A) Distribution of the relative frequency of 5′P end occurrences in a 50- nucleotide region upstream of the exon-exon junction. PARE data for Arabidopsis wild-type inflorescences (German et al., 2008) and yeast (Harigaya and Parker 2012), Degradome-Seq data for rice wild-type seedlings (Li et al., 2010) and adult stage worms (Park et al., 2013), and GMUCT data for HEK293 T cells (Willmann et al., 2014) were used in the metagene analyses. The illustration at the top shows the canonical position of an EJC and the size of a fragment protected by an EJC. The numbers of exons ≥ 50 nucleotides in size and with a 5′P end count ≥ 1 that were included in the analysis are shown in parentheses.

(B) Distribution of the maximum 5′P peaks within the transcripts of intron-containing genes in a 50-nucleotide region upstream of the exon-exon junction. The numbers of intron-containing genes that had the maximum 5′P peak with a count ≥ 3 and were included in the analyses are shown in parentheses.

In addition to the analysis of 5′P end distribution in a 50-nucleotide region upstream of the exon-exon junction, we also investigated the most abundant 5′P peaks (maximum 5′P peaks) occurring in the canonical EJC region (positions −26 to −30) within transcripts of intron-containing genes. Analysis of previously published 5′P end data revealed that 13.1 (2411 genes), 9.4 (1650 genes), 8.2 (710 genes), and 8.1% (1142 genes) of the maximum 5′P peaks within transcripts of Arabidopsis, rice, human cells, and worm, respectively, were located in the canonical EJC region (Figure 1B). These values were significantly higher than expected (P < 10⁻¹⁰⁰, χ² test), as the total length of the canonical EJC regions represents less than 2% of the total length of spliced transcripts. By contrast, the percentages of maximum 5′P peaks falling in the nearby 5-nucleotide regions in Arabidopsis, rice, and human cells were mostly lower than 2%, which is close to the expected value for a 5-nucleotide region calculated by dividing the total length of each 5-nucleotide region by that of spliced transcripts (Figure 1B). Hence, EJC footprints appear to account for a noticeable portion of major mRNA degradation intermediates.

Exoribonucleases Are Involved in in Vivo EJC Footprinting

To investigate whether EJC footprints in the RNA degradome are generated through 5′ to 3′ decay, as was previously reported for ribosome footprints (Pelechano et al., 2015; Yu et al., 2016), we profiled 5′P ends in the Arabidopsis wild type and two mutants defective in 5′ to 3′ RNA decay, xrn4-6 and fiery1-6 (fry1-6), and compared the abundances of 5′P ends mapping to canonical EJC regions. Although 5′P ends remained enriched in the canonical EJC region in the two mutants (Figure 2A), xrn4-6 had significantly fewer 5′P ends in the canonical EJC region than the wild type in a global analysis of Arabidopsis exons (P < 0.001, Kruskal-Wallis test followed by a pairwise Wilcoxon test with false discovery rate [FDR] correction; Figure 2B). The decrease in the number of 5′P ends in the canonical EJC region was more evident in fry1-6 (P < 0.001, Kruskal-Wallis test followed by a pairwise Wilcoxon test with FDR correction), wherein the activities of XRN2, XRN3, and XRN4 are likely all suppressed (Gy et al., 2007). Notably, the 5′P ends enriched around positions −20 to −15 are likely also fragments protected by RNA binding proteins, because the number of 5′P ends mapping to this region also appeared to be reduced in xrn4-6 and fry1-6. By contrast, the abundance of 5′P ends around the transcription start site (TSS) of genes harboring exons was largely increased in xrn4-6 and slightly but significantly elevated in fry1-6 compared with the wild type (P < 0.001, Kruskal-Wallis test followed by a pairwise Wilcoxon test with FDR correction; Figure 2C). This result strongly suggests that EJCs can sterically hinder XRN-mediated 5′ to 3′ decay, resulting in the EJC footprints observed in the RNA degradome. Moreover, fry1-6 exhibited a more pronounced reduction of EJC footprints than xrn4-6, implying that, in addition to the cytoplasm-localized XRN4 (Kastenmayer and Green, 2000), nuclear XRNs (XRN2 and XRN3) or other proteins with exoribonuclease activity may also contribute to in vivo EJC footprinting.

• Download figure
• Open in new tab
• Download powerpoint

Figure 2.

XRNs Are Involved in the Production of 5′P Ends Corresponding to the Canonical EJC Region in Arabidopsis.

(A) Distribution of 5′P ends in a 50-nucleotide region upstream of the exon-exon junction in wild-type, xrn4-6, and fry1-6 seedlings. The distributions for two independent biological replicates (R1 and R2) are shown. The numbers of exons ≥ 50 nucleotides in size with a 5′P end count ≥ 1 that were included in the analysis are shown in parentheses.

(B) and (C) Comparison of 5′P end abundance in the canonical EJC region (B) and at the TSS (C) of the genes harboring the selected exons in Arabidopsis wild-type, xrn4-6, and fry1-6 seedlings. The numbers of EJC regions (count > 10 TP40M) and TSSs (count > 5 TP40M) included in the analysis are shown in parentheses. Different letters above the box plots indicate significant differences between groups (P < 0.001) determined by the Kruskal-Wallis test followed by the pairwise Wilcoxon test with FDR correction. Box boundaries and whiskers indicate quartiles and the range of 5′P end abundance, respectively.

Prominent EJC Footprints Are Prevalent in NMD Targets

Given that ribosomal arrest by a PTC can trigger decapping, followed by 5′ to 3′ RNA decay, the EJCs deposited on the exons downstream of a PTC were predicted to hinder XRN-mediated degradation and lead to the accumulation of degradation intermediates with 5′P ends in the EJC region. Indeed, three Arabidopsis genes, TFIIIA, LPEAT2, and RS2Z33, which are known to produce PTC-containing mRNAs due to alternative splicing (AS; Yoine et al., 2006; Gloggnitzer et al., 2014), displayed prominent 5′P peaks in the EJC regions of the exons 3′ to the PTC (Figure 3; Supplemental Figure 2). Notably, in the wild type, the maximum 5′P peaks in these three genes were all located in the second exon 3′ end downstream of the PTC, although the maximum 5′P peak in LPEAT2 was not in the canonical EJC region but at a position 37 nucleotides upstream of the exon-exon junction. Moreover, the abundances of these maximum peaks, including the two in the canonical EJC region, were reduced in both xrn4-6 and fry1-6 compared with the wild type, but to a greater extent in fry1-6. However, a comparable or higher number of 5′P ends corresponding to the putative TSS and regions outside the canonical EJC site were observed in xrn4-6 and fry1-6 than in the wild type, which excluded the possibility that a change in gene expression accounted for the reduction of EJC footprints in the mutants.

• Download figure
• Open in new tab
• Download powerpoint

Figure 3.

The EJC Regions Downstream of PTCs Harbor SMG7-Dependent 5′P Peaks in Arabidopsis.

Positional distributions are shown for 5′P ends in three Arabidopsis genes producing PTC-containing mRNAs due to AS. The 5′P ends were identified in the PARE data generated by this study. Within the 5′P end plots, red peaks indicate peaks in the canonical EJC region, black arrows indicate maximum peaks in 10-d-old (10d) wild-type seedlings, and blue arrows indicate peaks corresponding to the putative TSS. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, thin lines indicate introns, orange boxes indicate regions altered by AS, and red asterisks indicate PTCs. Only the 5′P ends residing in exons are displayed in the plots. Data for replicate 2 are shown in Supplemental Figure 2.

A previous analysis of AS events suppressed in a double mutant of the NMD factors UPF1 and UPF3 of Arabidopsis revealed a large number of genes producing splicing variants with classical NMD features (Drechsel et al., 2013). We thus used this gene list to test the association between the NMD-suppressed AS events and the occurrence of maximum 5′P peaks in the canonical EJC region. Using a threshold applied in the previous report (upregulated in the double mutant of lba1 upf3-1 with FDR < 0.1; lba1 is also known as upf1-1; Drechsel et al., 2013), 1082 genes with NMD-suppressed AS events were identified. Among them, 37.2% (402 genes) and 32.8% (355 genes) possessed a maximum 5′P peak in the EJC region in replicates 1 and 2 of the Arabidopsis 10-d-old wild-type seedling PARE data we generated, respectively (Figures 4A and 4B). By contrast, the total incidence of maximum 5′P peaks in the EJC region for all intron-containing genes in the Arabidopsis genome was 15.3% (3379 genes) and 13.8% (3064 genes) in replicates 1 and 2, respectively, which is significantly lower than the incidence of genes with NMD-suppressed AS events (P < 10⁻¹⁰⁰, binomial test). This result implies that EJC footprints are the major and common degradation products of putative NMD targets. The predominance of EJC footprints in many putative NMD targets also validated their degradation during the pioneer round of translation. Because it has been shown that many intron-retention events are not sensitive to NMD (Kalyna et al., 2012; Drechsel et al., 2013), we also compared the incidence of maximum 5′P peaks in the EJC region for four types of AS events suppressed by NMD. Interestingly, a significantly higher incidence of maximum 5′P peaks in the EJC region was observed for alternative 5′ or 3′ splice site events than for intron-retention events (P < 10⁻⁵, binomial test), although these events were all suppressed by NMD (Figure 4B). However, for 2306 intron-containing genes harboring a maximum 5′P peak in the EJC region in both replicates of the PARE data we generated, only 293 genes (12.7%) had NMD-suppressed AS events (Figure 4A). A large number of maximum 5′P peaks in the EJC region did not overlap with NMD-suppressed AS events reported previously (Drechsel et al., 2013), implying that either the list was not complete or there are additional pathways involved in the degradation of EJC-bound RNA.

• Download figure
• Open in new tab
• Download powerpoint

Figure 4.

Genes Producing NMD-Suppressed AS Variants More Frequently Possess a Maximum 5′P Peak in the EJC Region.

(A) Venn diagram showing the overlap between genes producing NMD-suppressed AS variants (NMD⁺ AS) and genes possessing a maximum 5′P peak (M5′P) in the EJC region. The NMD-suppressed AS variants were extracted from the lba1 upf3-1 data reported by Drechsel et al. (2013). Two replicates (R1 and R2) of PARE data for 10-d-old (10d) seedlings generated by this study were used for the identification of genes possessing a maximum 5′P peak with a count ≥ 3 located in the EJC region.

(B) Relative frequency of genes possessing a maximum 5′P peak in the EJC region. Frequencies are shown for genes containing introns (intron⁺ gene) and those producing transcripts derived from the following AS events: NMD⁺ AS, intron retention, exon skipping, and alternative 5′ and 3′ splice sites (alt. 5′ss and alt. 3′ss). The number of genes in each category possessing a maximum 5′P peak in the EJC region is shown above the bar. Asterisks indicate significant differences between categories based on the binomial test (*, P < 10⁻⁵; **, P < 10⁻¹⁵; ***, P < 10⁻¹⁰⁰).

(C) Positional distribution of 5′P ends (from PARE data generated by this study) in three Arabidopsis genes producing NMD-suppressed AS variants. Within the 5′P end plots, red peaks indicate peaks in the canonical EJC region and black arrows indicate the maximum peaks in 10-d-old wild-type seedlings. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, thin lines indicate introns, red asterisks indicate PTCs, and orange boxes indicate regions altered by NMD-suppressed AS events according to Drechsel et al. (2013). Only the 5′P ends residing in exons are displayed in the plots. Data for replicate 2 are shown in Supplemental Figure 3.

AS events in CYP38, PPD6, and STR14 were suppressed by NMD (Drechsel et al., 2013) but were not annotated in The Arabidopsis Information Resource (TAIR) database. The maximum 5′P peaks in these three genes were located in the EJC region, and these peaks had a 5′P count ≥ 1000 tags per 40 million (TP40M; Figure 4C; Supplemental Figure 3). Although the 5′P end profiles of these three genes highly resembled those of validated plant miRNA targets reported previously (Addo-Quaye et al., 2008; German et al., 2008), the maximum 5′P peaks were dampened in xrn4-6 and fry1-6. This indicates that these extremely abundant 5′P peaks do not correspond to 5′P ends derived from endonucleolytic cleavage but instead the protected fragments of EJCs. Notably, these peaks were located in the first or third exon 3′ end downstream of the predicted PTCs caused by AS. This again supports the notion that the dominant EJC footprints in the exon 3′ ends downstream of predicted PTCs could be evidence of NMD triggered by those PTCs.

PTC-Triggered RNA Decay Leads to the Production of EJC Footprints

To verify that a PTC conditions the production of EJC footprints, we took advantage of existing nonsense mutants of PHYB and PHO2 (Reed et al., 1993; Delhaize and Randall, 1995), in which a PTC was introduced by ethyl methanesulfonate mutagenesis. The PTCs in the nonsense mutants should direct mRNA decay through the NMD pathway and are predicted to trigger the production of EJC footprints downstream of the PTC. We compared degradation intermediates of PHYB and PHO2 mRNAs in the wild type and the nonsense mutants phyB-9 and pho2 using a modified RNA ligase-mediated (RLM) 5′ RACE protocol. In both nonsense mutants, the 5′ ends of the most abundant amplification products of the mRNA degradation fragments were mainly mapped to positions 27 to 28 nucleotides upstream of the exon-exon junction downstream of the PTCs (Figures 5A and 5B). However, these degradation intermediates were not prominent in the wild type, although the miR159-guided 3′ cleavage product of MYB65, the positive control for modified RLM 5′ RACE, was detected in both the wild type and the nonsense mutants. The modified RLM 5′ RACE results for phyB-9 and pho2 mimic the PARE data of the three known NMD targets with PTCs introduced by AS (Figure 3), and together these results support the conclusion that PTCs elicit the accumulation of EJC footprints in downstream exons. Moreover, PTCs caused by AS or point mutations tended to result in the generation of 5′P ends mainly in the first or second EJC region downstream of PTCs (Figures 3, 4C, 5A, and 5B), suggesting that EJC footprints may serve as an immediate readout of PTC-mediated NMD irrespective of the cause of the PTC.

• Download figure
• Open in new tab
• Download powerpoint

Figure 5.

PTCs Trigger the Formation of EJC Footprints in Downstream Exons.

(A) and (B) Modified RLM 5′ RACE assays of RNA degradation intermediates generated from wild-type and nonsense mutant mRNAs of PHYB (A) and PHO2 (B). Arrowheads indicate the RACE products excised and cloned for Sanger sequencing (left panels). The positional distribution of the 5′P ends of the cloned RACE products from nonsense mutants is shown in the right panels. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, thin lines indicate introns, asterisks indicate PTCs, and half arrows indicate gene-specific primers used for modified RLM 5′ RACE. The 3′ remnant of MYB65 derived from miR159-guided cleavage was used as a positive control for modified RLM 5′ RACE.

(C) Comparison of EJC footprints from nonsense mutant mRNAs between phyB-9 and phyB-9 xrn4-6 and between pho2 and pho2 xrn4-6 using modified RLM 5′ RACE assays. Arrowheads indicate the RACE products sequenced in (A) and (B).

We further validated the contribution of XRN4 to EJC footprint production from nonsense mRNAs by generating the phyB-9 xrn4-6 and pho2 xrn4-6 double mutants. As expected, a lower number of putative EJC footprints was observed in the double mutants than in the nonsense mutants (Figure 5C). This result again confirms the role of XRN4 in the turnover of PTC-containing mRNAs and the generation of EJC footprints.

Analysis of EJC Footprints Reveals NMD Targets

Using the common features of EJC footprints in the NMD targets shown in Figures 3, 4, and 5, we globally validated potential targets of NMD triggered by AS or introns in the 3′ UTR that were annotated in the Arabidopsis and rice genome databases. For protein-coding transcripts with AS events, if the maximum 5′P peak was located in the first or second canonical EJC region downstream of an AS site, they were identified as putative targets of NMD triggered by AS, which frequently results in a PTC. In addition, transcripts possessing an intron-containing 3′ UTR were also identified as putative NMD targets if the maximum 5′P peak was located in the first or second canonical EJC region of the exon 3′ end at least 50 nucleotides downstream of the annotated termination codon. Analysis of the PARE data for Arabidopsis wild-type seedlings and flowers that we generated and the previously published Degradome-Seq data for rice wild-type seedlings (Li et al., 2010) revealed 254 Arabidopsis transcripts and 168 rice transcripts as putative NMD targets (Figure 6A; Supplemental Data Sets 1 and 2). Remarkably, 19.3 and 44% of the putative NMD targets identified from the analysis of Arabidopsis and rice 5′P end data, respectively, have an intron-containing 3′ UTR (Figure 6A). Compared with the large number of NMD-sensitive AS events reported by Drechsel et al. (2013), the number of Arabidopsis NMD targets we identified is small, because many AS events reported by Drechsel et al. (2013), including AS of CYP38, PPD6, and STR14, were not annotated in TAIR 10. However, 224 of the targets we identified in Arabidopsis did not show evidence of NMD-suppressed AS in the previous report (upregulated in the double mutant of lba1 upf3-1 with FDR < 0.1; Supplemental Data Set 1; Drechsel et al., 2013). Some of the novel NMD targets, such as Arabidopsis ARR8 and KCBP and rice transcripts encoding an SPA4-like protein (LOC_Os01g52640) and a transposon protein (LOC_Os04g03884), had high numbers of EJC footprints in the 3′ UTR (Figures 6B and 6C). For the Arabidopsis NMD targets we obtained, the analysis of gene Ontology (GO) terms revealed the highest enrichment for transcripts involved in “regulation of alternative mRNA splicing, via spliceosome” (22.27-fold enrichment, P < 0.001, Fisher’s exact test with Bonferroni correction; Supplemental Table 1). This enriched GO term is consistent with the previous finding that in Arabidopsis, NMD is an important mechanism regulating the splicing of genes encoding Ser/Arg proteins, which play key roles in splicing (Palusa and Reddy, 2010). Other significantly enriched GO terms unrelated to splicing included “organic cyclic compound metabolic process,” “nitrogen compound metabolic process,” and “chloroplast” (P < 0.05, Fisher’s exact test with Bonferroni correction). For the rice NMD targets we identified, however, no enriched GO terms were recovered. The difference in GO enrichment analysis results between the two species might be due to differences in the 5′P end data sets we analyzed, the frequency of introns in the 3′ UTR, or the annotations of splicing variants in the two genomes. It may also suggest that the role of the NMD pathway in gene regulation could be distinct among plant species.

• Download figure
• Open in new tab
• Download powerpoint

Figure 6.

Analysis of EJC Footprints Reveals Targets of NMD Triggered by AS and Intron-Containing 3′ UTRs.

(A) Fractions of putative NMD targets where NMD is triggered by AS or an intron-containing 3′ UTR (Intron⁺ 3′ UTR). The numbers of total putative NMD targets recovered from the analysis are in parentheses.

(B) and (C) Positional distribution of 5′P ends for targets of NMD triggered by intron-containing 3′ UTRs in Arabidopsis (B) and rice (C). 5′P ends in Arabidopsis were identified using the Arabidopsis PARE data (replicate 1) generated by this study, and those in rice were identified using the Degradome-Seq data for rice wild-type seedlings reported by Li et al. (2010). Within the 5′P end plots, red peaks indicate peaks in the canonical EJC region and black arrows indicate the maximum peaks. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, and thin lines indicate introns. Only the 5′P ends residing in exons are displayed in the plots.

Dysfunction of SMG7 Reduces the Number of EJC Footprints in Putative NMD Targets

To verify that EJC footprints in NMD targets are degradation products generated by the NMD pathway, we profiled 5′P ends in the Arabidopsis wild type and mutants of key NMD factors and compared the abundances of EJC footprints between them. In Arabidopsis, the crucial roles of UPF1 and SMG7 in the NMD pathway have been demonstrated (Arciga-Reyes et al., 2006; Riehs-Kearnan et al., 2012). Given that null mutants of UPF1 are seedling lethal (Arciga-Reyes et al., 2006), smg7-1, a severe but viable mutant of SMG7 with a T-DNA insertion in the CDS, was chosen for the comparison (Riehs et al., 2008). However, homozygous smg7-1 plants are sterile and morphologically different from wild-type plants when grown at 22°C due to an autoimmune response (Riehs et al., 2008; Riehs-Kearnan et al., 2012). To reduce the possibility that differences might be caused by the autoimmune response, we selected smg7-1 homozygous plants by genotyping and grew them at 28°C, under which the autoimmune response is repressed and wild-type and smg7-1 plants more closely resemble each other (Gloggnitzer et al., 2014). We then collected the inflorescences of wild-type and smg7-1 plants and isolated RNA for PARE library construction. Similar to the wild type, smg7-1 displayed an enrichment of 5′P ends in the canonical EJC region (Figure 7A). However, in smg7-1, the abundance of 5′P ends was significantly lower in the canonical EJC region but higher at the TSS compared with the wild type (P < 0.01, Wilcoxon rank-sum test; Figures 7B and 7C). Furthermore, the median of the log₂ fold change in maximum EJC footprint abundance between smg7-1 and the wild type was significantly lower for transcripts experiencing NMD-suppressed AS events than for other transcripts (P < 10⁻⁴, Kolmogorov-Smirnov two-sided test; Figure 7D), supporting the notion that SMG7 plays a crucial role in the production of EJC footprints from NMD targets. Consistent with the analysis of PARE data (Figures 3 and 4C), the analysis of modified RLM 5′ RACE data also confirmed the involvement of SMG7 in the production of EJC footprints for three NMD targets in both flowers and 14-d-old seedlings (Supplemental Figure 4). The reduction in the number of EJC footprints in smg7-1 seedlings supports the notion that the effect was due to the NMD function of SMG7 instead of its unique function in meiosis (Riehs et al., 2008). Like the 5′P peaks in the EJC regions of six NMD targets shown in Figures 3 and 4C, the maximum 5′P peak in LPEAT2, which was positioned outside the canonical EJC region, was also dampened in smg7-1 (Figure 3; Supplemental Figure 2), indicating that it might also be a degradation product of the SMG7-dependent NMD pathway. Moreover, for three previously reported NMD-sensitive transcripts that had xrn4-enhanced 5′P peaks in the CDS or 3′ UTR (Nagarajan et al., 2019), xrn4-enhanced 5′P peaks in eRF1-1 (Eukaryotic Release Factor1-1, NMD factor and target) were detected in the flowers of the wild type but were absent in those of smg7-1 (Supplemental Figure 5). This result further confirms the notion that the xrn4-enhanced 5′P peaks in eRF1-1 are also intermediates of degradation initiated by the NMD pathway (Nagarajan et al., 2019).

• Download figure
• Open in new tab
• Download powerpoint

Figure 7.

EJC Footprints in NMD Targets but Not miRNA Targets Are Less Abundant in smg7-1.

(A) Distribution of 5′P end counts in a 50-nucleotide region upstream of the exon-exon junction. Counts are shown for two replicates (R1 and R2) each from wild-type and smg7-1 flowers. The numbers of analyzed exons that are ≥ 50 nucleotides in size and with a 5′P end count ≥ 1 TP40M are shown in parentheses.

(B) and (C) Comparison of 5′P end counts in the canonical EJC region (B) and at the TSS (C) of the genes harboring the selected exons between the wild type and smg7-1. Asterisks indicate significant differences between groups (P < 0.01, Wilcoxon rank-sum test). Box boundaries and whiskers indicate quartiles and the range of 5′P end abundance, respectively. The numbers of analyzed EJC regions (count > 10 TP40M) or TSSs (count > 5 TP40M) are shown in parentheses.

(D) Scatter and density plots depicting log₂ fold change (smg7-1/wild type) and log₂ average 5′P end counts for the maximum EJC peaks in genes. Each data point represents a maximum EJC peak with an average 5′P end count ≥ 10. EJC peaks for genes producing NMD-suppressed AS variants (NMD⁺ AS), miRNA targets, and other genes (other) are shown. Arabidopsis genes producing NMD-suppressed AS variants and genes that are validated miRNA targets were extracted from data sets reported by Drechsel et al. (2013) and Zheng et al. (2012), respectively. The number of analyzed genes in each category is shown in parentheses. The density plots to the right cover all data points in each category, and straight lines indicate the median of each distribution. Asterisks indicate significant differences in the comparisons of the NMD⁺ AS group and the miRNA group with the group of other genes (Kolmogorov-Smirnov two-sided test). The colors and locations of asterisks indicate the associated category and the direction of the shift in the distribution.

Some miRNA Targets Possess EJC Footprints

Surprisingly, 5′P peaks corresponding to putative EJC footprints were also observed downstream of miRNA-guided cleavage sites in the targets of Arabidopsis miR159 (MYB33), miR160 (ARF10), and miR396 (GRF1; Figure 8A; Supplemental Figure 6). Likewise, rice miR159, miR160, and miR396 targets also harbored EJC footprints in the exon ends downstream of the miRNA target sites (Figure 8B). The accumulation of putative EJC footprints in the three Arabidopsis miRNA targets was reduced in xrn4-6 and nearly abolished in fry1-6; however, the amplitude of miRNA-guided cleavage peaks was increased in both xrn4-6 and fry1-6 compared with the wild type. We also confirmed the decreased number of EJC footprints and the increased number of miRNA-guided cleavage 3′ remnants for these three miRNA targets in xrn4-6 and fry1-6 by performing modified RLM 5′ RACE assays (Supplemental Figure 7). We also assayed XRN2 and XRN3 function in the production of EJC footprints from these three miRNA targets. The levels of 3′ remnants of miRNA-guided cleavage and EJC footprints were comparable between the wild type, xrn2-1, and xrn3-3 (Supplemental Figure 7). We thus further examined the level of EJC footprints in the XRN triple mutant xrn2-1 xrn3-3 xrn4-6. Like xrn4-6, the triple mutant also accumulated a higher number of miRNA-guided cleavage 3′ remnants but fewer EJC footprints compared with the wild type. However, the accumulation of these two types of degradation fragments in the triple mutant was neither different from that in xrn4-6 nor comparable to that in fry1-6. This result does not completely rule out a role for XRN3 in EJC footprinting because xrn3-3 is a knockdown mutant with a T-DNA insertion in the promoter region (Gy et al., 2007). These results confirm that XRN4 degrades the 3′ cleavage remnants of plant miRNA targets (Souret et al., 2004) and suggest that in addition to XRN4, other enzymes with 5′ to 3′ exonuclease activity also contribute to the turnover of miRNA cleavage remnants and EJC footprinting.

• Download figure
• Open in new tab
• Download powerpoint

Figure 8.

Some Arabidopsis and Rice miRNA Targets Possess Dominant EJC Footprints Downstream of miRNA-Guided Cleavage Sites.

Positional distribution is shown for 5′P ends of selected Arabidopsis (A) and rice (B) miRNA targets with prominent EJC footprints based on the Arabidopsis PARE data generated in this study and the rice Degradome-Seq data generated by Li et al. (2010). Within the 5′P end plots, red peaks indicate peaks in the canonical EJC region and red arrowheads indicate peaks corresponding to the miRNA-guided cleavage site. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, and thin lines indicate introns. Only 5′P ends residing in exons are displayed in the plots. Data for replicate 2 of Arabidopsis are shown in Supplemental Figure 6.

Notably, not all of the 5′P peaks evident at miRNA target sites were associated with prominent EJC footprints in downstream exon ends. For instance, while a high number of 5′P ends corresponding to the miR398 target site of Arabidopsis CSD1 and CSD2 were observed, no or poor 5′P peaks appeared in the canonical EJC regions 3′ to the miR398 target site (Figure 9A; Supplemental Figure 8). More intriguingly, the miR398 target site is close to the canonical EJC region in CSD1 and overlaps with the canonical EJC region in CSD2. Hence, EJCs might mask the miR398 target site on the CSD1 and CSD2 mRNAs that have not been translated and prevent miRNA regulation before the pioneer round of translation. Consistent with this result, rice miR398-regulated CSD1 and two miR444-regulated MADS genes, which possess a miRNA target site in or near the canonical EJC region, also lacked evident EJC footprints (Figure 9B). In addition to the miR398 and miR444 targets, a miR408 target encoding plantacyanin (ARPN) in Arabidopsis also possessed a prominent 5′P peak at the miRNA-guided cleavage site but poor peaks in the EJC region. However, unlike the miR398 and miR444 targets, the miRNA target site on ARPN is located at a position ∼200 nucleotides upstream of an exon-exon junction that is unlikely to be masked by EJCs. Taken together, these results suggest that there might be multiple factors controlling the interaction between plant miRNAs and EJC-bound targets.

• Download figure
• Open in new tab
• Download powerpoint

Figure 9.

Some Arabidopsis and Rice miRNA Targets Possess Poor EJC Footprints Downstream of miRNA-Guided Cleavage Sites.

Positional distribution is shown for 5′P ends of selected Arabidopsis (A) and rice (B) miRNA targets with poor EJC footprints based on the Arabidopsis PARE data generated in this study and the rice Degradome-Seq data generated by Li et al. (2010). Within the 5′P end plots, red peaks indicate peaks in the canonical EJC region and red arrowheads indicate peaks corresponding to the miRNA-guided cleavage site. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, and thin lines indicate introns. Only 5′P ends residing in exons are displayed in the plots. Data for replicate 2 of Arabidopsis are shown in Supplemental Figure 8.

Arabidopsis miRNAs Direct EJC-Bound mRNA Degradation

Unlike NMD targets, Arabidopsis miRNA targets did not exhibit changes in the number of EJC footprints in smg7-1 (Figure 8A; Supplemental Figure 6). Indeed, a global analysis of maximum EJC footprint abundance between smg7-1 and the wild type revealed a median log₂ fold change value close to 0 for miRNA targets. The log₂ fold change of miRNA targets is significantly higher than that of other transcripts (P < 0.01, Kolmogorov-Smirnov two-sided test; Figure 7D). This trend is opposite to that observed for NMD-suppressed AS events and suggests that most EJC footprints in miRNA targets are likely produced through an NMD-independent pathway.

Given that EJC footprints were evident downstream of miRNA target sites but rare in the upstream regions, we hypothesized that plant miRNAs are able to target EJC-bound mRNAs for degradation and trigger the production of EJC footprints. To test this hypothesis, we first assayed the role of AGO1, which mediates the function of most plant miRNAs, in the production of EJC footprints from the three miRNA targets shown in Figure 8A in the null ago1-3 mutant (Arribas-Hernández et al., 2016). In ago1-3, there were fewer miRNA-guided cleavage remnants and EJC footprints corresponding to the miR160 target ARF10 and the miR396 target GRF1 (Figure 10A). Interestingly, for the miR159 target MYB33, which exhibited no change in mRNA level in ago1-3 (Arribas-Hernández et al., 2016), the levels of miRNA-guided cleavage remnants and EJC footprints in ago1-3 were comparable to those in the wild type. Also, the mutation of AGO1 did not affect EJC footprint accumulation in the two NMD targets CYP38 and LPEAT2. Hence, these results confirm the link between the function of AGO1 and the production of EJC footprints in some miRNA targets.

• Download figure
• Open in new tab
• Download powerpoint

Figure 10.

AGO1 Function, Splicing, and miRNA Expression Determine EJC Footprint Production in miRNA Targets.

(A) RLM 5′ RACE assays of degradation intermediates generated from three miRNA targets and two NMD targets in the wild type and ago1-3. Open and solid arrowheads indicate miRNA-guided cleavage sites and EJC binding sites, respectively, in gene models and the corresponding RACE products separated by gel electrophoresis. The arrow indicates the maximum 5′ peak in the LPEAT2 gene model and in the corresponding RACE product on an electrophoretic gel. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, thin lines indicate introns, orange boxes indicate regions altered by AS, asterisks indicate PTCs, and half arrows indicate gene-specific primers. The steady state mRNA level of Arabidopsis ACT8 was used to control the input amount of total RNA and reverse transcription among samples. An equal amount of tobacco acid pyrophosphatase-treated total RNA from HeLa cells was spiked into each sample, and the TSS of human β-actin was used to control for differences in ligation efficiency among samples.

(B) Constructs expressing miR827 artificial targets: GUS with an intron (I-GUS827) and without an intron (GUS827). 35S, Cauliflower mosaic virus 35S promoter. Boxes indicate exons, thin lines indicate introns, and half arrows indicate gene-specific primers used for modified RLM 5′ RACE.

(C) RNA gel blot analysis of miR827 in Arabidopsis transgenic lines harboring miR827 artificial targets. Transgenic seedlings grown under phosphate-sufficient (+P) or phosphate-deficient (−P) conditions were used for the assays. U6 was used as a loading control. Numbers to the right of the blot show sizes in nucleotide.

(D) Modified RLM 5′ RACE assays of RNA degradation intermediates generated from miR827 artificial targets in Arabidopsis transgenic lines. The open arrowhead indicates the expected RACE product for the 3′ cleavage remnant of the artificial targets directed by miR827, and the solid arrowhead indicates the RACE product excised and cloned for Sanger sequencing. The miR159-guided 3′ cleavage remnant of MYB65 was used as a positive control for the modified RLM 5′ RACE assay.

(E) Positional distribution of the 5′P ends of cloned RACE products from I-GUS827 transgenic plants grown under phosphate-deficient conditions.

To further validate the ability of plant miRNAs to direct EJC footprint accumulation, we created an artificial miRNA target by fusing the 5′ UTR of Arabidopsis NITROGEN LIMITATION ADAPTATION (NLA), which contains a target site of miR827, to an intron-containing GUS gene and introduced this construct (I-GUS827) into Arabidopsis (Figure 10B). As miR827 is barely detectable when phosphate levels are sufficient but is highly induced upon phosphate starvation (Figure 10C; Hsieh et al., 2009), we thus were able to test the ability of miR827 to induce EJC footprint production by comparing RNA degradation intermediates produced from the I-GUS827 construct under phosphate-sufficient and -starvation conditions using modified RLM 5′ RACE. As expected, distinct amplification products of the RNA degradation intermediates with 5′P ends predominantly mapping to positions 26 to 27 nucleotides upstream of the exon-exon junction were detected in the plants grown under phosphate-starvation conditions but not in those grown with sufficient phosphate (Figures 10D and 10E). In contrast to the I-GUS827 construct, the control construct containing an intron-free GUS gene fused with the NLA 5′ UTR (GUS827) did not result in the accumulation of RNA degradation fragments truncated at the same site under phosphate-starvation conditions (Figure 10D). Nevertheless, the number of RNA degradation fragments with the 5′ end determined by miR827-guided cleavage was comparable between plants expressing the two constructs. This result validates the ability of plant miRNAs to target the mRNAs associated with EJCs for degradation and initiate the accumulation of transcripts with 5′P ends in the EJC region.

Unexpectedly, AS was observed in the I-GUS827 construct we created. In the AS form, the dominant sequence matching cloned 5′P ends was shifted farther upstream along with the alternative 5′ splice site compared with that in the normal splice form (Figure 10E). However, the distance between the predominant 5′P peak and the exon-exon junction remained 27 nucleotides in the AS form. As the deposition of EJC on mRNA is dependent on splicing, the shift of the 5′P peak in the AS form reinforces our hypothesis that the 5′P ends matching a region 26 to 30 nucleotides upstream of the exon-exon junction mostly represent in vivo EJC footprints.

Overexpression of an EJC Disassembly Factor Eliminates EJC Footprints

To provide direct evidence that the 5′P ends in the EJC region are EJC-protected fragments, we cloned an EJC disassembly factor, PYM, from Arabidopsis and transiently expressed a PYM overexpression construct in Arabidopsis protoplasts to disassemble EJCs (Figure 11). We predicted that disassembly of EJCs by PYM overexpression would promote the degradation of existing EJC-protected degradation intermediates and inhibit the production of new EJC footprints. For NMD targets, the production of intermediates with 5′P ends matching the canonical EJC regions of CYP38 and PPD6 was abolished when PYM was overexpressed (Figure 11). Degradation fragments corresponding to the maximum 5′P peak of LPEAT2, which was positioned outside the canonical EJC region and SMG7-dependent (Figure 3), were also absent in PYM-overexpressing protoplasts. This result indicates that this maximum 5′P peak is an EJC footprint and supports the notion that EJCs are also deposited on noncanonical sites in Arabidopsis. In addition to investigating the SMG7-dependent 5′P peaks in NMD targets (Figures 3 and 4C), we also examined the effect of PYM overexpression on the SMG7-independent 5′P peaks in the EJC region of miRNA targets shown in Figure 8A. Compared with the number of degradation intermediates in control protoplasts, the number of 3′ cleavage remnants directed by miRNAs in the three miRNA targets (MYB33, ARF10, and GRF1) was increased, whereas degradation intermediates with 5′P ends corresponding to the EJC region were barely detectable in the PYM-overexpressing protoplasts (Figure 11). Taken together, these results strongly indicate that the 5′P ends mapping to the canonical EJC regions of these transcripts are in vivo EJC footprints irrespective of their dependency on SMG7.

• Download figure
• Open in new tab
• Download powerpoint

Figure 11.

Overexpression of PYM Abolishes EJC Footprints in NMD and miRNA Targets.

Arabidopsis PYM was cloned and transiently overexpressed in Arabidopsis protoplasts. The steady state mRNA levels of PYM and a housekeeping gene, ACT8, in control (−) and PYM-transfected (+) protoplasts were measured using the cDNA from the modified RLM 5′ RACE assays. The amounts of RNA degradation intermediates from three NMD targets and three miRNA targets in control and PYM-overexpressing protoplasts were compared using modified RLM 5′ RACE assays. An equal amount of tobacco acid pyrophosphatase-treated total RNA from HeLa cells was spiked into each sample, and the TSS of human β-actin was used to infer ligation efficiency. Open and solid arrowheads indicate the miRNA-guided cleavage sites and EJC binding sites, respectively, in gene models and the corresponding RACE products separated by gel electrophoresis. The arrow indicates the maximum 5′ peak in the LPEAT2 gene model and the corresponding RACE product on an electrophoretic gel. The identity of all marked RACE products was confirmed by Sanger sequencing. Within the gene models, light gray boxes indicate UTRs, dark gray boxes indicate CDSs, thin lines indicate introns, orange boxes indicate regions altered by AS, asterisks indicate PTCs, and half arrows indicate gene-specific primers.