Transcriptomic analyses of rice (Oryza sativa) genes and non-coding RNAs under nitrogen starvation using multiple omics technologies
BMC Genomics volume 19, Article number: 532 (2018)
Nitrogen (N) is a key macronutrient essential for plant growth, and its availability has a strong influence on crop development. The application of synthetic N fertilizers on crops has increased substantially in recent decades; however, the applied N is not fully utilized due to the low N use efficiency of crops. To overcome this limitation, it is important to understand the genome-wide responses and functions of key genes and potential regulatory factors in N metabolism.
Here, we characterized changes in the rice (Oryza sativa) transcriptome, including genes, newly identified putative long non-coding RNAs (lncRNAs), and microRNAs (miRNAs) and their target mRNAs in response to N starvation using four different transcriptome approaches. Analysis of rice genes involved in N metabolism and/or transport using strand-specific RNA-Seq identified 2588 novel putative lncRNA encoding loci. Analysis of previously published RNA-Seq datasets revealed a group of N starvation-responsive lncRNAs showing differential expression under other abiotic stress conditions. Poly A-primed sequencing (2P-Seq) revealed alternatively polyadenylated isoforms of N starvation-responsive lncRNAs and provided precise 3′ end information on the transcript models of these lncRNAs. Analysis of small RNA-Seq data identified N starvation-responsive miRNAs and down-regulation of miR169 family members, causing de-repression of NF-YA, as confirmed by strand-specific RNA-Seq and qRT-PCR. Moreover, we profiled the N starvation-responsive down-regulation of root-specific miRNA, osa-miR444a.4-3p, and Degradome sequencing confirmed MADS25 as a novel target gene.
In this study, we used a combination of multiple RNA-Seq analyses to extensively profile the expression of genes, newly identified lncRNAs, and microRNAs in N-starved rice roots and shoots. Data generated in this study provide an in-depth understanding of the regulatory pathways modulated by N starvation-responsive miRNAs. The results of comprehensive, large-scale data analysis provide valuable information on multiple aspects of the rice transcriptome, which may be useful in understanding the responses of rice plants to changes in the N supply status of soil.
Nitrogen (N) is a key macronutrient for plants and has a strong influence on crop development and productivity. To increase crop yield, the application of synthetic N fertilizers to crops has increased substantially in recent decades. However, plants utilize less than half of the applied N because of low N use efficiency (NUE) and uptake saturation . The N fertilizers remaining in the soil poses several environmental problems, such as eutrophication. Moreover, increased total costs for N fertilizers lead to increases of product prices and a reduction in the farmer’s profitability. Therefore, improving crop NUE while maintaining crop productivity has several economic and environmental benefits.
Transgenic plants with improved NUE have been developed in which expression of protein-coding genes involved in N uptake, assimilation, and transport have been modulated by genetic engineering . Overexpression of nitrate transporters [3,4,5,6] or ammonium transporters [7,8,9,10] led to enhanced N source uptake ability and increased nitrate and ammonium contents in transgenic plants. Overexpression of N assimilation enzymes, including alanine aminotransferase (AlaAT) [11, 12] and glutamine synthetase (GS) increased total N content and plant dry biomass and produced yield increases. In parallel with these genetic approaches, a number of N-responsive genes were identified using high-throughput analysis tools such as microarrays and next-generation sequencing (NGS) platforms [13,14,15,16]. In these studies, the expression profiles of genes have been examined under N-deficient conditions, and novel candidate genes for improving the NUE of crop plants via genetic engineering have been identified. Of the identified N-responsive genes, overexpression of OsENOD93–1 (a rice early nodulin gene induced by N starvation) stimulated increases in shoot dry biomass and improvements in seed yield .
Recent studies of regulatory non-coding RNAs (ncRNAs), including microRNAs (miRNAs) and long non-coding RNAs (lncRNAs), identified their functional importance in many biological phenomena in plants including processes related to crop agricultural traits [17, 18]. Several large-scale investigations and functional studies of various stress-responsive ncRNAs in crops and other plants suggested that these regulatory ncRNAs had significant effects on physiological responses through regulation of gene [19,20,21,22,23,24,25,26,27,28]. Modulating expression of these regulatory ncRNAs could lead to improvements in some important agricultural traits such as productivity, male sterility, nutrient homeostasis, and floral organogenesis [19, 21, 29, 30]. Because of their functional significance and potential impact on agricultural traits, researchers are targeting ncRNAs as a novel resource for crop improvement [17, 31,32,33]. Researchers showed that changes in N supply status altered the expression of multiple miRNAs [34,35,36,37,38] and lncRNAs [25, 26] in several plant species. Many N-responsive miRNAs in plants were found to be involved in N metabolism [24, 36, 39], reprogramming root development , anthocyanin biosynthesis , and homeostasis of other nutrients . However, while species-wide studies on plant miRNAs have been performed, no studies examining the involvement of lncRNAs in regulating N responses or metabolism in plants have been reported, with the exception of research in Populus tomentosa .
Many previous studies investigated the rice transcriptome in N-starved conditions using large-scale datasets to identify N-responsive genes and miRNAs [13,14,15,16, 34]. However, these studies were limited to a single type of RNA molecules, such as mRNAs or small RNAs, and did not provide an overview of the transcriptome-wide responses to changes in N availability. In this study, we combined multiple RNA-Seq analyses to assess multiple aspects of transcriptomes in N-starved rice, and to provide detailed information on transcriptome-wide changes in response to N availability. We profiled and analyzed diverse aspects of the rice transcriptome including genes, putative lncRNAs, and miRNAs and their target mRNAs, using four different types of RNA-Seq methods in N-starved rice samples prepared over a time course. These large-scale datasets analyzed in this study revealed N-responsive characteristics and expanded our knowledge of regulatory activities in N-starved rice, and provided insights towards understanding the molecular mechanisms underpinning modulation of N homeostasis in rice.
Analysis of N-starved rice shoots and roots via various transcriptomic approaches
Four transcriptomic experimental tools were used to investigate the changes occurring in N-starved rice over a time course: strand-specific RNA-Seq, small RNA-Seq, 2P-Seq (poly A-primed sequencing), and Degradome sequencing (Fig. 1). N-starved samples were examined along with samples of 4-week-old rice plants grown under normal conditions (-N_0d) (Fig. 1a). All samples were analyzed by strand-specific RNA-Seq and small RNA-Seq. 2P-Seq and Degradome sequencing were applied to selected samples based on the analysis of the strand-specific RNA-Seq and small RNA-Seq datasets (Fig. 1b, Additional file 1: Figure. S1).
Approximately 513 million 101-bp paired-end reads generated from 12 strand-specific RNA-Seq datasets were mapped onto the rice genome (IRGSP-1.0). These reads were used for expression profiling of IRGSP-annotated genes (RAP-DB), and identifying previously unannotated putative lncRNAs based on the reference annotation-based transcript (RABT) assembly using Cufflinks . In addition, 117 million reads from three 2P-Seq analyses were examined to characterize the 3′ ends of IRGSP-annotated genes and novel putative lncRNAs. Together with strand-specific RNA-Seq reads, strand-specific 2P-Seq reads were helpful to confirm the expression of cis-natural antisense transcripts (cis-NATs). For rice miRNA profiling, 287 million genome-aligned small RNA-Seq reads were analyzed. After expression profiling of previously annotated rice miRNAs, Degradome data analysis was performed in 7 day N-starved roots to investigate targets of N starvation-responsive miRNAs. Altogether, these omics-based high-throughput analyses provided a spatio-temporal characterization of the rice transcriptome in response to N starvation.
Gene expression profiling of genes involved in N source transport and assimilation using strand-specific RNA-Seq data
First, the expression patterns of genes involved in N transport and assimilation were profiled via analysis of RNA-Seq data (Fig. 2a and b), since the transcript abundance of these genes were known to change dynamically in response to N availability [5, 43,44,45,46,47]. Transcript levels of genes encoding nitrate transporters (NRTs) and ammonium transporters (AMTs), varied more dynamically in roots than in shoots in response to N starvation (Fig. 2b). The transcription of NRT2.1, NRT2.2, NRT2.3a, NRT2.3b, and NAR2.2 increased in response to N starvation, whereas that of NRT/CHL (Os10g0554200, AtNRT1.1 homolog) and CLCa (Os12g0438600, AtCLCa homolog) decreased. The expression of NRT2.1 and NRT2.2 continued to increase gradually until 5 days after N starvation, and then reduced to un-induced levels at 7 days after N starvation. The two isoforms of NRT2.3 were gradually up-regulated until 7 days after N starvation (Fig. 2b and d), suggesting temporal-specific regulation of NRT expression in response to prolonged N starvation. Interestingly, one of the two down-regulated NRTs in both roots and shoots, CLCa, is a rice homolog of AtCLCa, which specifically localizes to the vacuole membrane and acts as an NO3−/H+ antiporter . Considering that leaf vacuoles function as major nitrate storage pools , the down-regulation of CLCa might induce the release of vacuole-reserved nitrate to the cytoplasm to maintain nitrate homeostasis in local tissue.
The AMT1.2 and AMT3.3 genes showed significant gradual increases in expression over time in N-starved rice roots, while AMT1.3 and AMT2.2 showed a gradual reduction in expression levels (Fig. 2b). These results are consistent with previous studies [45, 46]. These reciprocal expression changes in N-starved roots appeared to be part of a pre-programmed response to ammonium deficiency, because many isozymes involved in N assimilation showed similar responses. These results are consistent with previous expression analyses of rice GS, GOGAT, and asparagine synthetase (ASNS) [44, 50,51,52,53]. Taken together, these results suggest that genes that are down-regulated in response to N starvation (AMT1.3, GS1.2, and ASNS1) are mainly responsible for primary ammonium assimilation in N-sufficient conditions, whereas genes that are induced by N starvation (AMT1.2, GS1.1, and ASNS2) are active during N-deficient conditions in roots. This is supported by functional characterization studies of GS isozymes and two AMTs in rice . Overall, our profiling results indicate that plants use dynamic reprogramming of gene expression levels to adapt to N-deficient conditions.
Gene expression profiling of genes involved in amino acid/peptide transport
Effective distribution and localization of amino acids and peptides is necessary for rice growth and development. However, little is known about the genome-wide responses of rice genes encoding amino acid and peptide transporters in N-starved conditions. Here, genes encoding the amino acid transporter (AAT) family and peptide transporter (PTR) family of proteins were profiled in N-starved rice. A total of 191 genes with relevant annotations in Rice Annotation Project Database (RAP-DB; http://rapdb.dna.affrc.go.jp/). RAP-DB were identified; these included 82 out of 85 AAT genes , 82 out of 85 PTR genes , and 27 additional RAP-DB genes harboring the oligopeptide transporter domains (PF03169.8, E-value >5E-50) (Additional file 2: Table S3). Among these, 68 AAT, 68 PTR, and 13 OPT genes with expression levels > 1.0 FPKM in at least one out of ten N-starved rice RNA-Seq datasets were selected and further analyzed.
Of these genes, the expression of 39 AAT, 37 PTR, and 9 OPT genes was altered by more than 2-fold in at least one out of four N-starved conditions in roots and/or shoots. Hierarchical clustering results showed that AATs/PTRs/OPTs were generally up-regulated both in roots and/or in shoots in response to N starvation (Fig. 2c). Among the genes differentially expressed in 7 day N-starved roots and/or shoots, 31 and 24 genes were up-regulated, whereas only ten and two genes were down-regulated by > 2-fold in roots and shoots, respectively. These results indicated that amino acid/peptide transport systems were activated in roots and shoots in response to N starvation.
Among the amino acid/peptide transporter genes profiled, 64/30 genes were up−/down-regulated in at least one out of eight N-starved rice root and shoot datasets, respectively; 42/24 genes were up−/down-regulated in roots, and 33/9 were up−/down-regulated in shoots. These patterns suggested that most of the AATs/PTRs/OPTs were up-regulated in response to N starvation, and that the dynamic responses of these genes occurred in the roots rather than in the shoots; this was similar to the expression patterns of N metabolism-related genes. Additionally, many of the AATs/PTRs/OPTs showed tissue-specific regulation of expression, especially in roots. While only 15 of the genes exhibited expression changes in both roots and shoots at the same timing in N-starved conditions, 44 and 18 genes exhibited altered gene expression exclusively in roots and shoots, respectively. Tissue-dependent concordant/discordant regulation of some AATs/PTRs/OPTs was also observed. Seven genes showed concordant patterns of up-regulation in both roots and shoots: OsAAP13 (Os04g0470700), OsBAT1 (Os01g0607200), OsBAT3 (Os04g0470700), OsCAT3 (Os03g0641200), OsPTR12 (Os01g0872100), OsOPT09 (Os02g0695800), and OsOPT13 (Os04g0524500). Conversely, OsAAP5 (Os01g0878400), OsATL7 (Os01g0825800), OsPTR07 (Os01g0871500), and OsPTR08 (Os01g0871600) were discordantly regulated between two tissues.
Notably, one of the most abundantly expressed PTR genes, OsPTR07 was down-regulated in roots and up-regulated in shoots. OsPTR08, an evolutionary tandem-duplicate of OsPTR07, showed similar expression patterns to OsPTR07 both in N-starved roots and shoots; however, RNA-Seq data showed higher expression of OsPTR07 compared with OsPTR08 in RNA-Seq results (Additional file 2: Table S3). According to a previous report, many OsPTR genes are thought to be paralogous, and most (~ 67%) are tandemly duplicated; nevertheless, the paralogs appear to be differentially regulated in a tissue-specific manner, suggesting that these tandem duplicates of OsPTR genes are undergoing sub-functionalization and neo-functionalization . The different expression levels observed for OsPTR07 and OsPTR08 in shoots in response to N starvation may be an example of PTR sub-functionalization. This tissue-dependent concordant/discordant regulation of OsPTR genes, including OsPTR07 and OsPTR08, in response to N starvation might allow modulation of local amino acid/peptide transport in rice.
Profiling N starvation-responsive genes co-regulated by phosphate status
Along with N, inorganic phosphate (Pi) is a key macronutrient that is required for optimal growth and productivity of crop plants. Interactions between the regulatory systems governing these key nutrients are suspected at the physiological and gene expression levels [56,57,58,59]. To examine the regulatory effects of limiting N or Pi on homeostasis of the other nutrient at the transcriptome level in rice, we profiled genes responding to both N and Pi starvation using strand-specific RNA-Seq data collected in this study, and non-strand-specific RNA-Seq data obtained from a previous Pi starvation study  in which the time course and growth conditions were similar to those used in this study.
Expression profiling in 1 day, 3 day, and 7 day N- or Pi-starved rice indicated that some genes involved in the transport or metabolism of each nutrient were partly regulated in response to starvation conditions for the other nutrient. In 7 day N- or Pi-starved rice, 776 genes were up- or down- regulated by both N and Pi starvation (Fig. 3). Several of these genes were involved in N transport and metabolism such as AMTs, amino acid/peptide transporters, and N assimilating enzymes, and were also regulated by Pi starvation in roots (Table 1). Furthermore, Pi starvation-responsive genes like SPX2 and IPS1 were regulated in the roots of N-starved rice. Both AMT1.2 and AMT1.3 were down-regulated in Pi-starved roots, whereas only AMT1.2 was up-regulated in N-starved roots. IPS1, which is induced by Pi starvation in roots and in shoots of rice , exhibited discordant expression patterns between shoot and root in response to N starvation (Table 1). Collectively, these expression patterns suggested that, although the regulatory effects of N or Pi starvation on the metabolic pathways of the other nutrient appeared to be locally limited or partial, deficiency of N or Pi could regulate transcription of genes involved in their own and the other nutrient’s metabolic pathways.
The set of genes exhibiting responses to both N and Pi starvation was examined using GO analysis. GO results showed that genes up-regulated in 7 day N- and Pi-starved rice shoots were enriched for cellular homeostasis (GO:0019725) (Additional file 1: Figure. S2). Five of the six cellular homeostasis genes encoding thioredoxin fold domain-containing proteins (Os01g0194600, Os05g0198200, Os10g0482900) and two metallothioneins (Os03g0288000, Os050202800), were previously shown to be involved in antioxidant activity [61,62,63]. Oxidative stress is induced under nutrient-deficient conditions  and can cause deleterious damage to important metabolic enzymes such as GS . Increased expression of antioxidant genes in N- or Pi-starved conditions might provide protection from damage by scavenging harmful oxidants two cytokinin oxidase/dehydrogenases (OsCKX1, OsCKX4) involved in cytokinin degradation are up-regulated in 7 day N- and Pi-starved rice roots. The in planta status of cytokinin can affect the expression of genes involved in N and Pi pathways [65, 66]. N and Pi starvation may also induce cytokinin degradation, leading to modulation of root-to-shoot growth ratio and lateral root development [67, 68]. Induction of CKX genes in rice roots might be involved in regulating cytokinin status in response to N and Pi starvation.
Identification and profiling of putative lncRNAs under abiotic stress conditions, and characterization of the 3′-end of transcripts using 2P-Seq
After RABT assembly and a stepwise filtering, we identified a total of 2588 putative lncRNA-encoding loci (Fig. 1b, Additional file 2: Table S4). During the CPC value calculation step, we noticed that the rice protein sequences in the SwissProt non-redundant (NR) database included a number of short peptide sequences that were non-conserved and were not experimentally validated. This removed ~ 20% of the assembled transcripts, and we thus re-evaluated the CPC values of the filtered-out transcripts with SwissProt NR database that Oryza protein sequences were excluded, which resulted in the inclusion of 594 additional loci in our lncRNA dataset. To distinguish transcripts filtered by the different CPC evaluation processes, the re-evaluated lncRNAs were denoted as relatively high-coding (RHC) potential and relatively low-coding (RLC) potential lncRNAs. Among the lncRNA transcripts, 1255 and 1333 were expressed from intergenic regions and cis-NATs, respectively (Fig. 1b), and mostly comprised of single-exon transcripts (~ 74%). A large number of loci (918) were up- or down-regulated in response to N starvation in rice shoots and roots (Fig. 4a). Similar to the aforementioned profiling results for N metabolism and transport genes, many putative lncRNA loci showed more significant expression changes in N-starved roots than shoots (Fig. 4a), indicating that their responses to N starvation were tissue-dependent.
RNA-Seq data derived from rice samples subjected to Pi starvation and abiotic stress [60, 69] were analyzed to investigate stress-related characteristics of putative lncRNAs. Many of the putative lncRNAs that were unresponsive to N starvation were responsive to other stress conditions (Fig. 4b). A small number of putative lncRNAs were concordantly regulated in response to N and Pi starvation in 7 day rice samples (Fig. 4c). Analysis of N starvation and salt-, cold-, and drought-responsive transcription also showed a similar tendency (Fig. 4d–f). Some lncRNAs were specifically induced by cold stress, but the majority of these were down-regulated by salt and drought stress conditions (Fig. 4b). One of these cold stress-responsive putative lncRNAs, Chr03G0008, was expressed only in cold-treated seedlings (Additional file 1: Figure S3A). Chr03G0008 also showed N starvation- responsive induced expression patterns in N-starved roots, suggesting that Chr03G0008 responded to multiple stress conditions.
Expression of a well-known lncRNA, IPS1, changed in response to both N and Pi starvation (Table 1), suggesting that IPS1 might play a bridging role between the pathways regulating the two macronutrients. Similar putative lncRNAs that were co-regulated by both N and Pi deficiency were identified in our dataset. Expression of one such RHC-class lncRNA, Chr04G0017, was induced following both N and Pi starvation, especially in rice roots (Fig. 4g and h). Chr04G0017 comprised two exons with one distinct 2P-Seq signal at the 3′-end region, indicating strong polyadenylation. In the dataset analyzed by Secco et al., Chr04G0017 was ~ 10-fold induced after 21 days of Pi starvation . Polysomal fractionation analysis by Secco et al. suggested that Chr04G0017  was actively translated, which was consistent with the classification of Chr04G0017 as a RHC-class lncRNA in our filtering pipeline. Circumstantial and experimental evidence suggests that Chr04G0017 acts via its encoded peptide (116 a.a., 13.64 kDa) in nutrient-starved conditions. In a final example, expression of alternatively polyadenylated Chr07G0166 was induced after 7 days of both N and Pi starvation (Additional file 1: Figure S3B and S4E).
Although previous studies examined the expression of annotated genes and identified novel transcripts in rice, no studies have as yet defined or characterized the 3′-end of transcripts. Here, we used 2P-Seq analysis redefine models of various transcripts, including annotated genes, putative lncRNA loci, and primary miRNA (pri-miRNA) transcripts, with the aim to provide accurate models for future research. The 2P-Seq dataset examined sequences upstream of the poly A tail at the 3′ ends of transcripts. The results of 2P-Seq analysis revealed 647 (25%) putative lncRNAs with significant 2P-Seq signals on transcript models near the 3′ ends. Detailed scrutiny of N starvation-responsive putative lncRNAs revealed that a small number of these lncRNAs exhibited 3′ end-related characteristics, some of these were alternatively polyadenylated. One putative lncRNA, Chr04G0169, which was expressed exclusively in roots, showed two distinct groups of 2P-Seq signals at the 3′-end region: these signals were further validated by 3′ RACE (Fig. 5a). Another putative lncRNA, Chr07G0166, was present as three isoforms with different 3′ ends. Distances between 2P-Seq signals were ~ 1 kb, with similar distances between the polyadenylation sites of the short and long isoforms (Additional file 1: Figure S3B).
Integrated analysis of strand-specific RNA-Seq and 2P-Seq datasets also allowed the revision of transcript models, including those of previously annotated genes and pri-miRNAs. The expression of AMT1.1 and cis-NATAMT1.1 was detected by both strand-specific RNA-Seq and 2P-Seq (Fig. 5b). Read distribution and 5′ RACE identified 5′ extended regions of cis-NATAMT1.1 (Fig. 5b). The polyadenylation sites of AMT1.2 and cis-NATAMT1.2 transcripts were also detected by 2P-Seq, which showed longer 3′-ends of these transcripts than the IRGSP gene models (Additional file 1: Figure. S7A). In the case of pri-miR3979, two significant 3′ polyadenylation sites were observed for the single-exon transcript model of pri-miR3979 (Fig. 5c).
Analysis of N starvation-responsive miRNAs and their targets
The most substantial changes in miRNA expression were observed in 7 day N-starved samples (Fig. 6a and b). Shoots and roots of 7 day N-starved rice exhibited differential miRNA expression patterns. The majority of differentially expressed miRNAs were up-regulated in shoots of 7 day N-starved rice (Fig. 6a), and only five miRNAs were down-regulated. Almost three-fourths of differentially expressed miRNAs were down-regulated in roots of 7 day N-starved rice (Fig. 6b). Of these, members of the osa-miR169 family were significantly down-regulated in response to N starvation (Fig. 6c and d). Most of the osa-miR169 family exhibited increasing down-regulation over time in N-starved roots (Fig. 6d). Target prediction for osa-miR169 indicated that 8 of 11 rice genes encoding the nuclear transcription factor Y subunit alpha (NF-YA) were predicted to be targeted (Additional file 1: Figure S5A), and their expression levels were up-regulated in response to N starvation in rice roots (Fig. 6d and e). To confirm osa-miR169-mediated post-transcriptional regulation of NF-YAs, 5′ RACE was performed against all predicted NF-YA target genes. The 5′ RACE results showed that 7 of 8 NF-YAs were regulated by miR169-mediated cleavage (Additional file 1: Figure. S5C–I). These results suggested that up-regulation of NF-YAs in N-starved rice roots was the result of osa-miR169 de-repression in N starvation conditions.
We also observed the down-regulation of the osa-miR395 family and osa-miR399d-3p (Fig. 6b and c), which are responsible for modulating sulfur and Pi homeostasis, respectively, in response to changes in the sulfur/Pi supply status. Down-regulation of these miRNAs was consistent with previous reports from Arabidopsis and maize [35, 36, 38], and suggested that N starvation impacted the homeostasis of other nutrients regulated by nutrient-responsive miRNAs. Finally, root-specific expression of osa-miR444a.4-3p, one of the most abundant root miRNAs in our dataset, was also down-regulated in response to N starvation (Fig. 6b and c).
To examine the effect of changes in miRNA expression on target mRNAs, Degradome sequencing data from 7 day N-starved rice roots were analyzed using the CleaveLand pipeline. In total, 91 rice transcripts were identified as targets of 40 miRNA sequences (Additional file 2: Table S7). Several differentially expressed miRNAs were found to regulate their own target transcripts in 7 day N-starved rice roots (Table 2). One target gene identified by Degradome data analysis, OsMADS25, was shown to be targeted by osa-miR444a.4-3p (Fig. 7a). Osa-miR444 family sequences were derived from the NATs of four ANR1-like MADS-box genes (OsMADS23, OsMADS27a, OsMADS27b, OsMADS57) and were shown to target all four MADS-box genes [40, 70, 71]. However, osa-miR444a.4-3p-mediated regulation of OsMADS25 was not confirmed in a previous Degradome analysis of rice . In this study, Degradome analysis showed that only osa-miR444a.4-3p cleaved OsMADS25 (Fig. 7a), and osa-miR444a.4-3p only seemed to target OsMADS25. Expression levels of osa-miR444a.4-3p and OsMADS25 were negatively correlated in N-starved rice, suggesting that N starvation-mediated down-regulation of osa-miR444a.4-3p affected up-regulation of OsMADS25. In addition, osa-miR168a-5p.1-5p was involved in mediating the OsPTR29 cleavage pattern (Fig. 7b). OsPTR29, also known as OsNPF2.4, is a pH-dependent, low-affinity NRT that plays a role in nitrate acquisition and long-distance nitrate transport .
Some of the miRNAs that did not exhibit a significant change in the small RNA-Seq dataset were shown to regulate target genes through cleavage mechanisms (Additional file 1: Figure S6). The osa-miR1425-5p, which targeted genes encoding pentatricopeptide repeat (PPR) proteins (Additional file 1: Figure. S6C), was abundantly expressed in roots. Both osa-miR1425-5p and its primary transcript were slightly down-regulated in 7 day N-starved roots (Table 2). PPR is an organellar RNA-binding protein that participates in RNA editing in chloroplasts and mitochondria [73, 74]. Previously, functional analysis of mitochondrial PPRs in Arabidopsis has demonstrated that these proteins are involved in regulating root growth, energy metabolism, and abiotic stress responses [75,76,77,78]. These results suggest that osa-miR1425-5p has a potential role in modulating root architecture and energy metabolism in roots of rice plants.
Expression of cis-NATs and small RNAs at the genomic regions of AMT1 genes
Previous reports suggested that some cis-NATs could regulate target genes on the antisense genomic region corresponding to their protein-coding genes [79,80,81]. We examined whether genes involved in N metabolism and transport transcribed potential regulatory cis-NATs. Eight of nine AMT genes had IRGSP-annotated cis-NATs. In our strand-specific RNA-Seq dataset, most of these cis-NATAMTs were expressed at very low levels (< 0.1 FPKM) or were not detected. However, two cis-NATAMTs, cis-NATAMT1.1 (Os04g0509500) and cis-NATAMT1.2 (Os02g0620533), showed relatively high expression (> 1 FPKM), and their 3′-end poly A tails were supported by 2P-Seq (Fig. 5a and Additional file 1: Figure. S7A). These two polyadenylated cis-NATAMT1 transcripts and the corresponding AMT1 genes showed different expression patterns: cis-NATAMT1.2 and AMT1.2 were exclusively expressed in rice roots, and both showed increased expression levels in response to N starvation, whereas AMT1.1 and cis-NATAMT1.1 showed tissue-dependent discordant expression patterns. AMT1.1 was expressed more abundantly in roots, whereas cis-NATAMT1.1 was rarely expressed in roots and was expressed more abundantly in shoots (Fig. 5a).
During observation of RNA-Seq read distribution on the genomic region of AMT1.1 and AMT1.2, we observed that a few small RNA reads aligned to the genomic region of two AMT1 family genes and had expression patterns that generally correlated with the expression patterns of AMT genes and/or their cis-NATs described above (Fig. 5b and Additional file 1: Figure. S7A). First, we assumed that these small RNAs came from double-stranded RNAs (dsRNAs) consisting of cis-NAT pairs, since the genomic region of these genes overlapped extensively. However, expression patterns of AMT1 genes, their cis-NATs, and small RNAs suggest that these small RNAs are derived from single-stranded RNA of AMT1s and their cis -NATs. In the case of small RNAs in the AMT1.1 genomic region, expression levels of (+)-strand small RNAs correlated with expression levels of AMT1.1 rather than those of cis-NATAMT1.1 (Fig. 5a). Similarly, (+)-strand small RNAs sense to AMT1.2 were induced in concordance with AMT1.2 transcript levels in N-starved rice roots (Additional file 1: Figure. S7A). However, although expression patterns of cis-NATAMT1.1 and corresponding (−)-strand small RNAs were concordant in a tissue-dependent manner, cis-NATAMT1.1 and the corresponding (−)-strand small RNAs were not concordantly regulated in response to N starvation. Moreover, in silico RNA structure analysis of the homologous region of AMT1.1 and AMT1.2 sequences showed potent stem-loop structures that were sufficient to produce small RNAs (Additional file 1: Figure. S7B and S7C). Together, these results suggested that the small RNAs aligned to the genomic regions of AMT1.1 and AMT1.2 were derived from sense transcripts.
In this study, we profiled and characterized multiple aspects of the rice transcriptome in N-starved rice plants in a spatio-temporal manner, via four different transcriptomic approaches. Expression patterns of genes involved in N metabolism and transport were characterized, and a set of putative lncRNA loci was identified. In addition, entities and characteristics of 3′-ends of rice transcripts were identified and analyzed by 2P-Seq. We also profiled annotated rice genes and putative lncRNAs in Pi starvation-treated rice samples to gain insights into the relationships between gene regulatory pathways involved in the metabolism of N and other nutrients. Comparative analysis of RNA-Seq datasets from N- and Pi-starved rice identified co-regulated genes and putative lncRNAs, and one of those lncRNAs (Chr04G0017) was identified to be regulated by both N and Pi starvation in a root-specific manner. Additionally, analysis of putative lncRNAs in Pi-starved samples as well as salt-, cold-, and drought-stressed rice samples suggested that a set of putative lncRNAs was regulated in a stress-specific manner. Lastly, we profiled miRNA expression in N-starved rice, and identified their target genes in 7 day N-starved roots. These unprecedented large-scale, integrative analyses of the N-starved rice transcriptome will provide valuable resources for researchers developing novel gene resources for NUE-improved crop development.
We attempted to predict the potential functions of the N and Pi starvation-responsive putative RHC-class lncRNA Chr04G0017. Previous reports have indicated that short peptide-coding RNAs function as both RNAs and peptides [82,83,84]. To assess the potential functions of Chr04G0017 as a 116 amino acid peptide, we performed protein BLAST and peptide structure prediction. No homologous plant proteins were identified other than those in Oryza subspecies via protein BLAST analysis; however, weak similarities were observed with some cytoskeleton-like proteins, such as Type-II keratin, of other species. Peptide structure prediction using PEP-FOLD3  showed the presence of plausible alpha-helical structure models (coiled-coil) between amino acid positions 31 and 100 (Additional file 1: Figure. S10), supporting the existence of cytoskeletal protein-like structures. Because protein BLAST results of peptide sequences from Nipponbare cultivar showed high similarities with those of other Oryza subspecies, we analyzed the sequence alignment of Chr04G0017 genomic region between Oryza subspecies. Results showed near-perfect similarities in both intronic and exonic regions of Chr04G0017, including ORF sequences (blue-colored in Additional file 1: Figure. S9), suggesting that Chr04G0017 is well-conserved among Oryza subspecies. Although sufficient clues could not be identified for understanding the function of Chr04G0017, sequence alignment indicated that this putative RHC-class lncRNA may play a conserved role in nutrient-deficient conditions, such as during N and Pi starvation. Further investigation is needed to understand the potential functions of both peptide and RNA forms of Chr04G0017.
Information on strand-specificity from strand-specific RNA-Seq and 2P-Seq data was helpful in defining transcriptional direction of putative lncRNAs, for identifying exact transcript models, and for measuring the expression of each transcript of cis-NAT pairs. This allowed us to confirm the expression of two cis-NAT pairs, AMT1.1–cis-NATAMT1.1 and AMT1.2–cis-NATAMT1.2. Moreover, we observed small RNA reads covering the overlapped region of cis-NATAMT1 pairs, whereas small RNAs were not significantly detected in the overlapped region of AMT2–cis-NATAMT2 or AMT3–cis-NATAMT3 pairs (Additional file 1: Figure. S8). Our data suggested that these small RNAs were generated from corresponding sense transcripts rather than dsRNAs formed by complementary binding of cis-NAT pairs. Previous studies of the molecular functions of cis-NATs reported RNA-RNA interaction-based regulatory mechanisms through complementary sequence binding between cis-NAT pairs, leading to dsRNA-dependent RNA interference, targeted RNA protection by RNA masking and translational enhancements [21, 80]. Because these two cis-NATAMT1 pairs overlapped significantly, and the nucleotide sequences of the three rice AMT1 genes were highly similar to one other, the possibilities remain of in cis or in trans interactions between cis-NATAMT1 and AMT1 genes. Similar cis-NATAMT1 transcripts also appeared to be expressed in Arabidopsis (At4G13505 for AtAMT1.1 cis-NAT pair), Hordeum vulgare (PUT-169a-Hordeum_vulgare-46,748), and maize (GRMZM2G474905 for cis-NATZmAMT1, GRMZM2G332891 for cis-NATZmAMT2), suggesting that cis-NATAMT expression and function might be conserved among plants. Further research is needed to determine the molecular functions of these cis-NATs in N-starved rice.
By profiling miRNA expression and analyzing Degradome sequencing datasets, we discovered a novel target of osa-miR444a.4-3p, OsMADS25. OsMADS25 is one out of five ANR1-like rice MADS-box proteins [40, 86], all of which, except OsMADS25, are regulated by members of the osa-miR444 family derived from their NATs [70, 72]. Evidence from recent studies demonstrated the importance of miR444-mediated regulation of ANR1-like MADS-box genes in modulating N homeostasis in rice . Furthermore, previous research showed that overexpression of OsMADS25 in rice promoted primary/lateral root growth, increased shoot fresh weight, and increased nitrate accumulation in the presence of nitrate . However, the relationship between the osa-miR444 family and OsMADS25 was not noted previously as the genomic region of OsMADS25 did not transcribe miR444-encoding cis-NATs. Our discoveries of N starvation-responsive characteristics of osa-miR444a.4-3p and its in trans regulation of OsMADS25 expanded the regulatory pathways of the osa-miR444 family and reinforced the significance of osa-miR444 in modulating N homeostasis in rice. As a significant proportion of osa-miR444 family expression in roots (~ 47%) is attributable to osa-miR444a.4, and osa-miR444a.4 is only minimally expressed in shoots (< 0.01%), the roles of osa-miR444a.4 in N acquisition and homeostasis in roots merit further investigation.
In this study, multiple transcriptomic investigations on N-starved rice plants improved our understanding of the transcriptomic responses to N starvation by providing detailed and intricate information on changes in the rice transcriptome. Strand-specific RNA-Seq datasets not only provided information on the responses of genes involved in N metabolism- and transport-involved genes, but also newly identified 2588 novel putative lncRNA encoding genomic loci. Information on the transcript models of these lncRNAs combined with the results of 2P-Seq analysis showed examples of alternatively polyadenylated isoforms of N starvation-responsive lncRNAs, providing precise information on transcript models of these lncRNAs. Analyses of lncRNAs using previously published RNA-Seq datasets revealed lncRNAs that not only responded to N starvation but also showed differential expression in response to various kinds of abiotic stresses. We also reported N-responsive characteristics of the root-specific osa-miR444a.4-3p and its novel target gene, OsMADS25. Overall, these large-scale datasets provide valuable information for the generation of new rice cultivars with higher NUE or greater resistance to N starvation in future breeding programs.
Plant material and growth conditions
Rice (Oryza sativa cv. Nipponbare) seeds were germinated in MS media for 4 d, and then transferred to tap water for 3d before being transferred into the hydroponic solution. Rice seedlings were grown in the modified Yoshida solution for 10 d . The solution was renewed every 3 d. For preparing N-starved rice samples, seedlings were transferred to solution lacking N (0 mM of NH4NO3), and roots and shoots were harvested separately at 1, 3, 5, 7 d of N starvation (Fig. 1a). All samples were harvested at the same time of the day (i.e., 2 h after the onset of subjective day) to minimize potential circadian effects.
Total RNA isolation and library preparation for high-throughput sequencing
Root and shoot samples were ground in liquid N, separately. Total RNA was extracted from the samples using TRIzol Reagent (Invitrogen), according to the manufacturer’s instructions, and the integrity and quality of RNA samples was analyzed. Strand-specific RNA-Seq libraries were constructed using 5 μg total RNA, according to the modified protocol previously described . For poly A-primed sequencing, 2P-Seq library was constructed using 60 μg total RNA extracted from rice shoots and roots under normal and N-starvation conditions . Using Illumina HiSeq 2500, strand-specific RNA-Seq libraries were analyzed with 101-bp paired-end sequencing, and 2P-Seq libraries were analyzed with 101-bp single-end sequencing. The construction and sequencing of small RNA-Seq libraries were performed according to protocols described previously . The Degradome library of 7 d N-starved rice roots was constructed as described previously . Degradome library were analyzed using Illumina HiSeq 2500 51-bp single-end sequencing.
Analysis of RNA-Seq data using bioinformatics
Strand-specific RNA-Seq reads were aligned to the rice IRGSP-1.0 genome  using TopHat in the Cufflinks package , and RABT assembly was performed using Cufflinks package using the rice gene model annotation from RAP-DB (http:// http://rapdb.dna.affrc.go.jp/). Identification of putative lncRNAs was performed with analysis pipeline described in Fig. 1b. Transcript abundance of genes and putative lncRNAs was then estimated as fragments per kilobase of exon (FPKM). RNA-Seq datasets of Pi-starved rice plants were obtained from Gene Expression Omnibus (GEO), with the accession numbers SRR1005258, SRR1005300, SRR1005306, SRR1005318, SRR1005321, SRR1005363, SRR1005369, and SRR1005381 (https://www.ncbi.nlm.nih.gov/bioproject/PRJNA215013) , and RNA-Seq datasets of abiotic stress-treated rice plants were obtained with the accession numbers of ERR037679, ERR037681, ERR037683 and ERR037687 (https://trace.ncbi.nlm.nih.gov/Traces/study/?acc=ERP000760) .
For 2P-Seqdata analysis, sequence reads ranging from 20 to 95-nt were aligned to the rice genome using TopHat for further analysis. The 5′-ends of sequence reads were used to represent the genomic position of aligned 2P-Seq data, and peak signals were calculated. Peak signals up to 1000-bp distant from the 3′ end of assembled transcript models were selected as candidate polyadenylation sites.
Raw sequence reads generated from small RNA-Seq and Degradome sequencing were cleaned by removing adapter sequences. Rice miRNA prediction was performed as described previously . Expression levels of miRNAs were normalized and estimated as reads per 40 million of sequence reads (RP40M). Degradome data analysis was performed to identify miRNA target genes. High-quality Degradome sequence reads were obtained from rawdata by filtering out poor quality reads and removing adapter sequences using FASTX toolkit (http://hannonlab.cshl.edu/fastx_toolkit/). Reads corresponding to structural non-coding RNAs and repeat sequences in the rice genome database  were filtered-out. The processed Degradome reads were analyzed using CleaveLand4 pipeline  with rice transcript models and miRNA sequences.
Gene ontology (GO) analysis
GO analysis was performed using AgriGO web server with default options  (http://bioinfo.cau.edu.cn/agriGO/). MSU (v7.0) gene ID corresponding to differentially expressed (up- or down-regulated by more than 2-fold) RAP-DB gene ID was used in GO analysis.
Identification of 5′ and 3′ends of transcripts and miRNA cleavage sites, and quantitative real-time-PCR (qRT-PCR)
To identify the 5′ and 3′ ends of transcripts, 5′ RNA ligase-mediated rapid amplification of cDNA ends (RLM-RACE) and 3′ RACE were performed respectively, using 1 μg total RNA and GeneRacer Kit (Invitrogen), following the manufacturer’s instructions. To identify miRNA cleavage sites in target mRNAs, 5′ RLM-RACE without serial dephosphorylation-decapping treatment was performed . For quantitative RT-PCR, cDNA was synthesized from 1 μg of total RNA using oligo(dT)18 primers and random hexamers. Results were normalized against UBQ5 (Os01g0328400) for shoots and eEF-1α (Os03g0177500) for roots. To measure the expression level of Chr04G0017, results from both shoots and roots were normalized relative to eEF-1α.
Small RNA northern blot analysis
To examine the expression level of miRNAs in N-starved rice roots and shoots, 10 μg of total RNA from roots and shoots of N-starved rice samples were resolved by 15% urea-PAGE, and transferred to Hybond-NX membrane (GE Healthcare), followed by UV-crosslinking. 10 nmole of each probe corresponding to rice miRNA sequences was radiolabeled by standard 5′ end-labeling reaction using T4 polynucleotide kinase (Takara). Probe sequences used in this study were shown in Additional file 2: Table S3.
Poly A-primed sequencing
Amino acid permease
Amino acid transporter
International Rice Genome Sequencing Project
Long non-coding RNA
Natural antisense transcript
Nuclear transcription factor Y subunit alpha
Nitrogen Use Efficiency
quantitative real-time PCR
Reference annotation-based transcript
Rice Annotation Project Database
Raun WR, Johnson GV. Improving nitrogen use efficiency for cereal production. Agron J. 1999;91(3):357–63.
Xu GH, Fan XR, Miller AJ. Plant nitrogen assimilation and use efficiency. Annu Rev Plant Biol, Vol 63. 2012;63:153–82.
Chen JG, Zhang Y, Tan YW, Zhang M, Zhu LL, Xu GH, Fan XR. Agronomic nitrogen-use efficiency of rice can be increased by driving OsNRT2.1 expression with the OsNAR2.1 promoter. Plant Biotechnol J. 2016;14(8):1705–15.
Fang ZM, Xia KF, Yang X, Grotemeyer MS, Meier S, Rentsch D, Xu XL, Zhang MY. Altered expression of the PTR/NRT1 homologue OsPTR9 affects nitrogen utilization efficiency, growth and grain yield in rice. Plant Biotechnol J. 2013;11(4):446–58.
Feng HM, Yan M, Fan XR, Li BZ, Shen QR, Miller AJ, Xu GH. Spatial expression and regulation of rice high-affinity nitrate transporters by nitrogen and carbon status. J Exp Bot. 2011;62(7):2319–32.
Xia X, Fan X, Wei J, Feng H, Qu H, Xie D, Miller AJ, Xu G. Rice nitrate transporter OsNPF2.4 functions in low-affinity acquisition and long-distance transport. J Exp Bot. 2015;66(1):317–31.
Bao AL, Liang ZJ, Zhao ZQ, Cai HM. Overexpressing of OsAMT1-3, a high affinity ammonium transporter gene, modifies Rice growth and carbon-nitrogen metabolic status. Int J Mol Sci. 2015;16(5):9037–63.
Hoque MS, Masle J, Udvardi MK, Ryan PR, Upadhyaya NM. Over-expression of the rice OsAMT1-1 gene increases ammonium uptake and content, but impairs growth and development of plants under high ammonium nutrition. Funct Plant Biol. 2006;33(2):153–63.
Kumar A, Kaiser BN, Siddiqi MY, Glass ADM. Functional characterisation of OsAMT1.1 overexpression lines of rice, Oryza sativa. Funct Plant Biol. 2006;33(4):339–46.
Ranathunge K, El-kereamy A, Gidda S, Bi YM, Rothstein SJ. AMT1;1 transgenic rice plants with enhanced NH4 permeability show superior growth and higher yield under optimal and suboptimal NH4 conditions. J Exp Bot. 2014;65(4):965–79.
Beatty PH, Carroll RT, Shrawat AK, Guevara D, Good AG. Physiological analysis of nitrogen-efficient rice overexpressing alanine aminotransferase under different N regimes. Botany. 2013;91(12):866–83.
Shrawat AK, Carroll RT, DePauw M, Taylor GJ, Good AG. Genetic engineering of improved nitrogen use efficiency in rice by the tissue-specific expression of alanine aminotransferase. Plant Biotechnol J. 2008;6(7):722–32.
Bi YM, Kant S, Clark J, Gidda S, Ming F, Xu JY, Rochon A, Shelp BJ, Hao LX, Zhao R, et al. Increased nitrogen-use efficiency in transgenic rice plants over-expressing a nitrogen-responsive early nodulin gene identified from rice expression profiling. Plant Cell Environ. 2009;32(12):1749–60.
Cai HM, Lu YG, Xie WB, Zhu T, Lian XM. Transcriptome response to nitrogen starvation in rice. J Biosci. 2012;37(4):731–47.
Lian XM, Wang SP, Zhang JW, Feng Q, Zhang LD, Fan DL, Li XH, Yuan DJ, Han B, Zhang QF. Expression profiles of 10,422 genes at early stage of low nitrogen stress in rice assayed using a cDNA microarray. Plant Mol Biol. 2006;60(5):617–31.
Yang W, Yoon J, Choi H, Fan Y, Chen R, An G. Transcriptome analysis of nitrogen-starvation-responsive genes in rice. BMC Plant Biol. 2015;15:31.
Fischer JJ, Beatty PH, Good AG, Muench DG. Manipulation of microRNA expression to improve nitrogen use efficiency. Plant Sci. 2013;210:70–81.
Nguyen GN, Rothstein SJ, Spangenberg G, Kant S. Role of microRNAs involved in plant response to nitrogen and phosphorous limiting conditions. Front Plant Sci. 2015;6:629.
Ding J, Lu Q, Ouyang Y, Mao H, Zhang P, Yao J, Xu C, Li X, Xiao J, Zhang Q. A long noncoding RNA regulates photoperiod-sensitive male sterility, an essential component of hybrid rice. Proc Natl Acad Sci U S A. 2012;109(7):2654–9.
Franco-Zorrilla JM, Valli A, Todesco M, Mateos I, Puga MI, Rubio-Somoza I, Leyva A, Weigel D, Garcia JA, Paz-Ares J. Target mimicry provides a new mechanism for regulation of microRNA activity. Nat Genet. 2007;39(8):1033–7.
Jabnoune M, Secco D, Lecampion C, Robaglia C, Shu QY, Poirier Y. A Rice cis-natural antisense RNA acts as a translational enhancer for its cognate mRNA and contributes to phosphate homeostasis and plant fitness. Plant Cell. 2013;25(10):4166–82.
Liang G, Yang F, Yu D. MicroRNA395 mediates regulation of sulfate accumulation and allocation in Arabidopsis thaliana. Plant J. 2010;62(6):1046–57.
Zhao B, Ge L, Liang R, Li W, Ruan K, Lin H, Jin Y. Members of miR-169 family are induced by high salinity and transiently inhibit the NF-YA transcription factor. BMC Mol Biol. 2009;10:29.
Zhao M, Ding H, Zhu JK, Zhang F, Li WX. Involvement of miR169 in the nitrogen-starvation responses in Arabidopsis. The New phytologist. 2011;190(4):906–15.
Chen M, Wang CL, Bao H, Chen H, Wang YW. Genome-wide identification and characterization of novel lncRNAs in Populus under nitrogen deficiency. Mol Gen Genomics. 2016;291(4):1663–80.
Lv YD, Liang ZK, Ge M, Qi WC, Zhang TF, Lin F, Peng ZH, Zhao H. Genome-wide identification and functional prediction of nitrogen-responsive intergenic and intronic long non-coding RNAs in maize (Zea mays L.). BMC Genomics. 2016;17:350.
Shuai P, Liang D, Tang S, Zhang ZJ, Ye CY, Su YY, Xia XL, Yin WL. Genome-wide identification and functional prediction of novel and drought-responsive lincRNAs in Populus trichocarpa. J Exp Bot. 2014;65(17):4975–83.
Zhang H, Hu WG, Hao JL, Lv SK, Wang CY, Tong W, Wang YJ, Wang YZ, Liu XL, Ji WQ. Genome-wide identification and functional prediction of novel and fungi-responsive lincRNAs in Triticum aestivum. BMC Genomics. 2016;17:238.
Zhang YC, Liao JY, Li ZY, Yu Y, Zhang JP, Li QF, Qu LH, Shu WS, Chen YQ. Genome-wide screening and functional analysis identify a large number of long noncoding RNAs involved in the sexual reproduction of rice. Genome Biol. 2014;15(12):512.
Zhang YC, Yu Y, Wang CY, Li ZY, Liu Q, Xu J, Liao JY, Wang XJ, Qu LH, Chen F, et al. Overexpression of microRNA OsmiR397 improves rice yield by increasing grain size and promoting panicle branching. Nat Biotechnol. 2013;31(9):848.
Shin SY, Shin C. Regulatory non-coding RNAs in plants: potential gene resources for the improvement of agricultural traits. Plant Biotechnol Rep. 2016;10(2):35–47.
Zheng LL, Qu LH. Application of microRNA gene resources in the improvement of agronomic traits in rice. Plant Biotechnol J. 2015;13(3):329–36.
Zhou M, Luo H. MicroRNA-mediated gene regulation: potential applications for plant genetic engineering. Plant Mol Biol. 2013;83:59–75.
Jeong DH, Park S, Zhai J, Gurazada SG, De Paoli E, Meyers BC, Green PJ. Massive analysis of rice small RNAs: mechanistic implications of regulated microRNAs and variants for differential target RNA cleavage. Plant Cell. 2011;23(12):4185–207.
Liang G, Ai Q, Yu D. Uncovering miRNAs involved in crosstalk between nutrient deficiencies in Arabidopsis. Sci Rep. 2015;5:11813.
Liang G, He H, Yu D. Identification of nitrogen starvation-responsive microRNAs in Arabidopsis thaliana. PLoS One. 2012;7(11):e48951.
Wang YJ, Zhang CJ, Hao QN, Sha AH, Zhou R, Zhou XA, Yuan LP. Elucidation of miRNAs-mediated responses to low nitrogen stress by deep sequencing of two soybean genotypes. PLoS One. 2013;8(7):e67423.
Zhao M, Tai HH, Sun SZ, Zhang FS, Xu YB, Li WX. Cloning and characterization of maize miRNAs involved in responses to nitrogen deficiency. PLoS One. 2012;7(1):e29669.
Lin WY, Huang TK, Chiou TJ. NITROGEN LIMITATION ADAPTATION, a target of MicroRNA827, mediates degradation of plasma membrane-localized phosphate transporters to maintain phosphate homeostasis in Arabidopsis. Plant Cell. 2013;25(10):4061–74.
Yan Y, Wang H, Hamera S, Chen X, Fang R. miR444a has multiple functions in the rice nitrate-signaling pathway. Plant J. 2014;78(1):44–55.
Hsieh LC, Lin SI, Shih ACC, Chen JW, Lin WY, Tseng CY, Li WH, Chiou TJ. Uncovering small RNA-mediated responses to phosphate deficiency in Arabidopsis by deep sequencing. Plant Physiol. 2009;151(4):2120–32.
Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and cufflinks. Nat Protoc. 2012;7(3):562–78.
Araki R, Hasegawa H. Expression of rice (Oryza sativa L.) genes involved in high-affinity nitrate transport during the period of nitrate induction. Breed Sci. 2006;56(3):295–302.
Ishiyama K, Inoue E, Tabuchi M, Yamaya T, Takahashi H. Biochemical background and compartmentalized functions of cytosolic glutamine synthetase for active ammonium assimilation in rice roots. Plant Cell Physiol. 2004;45(11):1640–7.
Li SM, Li BZ, Shi WM. Expression patterns of nine ammonium transporters in Rice in response to N status. Pedosphere. 2012;22(6):860–9.
Sonoda Y, Ikeda A, Saiki S, von Wiren N, Yamaya T, Yamaguchi J. Distinct expression and function of three ammonium transporter genes (OsAMT1;1-1;3) in rice. Plant Cell Physiol. 2003;44(7):726–34.
Sonoda Y, Ikeda A, Saiki S, Yamaya T, Yamaguchi J. Feedback regulation of the ammonium transporter gene family AMT1 by glutamine in rice. Plant Cell Physiol. 2003;44(12):1396–402.
De Angeli A, Monachello D, Ephritikhine G, Frachisse JM, Thomine S, Gambale F, Barbier-Brygoo H. The nitrate/proton antiporter AtCLCa mediates nitrate accumulation in plant vacuoles. Nature. 2006;442(7105):939–42.
Granstedt RC, Huffaker RC. Identification of the leaf vacuole as a major nitrate storage pool. Plant Physiol. 1982;70(2):410–3.
Funayama K, Kojima S, Tabuchi-Kobayashi M, Sawa Y, Nakayama Y, Hayakawa T, Yamaya T. Cytosolic glutamine synthetase1;2 is responsible for the primary assimilation of ammonium in rice roots. Plant Cell Physiol. 2013;54(6):934–43.
Kawachi T, Sueyoshi K, Nakajima A, Yamagata H, Sugimoto T, Oji Y. Expression of asparagine synthetase in rice (Oryza sativa) roots in response to nitrogen. Physiol Plantarum. 2002;114(1):41–6.
Ohashi M, Ishiyama K, Kojima S, Konishi N, Nakano K, Kanno K, Hayakawa T, Yamaya T. Asparagine synthetase1, but not asparagine synthetase2, is responsible for the biosynthesis of asparagine following the supply of ammonium to rice roots. Plant & cell physiol. 2015;56(4):769–78.
Yamaya T, Kusano M. Evidence supporting distinct functions of three cytosolic glutamine synthetases and two NADH-glutamate synthases in rice. J Exp Bot. 2014;65(19):5519–25.
Zhao H, Ma H, Yu L, Wang X, Zhao J. Genome-wide survey and expression analysis of amino acid transporter gene family in rice (Oryza sativa L.). PLoS One. 2012;7(11):e49210.
Zhao X, Huang J, Yu H, Wang L, Xie W. Genomic survey, characterization and expression profile analysis of the peptide transporter family in rice (Oryza sativa L.). BMC Plant Biol. 2010;10:92.
Cai HM, Xie WB, Lian XM. Comparative analysis of differentially expressed genes in Rice under nitrogen and phosphorus starvation stress conditions. Plant Mol Biol Rep. 2013;31(1):160–73.
Cerutti T, Delatorre CA. Nitrogen and phosphorus interaction and cytokinin: responses of the primary root of Arabidopsis thaliana and the pdr1 mutant. Plant Sci. 2013;198:91–7.
Wu P, Ma LG, Hou XL, Wang MY, Wu YR, Liu FY, Deng XW. Phosphate starvation triggers distinct alterations of genome expression in Arabidopsis roots and leaves. Plant Physiol. 2003;132(3):1260–71.
Ziadi N, Belanger G, Cambouris AN, Tremblay N, Nolin MC, Claessens A. Relationship between phosphorus and nitrogen concentrations in spring wheat. Agron J. 2008;100(1):80–6.
Secco D, Jabnoune M, Walker H, Shou H, Wu P, Poirier Y, Whelan J. Spatio-temporal transcript profiling of rice roots and shoots in response to phosphate starvation and recovery. Plant Cell. 2013;25(11):4285–304.
Chiaverini N, De Ley M. Protective effect of metallothionein on oxidative stress-induced DNA damage. Free Radic Res. 2010;44(6):605–13.
Lee KO, Lee JR, Yoo JY, Jang HH, Moon JC, Jung BG, Chi YH, Park SK, Lee SS, Lim CO, et al. GSH-dependent peroxidase activity of the rice (Oryza sativa) glutaredoxin, a thioltransferase. Biochem Biophys Res Commun. 2002;296(5):1152–6.
Ruttkay-Nedecky B, Nejdl L, Gumulec J, Zitka O, Masarik M, Eckschlager T, Stiborova M, Adam V, Kizek R. The role of metallothionein in oxidative stress. Int J Mol Sci. 2013;14(3):6044–66.
Kandlbinder A, Finkemeier I, Wormuth D, Hanitzsch M, Dietz KJ. The antioxidant status of photosynthesizing leaves under nutrient deficiency: redox regulation, gene expression and antioxidant activity in Arabidopsis thaliana. Physiol Plantarum. 2004;120(1):63–73.
Hirose N, Makita N, Kojima M, Kamada-Nobusada T, Sakakibara H. Overexpression of a type-a response regulator alters rice morphology and cytokinin metabolism. Plant Cell Physiol. 2007;48(3):523–39.
Kiba T, Kudo T, Kojima M, Sakakibara H. Hormonal control of nitrogen acquisition: roles of auxin, abscisic acid, and cytokinin. J Exp Bot. 2011;62(4):1399–409.
Horgan JM, Wareing PF. Cytokinins and the growth responses of seedlings of Betula pendula Roth. And Acer pseudoplatanus L. to nitrogen and phosphorus deficiency. J Exp Bot. 1980;31(2):525–32.
Kamada-Nobusada T, Makita N, Kojima M, Sakakibara H. Nitrogen-dependent regulation of De novo Cytokinin biosynthesis in Rice: the role of glutamine metabolism as an additional signal. Plant Cell Physiol. 2013;54(11):1881–93.
Lu TT, Zhu CR, Lu GJ, Guo YL, Zhou Y, Zhang ZY, Zhao Y, Li WJ, Lu Y, Tang WH, et al. Strand-specific RNA-seq reveals widespread occurrence of novel cis-natural antisense transcripts in rice. BMC Genomics. 2012;13:721.
Lu C, Jeong DH, Kulkarni K, Pillay M, Nobuta K, German R, Thatcher SR, Maher C, Zhang L, Ware D, et al. Genome-wide analysis for discovery of rice microRNAs reveals natural antisense microRNAs (nat-miRNAs). Proc Natl Acad Sci U S A. 2008;105(12):4951–6.
Wu L, Zhang Q, Zhou H, Ni F, Wu X, Qi Y. Rice MicroRNA effector complexes and targets. Plant Cell. 2009;21(11):3421–35.
Li YF, 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(5):742–59.
Barkan A, Small I. Pentatricopeptide repeat proteins in plants. Annu Rev Plant Biol. 2014;65:415–42.
Yagi Y, Tachikawa M, Noguchi H, Satoh S, Obokata J, Nakamura T. Pentatricopeptide repeat proteins involved in plant organellar RNA editing. RNA Biol. 2013;10(9):1419–25.
Sung TY, Tseng CC, Hsieh MH. The SLO1 PPR protein is required for RNA editing at multiple sites with similar upstream sequences in Arabidopsis mitochondria. Plant J. 2010;63(3):499–511.
Xie T, Chen D, Wu J, Huang X, Wang Y, Tang K, Li J, Sun M, Peng X. Growing slowly 1 locus encodes a PLS-type PPR protein required for RNA editing and plant development in Arabidopsis. J Exp Bot. 2016;67(19):5687–98.
Yuan H, Liu D. Functional disruption of the pentatricopeptide protein SLG1 affects mitochondrial RNA editing, plant development, and responses to abiotic stresses in Arabidopsis. Plant J. 2012;70(3):432–44.
Zhu Q, Dugardeyn J, Zhang CY, Takenaka M, Kuhn K, Craddock C, Smalle J, Karampelias M, Denecke J, Peters J, et al. SLO2, a mitochondrial pentatricopeptide repeat protein affecting several RNA editing sites, is required for energy metabolism. Plant J. 2012;71(5):836–49.
Faghihi MA, Wahlestedt C. Regulatory roles of natural antisense transcripts. Nat Rev Mol Cell Bio. 2009;10(9):637–43.
Rosikiewicz W, Makalowska I. Biological functions of natural antisense transcripts. Acta Biochim Pol. 2016;63(4):665–73.
Werner A. Biological functions of natural antisense transcripts. BMC Biol. 2013;11:31.
Campalans A, Kondorosi A, Crespi M. Enod40, a short open reading frame-containing mRNA, induces cytoplasmic localization of a nuclear RNA binding protein in Medicago truncatula. Plant Cell. 2004;16(4):1047–59.
Laporte P, Satiat-Jeunemaitre B, Velasco I, Csorba T, Van de Velde W, Campalans A, Burgyan J, Arevalo-Rodriguez M, Crespi M. A novel RNA-binding peptide regulates the establishment of the Medicago truncatula-Sinorhizobium meliloti nitrogen-fixing symbiosis. Plant J. 2010;62(1):24–38.
Sousa C, Johansson C, Charon C, Manyani H, Sautter C, Kondorosi A, Crespi M. Translational and structural requirements of the early nodulin gene enod40, a short-open reading frame-containing RNA, for elicitation of a cell-specific growth response in the alfalfa root cortex. Mol Cell Biol. 2001;21(1):354–66.
Shen Y, Maupetit J, Derreumaux P, Tuffery P. Improved PEP-FOLD approach for peptide and Miniprotein structure prediction. J Chem Theory Comput. 2014;10(10):4745–58.
Yu C, Liu Y, Zhang A, Su S, Yan A, Huang L, Ali I, Liu Y, Forde BG, Gan Y. MADS-box transcription factor OsMADS25 regulates root development through affection of nitrate accumulation in rice. PLoS One. 2015;10(8):e0135196.
Yoshida S, Forno DA, Cock JH, Gomez KA. Laboratory manual for physiological studies of rice, Ed 3. Manila: International Rice Research Institute; 1976.
Zhong S, Joung JG, Zheng Y, Chen YR, Liu B, Shao Y, Xiang JZ, Fei Z, Giovannoni JJ. High-throughput illumina strand-specific RNA sequencing library preparation. Cold Spring Harb Protoc. 2011;2011(8):940–9.
Spies N, Burge CB, Bartel DP. 3 ' UTR-isoform choice has limited influence on the stability and translational efficiency of most mRNAs in mouse fibroblasts. Genome Res. 2013;23(12):2078–90.
Kim T, Park JH, Lee SG, Kim S, Kim J, Lee J, Shin C. Small RNA transcriptome of hibiscus Syriacus provides insights into the potential influence of microRNAs in flower development and terpene synthesis. Mol cells. 2017;40(8):587–97.
German MA, Luo SJ, 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(3):356–62.
Kawahara Y, de la Bastide M, Hamilton JP, Kanamori H, McCombie WR, Ouyang S, Schwartz DC, Tanaka T, Wu JZ, Zhou SG, et al. Improvement of the Oryza sativa Nipponbare reference genome using next generation sequence and optical map data. Rice. 2013;6:4.
Hwang DG, Park JH, Lim JY, Kim D, Choi Y, Kim S, Reeves G, Yeom SI, Lee JS, Park M, et al. The hot pepper (Capsicum annuum) MicroRNA transcriptome reveals novel and conserved targets: a Foundation for Understanding MicroRNA functional roles in hot pepper. PLoS One. 2013;8(5):e64238.
Brousse C, Liu Q, Beauclair L, Deremetz A, Axtell MJ, Bouché N. A non-canonical plant microRNA target site. Nucleic Acids Res. 2014;42(8):5270–9.
Du Z, Zhou X, Ling Y, Zhang ZH, Su Z. agriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010;38:W64–70.
Park JH, Shin C. MicroRNA-directed cleavage of targets: mechanism and experimental approaches. BMB Rep. 2014;47(8):417–23.
We are grateful for helpful discussions with members of the Shin laboratory. We also thank Dr. Hyun-Ah Lee and Prof. Do Il Choi for providing strand-specific RNA-Seq protocols and materials. We are also thankful to Soyoung Kim and Dooyoung Lee for helping small RNA northern blot analysis.
This work was supported by the Next-Generation BioGreen 21 Program (No. PJ01332501), Rural Development Administration, Republic of Korea.
Availability of data and materials
The datasets generated and analyzed during the current study are available in the Sequence Read Archive repository (http://www.ncbi.nlm.nih.gov/sra/) under the accession number SRP109719 (strand-specific RNA-Seq: SRR571384, SRR5713898, SRR5713899, SRR5713900, SRR5713901, SRR5713902, SRR5713903, SRR5713904, SRR5713905, SRR5713906, SRR5713907; Small RNA-Seq: SRR5713880, SRR5713881, SRR5713882, SRR5713883, SRR5713885, SRR5713886, SRR5713887, SRR5713888, SRR5713889, SRR5713890; 2P-Seq: SRR5713895, SRR5713896, SRR5713897; Degradome: SRR5713893). Additional supporting tables are included as Additional files.
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Figure S1. Correlation analysis of strand-specific RNA-Seq and small RNA-Seq. Figure S2. GO analysis of up-regulated and down-regulated genes under nitrogen and phosphate starvation in rice. Figure S3. RNA-Seq read distribution of putative lncRNAs responsive to nitrogen starvation and other stressors. Figure S4. Expression patterns of N-responsive putative lncRNAs and quantitative PCR validation results. Figure S5. Predictions of miR169-targeted rice NF-YAs and 5′ RACE results. Figure S6. Degradome sequencing analysis of genes targeted by rice microRNAs. Figure S7. RNA-Seq read distribution on the genomic region of ammonium transporters, and predicted secondary structures. Figure S8. RNA-Seq and Small RNA-Seq read distribution on the genomic region of AMT2.1 and AMT3.3. Figure S9. Sequence alignment of Chr04G0017 genomic region in rice subspecies. Blue characters indicate predicted open reading frame. Figure S10. Peptide modeling results of Chr04G0017-encoding open reading frame. (PDF 33992 kb)
Table S1. Read statistics of the four types of transcriptome data used in this study. Table S2. List of genes involved in nitrogen source uptake, assimilation and transport and their expression level in rice. Table S3. List of amino acid/peptide transporter genes and their expression level in rice. Table S4. Pipeline for identification of putative long non-coding RNAs in rice. Table S5. List of putative long non-coding RNAs identified in this study and their expression level in rice. Table S6. List of miRNAs analyzed in this study and their expression level in rice. Table S7. List of microRNA-targeted genes identified by Degradome-Seq. Table S8. Oligonucleotide sequences used in this study. (XLSX 1030 kb)
About this article
Cite this article
Shin, SY., Jeong, J.S., Lim, J.Y. et al. Transcriptomic analyses of rice (Oryza sativa) genes and non-coding RNAs under nitrogen starvation using multiple omics technologies. BMC Genomics 19, 532 (2018). https://doi.org/10.1186/s12864-018-4897-1