- Open Access
A global profiling of uncapped mRNAs under cold stress reveals specific decay patterns and endonucleolytic cleavages in Brachypodium distachyon
Genome Biology volume 14, Article number: R92 (2013)
mRNA degradation is a critical factor in determining mRNA abundance and enables rapid adjustment of gene expression in response to environmental stress. The involvement of processing bodies in stress response suggests a role for decapping-mediated mRNA degradation. However, little is known about the role of mRNA degradation under stressful environmental conditions.
Here, we perform a global study of uncapped mRNAs, via parallel analysis of RNA ends (PARE), under cold stress in Brachypodium distachyon. Enrichment analysis indicates that degradation products detected by PARE are mainly generated by the decapping pathway. Endonucleolytic cleavages are detected, uncovering another way of modulating gene expression. PARE and RNA-Seq analyses identify four types of mRNA decay patterns. Type II genes, for which light-harvesting processes are over-represented in gene ontology analyses, show unchanged transcript abundance and altered uncapped transcript abundance. Uncapping-mediated transcript stability of light harvesting-related genes changes significantly in response to cold stress, which may allow rapid adjustments in photosynthetic activity in response to cold stress. Transcript abundance and uncapped transcript abundance for type III genes changes in opposite directions in response to cold stress, indicating that uncapping-mediated mRNA degradation plays a role in regulating gene expression.
To our knowledge, this is the first global analysis of mRNA degradation under environmental stress conditions in Brachypodium distachyon. We uncover specific degradation and endonucleolytic cleavage patterns under cold stress, which will deepen our understanding of mRNA degradation under stressful environmental conditions, as well as the cold stress response mechanism in monocots.
Degradation of messenger RNAs (mRNAs) is one of the key factors determining mRNA abundance in a cell. Proper degradation of mRNA is crucial for the maintenance of cellular homeostasis and enables rapid and precise adjustment of gene expression in response to various environmental stimulus[1, 2]. The study of mRNA degradation in plants has fallen behind compared with what has been done in other organisms and recent advances in this area provide us with only preliminary knowledge about plant RNA decay pathways[3–5]. For normal mRNAs, degradation is believed to initiate with deadenylation, the deadenylated mRNA can then either be degraded in the 3'-5' direction in a process mediated by the exosome or be decapped and then degraded in the 5'-3' direction. Alternatively, mRNA decay in plants can be initiated with internal cleavage(s), for example, by miRNA- or siRNA-directed RISC degradation[6–9]. The mechanism of aberrant mRNAs degradation in plants is largely unknown, except for the identification of some key enzymes in the plant nonsense-mediated decay (NMD) pathway[10, 11]. Obviously, our understanding of the degradation pathways in plants is still far from complete.
In eukaryotic cells, decapping is a critical step in mRNA turnover as it requires the stop of translation and the change of mRNAs into a non-translating state, which makes them accessible for further degradation. Decapping, to a large extent, competes with translation for the access to the 5' cap structure, which determines whether the mRNAs will be degraded or translated[12–15]. Recently, ribosomes are reported to be allowed to bind to mRNAs and complete translation before the transcript is ultimately degraded. Messenger ribonucleoproteins (mRNPs) that are translationally repressed and associated with the decapping machinery can assemble into granules, named as processing bodies (PBs) [17, 18]. The assembly of PBs occurs primarily in response to stress stimuli, and the size of PBs has been found to increase significantly under stress conditions, implying that PBs play an important role during plant stress responses. Furthermore, PBs co-localize with stress granules (SGs), suggesting that they probably interact with each other and function together in stress responses[19–21]. The significant increase of PBs in response to stress and their close connection with SGs indicate that mRNA decay, especially the decapping pathway, plays an important role in cellular stress responses[22, 23]. However, little information is available for a complete view of mRNA degradation under stress conditions.
In this study, uncapped transcriptome has been analyzed with a deep sequencing approach, parallel analysis of RNA ends (PARE). Taking advantage of the modified 5' rapid amplification of cDNA ends (RACE) and high-throughput deep sequencing, PARE was developed to identify 5'-phosphorylated, polyadenylated RNAs. In most cases, it was used to identify potential targets for miRNA-directed cleavages[25–28]. In mouse, this approach was employed to perform a global analysis of the RNA degradome[29, 30]. This kind of application, as far as we know, was seldom reported in plant, except for a microarray study with similar experimental principle. Here, this approach was employed to study mRNA degradation under cold stress condition in Brachypodium distachyon, a new monocot model plant . Our studies will not only improve our understanding of mRNA degradation mechanism under stress conditions, but also provide deeper insight into how monocot plants respond to low temperature stress.
Global analysis of the uncapped mRNAs under cold stress
To obtain a complete view of the uncapped transcriptome and understand its role in plant cold stress responses, PARE was performed to identify 5'-phosphorylated RNAs using high-throughput sequencing. The Brachypodium diploid ecotype ABR5, with obviously higher cold tolerance and stronger vernalization requirement than the sequenced ecotype BD21 (Figure S1 in Additional file 1), was chosen for our analysis. Four PARE libraries were prepared. DC and DW were for Brachypodium seedlings with and without cold treatment, respectively. DC replicate and DW replicate were biological replicates for DC and DW (Table S1 in Additional file 2). For direct comparison, samples used for PARE analysis were prepared in the same manner and run side by side on Illumina GAIIx platform for high-throughput sequencing.
For all the four PARE libraries, over 10 million reads were generated through sequencing, respectively (Table S1 in Additional file 2). After removing sequences corresponding to known non-coding RNAs (rRNAs, tRNAs, small nuclear RNAs, and small nucleolar RNAs) and repeats/transposons, the remaining unique signatures were mapped to Brachypodium transcriptome. For DW and DC libraries, 65.48% and 67.68% of the generated sequences could be mapped to annotated transcripts. For DW replicate and DC replicate libraries, 56.56% and 57.89% of the generated sequences could be mapped to annotated transcripts (Table S1 in Additional file 2). For direct comparison, the DC library was normalized to the DW library according to the population size of mapped transcripts (Table S1 in Additional file 2). Biological replicates libraries for DW and DC exhibited good reproducibility (Figure S2 in Additional file 1). DW and DC libraries were used for further analysis, and the obtained results were verified in their biological replicates.
About 12 thousand genes were found to be represented by PARE reads in the uncapped transcriptome. The most abundant uncapped transcript in the DW library was from a gene encoding chlorophyll a-b binding protein, whereas in the DC library, except for a protein with unknown function, it was from agene encoding cold acclimation protein (Table S2 in Additional file 2). When genes with the top 20 highest PARE reads were considered, the DW library was composed mainly of degradation products for chlorophyll-related genes (Table S2 in Additional file 2). After cold treatment, mRNA degradation products for cold-stress-related genes appeared in the DC library and accounted for a proportion of sequences with high sequencing reads, while the uncapped transcript abundance for chlorophyll-related genes was still high (Table S2 in Additional file 2). Cold stress was observed to exert a clear effect on the uncapped transcriptome. Overall, approximately half of uncapped transcripts showed obviously altered abundance (>5-fold) in the cold-stressed Brachypodium seedlings. Interestingly, several phytohormone-related genes showed significant up-regulation after cold treatment (Table S3 in Additional file 2).
It is interesting to note that a considerable portion of the unique reads for both DW and DC libraries were mapped to non-coding RNAs (ncRNAs), which has also been observed for PARE libraries generated from Arabidopsis, mouse and human[29, 34, 35]. No further analysis was performed on ncRNAs in this paper, which mainly focused on the degradation of mRNAs.
According to the principle of PARE analysis, only polyadenylated RNAs will be detected. Deadenylation precedes decapping in mRNA degradation pathways [3–5], which means some mRNAs, whose poly(A) tails length are too short (<18nt), can not be represented in the PARE library. Although the information about poly(A) tail length in plants is hardly available, mRNAs are generally believed to possess long poly(A) tails (about 250 nucleotides), which is shortened to dozens of nucleotides in deadenylation[36, 37]. Only 18 adenosine residues of poly(A) tail are needed for PARE analysis , suggesting that the negative effect caused by deadenylation may be limited to a small range.
Relating uncapping-mediated mRNA decay profiles to transcriptome
To globally analyze the contribution of mRNA decay to gene expression regulation, mRNA abundance was analyzed by RNA-Seq. Two cDNA libraries, with (RC) and without (RW) cold treatment, were prepared from Brachypodium seedlings. High-throughput sequencing was performed and the obtained results were verified by real-time RT-PCR (Figure S3 in Additional file 1). Cold stress had obvious effect on Brachypodium transcriptome. According to the results of RNA-Seq analysis, approximately one-third of genes showed significant changes (>5-fold) in their expression levels after cold treatment.
About eight thousand genes could be found in both PARE and RNA-Seq libraries. The PARE dataset was analyzed for degradation patterns, in conjunction with them RNA abundance data provided by RNA-Seq, with log2 fold change of ±1.5 as the threshold. Four major change patterns were identified (Table 1). For type I genes, the abundance of their transcripts and uncapped transcripts changed in the same direction after cold treatment. For genes in type II, uncapped transcripts showed clear changes after cold treatment, but their transcript abundance remained unchanged. For genes belonging to type III, the abundance of their transcripts and uncapped transcripts changed in the opposite direction after cold treatment. Transcripts of class IV genes showed clear changes in abundance after cold treatment, but their uncapped transcript abundance remained unchanged. About three-eighths of the genes belonged to type II group, while only one-eighth of genes were grouped as type III. The rest of the genes distributed almost equally between I and type IV (Table 1). Similar results were obtained for biological replicate libraries of DW and DC (data not shown).
Relationship between uncapping-mediated mRNA decay and gene function
Gene Ontology (GO) analysis was performed for type I-IV genes using the GO database, which organizes information based on molecular function, biological process, and cellular component categories . In the biological process category, obvious over-representation against background was detected for the light-harvesting process of photosynthesis and nitrogen compound metabolic process for the type II population (Figure 1). When just specific biological processes were considered, only light-harvesting-related genes involved in photosynthesis were significantly enriched in type II. Obvious enrichment was observed for ATP-dependent enzymes in the molecular function category (Figure 1). These kinds of enzyme activity are involved in multiple biological processes and it is hard to deduce the biological significance of their enrichment in GO analysis.
For type IV genes, translation-related genes were obviously over-represented in the biological process category. Consistently, genes encoding ribosome-localized proteins were clearly over-represented in the cellular component category, and genes encoding ribosome-related proteins were clearly over-represented in the gene function category (Figure 1).
For type I genes, enrichment was observed for oxidoreductase and transferase activity in the molecular function category, but no corresponding enrichment was detected in both biological process category and cellular component category. Some of the stress-responsive genes were identified in type I, but no obviously significant enrichment was detected (Figure S4 in Additional file 3). For type III, genes involved in several important biological processes were represented at similar levels in GO analysis, with no obvious enrichment (Figure S4 in Additional file 3).
The overall mRNA stabilities for transcripts with different change tendencies during cold stress responses
To study the effect of cold stress on the overall mRNA stability of type I-IV genes, mRNA half-lives before and after cold treatment were analyzed with cordycepin, a global transcription inhibitor acting as a chain-terminating adenosine analog . Genes randomly selected from type I, III, and IV were analyzed. For type II genes, substantial enrichment was detected for light-harvesting-related genes in the GO analysis. This finding, along with the fundamental role of photosynthesis in plants[40–42], led us to further study this group of genes. Both light-harvesting-related and light-harvesting-unrelated type II genes were analyzed for their overall mRNA stabilities. Interestingly, for light-harvesting-related type II genes, all of the mRNAs with decreased uncapping-mediated degradation products (detected by PARE) showed increased overall stability, and all of the mRNAs with increased uncapping-mediated degradation products (detected by PARE) showed decreased overall stability, implying that uncapping-mediated mRNA degradation may be the key factor in determining the overall mRNA stability for this group of genes (Figure 2, Figure S5 in Additional file 1). For the other groups of genes, their overall mRNA stabilities were either up- or downregulated, or remained unchanged, no clear relationship was detected between the overall mRNA stability and uncapping-mediated mRNA degradation (Figure 2, Figure S5 in Additional file 1).
Sequence features associated with different mRNA decay patterns
To explore the features of mRNAs with different change patterns, mRNA length, untranslated region (UTR) length, GC content as well as the number of introns were analyzed for type I-IV mRNAs. The connection between introns and mRNA stability has been uncovered recently[43–45]. Interestingly, the average intron number for type II genes was higher than for the other genes, but the percentage of intronless genes was lower in this population (Figure S6 in Additional file 1). These data indicated that a mechanism, more complicated than expected, may exist for uncapping-mediated mRNA stability under stress conditions.
Regulatory sequence motifs in mRNA UTR regions have a close relationship with mRNA stability[46, 47]. Genes with the same uncapped/total transcript change patterns, therefore, may share common regulating elements in their mRNA UTR regions. Type I-IV genes were analyzed with an integrated motif-discovery program to search for possible motifs located in 5' or 3'UTRs of corresponding transcripts. Conserved sequence patterns in the 5'UTRs of mRNAs were identified for all four types of genes (Figure S6 in Additional file 1), implying the existing of a specific regulator or regulation mechanism. Interestingly, the identified motifs for type II and type IV genes are similar, providing a clue that similar regulation factor may exist for these two types of genes. The motif-discovery analysis was also performed for the subgroups of type I, II, III, and IV genes, according to the trend of change in PARE and RNA-Seq analysis. More conserved motifs were identified, including several GA-rich motifs (Figure S7 in Additional file 1). It is interesting to note that similar motifs were also found in the 3' UTRs for genes with decreased transcript abundance in both type I and type IV, implying similar expressional regulator for these genes (Figure S7 in Additional file 1).
Decapping pathway, with clear response to cold stress, makes the most significant contribution to the uncapped transcriptome
A high degree of conservation has been observed for degradation complexes or their core components in plant decay pathways, including the decapping complex, the deadenylation complex, the NMD complex and the exosome[4, 5]. Although target sets for the core components of these degradation complexes have been reported in Arabidopsis[49–52], such kind of information is hardly available for monocotyledonous plants. To investigate the contribution of different decay pathways to the uncapped transcriptome, enrichment analysis was carried out with target sets in Arabidopsis, based on the assumption that degradation complexes have conserved target sets in dicots and monocots (Figure 3, Table S4 in Additional file 2).
Among all the tested mRNA degradation complexes, substrates for the decapping complex showed significant enrichment before and after cold treatment (Figure 3), consistent with the origin of PARE libraries. Significant enrichment was also detected for the substrates of deadenylation complex after cold treatment (Figure 3), suggesting that deadenylation may make more significant contribution to PARE libraries under cold stress condition.
To gain a complete view of the mRNA degradation network's response to cold stress, the targets of the exosome and NMD pathway were also analyzed. No statistically significant enrichment was detected for exosome pathway (Figure 3), which is not unexpected because the exosome mediates the 3'end decay pathway and its degradation products could not be detected by PARE [53, 54]. No statistically significant enrichment was detected for targets of the RNA helicase UPF1, the key component of the NMD pathway (Figure 3), implying no direct connection between NMD- and uncapping-mediated degradation.
miRNA-mediated RNA degradation also makes contribution to the uncapped transcriptome
Endonucleolytic cleavage(s) mediated by miRNA-programmed RNA-induced silencing complexes (RISCs) could also produce 5'-phosphorylated RNAs, which could be detected by PARE [55, 56]. In Brachypodium, approximately one-fourth of miRNAs showed altered expression (>3-fold) in cold stress. Among them, three conserved (miR169e, miR172b, miR397) and 25 predicted miRNAs respond significantly (≥5-fold) to cold treatment. The target genes for most conserved miRNAs were identified by PARE analysis and most of them were conserved targets (Table 2). The obtained targets were confirmed with biological replicate libraries for DW and DC (Table 2). For the 25 cold-responsive predicted miRNAs, corresponding cleavage remnants for some of them could be found in PARE libraries (Table S5 in Additional file 2).
The miRNA libraries, with (WC) and without (NC) cold-treatment, were prepared in the same way as PARE and RNA-Seq libraries, enabling direct comparison between different libraries for changes in miRNAs, their target genes, and corresponding degradation products. In Brachypodium, miRNAs play different kinds of roles in gene expression regulation. For example, miR393 was a cold-induced miRNA and its target, an F-box gene, was downregulated by cold stress (Table 2). Consistently, more miR393 cleavage products were detected by PARE after cold treatment (Table 2), revealing the dominant role of miRNA-mediated regulation. Similar miRNA-target interactions were also observed between miR396 and some of its targets (Table 2). For the other cold-induced conserved miRNAs, miR172, some of their target genes showed obvious decreases in transcript abundance after cold treatment. However, for these target genes, their corresponding degradation products were unchanged in cold stress response (Table 2), suggesting that miRNA-mediated degradation is not the only influencing factor involved in expressional regulation of these genes.
Endonucleolytic cleavages, another contributing factor to the uncapped transcriptome, were identified and analyzed for its functional mechanism
In addition to the miRNA-directed internal cleavage pathway, is there any other endogenous cleavage pathway existing in plant? In our PARE analysis, several specific endogenous cleavages were mapped to sites unrecognized by any conserved or predicted miRNAs (Figure 4, Figure S8 in Additional file 1). For all of the transcripts involved in these endogenous cleavage events, PARE reads were identified only at these specific cleavage sites, no reads or reads with obviously lower abundance were identified at other positions (Figure 4b). Moreover, conserved motifs were found at the cleavage sites (Figure 4a, Figure S8 in Additional file 1).The way mRNAs cut in these endogenous cleavage events was quite different from that for miRNA- or siRNA-mediated degradation [55, 58], excluding the possibility that they are guided by miRNAs or siRNAs. Plant miRNAs/siRNAs have perfect or near-perfect complementarity to their targets and the cleavage site is located in the central region of the sequence spanned by the guiding small RNAs[55, 58]. However, it is not the case for the identified endogenous cleavages, in which the cleavage sites are located at the 5' end of the conserved motifs. For example, 12 mRNAs in the DC library shared a consensus motif, in which '↓' stands for the cleavage site (Figure 4a).
In PARE libraries, approximately 90 transcripts were found to be involved in the identified endogenous cleavage events, when only prominent endogenous cleavages were considered. Ten endogenous cleavage groups were found based on conserved motifs identified at the cleavage sites. Endogenous cleavage groups with fewer than six group members were omitted from further analysis. Interestingly, about four-tenths of the identified endogenous cleavages were located in exons, and approximately six-tenths of them were located in the 3' UTR. Only a few of them were found in the 5' UTR region. These endogenous cleavages distributed randomly in PARE libraries, with no distribution preference identified for any type of genes, implying that this degradation mechanism may be independent of uncapping.
A significant portion of the identified endogenous cleavage groups were cold-stress-responsive. Four groups were identified only in the DW library and four groups were found only in the DC library. These endogenous cleavage groups were also identified in corresponding biological replicate PARE libraries, except for a few genes with relatively low reads (Figure S8 in Additional file 1). Conservation analysis was performed for the conserved motifs identified in these endogenous cleavage groups. The identified motifs expressed different levels of conservation in barley, wheat and rice, but not in Arabidopsis (Table S6 in Additional file 2).
No clear functional connection was detected among group members, but a significant portion of them encoded enzymes catalyzing various kinds of biological reactions (Figure S8 in Additional file 1). Of the 12 mRNAs shown in Figure 4, about half of them were enzyme genes (Figure 4). For the DW- or DC-specific endogenous cleavage groups, group members always encoded proteins involved in stress responses. For example, the DC-specific endogenous cleavage group was shown in Figure 4a. Among them, sterile (STE) kinases, serine/threonine kinase, erine/threonine phosphatase, and shaggy-related protein kinase are well known for their roles in stress signal transduction[59–62]. Glutaredoxins belong to glutathione-dependent disulfide oxidoreductases that are involved in oxidative stress responses . Annexins, WEE1-like protein kinase and importin-α have also been reported to be involved in plant stress responses[64–67]. The premRNA splicing factor, ribosomal protein as well as zinc finger protein are implicated in various kinds of cellular processes, including stress responses[68–71].
These observed endogenous cleavages could be made by endonucleolytic enzymes, or be the results of stalled exonuclease activity. If the second case is true, the 5' cap structure and 3' poly(A) tail of their mRNAs would probably be destroyed. To explore the mechanism responsible for these endogenous cleavages, quantitative splinted-ligation reverse transcriptase polymerase chain reaction assay (qSL-RT-PCR) and ligation-mediated poly(A) test (LM-PAT) assay were employed to determine whether these endogenously cleaved mRNA shad intact 5' and 3' end structures, respectively[72, 73]. Genes involved in endogenous cleavage events were randomly selected for further analyses. To obtain clear and stable results, genes with low expression levels were not included.
The qSL-RT-PCR assay showed that no RT-PCR products were identified for Bradi2g24790, Bradi2g35600, Bradi3g08917, and Bradi4g28400 (Figure 5), indicating that the anchor RNA could not be ligated to 5'end of these mRNAs because of their intact 5' cap structures. After mRNAs were decapped with Tobacco Acid Pyrophosphatase (TAP), RT-PCR products were detected, suggesting that the qSL-RT-PCR system works well (Figure 5). The results of the qSL-RT-PCR indicated that all the analyzed transcripts had intact 5'cap structures.
In the LM-PAT assay, poly(A) tail lengths were measured for Bradi2g24790 Bradi2g35600 and Bradi4g28400, which were involved in cold-responsive endogenous cleavage events (Figure S8 in Additional file 1). An endogenous-cleavage-unrelated type II gene, Bradi1g56580, was used as a positive control. As shown in Figure 6, the size of the PCR product for Bradi1g56580 decreased after cold treatment, indicating that the LM-PAT analysis system works. No obvious change was observed in the LM-PAT assay for Bradi2g24790 Bradi2g35600 and Bradi4g28400 (Figure 6), indicating that their mRNA poly(A) tail length did not change after cold treatment. The results of the LM-PAT indicated that although the analyzed transcripts were involved cold-responsive endogenous cleavages, their mRNA poly(A) tail lengths did not respond to cold treatment.
In response to cold stress, eukaryotic cells effectively regulate gene expression through rapidly adjusting their transcriptome[74–76]. Transcript abundance is the equilibrium between the rates of mRNA synthesis and degradation. Up to now, relatively little is known about mRNA degradation. Although the involvement of PBs in stress responses suggested that mRNA decay, especially the decapping pathway, may play a role under stress condition, further experimental evidence is needed in this area. To a certain extent, mRNA degradation is only assumed to be a factor determining the speed by which the transcript abundance is adjusted. Our results indicated that mRNA degradation did play a role in the regulation of steady-state transcript abundance. For example, the opposite change of transcript and uncapped transcript abundance for type III genes uncovered the obvious contribution of uncapping to gene expression regulation (Table 1). Moreover, genome-wide analysis of mRNA degradation under cold stress will deepen our understanding of mRNA degradation and stress response mechanism in plants.
Different mRNA degradation patterns were uncovered during cold stress responses
Cold stress exerted a clear effect on the uncapped transcriptome. When both gene transcription and transcript degradation were taken into consideration, genes were classified into four groups according to their alternation patterns in response to cold stress (Table 1).
For type I genes, the abundance of transcripts and uncapped transcripts changed in the same direction after cold treatment (Table 1). The mRNA half-life analysis indicated that the overall mRNA stability of type I genes changed during cold stress (Figure 2, Figure S5 in Additional file 1), suggesting that mRNA stability, in addition to transcription, also plays an important role in determining steady-state mRNA abundance. Base on present data, it is not easy to deduce how the uncapping-dependent stability changed for type I genes under cold stress condition.
For type II genes, uncapped transcripts changed significantly after cold treatment, but their transcript abundance remained unchanged (Table 1), suggesting a balance between transcription and transcript degradation. Specific sequence features were identified for type II genes (Figure S6 in Additional file 1), implying a specific regulation mechanism for this group of genes.
For genes belonging to type III, the abundance of their transcripts and uncapped transcripts changed in opposite directions after cold treatment (Table 1). Because the uncapped transcript level is affected by both transcript abundance and its degradation rate, the opposite changes in transcript and uncapped transcript abundance indicate that uncapping did play a role in the regulation of gene expression. For type III genes, no clear relationship was detected between the overall mRNA stability measured with cordycepin-treatment and uncapping-mediated mRNA degradation, implying that uncapping may not be the key factor determining the overall mRNA stability for this group of genes (Figure 2, Figure S5 in Additional file 1).
For type IV genes, their uncapped transcripts remained unchanged after cold treatment although their transcript levels showed obvious alterations (Table 1). For these genes, the opposite changes in transcription and uncapping-mediated degradation reached a balance to ensure that their uncapped transcripts remained stable during cold stress. Now it is not easy to deduce the biological significance of this balance. Translation-related genes were obviously over-represented for type IV genes in GO analysis (Figure 1), implying that some kind of regulation mechanism, preferentially controlling translation, may exist for this group of genes.
The mRNA degradation products from cold-responsive genes account for a proportion of the sequences with high sequencing reads in the DC library (Table S2 in Additional file 2). However, only slight enrichment was observed for stress-related genes in GO analysis (Figure S4 in Additional file 3), implying that the effect of uncapping-mediated degradation cold-stress-related transcripts may be limited to a small group of genes.
Uncapping plays an important role in gene expression regulation for the light-harvesting process of photosynthesis under cold stress
During cold stress response, the steady-state transcript abundance for type II genes was stable during cold stress, but their uncapping-dependent transcript stability showed significant changes (Table 1). According to GO analysis on type II population, clear over-representation was detected for light harvesting process of photosynthesis (Figure 1). This result, together with the widely accepted important roles of photosynthesis[40–42], led us to further dissect this group of genes. The overall mRNA stability was analyzed for all types of genes, and only light-harvesting-related genes in type II showed opposite alternation trends for uncapping-mediated degradation products and the overall mRNA stability during cold stress response (Figure 2, Figure S5 in Additional file 1), implying that uncapping might be the key factor in determining the overall mRNA stability for this group of genes.
What is the functional significance of this uncapping-mediated regulation mechanism for light-harvesting-related genes during cold stress? Photosynthesis is involved in the first group of biological processes influenced by temperature fluctuation [78, 79]. Low temperature causes an energy imbalance in photosynthesis, which results in photoinhibition[78, 80]. Annual cold tolerant plants show a decrease in light-saturated rates of CO2 uptake in response to a sudden low temperatures stimulus, which is followed by a strong recovery of photosynthesis once their leaves are acclimated [81, 82]. For these plants, the photosynthesis activity shows reverse trends of change in short-term and long-term cold responses. One possible explanation of our results is that the expression levels of light-harvesting-related genes remained unchanged to save energy and materials at the middle stage (the stage during which our experiment was performed) of cold response, but their uncapping-dependent transcript stability altered markedly, to prepare for a quick and significant change at gene expression level in case the energy balance in photosynthesis broken and the expression of light-harvesting-related genes needs to be altered instantly under cold stress.
Decapping, among all the degradation pathways, contributes most significantly to the uncapped transcriptome
In this study, a complete view of the uncapped transcriptome was obtained through PARE analysis. There are different degradation pathways in plants. Up to now, our understanding of them is fairly rudimentary, with several issues remaining to be addressed[3–5]. Based on currently available information, enrichment analysis was carried out and the results indicated that decapping pathway, whose products are 5'-phosphorylated RNAs, is the main contributor to the uncapped polyadenylated degradome (Figure 3).
In addition to the decapping pathway, the NMD and miRNA-directed RISC degradation can also produce 5'-phosphorylated RNAs, which could be identified by PARE. The degradation products for thousands of genes were identified in PARE libraries, and only dozens of them were found to be miRNA targets in our study (Table 2, Table S5 in Additional file 2). To exclude the confusion caused by false miRNA targets, targets with relatively low credibility were not included in our analysis (Table 2). The total number of miRNA targets, therefore, may be underestimated in our study. Up to now, the number of experimental identified miRNAs is limited for monocots[28, 83]. However, hundreds or even thousands of miRNA targets have been predicted by software according to the information provided by presently available database for plant miRNA targets, suggesting the existence of unidentified miRNA targets. Thus, miRNA-mediated degradation did make contribution to PARE libraries, but presently it is hard to tell its extent when taking the limitation of tools used for miRNA target identification into consideration.
Only limited information is available for NMD pathway in plants . Recent studies showed that both long 3'UTRs and 3' UTR-located introns can efficiently induce NMD in plants and 3'UTR-bound exon junction complex (EJC) is required for intron-based plant NMD. However, the position of NMD in the whole mRNA decay system is still unclear in plant. In the enrichment analysis of PARE library, no statistically significant enrichment was detected for targets of the key component of NMD (Figure 3), providing a clue that NMD may not have close connection with uncapping in plant.
Although the degradation complexes are conserved, they may have different sets of targets in dicots and monocots. Also, it is possible that other undetected degradation complexes may exist in plants. Up to now, relatively little is known about plant degradation pathways. More information is needed for us to better understand their contribution to uncapped transcriptome.
Endonucleolytic cleavage, identified by PARE, plays a role in the regulation of gene expression, especially under cold stress
The decay of RNA in eukaryotes has long been assumed to be executed mainly by exoribonucleases[86–89]. Recently, however, there are more and more reports demonstrating the involvement of endoribonucleases in RNA turnover [29, 90]. In addition to the enzymes displaying endoribonucleolytic activity in general processes associated with RNA metabolism, several endoribonucleases have close connections with specific signals, including stress stimuli. For example, the angiogenin, an RNase A protein with endonucleolytic activity, and endonuclease PMR1 have been proved to be involved in stress response [92–95]. As far as we know, there is no such kind of report in plants. In this study, a series of cold-stress-responsive endogenous cleavage events were detected in PARE libraries, demonstrating the existence of this mechanism in plants. Similar endogenous cleavage events have been reported in mouse. Although PARE analyses have also been performed in plants[25–28], the purpose of those studies is to identify potential miRNA targets, with almost no attention paid on endogenous cleavages.
The mRNAs involved in cold-responsive endogenous cleavage events had intact 5'cap structures and 3'poly(A) tails with unchanged length after cold treatment (Figure 5, Figure 6), suggesting that these endogenous cleavages have no direct connection with decapping and deadenylation. The motifs identified for endogenous cleavage events were not highly conserved (Table S6 in Additional file 2), suggesting that endogenous cleavages may choose different groups of genes as targets in different plant species. The endogenous cleavages distributed randomly among different types of genes. Most of these endonucleolytic cleavages were either cold-induced or cold-suppressed, and their targets were mostly genes involved in cold stress responses (Figure 4, Figure S8 in Additional file 1), suggesting that endonucleolytic cleavage may exert its function in cold stress response. After endonucleolytic cleavage, mRNA degradation can proceed in both 3'-5' and 5'-3' directions, which enabled prompt degradation of the entire transcript. Thus, it is quite possible that the identification of endonucleolytic cleavages uncovers an efficient way to modulate gene expression during cold stress responses.
The mRNA degradation, one of the key parameters determining steady-state gene expression levels, is a complicated and well-ordered process, controlled by multiple mechanisms. Our understanding of mRNA decay system is far from complete, with several pathways remaining unknown. The identification of endonucleolytic cleavages revealed another way of degradation, especially under cold stress condition. With current knowledge of mRNA degradation, the biological significance of this mechanism is still not clear. More efforts are needed in this area to further understand this mechanism and its connection with other degradation pathways.
Our global analysis provides a complete view of uncapping-mediated mRNA degradation and its related degradation pathways under cold stress. Specific degradation patterns have been uncovered. Endonucleolytic cleavages, another way to modulate gene expression in cold stress response, were revealed. Our results will not only deepen our understanding of mRNA degradation under stress condition, but also help us to gain deeper insight into the cold stress response mechanism of economically important winter-habit crops and biofuel grasses, which has close evolutionary relationship with Brachypodium.
Materials and methods
Plant material and growth conditions
Seeds of the Brachypodium distachyon (L.) Beauv. diploid line ABR5 were placed in petri dishes containing 2 layers of damp sterile filter paper. The seeds were first stratified at 4°C for 1 week to promote synchronous germination, then grown in a growth chamber at 24°C under a 16 h/8 h (light/dark) photoperiod with light intensity of approximately 5,000 lux.
Total RNA was extracted from 12-day-old ABR5 seedlings using the mirVana miRNA Isolation Kit (Ambion, Austin, USA) following the manufacturer's protocol. The aerial part of seedlings treated or untreated with cold stress (4°C for 24 h) were pooled and used for construction of the DC or DW library, respectively. PARE libraries were prepared according to  in principle. Modification was made according to . Sequencing was performed on Illumina GAIIx by LC Sciences (Houston, TX, USA).
Raw sequencing reads were obtained using Illumina's Pipeline V1.5 software. After removing sequences corresponding to known rRNAs, tRNAs, small nuclear RNAs, small nucleolar RNAs, and repeats/transposons[98, 99], the remaining sequences were mapped to the Brachypodium genome using bowtie, allowing one mismatch. Low quality reads and low frequency reads (≤5) were omitted from further analysis. For PARE reads mapped to multiple positions, reads were divided among different positions. Brachypodium draft genome sequences and V1.2 annotation was obtained from . Custom Perl scripts were used to identify the total reads for specific transcripts. For comparison, the DC library was normalized to the DW library based on the population size of the genome-mapped reads (Table S1 in Additional file 2). Genes with log2 fold of changes ≥1.5 were considered as upregulated genes and genes with log2 fold of changes ≤-1.5 were considered as downregulated genes.
Total RNA was extracted with the mirVana miRNA Isolation Kit (Ambion, Austin, TX, USA) from 12-day-old ABR5 seedlings, which were grown and cold-treated in the same manner as for the PARE libraries. The mRNAs were enriched by using oligo(dT) magnetic beads and broken into short fragments (about 200 bp). The first strand cDNA was synthesized using random hexamer primer with mRNA fragments as templates. DNA polymerase I (Epicentre® Biotechnologies, Madison, WI, USA) was used to synthesize the second strand. The double-strand cDNA was purified with the QiaQuick PCR extraction kit (Qiagen). Then, sequencing adaptors were added and the obtained fragments were purified using agarose gels and enriched by PCR amplification. The resultant products were used for sequencing analysis via Illumina HiSeq™ 2000.
After low quality reads were removed, Illumina sequencing reads were mapped to the Brachypodium genome using SOAP, allowing two mismatch. Gene expression profiles were obtained with reference to the Brachypodium V1.2 annotation. The population size of the genome-mapped unique reads for RC library and RW library was similar, enabling direct comparison. Genes with log2 fold of changes ≥1.5 were considered as upregulated genes and genes with log2 fold of changes ≤-1.5 were considered as downregulated genes.
All the sequences used for mRNA feature analysis were downloaded from the Brachypodium website. As described above, the Brachypodium V1.2 annotation was employed in this analysis. Conserved motifs were identified by MEME software . Gene Ontology analysis was performed using agriGO. All other sequence feature analyses were performed using custom scripts in Perl.
For miRNA target identification, PAREsnip was used to detect potentially cleaved targets based on degradome sequences, with Brachypodium transcripts (V1.2) and miRNA sequence as input. Conserved miRNA sequences were downloaded from miRBase and Brachypodium-specific miRNA sequences were obtained from previous publication .
Enrichment analyses were performed with Fisher's exact test. For the decapping pathway, target sets for TDT (88 genes), the mRNA-decapping enzyme, were analyzed. For the deadenylation pathway, target sets of PARN (135 genes), the poly(A)-specific ribonuclease, were analyzed. For the nonsense-mediated decay (NMD) pathway, target sets of UPF (58 genes), an RNA helicase in the NMD, were analyzed. Target sets of RRP4 (119 genes), a core component of the exosome, were also analyzed. Their presence in the DW and DC libraries was analyzed using the whole Brachypodium genome as background. The two-tail P value was calculated using online software . Only genes represented in the Brachypodium V1.2 annotation  were considered in the analysis.
Total mRNA stability measurement
Total mRNA stability was measured through inhibiting transcription with cordycepin treatment[39, 106–108]. The middle part of the primary leaves from 12-day-old Brachypodium seedlings with (CD) and without (WT) cold treatment (4°C, 24h) were used. Cordycepin treatment was performed as described. The 3-cm leaf segments were incubated with 600 µM cordycepin solution with regular shaking at 22°C for 7h. Then tissue samples were harvested at regular intervals (1h, 2h, and 4h) and quickly frozen in liquid nitrogen. Total RNA was isolated as described above. Two microgram of total RNA was used for reverse transcription with SuperScript II reverse transcriptase (Invitrogen). The obtained cDNAs were quantified and used as template for quantitative PCR using SYBR Green PCR Master mix (Applied Biosystems, Foster City, CA, USA). For both WT and CD samples, mRNA abundance for the time zero was arbitrarily set to 1 and mRNA abundance for the other time points is shown as the relative value of time zero. Quantitative PCR reactions were repeated three times.
RLM 5'-RACE was conducted using the GeneRacer kit (Invitrogen). Total RNA was isolated as described above. Total RNA was ligated to the GeneRacer RNA Oligo adaptor (5'-CGA CUG GAG CAC GAG GAC ACU GAC AUG GAC UGA AGG AGU AGA AA-3'), then reverse transcription was performed with SuperScript II reverse transcriptase (Invitrogen) using oligo(dT) or gene-specific primers. The obtained cDNAs were subjected to two rounds of PCR amplification (30 to 35 cycles). Primers specific to the RNA Oligo adapter and primers specific to the target genes were designed. For Bradi1g31450,5'-CGA GGA CAC TGA CAT GGA CTG AAG GAG TAG AA-3' (BRADI1G31450-1-F) and 5'-TTG TAG ACA GTG CTA AGC TTC AAG GAC AGG CG-3' (BRADI1G31450-1-R) were used for the first round PCR; BRADI1G31450-1-F and 5'-CGA GCT ATG TCT GAT AGG ACA GAA CAA ATT TCA-3' (BRADI1G31450-2-R) were used for the second round PCR. For Bradi2g47100, 5'-GAG GAC ACT GAC ATG GAC TGA AGG AGT AGA-3' (BRADI2G47100-1-F) and 5'-AGC TGC AAC TCA GAA CGA AAG CTA AAG AG-3' (BRADI2G47100-1-R) were used for the first round PCR; BRADI2G47100-1-F and 5'-AAC AGA GCC AAG CCT GTT TGT GTT GC-3' (BRADI2G47100-2-R) were used for the second round PCR. For Bradi4g29680, 5'-ACT GGA GCA CGA GGA CAC TGA CAT-3' (BRADI4G29680-1-F) and 5'-AAC AGC GTG CAA CCA TTG AAA CCA GGA CAC-3' (BRADI4G29680-1-R) were used for the first round PCR; BRADI4G29680-1-F and 5'-CCG ATT ACA GTC ACC GAA CGA AAT GAC-3' (BRADI4G29680-2-R) were used for the second round PCR. Amplification products from a single band were cloned using the TOPO TA cloning kit (Invitrogen) and sequenced.
qSL-RT-PCR was performed as described previously. Total RNA was extracted as described above and mixed with RNA Oligo adaptor (sequence was provided in RLM 5'-RACE) and complementary DNA Splint oligonucleotides, catalyzed by T4 DNA ligase (Promega, Madison, WI, USA). The mRNAs containing a 5' cap could not be ligated, whereas decapped RNAs with a 5' monophosphate could be ligated with RNA Oligo adaptor. After ligation, the DNA splints were destroyed by RQ1 DNase (Promega).The ligated RNAs were converted to cDNAs by reverse transcription with a reverse gene-specific primer (GSP-R). The resulting cDNAs were then detected by PCR using GSP-R and a forward primer annealing to the RNA Oligo adaptor (adaptor Primer). The reverse transcription and PCR amplification was performed as described above. PCR amplification was run for 30 cycles (+TAP) and 37 cycles (-TAP) for target genes and 22 cycles for UBQ10. Primers and DNA Splints used were: for Bradi2g35600, 5'-CGC TGC TAG GGT TTG ACG GTT TTG GTT GGT TTT CTA CTC CTT CAG TCC ATG TCA GTG TCC TCG TGC TCC AGT CG-3' (Bradi2g35600DNAsplint), 5'-GAC TGG AGC ACG AGG ACA CTG ACA T-3' (Bradi2g35600anchor) and5'-GTC GTT GCT GGA CTG GGC CAT CTC TT-3' (Bradi2g35600reverse); for Bradi2g24790, 5'-GCG AGA AGG TGG GGA GGT TTT GGA GGA AGG TTT CTA CTC CTT CAG TCC ATG TCA GTG TCC TCG TGC TCC AGT CG-3' (Bradi2g24790DNAsplint), 5'-GAG GAC ACT GAC ATG GAC TGA AGG AGT AGA-3' (Bradi2g24790anchor) and 5'-ACG AGA GGA AAT AAG GAG ATC CCT CGC T-3' (Bradi2g24790Reverse); for Bradi3g08917, 5'- GAC GGA GCA AGA ACC CGC CGCCGC CAA CCA TTT TCT ACT CCT TCA GTC CAT GTC AGT GTC CTC GTG CTC CAG TCG-3' (Bradi3g08917DNAsplint), 5'-AGG ACA CTG ACA TGG ACT GAA GGA GT-3' (Bradi3g08917anchor) and 5'-TCC AGA GCC AAT GCA TGG ATC AAC AC-3' (Bradi3g08917Reverse); for Bradi4g28400, 5'-GGT TCA GCT CCA ACA ACC TCC ATT CCC ATG TCC CTT TCT ACT CCT TCA GTC CAT GTC AGT GTC CTC GTG CTC CAG TCG-3' (Bradi4g28400DNAsplint), 5'-CTG ACA TGG ACT GAA GGA GTA GAA AGG G-3' (Bradi4g28400anchor) and 5'-CTT AGG GTC CTC AAA TGA TCG GAC CTT G-3'(Bradi4g28400Reverse); for ubiquitin10 (UBQ10),5'-TCC TCT GAC ACA ATC GAC AAC-3' (UBQ10For) and 5'-TCC TGG ATC TTT GCC TTC AC-3' (UBQ10Rev).
LM-PAT assays were performed as previously described . The total RNA (250 ng) was incubated with 20 ng phosphorylated oligo(dT)12-18 [p(dT)12-18] at 65°C for 10 min, then 42°C for 30 min in the presence of T4 DNA ligase (Promega). Then, 200 ng oligo(dT) anchor (5'-GCG AGC TCC GCG GCC GCG TTT TTTTTT TTT-3') was added and the temperature was lowered to 12°C for 2 h so that the oligo(dT) anchor could be annealed to the unpaired ends. This ligated primer was then used for reverse transcription, catalyzed by AMV reverse transcriptase (Promega) at 42°C for 1 h. The first-strand cDNAs were synthesized and used as template for PCR with oligo(dT) anchor primer and gene-specific primer. PCR was performed as described above and primer used was: for Bradi1g56580, 5'-GGC AAT GCG CCA TGT GTG TAT GTA-3'; for Bradi2g24790, 5'-TGG CTC TGC CTT ATT TCA ATT ACT TGA TGC TT-3'; for Bradi2g35600, 5'-CGA GTG TTC GAA ATA AGT AAG CCG TCA GG-3'; forBradi4g28400, 5'-AGT CAA TTC AGC TTG GAT CAG CCC-3'.
Raw data files have been deposited at the National Center for Biotechnology Information Gene Expression Omnibus (NCBI GEO) with accession number GSE48040 and GSE48234.
- COR :
cold regulated gene
PARE library to analyze the RNA degradome of Brachypodium seedlings with cold treatment
PARE library to analyze the RNA degradome of Brachypodium seedlings without cold treatment
ligation-mediated poly(A) test
messenger ribonucleo proteins
parallel analysis of RNA ends
pyruvate phosphate dikinase
quantitative splinted-ligation reverse transcriptase polymerase chain reaction assay
RNA-high-throughput deep sequencing
RNA-Seq library for Brachypodium seedlings with cold treatment
RNA-Seq library for Brachypodium seedlings without cold treatment
RNA-induced silencing complexes
- STE kinases:
Shim J, Karin M: The control of mRNA stability in response to extracellular stimuli. Mol Cells. 2002, 14: 323-331.
Nakaminami K, Matsui A, Shinozaki K, Seki M: RNA regulation in plant abiotic stress responses. Biochim Biophys Acta. 2012, 1819: 149-153.
Belostotsky DA: State of decay: an update on plant mRNA turnover. Curr Top Microbiol Immunol. 2008, 326: 179-199.
Belostotsky DA, Sieburth LE: Kill the messenger: mRNA decay and plant development. Curr Opin Plant Biol. 2009, 12: 96-102.
Xu J, Chua N-H: Processing bodies and plant development. Curr Opin Plant Biol. 2011, 14: 88-93.
Chen X: MicroRNA metabolism in plants. Rna Interference. Edited by: Paddison PJVPK. 2008, 320: 117-136. Curr Top Microbiol Immunol
Ramachandran V, Chen X: Small RNA metabolism in Arabidopsis. Trends Plant Sci. 2008, 13: 368-374.
Jones-Rhoades MW, Bartel DP, Bartel B: MicroRNAs and their regulatory roles in plants. Annu Rev Plant Biol. 2006, 57: 19-53.
Guleria P, Mahajan M, Bhardwaj J, Yadav SK: Plant small RNAs: biogenesis, mode of action and their roles in abiotic stresses. Genomics Proteomics Bioinformatics. 2011, 9: 183-199.
Hori K, Watanabe Y: UPF3 suppresses aberrant spliced mRNA in Arabidopsis. Plant J. 2005, 43: 530-540.
Arciga-Reyes L, Wootton L, Kieffer M, Davies B: UPF1 is required for nonsense-mediated mRNA decay (NMD) and RNAi in Arabidopsis. Plant J. 2006, 47: 480-489.
Li Y, Kiledjian M: Regulation of mRNA decapping. Wiley Interdiscip RevRNA. 2010, 1: 253-265.
Ling SHM, Qamra R, Song H: Structural and functional insights into eukaryotic mRNA decapping. Wiley Interdiscip RevRNA. 2011, 2: 193-208.
Franks TM, Lykke-Andersen J: The control of mRNA decapping and P-Body formation. Mol Cell. 2008, 32: 605-615.
Shyu A-B, Wilkinson MF, van Hoof A: Messenger RNA regulation: to translate or to degrade. EMBO J. 2008, 27: 471-481.
Hu W, Sweet TJ, Chamnongpol S, Baker KE, Coller J: Co-translational mRNA decay in Saccharomyces cerevisiae. Nature. 2009, 461: 225-U296.
Eulalio A, Behm-Ansmant I, Izaurralde E: P bodies: at the crossroads of post-transcriptional pathways. Nat Rev Mol Cell Biol. 2007, 8: 9-22.
Parker R, Sheth U: P bodies and the control of mRNA translation and degradation. Mol Cell. 2007, 25: 635-646.
Weber C, Nover L, Fauth M: Plant stress granules and mRNA processing bodies are distinct from heat stress granules. Plant J. 2008, 56: 517-530.
Anderson P, Kedersha N: Stress granules. Curr Biol. 2009, 19: R397-R398.
Buchan JR, Muhlrad D, Parker R: P bodies promote stress granule assembly in Saccharomyces cerevisiae. J Cell Biol. 2008, 183: 441-455.
Balagopal V, Parker R: Polysomes, P bodies and stress granules: states and fates of eukaryotic mRNAs. Curr Opin Cell Biol. 2009, 21: 403-408.
Bruno I, Wilkinson MF: P-bodies react to stress and nonsense. Cell. 2006, 125: 1036-1038.
German MA, Luo S, Schroth G, Meyers BC, Green PJ: Construction of Parallel Analysis of RNA Ends (PARE) libraries for the study of cleaved miRNA targets and the RNA degradome. Nat Protoc. 2009, 4: 356-362.
Tang Z, Zhang L, Xu C, Yuan S, Zhang F, Zheng Y, Zhao C: Uncovering small RNA-mediated responses to cold stress in a wheat thermosensitive genic male-sterile line by deep sequencing. Plant Physiol. 2012, 159: 721-738.
Song Q-X, Liu Y-F, Hu X-Y, Zhang W-K, Ma B, Chen S-Y, Zhang J-S: Identification of miRNAs and their target genes in developing soybean seeds by deep sequencing. BMC Plant Biology. 2011, 11: 5-
Pantaleo V, Szittya G, Moxon S, Miozzi L, Moulton V, Dalmay T, Burgyan J: Identification of grapevine microRNAs and their targets using high-throughput sequencing and degradome analysis. Plant J. 2010, 62: 960-976.
Li Y-F, Zheng Y, Addo-Quaye C, Zhang L, Saini A, Jagadeeswaran G, Axtell MJ, Zhang W, Sunkar R: Transcriptome-wide identification of microRNA targets in rice. Plant J. 2010, 62: 742-759.
Bracken CP, Szubert JM, Mercer TR, Dinger ME, Thomson DW, Mattick JS, Michael MZ, Goodall GJ: Global analysis of the mammalian RNA degradome reveals widespread miRNA-dependent and miRNA-independent endonucleolytic cleavage. Nucleic Acids Res. 2011, 39: 5658-5668.
Mercer TR, Dinger ME, Bracken CP, Kolle G, Szubert JM, Korbie DJ, Askarian-Amiri ME, Gardiner BB, Goodall GJ, Grimmond SM, Mattick JS: Regulated post-transcriptional RNA cleavage diversifies the eukaryotic transcriptome. Genome Res. 2010, 20: 1639-1650.
Jiao Y, Riechmann JL, Meyerowitz EM: Transcriptome-wide analysis of uncapped mRNAs in Arabidopsis reveals regulation of mRNA degradation. Plant Cell. 2008, 20: 2571-2585.
Opanowicz M, Vain P, Draper J, Parker D, Doonan JH: Brachypodium distachyon: making hay with a wild grass. Trends Plant Sci. 2008, 13: 172-177.
The Resource for Brachypodium Genomics. [http://www.brachypodium.org]
Shin C, Nam J-W, Farh KK-H, Chiang HR, Shkumatava A, Bartel DP: Expanding the microRNA targeting code: functional sites with centered pairing. Mol Cell. 2010, 38: 789-802.
German MA, Pillay M, Jeong D-H, Hetawal A, Luo S, Janardhanan P, Kannan V, Rymarquis LA, Nobuta K, German R, De Paoli E, Lu C, Schroth G, Meyers BC, Green PJ: Global identification of microRNA-target RNA pairs by parallel analysis of RNA ends. Nature Biotechnol. 2008, 26: 941-946.
Eckmann CR, Rammelt C, Wahle E: Control of poly(A) tail length. Wiley Interdiscip RevRNA. 2011, 2: 348-361.
Chen CY, Shyu AB: Deadenylation and p-bodies. AdvExpMedBiol. 2013, 768: 183-195.
Du Z, Zhou X, Ling Y, Zhang Z, Su Z: agriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010, 38: W64-W70.
Schmidt M, Dehne S, Feierabend J: Post-transcriptional mechanisms control catalase synthesis during its light-induced turnover in rye leaves through the availability of the hemin cofactor and reversible changes of the translation efficiency of mRNA. Plant J. 2002, 31: 601-613.
Nevo R, Charuvi D, Tsabari O, Reich Z: Composition, architecture and dynamics of the photosynthetic apparatus in higher plants. Plant J. 2012, 70: 157-176.
Chaves MM, Flexas J, Pinheiro C: Photosynthesis under drought and salt stress: regulation mechanisms from whole plant to cell. Ann Bot. 2009, 103: 551-560.
Makino A: Photosynthesis, grain yield, and nitrogen utilization in rice and wheat. Plant Physiol. 2011, 155: 125-129.
Wang H-F, Feng L, Niu D-K: Relationship between mRNA stability and intron presence. Biochem Biophys Res Commun. 2007, 354: 203-208.
Zhao C, Hamilton T: Introns regulate the rate of unstable mRNA decay. J Biol Chem. 2007, 282: 20230-20237.
Narsai R, Howell KA, Millar AH, O'Toole N, Small I, Whelan J: Genome-wide analysis of mRNA decay rates and their determinants in Arabidopsis thaliana. Plant Cell. 2007, 19: 3418-3436.
Chakkalakal JV, Miura P, Belanger G, Michel RN, Jasmin BJ: Modulation of utrophin A mRNA stability in fast versus slow muscles via an AU-rich element and calcineurin signaling. Nucleic Acids Res. 2008, 36: 826-838.
Masaru OT: Control of mRNA stability in plants cis-acting elements which affect mRNA stability. Tanpakushitsu Kakusan Koso. 1993, 38: 2367-2371.
Bailey TL, Williams N, Misleh C, Li WW: MEME: discovering and analyzing DNA and protein sequence motifs. Nucleic Acids Res. 2006, 34: W369-W373.
Goeres DC, Van Norman JM, Zhang W, Fauver NA, Spencer ML, Sieburth LE: Components of the Arabidopsis mRNA decapping complex are required for early seedling development. Plant Cell. 2007, 19: 1549-1564.
Nishimura N, Kitahata N, Seki M, Narusaka Y, Narusaka M, Kuromori T, Asami T, Shinozaki K, Hirayama T: Analysis of ABA hypersensitive germination2 revealed the pivotal functions of PARN in stress response in Arabidopsis. Plant J. 2005, 44: 972-984.
Yoine M, Ohto M-a, Onai K, Mita S, Nakamura K: The lba1 mutation of UPF1 RNA helicase involved in nonsense-mediated mRNA decay causes pleiotropic phenotypic changes and altered sugar signalling in Arabidopsis. Plant J. 2006, 47: 49-62.
Chekanova JA, Gregory BD, Reverdatto SV, Chen H, Kumar R, Hooker T, Yazaki J, Li P, Skiba N, Peng Q, Alonso J, Brukhin V, Grossniklaus U, Ecker JR, Belostotsky DA: Genome-wide high-resolution mapping of exosome substrates reveals hidden features in the Arabidopsis transcriptome. Cell. 2007, 131: 1340-1353.
Estevez AM, Lehner B, Sanderson CM, Ruppert T, Clayton C: The roles of intersubunit interactions in exosome stability. J Biol Chem. 2003, 278: 34943-34951.
Mitchell P, Petfalski E, Shevchenko A, Mann M, Tollervey D: The exosome: a conserved eukaryotic RNA processing complex containing multiple 3'->5' exoribonucleases. Cell. 1997, 91: 457-466.
Carthew RW, Sontheimer EJ: Origins and mechanisms of miRNAs and siRNAs. Cell. 2009, 136: 642-655.
Bartel DP: MicroRNAs: Genomics, biogenesis, mechanism, and function. Cell. 2004, 116: 281-297.
Zhang J, Xu Y, Huan Q, Chong K: Deep sequencing of Brachypodium small RNAs at the global genome level identifies microRNAs involved in cold stress response. BMC Genomics. 2009, 10:
Huntzinger E, Izaurralde E: Gene silencing by microRNAs: contributions of translational repression and mRNA decay. Nat Rev Genet. 2011, 12: 99-110.
Diedhiou CJ, Popova OV, Dietz K-J, Golldack D: The SNF1-type serine-threonine protein kinase SAPK4 regulates stress-responsive gene expression in rice. BMC Plant Biol. 2008, 8: 49-
Pais SM, Tellez-Inon MT, Capiati DA: Serine/threonine protein phosphatases type 2A and their roles in stress signaling. Plant SignalBehav. 2009, 4: 1013-1015.
Kido EA, de Araujo Barbosa PK, Costa Ferreira Neto JR, Pandolfi V, Houllou-Kido LM, Crovella S, Molina C, Kahl G, Benko-Iseppon AM: Identification of plant protein kinases in response to abiotic and biotic stresses using SuperSAGE. Curr Protein Pept Sci. 2011, 12: 643-656.
Saidi Y, Hearn TJ, Coates JC: Function and evolution of 'green' GSK3/Shaggy-like kinases. Trends Plant Sci. 2012, 17: 39-46.
Rouhier N, Lemaire SD, Jacquot J-P: The role of glutathione in photosynthetic organisms: Emerging functions for glutaredoxins and glutathionylation. Ann Rev Plant Biol. 2008, 59: 143-166.
Clark G, Konopka-Postupolska D, Hennig J, Roux S: Is annexin 1 a multifunctional protein during stress responses?. Plant SignalBehav. 2010, 5: 303-307.
Lee S, Lee EJ, Yang EJ, Lee JE, Park AR, Song WH, Park OK: Proteomic identification of annexins, calcium-dependent membrane binding proteins that mediate osmotic stress and abscisic acid signal transduction in Arabidopsis. Plant Cell. 2004, 16: 1378-1391.
Vanderauwera S, Suzuki N, Miller G, van de Cotte B, Morsa S, Ravanat J-L, Hegie A, Triantaphylides C, Shulaev V, Van Montagu MCE, et al: Extranuclear protection of chromosomal DNA from oxidative stress. Proc Natl Acad Sci USA. 2011, 108: 1711-1716.
Yasuda Y, Miyamoto Y, Yamashiro T, Asally M, Masui A, Wong C, Loveland KL, Yoneda Y: Nuclear retention of importin alpha coordinates cell fate through changes in gene expression. EMBO J. 2012, 31: 83-94.
Lee B-h, Kapoor A, Zhu J, Zhu J-K: STABILIZED1, a stress-upregulated nuclear protein, is required for pre-mRNA splicing, mRNA turnover, and stress tolerance in Arabidopsis. Plant Cell. 2006, 18: 1736-1749.
Lorkovic ZJ, Kirk DAW, Lambermon MHL, Filipowicz W: Pre-mRNA splicing in higher plants. Trends Plant Sci. 2000, 5: 160-167.
Floris M, Mahgoub H, Lanet E, Robaglia C, Menand B: Post-transcriptional regulation of gene expression in plants during abiotic stress. Int J Mol Sci. 2009, 10: 3168-3185.
Pomeranz M, Finer J, Jang J-C: Putative molecular mechanisms underlying tandem CCCH zinc finger protein mediated plant growth, stress, and gene expression responses. Plant SignalBehav. 2011, 6: 647-651.
Blewett N, Coller J, Goldstrohm A: A quantitative assay for measuring mRNA decapping by splinted ligation reverse transcription polymerase chain reaction: qSL-RT-PCR. RNA. 2011, 17: 535-543.
Murray EL, Schoenberg DR: Assays for determining poly(A) tail length and the polarity of mRNA decay in mammalian cells. Methods Enzymol. 2008, 448: 483-504.
Winfield MO, Lu C, Wilson ID, Coghill JA, Edwards KJ: Plant responses to cold: transcriptome analysis of wheat. Plant Biotechnol J. 2010, 8: 749-771.
Deyholos MK: Making the most of drought and salinity transcriptomics. Plant Cell Environ. 2010, 33: 648-654.
Hazen SP, Wu Y, Kreps JA: Gene expression profiling of plant responses to abiotic stress. Funct Integr Genomics. 2003, 3: 105-111.
Munchel SE, Shultzaberger RK, Takizawa N, Weis K: Dynamic profiling of mRNA turnover reveals gene-specific and system-wide regulation of mRNA decay. Mol Biol Cell. 2011, 22: 2787-2795.
Ensminger I, Busch F, Huner NPA: Photostasis and cold acclimation: sensing low temperature through photosynthesis. Physiologia Plantarum. 2006, 126: 28-44.
Kocova M, Hola D, Wilhelmova N, Rothova O: The influence of low-temperature on the photochemical activity of chloroplasts and activity of antioxidant enzymes in maize leaves. Biologia Plantarum. 2009, 53: 475-483.
Huner NPA, Oquist G, Hurry VM, Krol M, Falk S, Griffith M: Photosynthesis, photoinhibition and low-temprature acclimation in cold tolerant plants. Photosynthesis Research. 1993, 37: 19-39.
Hurry V, Strand A, Furbank R, Stitt M: The role of inorganic phosphate in the development of freezing tolerance and the acclimatization of photosynthesis to low temperature is revealed by the pho mutants of Arabidopsis thaliana. Plant J. 2000, 24: 383-396.
Oquist G, Huner NPA: Cold-hardening-induced resistance to photoinhibition of photosynthesis in winter rye is dependent upon an increased capacity for photosynthesis. Planta. 1993, 189: 150-156.
Ming Zhou, Lianfeng Gu, Pingchuan Li, Xianwei Song, Liya Wei, Zhiyu CHen, Cao X: Degradome sequencing reveals endogenous small RNA targets in rice (Oryza sativa L. ssp. indica). Frontiers in Biology. 2010, 5: 67-90.
Zhang Z, Yu J, Li D, Zhang Z, Liu F, Zhou X, Wang T, Ling Y, Su Z: PMRD: plant microRNA database. Nucleic Acids Res. 2010, 38: D806-D813.
Nyiko T, Kerenyi F, Szabadkai L, Benkovics AH, Major P, Sonkoly B, Merai Z, Barta E, Niemiec E, Kufel J, Silhavy D: Plant nonsense-mediated mRNA decay is controlled by different autoregulatory circuits and can be induced by an EJC-like complex. Nucleic Acids Res. 2013, 41: 6715-6728.
Houseley J, Tollervey D: The many pathways of RNA degradation. Cell. 2009, 136: 763-776.
Deutscher MP: Degradation of stable RNA in bacteria. J Biol Chem. 2003, 278: 45041-45044.
Deutscher MP: Degradation of RNA in bacteria: comparison of mRNA and stable RNA. Nucleic Acids Res. 2006, 34: 659-666.
Kushner SR: mRNA decay in prokaryotes and eukaryotes: different approaches to a similar problem. Iubmb Life. 2004, 56: 585-594.
Liu H, Kiledjian M: An erythroid-enriched endoribonuclease (ErEN) involved in alpha-globin mRNA turnover. Protein Pept Lett. 2007, 14: 131-136.
Tomecki R, Dziembowski A: Novel endoribonucleases as central players in various pathways of eukaryotic RNA metabolism. RNA. 2010, 16: 1692-1724.
Shapiro R, Riordan JF, Vallee BL: Characteristic ribonucleolytic activity of human angiogenin. Biochemistry. 1986, 25: 3527-3532.
Fu H, Feng J, Liu Q, Sun F, Tie Y, Zhu J, Xing R, Sun Z, Zheng X: Stress induces tRNA cleavage by angiogenin in mammalian cells. Febs Letters. 2009, 583: 437-442.
Gilks N, Kedersha N, Ayodele M, Shen L, Stoecklin G, Dember LM, Anderson P: Stress granule assembly is mediated by prion-like aggregation of TIA-1. Mol Biol Cell. 2004, 15: 5383-5398.
Yang F, Peng Y, Murray EL, Otsuka Y, Kedersha N, Schoenberg DR: Polysome-bound endonuclease PMR1 is targeted to stress granules via stress-specific binding to TIA-1. Mol Cell Biol. 2006, 26: 8803-8813.
Jackowiak P, Nowacka M, Strozycki PM, Figlerowicz M: RNA degradome-its biogenesis and functions. Nucleic Acids Res. 2011, 39: 7361-7370.
Ma Z, Coruh C, Axtell MJ: Arabidopsis lyrata Small RNAs: Transient MIRNA and Small Interfering RNA Loci within the Arabidopsis Genus. Plant Cell. 2010, 22: 1090-1103.
The Rfam database. [http://www.sanger.ac.uk/Software/Rfam]
Repbase Update. [http://www.girinst.org/server/RepBase]
Langmead B: Aligning short sequencing reads with Bowtie. CurrProtocBioinformatics. 2010, Chapter 11: Unit 11.17
International Brachypodium I: Genome sequencing and analysis of the model grass Brachypodium distachyon. Nature. 2010, 463: 763-768.
Li R, Yu C, Li Y, Lam T-W, Yiu S-M, Kristiansen K, Wang J: SOAP2: an improved ultrafast tool for short read alignment. Bioinformatics. 2009, 25: 1966-1967.
Folkes L, Moxon S, Woolfenden HC, Stocks MB, Szittya G, Dalmay T, Moulton V: PAREsnip: a tool for rapid genome-wide discovery of small RNA/target interactions evidenced through degradome sequencing. Nucleic Acids Res. 2012, 40: e103-
The microRNA database. [http://www.mirbase.org]
Fisher's Exact Test. [http://www.langsrud.com/fisher.htm]
Gutierrez RA, Ewing RM, Cherry JM, Green PJ: Identification of unstable transcripts in Arabidopsis by cDNA microarray analysis: Rapid decay is associated with a group of touch- and specific clock-controlled genes. Proc Natl Acad Sci USA. 2002, 99: 11513-11518.
Noel GM, Tognetti JA, Pontis HG: Protein kinase and phosphatase activities are involved in fructan synthesis initiation mediated by sugars. Planta. 2001, 213: 640-646.
Seeley KA, Byrne DH, Colbert JT: Red light-independent instability of oat phytochrome messenger-RNA invivo. Plant Cell. 1992, 4: 29-38.
This study was supported by the National Natural Science Foundation of China (NSFC: 30970241), the Main Direction Program of Knowledge Innovation of Chinese Academy of Sciences (KSCX3-EW-N-07 and KSCX2-EW-J-1).
KC and JZ conceived and designed the experiments; JZ performed the experiments; ZM conducted most of the bioinformatic analysis; JZ and KC analyzed the data; JZ and KC wrote the paper. All authors have read and approved the manuscript for publication.
Electronic supplementary material
Additional file 1: . Figure S1. Cold tolerance and vernalization responsiveness assessment of two Brachypodium diploid genotypes, ABR5 and BD21. (a) Seeds for ABR5 and BD21 were germinated and grown for 2 weeks under the growth condition described above. Then the seedlings were treated with -5°C for 6h, 12h, 18h, 24h, and recovered for 10 days. Survival rate was investigated after cold treatment and recovery. Error bars represent standard error. Each time point consisted of 10 to 15 plants, and the experiment was repeated three times.(b) Seeds for ABR5 and BD21 were treated with °C for 1, 2, 3, 4, 5, and 6 weeks as vernalization treatment and grown in a growth chamber at 24°C, with 6 h/8 h (light/dark) photoperiod and approximately 5,000 lux light intensity. Flowering time was measured by the number of days after germination until the first flower bud appears. Error bars represent standard error. Each time point consisted of 10 to 15 plants, and the experiment was repeated three times. Figure S2. Reproducibility of PARE libraries. Scatter plot showing the reproducibility between PARE libraries and their biological repeats. The reads from DW (a) and DC (b) library was mapped to each transcript and the sum of total reads was plotted against corresponding reads from their biological replicate libraries, DW replicate and DC replicate, respectively. Figure S3. Real-time RT-PCR validation of gene expression data provided by RNA-Seq analysis. Total RNA was isolated from 12-day-old Brachypodium seedlings with (WC) and without (WO) cold treatment. Real-time PCR analysis was performed as described by Dai et al. (Plant Physiol 2007, 143:1739). The SuperScripts II reverse transcriptase (Invitrogen, Carlsbad, CA, USA) was used for reverse transcription and SYBR Green PCR Master mix (Applied Biosystems, Foster City, CA, USA) was used for quantitative PCR. The amplification of a Brachypodium Tubulin gene was used as an internal control to normalize all data. The normalized mRNA levels in the WO samples were arbitrarily set to 1. Quantitative PCR reactions were repeated three times. Error bars represent the standard error. Figure S5. Overall mRNA stability estimation through quantification of mRNA abundance before and after cordycepin (transcription inhibitor) treatment. The middle part of the primary leaves from 12-day-old Brachypodium seedlings with (WC) and without (WO) cold treatment (4°C, 24h) were used. Cordycepin treatment was performed as described (Plant J 2002, 31:601). The 3-cm leaf segments were incubated with 600µM cordycepin solution with regular shaking at 22°C for 7h. Then tissue samples were harvested at regular intervals (1h, 2h, and 4h) and quickly frozen in liquid nitrogen. Total RNA was isolated with mirVana miRNA Isolation Kit (Ambion, Austin, TX, USA). Two microgram of total RNA was used for reverse transcription with SuperScript II reverse transcriptase (Invitrogen). The obtained cDNAs were quantified and used as template for quantitative PCR using SYBR Green PCR Master mix (Applied Biosystems, Foster City, USA). For both WO and WC samples, mRNA abundance for the time zero was arbitrarily set to 1 and mRNA abundance for the other time points were shown as the relative value of time zero. Quantitative PCR reactions were repeated three times. Error bars indicate standard deviation. Figure S6. mRNA features correlated with different degradation patterns. The mRNA length, UTR length, GC content, and the number of introns were analyzed for type I-IV mRNAs. N: no conserved motif (E-value <0.001) was identified; Y: conserved motifs were identified and the motifs for Y1, Y2, Y3, and Y4 were indicated below the table. The x-axis represented nucleotide position. The y-axis represented the information content measured in bits. The overall height of each stack indicated the sequence conservation at that position, and the height of a letter within the stack indicated the relative frequency of corresponding nucleotide in the motif. Figure S7. Conserved motifs identified in the 5' and 3' UTRs for subgroups of type I, II, III, and IV genes. I, II, III, and IV: different mRNA decay patterns classified according to the variation tend of uncapped transcript/transcript abundance in cold stress. D: uncapped transcript abundance indicated by PARE reads; R: transcript abundance indicated by RNA-Seq reads; N: no conserved motif (E-value <0.001) was identified. Conserved motifs were analyzed with MEME software (Nucleic Acids Res. 2006, 34:W369) and identified conserved motifs with E-value <0.001 were shown in the table. Figure S8. Endonucleolytic cleavages identified in the DW and DC library. Transcripts with endonucleolytic cleavages shared different conserved motifs. Conserved motifs, at the endonucleolytic cleavage sites, were marked with gray background. Cleavage sites were indicated by the asterisk. The putative functions of all the transcripts with endonucleolytic cleavages were shown. The majority of the detected endonucleolytic cleavages were also found in corresponding biological replicate PARE libraries. Underlined gene names denoted endonucleolytic cleavages that were not detected in corresponding biological replicate PARE libraries.EW: endonucleolytic cleavage groups only found in the DW library;EC: endonucleolytic cleavage groups only found in the DC library;GONG: endonucleolytic cleavage groups found in both DW and DC library. (PDF 1 MB)
Additional file 2: . Table S1. Summary of PARE analysis. Table S2. Top 20 genes in DW and DC library. Table S3. Top 10 upregulated and downregulated genes in PARE analysis. Table S4. Enrichment analysis with Fisher's exact test. Table S5. Predicted cold-responsive Brachypodium miRNAs and their putative targets. Table S6. The conservation of motifs identified at the endonucleolytic cleavage sites. (XLSX 20 KB)
Additional file 3: . Figures S4. Relationship between mRNA decay pattern and gene function. Gene Ontology (GO) analysis was performed for the type I, II, III, and IV genes using agriGO which organizes information for molecular function, biological process, and cellular component categories.type I: (a), molecular function; (b),biological process; (c) cellular component;type II: (d), molecular function; (e),biological process; (f) cellular component;type III: (g), molecular function; (h),biological process; (i) cellular component;type IV: (j), molecular function; (k),biological process; (l) cellular component. (RAR 4 MB)
Authors’ original submitted files for images
About this article
Cite this article
Zhang, J., Mao, Z. & Chong, K. A global profiling of uncapped mRNAs under cold stress reveals specific decay patterns and endonucleolytic cleavages in Brachypodium distachyon. Genome Biol 14, R92 (2013). https://doi.org/10.1186/gb-2013-14-8-r92
- mRNA stability
- cold stress
- Brachypodium distachyon