Skip to main content

Advertisement

Transcriptomic analysis reveals adaptive strategies to chronic low nitrogen in Tibetan wild barley

Article metrics

  • 740 Accesses

Abstract

Background

Development of crop cultivars with high low nitrogen (LN) tolerance or nitrogen use efficiency (NUE) is imperative for sustainable agriculture development. Tibetan wild barley is rich in genetic diversity and may provide elite genes for LN tolerance improvement. Little has been known about transcriptional responses of the wild barley to chronic LN stress.

Results

In this study, transcriptomic profiling of two Tibetan wild barley genotypes, LN- tolerant XZ149 and LN-sensitive XZ56 has been conducted using RNA-Seq to reveal the genotypic difference in response to chronic LN stress. A total of 520 differentially expressed genes (DEGs) were identified in the two genotypes at 12 d after LN stress, and these DEGs could be mainly mapped to 49 metabolism pathways. Chronic LN stress lead to genotype-dependent responses, and the responsive pattern in favor of root growth and stress tolerance may be the possible mechanisms of the higher chronic LN tolerance in XZ149.

Conclusion

There was a distinct difference in transcriptional profiling between the two wild barley genotypes in response to chronic LN stress. The identified new candidate genes related to LN tolerance may cast a light on the development of cultivars with LN tolerance.

Background

As an essential component of key macromolecules, nitrogen (N) is quantitatively the most important mineral nutrient for plants [1]. In soils, N is often the most important factor limiting plant growth, and plants frequently encounter N deficiency in their natural habitats. In the past several decades, the increasing use of N fertilizers in crop production has played a major role in increasing yields [2]. However, the main problem is the fact that crop plants only use less than half of the applied N [3], with the remaining N causing severe environmental pollution. There is thus an impending need to realize high productivity while decreasing the rate of N application. This necessitates a comprehensive understanding of molecular mechanisms underlying morphological and physiological adaptation to low N (LN) stress in crops.

Under LN stress, the sessile plants have evolved many adaptive responses. Plants improves the efficiency of N uptake by modifying root architecture, enlarging root system and enhancing the expression of high-affinity transport systems for nitrate and ammonium [4,5,6,7,8,9]. Meanwhile, N utilization efficiency could be also improved in the plants [10]. In addition, the remobilization of N from source organs might be stimulated when plants are subjected to N limitation [11], resulting in the enhancement of N re-assimilation, and maintenance of N economy in plants [12, 13]. Furthermore, the expression levels of a number of genes associated with N metabolism in plants were altered to ensure the survival or complete their life cycle [14, 15].

Tibetan annual wild barley is regarded as one of the progenitors of cultivated barley [16, 17], possessing wider genetic diversity and generally better adaption to N deficiency in comparison with cultivated barley [18]. Some wild barley genotypes with high LN tolerance were identified, providing the elite genetic materials for improving LN tolerance of barley as well as other cereal crops [19]. Two wild barley genotypes differing dramatically in LN tolerance were used for transcriptome analysis at early stage of LN stress (6 h and 48 h) in our previous study [20]. However, plants respond to nutrient deficiency by inducing or repressing different sets of genes at special time [21]. Obviously, the mechanisms of LN tolerance in wild barley still remain to be revealed.

Although understanding of immediate responses is necessary for revealing molecular mechanisms of LN tolerance in plants, it still needs knowing the responses to long-term LN stress. In this study, we investigated genotypic differences in transcriptomic responses to chronic LN stress using the two wild barley genotypes (XZ149 and XZ56) differing in LN tolerance. The major objective of this study was to identify genes associated with the LN tolerance and to understand the molecular mechanisms underlying LN tolerance in wild barley.

Results

Effect of LN stress on growth performances of two wild barley genotypes

XZ149 and XZ56 used in this study exhibited a distinct difference in plant growth at 12 d under LN stress (Fig. 1). Shoot dry weight of the two genotypes was significantly reduced under LN stress, with 15 and 37% reduction in comparison with control for XZ149 and XZ56, respectively (Fig. 1a). On the other hand, root dry weight was increased by 37 and 22% under LN stress for XZ149 and XZ56, respectively (Fig. 1b). As a result, total dry weight per plant under LN stress was only reduced by 7% in XZ149, but as much as 29% in XZ56 (Fig. 1c). LN stress also caused a significant reduction in SPAD value (chlorophyll content), N concentration and accumulation of the two genotypes, with XZ149 being reduction of 13, 11, 25% and XZ56 of 27, 19, 49%, respectively (Fig. 1d, e and f). The less reduction of shoot dry weight, N concentration and accumulation, and SPAD value, as well as the more increase in root dry weight in XZ149 than in XZ56 proved previous finding that XZ149 is more LN tolerant than XZ56 [18].

Fig. 1
figure1

Growth performances of the two wild barley genotypes at 12 d after low N stress. a Shoot dry weight (n = 6); b Root dry weight (n = 6); c Total dry weight (n = 6); (d) SPAD value (n = 10); e N concentration (n = 3); f N accumulation (n = 6). CK: Normal N level (2 mM N); LN: Low N level (0.2 mM N); DW: dry weight. The different letters mean significant difference among treatments and genotypes according to the Duncan’s multiple range, P < 0.05

DEGs at the two N levels

The responding kinetics of HvNRT2.1 to LN stress was investigated in both XZ149 and XZ56. The relative expression was obtained using pair-wise comparison between normal and LN conditions at 1 d, 2 d, 4 d, 8 d and 12 d after LN treatment, respectively (Additional file 1: Figure S1). The transcript level of HvNRT2.1 was the maximum at 12 d after LN stress in both XZ149 and XZ56 (Additional file 1: Figure S1). In our previous study, 12 d was the time for sampling to determine the long-term changes of metabolites under LN stress [22]. Accordingly, in this study we took the samples at 12 d after LN stress for RNA-Seq analysis, to reveal the long-term responses to LN treatment.

In this study, eight cDNA libraries were constructed, and a total of 132,152,200 clean reads were obtained in the all samples. For most samples, 70% of the sequenced reads could be uniquely mapped (Additional file 2: Table S2). The results from real time PCR of 15 responsive genes highly validated the RNA-Seq data (Fig. 2).

Fig. 2
figure2

Quantitative real-time PCR validation of 15 differentially expressed genes (DEGs). a Transcript levels of 15 DEGs and the corresponding expression data of RNA-Seq. The grey columns represent relative expression obtained by qRT-PCR, and the black columns represent relative expression obtained by RNA-Seq. b Comparison between the relative expression obtained from RNA-Seq data and qRT-PCR. The RNA-Seq log2 value of the relative expression (y-axis) has been plotted against qRT-PCR log2 value of the relative expression (x-axis)

The transcriptional levels were analyzed by calculating the FPKM. Meanwhile, FPKM ≥1 at least in one of the samples and FDR < 0.05 was set as screening thresholds for DEGs. Consequently, 520 DEGs were identified using pair-wise comparison for each accession (Additional file 3: Table S3, all the sequences of the DEGs were listed in Additional file 4: Table S4), including both up-regulated (229) and down-regulated (296) genes (Fig. 3). Notably, DEGs in XZ149 (278) were almost equal to those in XZ56 (290) (Fig. 3). Meanwhile, there were 48 DEGs shared commonly in XZ149 and XZ56 (Additional file 5: Table S5). The two accessions showed the different expression patterns, with XZ149 being basically equal in up-regulated and down-regulated DEGs and XZ56 having less up-regulated DEGs than down-regulated ones (Fig. 3).

Fig. 3
figure3

A Venn diagram describing overlaps among differentially expressed genes (DEGs) in XZ149 and XZ56

DEGs over the different time under LN stress

Combined with the previous data on the samples taken in 6 h and 48 h after treatment, the expression patterns of different DEGs were clustered into 13 classes (Additional file 6: Figure S2). Across three time points examined, most of DGEs were only transiently down- or up- expressed at a specific time point in each genotype (Additional file 6: Figure S2). Only 29 DEGs in XZ149 and 11 in XZ56 were consistently expressed at three time points, with only one DEG being consistently up-regulated in XZ149, which encoded aminocyclopropane-1-carboxylate oxidase (ACO) (Additional file 7: Table S6). No gene was consistently down-regulated and only two DEGs, encoding nitrate transporter, were consistently up-regulated in both XZ149 and XZ56 at three time points (Additional file 7: Table S6). The results suggest that plants respond to LN stress by taking different countermeasures, inducing or repressing different sets of genes at certain time. In other words, the molecular mechanisms of LN tolerance could vary with growth stages.

Functional classification of DEGs

By hierarchical clustering analysis, 520 DEGs could be clarified into four classes (Fig. 4). Furthermore, these DEGs were enriched in 10 functional groups through GO enrichment analysis (Fig. 4b, Additional file 8: Figure S3). Molecular function Ontology included four GO terms ‘catalytic activity’, ‘antioxidant activity’, ‘molecular function regulator’ and ‘binding’, while cellular component Ontology only consisted of term ‘extracellular region’. The DEGs associated with ‘metabolic process’ were the most enriched and accounted for 35% of the biological process Ontology (Fig. 4b).

Fig. 4
figure4

Hierarchical cluster, gene ontology (GO) enrichment and KEGG analysis of 520 DEGs at 12 d after low N stress. A total of 520 DEGs were performed on a Hierarchical cluster analysis. The samples and treatments are displayed above each column. Genes are displayed by different colors and relative levels of expression are showed by a color gradient from low (blue) to high (red). b GO enrichment. c KEGG overview. X- axis depicts the different pathway; Y- axis represents the number of DEGs involving in each pathway

KEGG pathway analysis showed that 53 encoded enzymes by the 520 DEGs were mapped to 49 KEGG pathways (Fig. 4c), such as amino acid, carbohydrate, energy and secondary metabolite metabolisms. Among these pathways, 11 DEGs were involved in sucrose and starch metabolism, and three of them were up-regulated only in XZ149 and one down-regulated only in XZ56 (Additional file 9: Table S7).

DEGs involved in nitrogen/carbon metabolism and other nutrient uptake

Many genes involved in N absorption were differentially expressed, with ten DEGs encoding nitrate transporters, which were found at 12 d after LN stress (Table 1). Most of these DEGs were up-regulated and their relative expression was higher in XZ149 than in XZ56 except MLOC_65110 and MLOC_70747 (Table 1). In addition, asparagine synthetase (ASN) DEG (MLOC_63089) involved in N assimilation was up-regulated in XZ149 but not changed in XZ56 (Table 1), and relative content of asparagine (Asn) was also higher in XZ149 than in XZ56 in both leaves and roots (Additional file 10: Table S8). Furthermore, in comparison with XZ56, XZ149 showed higher relative content of soluble protein and activities of the two key assimilating enzymes in leaves, namely NR and GS (Additional file 10: Table S8).

Table 1 Genes encoding protein transporters and enzymes involve in C/N metabolism showing genotypic difference expression in response to chronic low N stress

Carbon and N metabolisms are highly interconnected [23], and LN stress caused alteration in the expression of the genes related to carbon metabolism in barley roots. In the current study, the DEG encoding soluble acid invertase (MLOC_60412) was up-regulated in XZ149, but not changed in XZ56 under LN stress (Table 1). Two trehalose-6-phosphate phosphatase (TPP) DEGs were identified, with one (MLOC_59991) being up-regulated in XZ149 and unchanged in XZ56, and the other one (MLOC_4115) being down-regulated in XZ56 and not changed in XZ149 (Table 1). Moreover, XZ149 had higher relative contents of sucrose and trehalose in comparison with XZ56 (Additional file 10: Table S8). Futhermore, the DEGs related to cell wall remodeling, including pectinesterase (PE, MLOC_65390), pectate lyase (PEL, MLOC_59629), and pectin acetylesterase (PAE, MLOC_55102), were significantly changed under LN stress. MLOC_55102 showed elevated expression in XZ149 but not changed in XZ56, and MLOC_65390 and MLOC_59629 showed the reduced expression in XZ56 but not changed in XZ149 (Table 1).

In addition to nitrate transporter, many DEGs encoding other transporters were also identified, including phosphate transporters, potassium transporters, amino acid transporters, sugar transporters and ABC transporters (Table 1), indicating that N metabolism could affect uptake of many other nutrients in barley through interactive regulation under long-term LN stress. It was also noted that all of the four DEGs encoding sugar transporters were up-regulated only in XZ149 or down-regulated only in XZ56 (Table 1).

Transcription factors and protein kinases

Totally 37 DEGs encoding transcription factors (TFs) were identified in XZ149 and XZ56 under long-term LN stress, less than half of the early responsive TFs (89). These TFs belonged to different families, such as Zinc finger (8), ERF (8), MYB (6), HSF (3), CBF/DREB (3), bZIP (3), NAC (2), WRKY (2), MADS (1) and ABF/AREB (1) (Fig. 5a). Among the 37 TFs, only one DEG was shared by the two genotypes. All other TFs showed different responses to LN stress between the two genotypes, while most DEGs were up-regulated in XZ149 and down-regulated in XZ56. It is noteworthy that only 13 TFs were found in response to long-term LN stress in XZ149, with nine TFs being up-regulated (Fig. 5a).

Fig. 5
figure5

Average linkage hierarchical cluster analysis of transcription factors (TFs) and protein kinases (PKs) identified in differentially expressed genes (DEGs). a Transcription factors (TFs); b Protein kinases (PKs). The samples and treatments are displayed above each column. Genes are displayed by different colors. Relative levels of expression are showed by a color gradient from low (blue) to high (red)

In the present study, 29 DEGs encoding different groups of protein kinases (PKs) were identified under LN stress (Fig. 5b). The PKs belonging to receptor-like kinase (RLK) family was the most enriched among the PKs. Six of the fifteen RLKs carried leucine-rich repeat (LRR) domains, and three contained lectin domains (Fig. 5b). Among RLKs, a wall associated kinase 3 homolog (MLOC_11098) was up-regulated in XZ149 but not changed in XZ56 (Fig. 5b). Besides, a N-acetyl-l-glutamate kinase (NAGK, MLOC_77291) was up-regulated in XZ149 and down-regulated in XZ56, while one NSP-interacting kinase 3-like DEG (MLOC_52843) was down-regulated in XZ56 but not changed in XZ149 (Fig. 5b).

Discussion

Nitrogen deficiency in soil is a major factor limiting crop production worldwide. It is imperative to understand the molecular mechanisms of LN tolerance in plants in order to improve LN tolerance through breeding. It is generally believed that the early responsive DEGs play more important roles for plant survival. Indeed, more DEGs were identified at 6 h (922) and 48 h (842) than at 12 d (520) after LN stress, suggesting the more genes are involved in early response to LN tolerance in plants. Moreover, plants respond to nutrient deficiency by regulating different sets of genes at special time [21]. Actually, the DEGs identified at 12 d of LN stress were different from those identified at the earlier time (Additional file 6: Figure S2). In the present study, we analyzed the transcriptomic responses of the two wild barley genotypes to LN stress by the RNA-Seq. Clearly, there is a distinct difference at the transcriptional level between the two genotypes in their responses to chronic LN stress.

The DEGs associated directly with N metabolism

Nitrate metabolism starts from N acquisition from soil by nitrate transporters. N deficiency enhances the expression of nitrate transporters, and high affinity nitrate transporters played critical roles in N acquisition and remobilization in N-starved plants [5, 6]. The higher relative expression of most nitrate transporter DEGs in XZ149 shows its higher ability of nitrate uptake and N reallocation compared with XZ56 under LN stress, which is also confirmed by higher relative soluble protein and activity of NR and GS in the leaves of XZ149. Consequently, XZ149 had higher N concentration and more N accumulation in the plant tissues in comparison with XZ56 (Additional file 10: Table S8, Fig. 1e and f).

It is well documented that N limitation acts as an exogenous senescence-triggering factor and enhances N remobilization [11]. Amino acids from protein degradation or assimilation will be translocated and reallocated as glutamine (Gln) and Asn, the major N carriers in the phloem of many plants [24,25,26]. Asn synthetase gene 1 (ASN1), controlling Asn synthesis will be up-expressed in the plants exposed to LN stress. In addition, over-expression of ASN1 in Arabidopsis enhanced its tolerance to limited N condition [27]. In the current study, an Arabidopsis ASN1 homolog DEG (MLOC_63089) was only up-regulated in XZ149, but not changed in XZ56 at 12 d under LN stress. In view of higher relative content of Asn in XZ149 than in XZ56 in both leaves and roots (Additional file 10: Table S8), it may be suggested that higher LN tolerance in XZ149 is partially contributed by its better reallocation of N.

Chronic LN stress alters carbon metabolism

It is well known that carbon and N metabolism are highly interconnected [23, 28, 29]. Indeed, 11 DEGs involved in sucrose and starch metabolism were found in barley roots under LN stress (Additional file 9: Table S7). Moreover, four DEGs encoding sugar transporters were up-regulated in XZ149 but down-regulated in XZ56 (Table 1). Meanwhile, more accumulation of sucrose was also observed in the roots of LN-treated plants relative to control (Additional file 10: Table S8), and moreover, XZ149 accumulated more sucrose in its roots than XZ56 under LN stress. As sucrose plays a role as a stimulant for lateral root formation [30], thus it could be understood that XZ149 showed better growth than XZ56 under LN stress (Fig. 1).

The gene encoding soluble acid invertase (MLOC_60412) was up-regulated in XZ149 but not changed in XZ56 after 12 d under LN stress, and the results were consistent with those in the previous study, where the sampling were taken at 6 h and 48 h after LN treatment [20]. As invertase could irreversibly catalyze the hydrolysis of sucrose [31], it may be suggested that sucrose degradation was continuously enhanced in XZ149. In fact, it is well known that invertase is crucial in the carbon supply from sucrose to the non-photosynthetic cells of plants [32]. Thus, the consistently enhanced expression of invertase in the roots of XZ149 may provide more carbon source and energy, accounting for its better root growth under LN stress.

Trehalose-6-phosphate phosphatase (TPP) involves in trehalose biosynthesis, and trehalose plays a protective role against stress in plants [33]. Currently, one DEG encoding TPP was up-regulated in XZ149 but unchanged in XZ56, and another one was down-regulated in XZ56 but not changed in XZ149. The higher relative expression of DEGs encoding TPP in XZ149 may enhance its higher biosynthesis of stress protectant trehalose (Additional file 10: Table S8). It was reported that OsTPP1 was induced by stress [34], and its over-expression conferred stress tolerance in rice [35]. Therefore, it may be assumed that the higher relative expression of TPP and biosynthesis of trehalose in XZ149 may be a mechanism for protecting the plants against LN stress.

The DEGs associated with cell wall remodeling

The dynamic and active structure of cell wall is necessary for plants to remodel when responding to external stimuli [36, 37]. Cell wall would be remodeled when plants are subjected to nutrition deficiency, in particular N deficiency [36]. Cell wall related genes play a role in inducing cell wall relaxation and maintaining N assimilation under N starvation [38]. In this study, we found that cell wall remodeling DEGs, including PE, PEL, and PAE were more highly expressed in XZ149 than in XZ56 under LN stress. In addition, it should be noted that above-mentioned three enzymes are involved in pectin degradation. Hence, in comparison with XZ56, more pectin breakdown in XZ149 may provide it more materials for use in other biological processes, compensating for the depressed photosynthetic carbon assimilation under LN stress [15]. Thus, higher relative expression of cell wall related DEGs in XZ149 may be beneficial for its coping with LN stress.

Transcription factors and hormone signaling

Expression of plant-specific Dof1 TFs has been proved to improve plant growth under LN condition [39, 40]. Several R2R3-type MYB TFs are involved in plant stress responses [41], and over-expression of OsMYB48–1 enhanced drought and salinity tolerance [42]. In this study, we found that one Dof and two MYB were down-regulated in XZ56, while three MYB were up-regulated in XZ149 (Fig. 5a). It may be assumed that these IFs may contribute to the development of complex signaling webs and should be potential candidates for improving LN stress tolerance.

The synthesis of plant hormone ethylene is tightly regulated in response to stresses [43]. Genes encoding ACO, a rate-limiting enzyme in ethylene synthesis, have been reported to be response to N deficiency [44, 45]. In our previous study, three ACO homologs were up-regulated only in XZ149 under LN stress, however all of them were transiently expressed at either 6 h or 48 h [20]. Here, another ACO DEG (MLOC_34709) was detected, which was up-regulated in XZ149 at all the three time points, but unchanged in XZ56 at the later two time points (Additional file 7: Table S6). Considering the fine tuning between NRT2.1 expression and ethylene biosynthesis under LN level [46], it may be suggested that the up-regulation of MLOC_34709 in XZ149 contributes to its early and long-term LN tolerance.

Protein kinases

Protein kinases are important for adaptation to abiotic stress in plants [47]. The NSP-interacting kinase (NIK), belonging to the five leucine-rich repeats-containing receptor-like kinase subfamily, is a transducer of plant defense signaling [48]. In this study, we found that NSP-interacting kinase 3-like (MLOC_52843) was down-regulated only in XZ56 (Fig. 5b), suggesting that XZ149 may enhance its LN stress tolerance by maintaining its defense capability. In addition, a wall associated kinase 3 (WAK) homologous DEG (MLOC_11098) was up-regulated in XZ149 but unchanged in XZ56 (Fig. 5b). WAK is involved in root growth under LN stress [49], thus its higher expression in XZ149 could partially account for its better adaptation to LN stress.

N-acetyl-l-glutamate kinase (NAGK), the prototype of the amino acid kinase family, can catalyzes a key step of Arginine (Arg) synthesis. Arg content was increased in plants under N-deficient condition [50]. It was reported that some amino acids play roles in alleviating abiotic stress [51]. Currently, a NAGK homolog (MLOC_77291) was up-regulated in XZ149 and down-regulated in XZ56 (Fig. 5b), indicating that the increase of NAGK in XZ149 could confer its LN tolerance.

Conclusion

The current study identified the DEGs in barley roots, which respond to chronic (12 days) LN stress. There was a dramatic difference between the two Tibetan wild barley genotypes in transcriptomic response to chronic LN stress. The pattern of genetic responses in favor of better root growth and higher stress tolerance may be the possible mechanisms of the higher LN tolerance in XZ149. In addition, some new candidate genes related to LN tolerance were identified, which could be useful for developing barley cultivars with LN tolerance.

Methods

Plant materials and LN treatment

The experiment using two Tibetan wild barley accessions XZ149 and XZ56 (LN- tolerant and sensitive genotypes, respectively) was carried out in black plastic pots (5 L) in a greenhouse with natural light. The wild barley accessions were collected from Tibetan area in last century and kindly presented by professor Sun of Huazhong Agricultural University, China. The nutrition solution was used according to Quan et al. [20]. The solution was renewed every five days, continuously aerated with pumps. Treatments were conducted on three-leaf-stage seedlings with two N levels (0.2 mM N as LN treatment, 2 mM N as control).

For biomass determination, the plants were harvested at 12 d after LN treatment and separated into shoots and roots. Dry weight was recorded after the samples were dried at 105 °C for 30 min and to constant weight at 80 °C. Meanwhile fresh plant tissues were taken for use in determining nitrate reductase (NR) activity, glutamine synthetase (GS) activity and soluble protein content with three biological replications, and the content of metabolites with four biological replications. The roots of the two barley accessions under both N treatments were sampled with three biological replicates at 1 d, 2 d, 4 d, 8 d and 12 d after treatments for the time course analysis of gene HvNRT2.1 expression. For RNA-Seq analysis, the samples were taken at 12 d after treatments. Roots of four seedlings for each treatment were pooled as one biological replication. Totally eight samples [2 genotypes × 2 treatments × 2 biological replications] were taken for analysis.

Physiological measurement

N concentration in plant tissue was determined using Foss Kjeltec 8400. Soluble protein content was measured as described by Andrews et al. [52]. NR and GS activity were determined according to Kaiser et al. [53] and Masclaux-Daubresse et al. [54], respectively. The metabolites were extracted and analyzed using 7890A/5975C GC–MS system (Agilent, USA) and AMDIS 32 software according to Quan et al. [22].

cDNA library construction and sequencing for RNA-seq

The total RNA was extracted using miRNeasy mini kit (QIAGEN, Germany) following the manufacturer’s specification. RNA degradation, integrity, abundances and purity were assessed for meeting the requirements [20]. cDNA libraries were constructed using the Illumina TruSeq™ RNA Sample Preparation Kit (Illumina, San Diego, CA, USA) according to the manufacturer’s instructions. Briefly, mRNA was obtained from the total RNA using magnetic beads with poly-T oligonucleotide. Then the random fragmentation of the purified mRNA was reversely transcribed into cDNA. After ligated with the adapters on both ends, DNA fragments were selectively amplified and enriched. Subsequently, the purified PCR products were quantified using Agilent Bioanalyzer 2100 system. After cluster generation, the final cDNA library was sequenced on an Illumina NextSeq 500 platform.

Raw reads with 75 bp single-end were initially processed to remove adapters sequences, empty sequences and low-quality bases, and then the Q20, Q30, GC contents, and sequence duplication level of the clean data were analyzed. Then the clean reads were mapped against the barley reference genomes using TopHat (http://tophat.cbcb.umd.edu/), and finally the mapping results were analyzed to identify splice junctions.

Identification of the DEGs and validation of RNA-Seq by real time PCR

Gene expression levels were calculated by the FPKM values (fragments per kilobase of exon per million fragments mapped reads) [55]. Fold-changes were defined as normalized read count abundance for the LN-stressed samples divided by that of the control samples. To identify differentially expression genes (DEGs), the difference in expression between control and LN treatment was analyzed using the DESeq R package (1.10.1) [56]. An FDR (false discovery rate) of 0.05 was used for determining significant DEGs [57].

To validate the reliability of the RNA-Seq results, the expression of candidate genes was determined by real time PCR assay using the RNA for RNA-Seq. The first strand cDNA was synthesized using PrimerScript™ RT reagent Kit with gDNA Eraser (Takara, Japan). The gene-specific primers, designed by primer-blast (http:/www.ncbi.nlm.nih. gov/tools/primer-blast/), were presented in Additional file 11: Table S1. All real time PCR analyses were performed on a CFX96 system (Bio-Rad USA) with two biological replicates and three technical replicates. HvGAPDH was used as an internal control. The relative expression was calculated by the comparative CT method and expressed as the fold change referred to the expression in the control plants [58].

Statistical analysis

Gene Ontology (GO) annotation and KEGG (Kyoto Encyclopedia of Genes and Genomes) analysis for the DEGs were conducted using the Blast2GO program [59] according to Zeng et al. [60]. GO terms were tested by applying tools for GO enrichment (http://systemsbiology.cau.edu.cn/agriGOv2/) at p-values ≤0.05 [61]. Venn diagram was made on jvenn (http://jvenn.toulouse.inra.fr/app/example.html) [62]. Heatmaps and hierarchical clustering were generated with genesis 1.8.1. Significant differences of gene expression between treatments were tested using a DPS statistical software, and the difference at P < 0.05 was considered as significant.

Abbreviations

ACO:

Aminocyclopropane-1-carboxylate oxidase

Arg:

Arginine

Asn:

Asparagine

ASN:

Asparagine synthetase

DEGs:

Differentially expressed genes

FDR:

False discovery rate

FPKM:

Fragments per kilobase of exon per million fragments mapped reads

Gln:

Glutamine

GO:

Gene Ontology

GS:

Glutamine synthetase

KEGG:

Kyoto Encyclopedia of Genes and Genomes

LN:

Low nitrogen

LRR:

Leucine-rich repeat domain containing kinase

N:

Nitrogen

NAGK:

N-acetyl-l-glutamate kinase

NR:

Nitrate reductase

PAE:

Pectin acetylesterase

PE:

Pectinesterase

PEL:

Pectate lyase

PKs:

Protein kinases

RLK:

Receptor-like kinase

TFs:

Transcription factors

TPP:

Trehalose-6-phosphate phosphatase

WAK:

Wall associated kinase

References

  1. 1.

    Crawford NM, Forde BG. Molecular and developmental biology of inorganic nitrogen nutrition. In: The Arabidopsis Book, vol. 1; 2002. p. e0011.

  2. 2.

    Tilman D, Reich PB, Knops J, Wedin D, Mielke T, Lehman C. Diversity and productivity in a long-term grassland experiment. Science. 2001;294:843–5.

  3. 3.

    Socolow RH. Nitrogen management and the future of food: lessons from the management of energy and carbon. Proc Natl Acad Sci. 1999;96:6001–8.

  4. 4.

    Engineer CB, Kranz RG. Reciprocal leaf and root expression of AtAmt1.1 and root architectural changes in response to nitrogen starvation. Plant Physiol. 2007;143:236–50.

  5. 5.

    Kiba T, Feria-Bourrellier AB, Lafouge F, Lezhneva L, Boutet-Mercey S, Orsel M, et al. The Arabidopsis nitrate transporter NRT2.4 plays a double role in roots and shoots of nitrogen-starved plants. Plant Cell. 2012;24:245–58.

  6. 6.

    Lezhneva L, Kiba T, Feria-Bourrellier AB, Lafouge F, Boutet-Mercey S, Zoufan P, et al. The Arabidopsis nitrate transporter NRT2.5 plays a role in nitrate acquisition and remobilization in nitrogen-starved plants. Plant J. 2014;80:230–41.

  7. 7.

    López-Bucio J, Cruz-Ramırez A, Herrera-Estrella L. The role of nutrient availability in regulating root architecture. Curr Opin Plant Biol. 2003;6:280–7.

  8. 8.

    Remans T, Nacry P, Pervent M, Girin T, Tillard P, Lepetit M, et al. A central role for the nitrate transporter NRT2.1 in the integrated morphological and physiological responses of the root system to nitrogen limitation in Arabidopsis. Plant Physiol. 2006;140:909–21.

  9. 9.

    von Wire’n N, Lauter FR, Ninnemann O, Gillissen B, Walch-Liu P, Engels C, et al. Differential regulation of three functional ammonium transporter genes by nitrogen in root hairs and by light in leaves of tomato. Plant J. 2000;21:167–75.

  10. 10.

    Hirel B, Le Gouis J, Ney B, Gallais A. The challenge of improving nitrogen use efficiency in crop plants: towards a more central role for genetic variability and quantitative genetics within integrated approaches. J Exp Bot. 2007;58:2369–87.

  11. 11.

    Wingler A, Marès M, Pourtau N. Spatial patterns and metabolic regulation of photosynthetic parameters during leaf senescence. New Phytol. 2004;161:781–9.

  12. 12.

    Bernard SM, Habash DZ. The importance of cytosolic glutamine synthetase in nitrogen assimilation and recycling. New Phytol. 2009;182:608–20.

  13. 13.

    Miflin BJ, Habash DZ. The role of glutamine synthetase and glutamate dehydrogenase in nitrogen assimilation and possibilities for improvement in the nitrogen utilization of crops. J Exp Bot. 2002;53:979–87.

  14. 14.

    Gelli M, Duo Y, Konda AR, Zhang C, Holding D, Dweikat I. Identification of differentially expressed genes between sorghum genotypes with contrasting nitrogen stress tolerance by genome-wide transcriptional profiling. BMC Genomics. 2014;15:179.

  15. 15.

    Zhao W, Yang X, Yu H, Jiang W, Sun N, Liu X, et al. RNA-Seq-based transcriptome profiling of early nitrogen deficiency response in cucumber seedlings provides new insight into the putative nitrogen regulatory network. Plant Cell Physiol. 2015;56:455–67.

  16. 16.

    Dai F, Nevo E, Wu DZ, Comadran J, Zhou MX, Qiu L, et al. Tibet is one of the centers of domestication of cultivated barley. Proc Natl Acad Sci U S A. 2012;109:16969–73.

  17. 17.

    Dai F, Wang X, Zhang XQ, Chen Z, Nevo E, Jin G, et al. Assembly and analysis of a qingke reference genome demonstrate its close genetic relation to modern cultivated barley. Plant Biotechnol J. 2017. https://doi.org/10.1111/pbi.12826.

  18. 18.

    Quan X, Zeng J, Han Z, Zhang G. Ionomic and physiological responses to low nitrogen stress in Tibetan wild and cultivated barley. Plant Physiol Bioch. 2017;111:257–65.

  19. 19.

    Yang L, Hu H, Zhu B, Jin X, Wu F, Zhang G. Genotypic variations of N use efficiency in Tibetan wild and cultivated barleys. J Zhejiang Univ (Agric & Life Sci). 2014;40:155–64.

  20. 20.

    Quan X, Zeng J, Ye L, Chen G, Han Z, Shah JM, et al. Transcriptome profiling analysis for two Tibetan wild barley genotypes in responses to low nitrogen. BMC Plant Biol. 2016;16:1.

  21. 21.

    Cai H, Lu Y, Xie W, Zhu T, Lian X. Transcriptome response to nitrogen starvation in rice. J Biosci. 2012;37:731–47.

  22. 22.

    Quan X, Qian Q, Ye Z, Zeng J, Zhang G. Metabolic analysis of two contrasting wild barley genotypes grown hydroponically reveals adaptive strategies in response to low nitrogen stress. J Plant Physiol. 2016b;206:59–67.

  23. 23.

    Schluter U, Mascher M, Colmsee C, Scholz U, Bräutigam A, Fahnenstich H, et al. Maize source leaf adaptation to nitrogen deficiency affects not only nitrogen and carbon metabolism but also control of phosphate homeostasis. Plant Physiol. 2012;160:1384–406.

  24. 24.

    Gaufichon L, Reisdorf-Cren M, Rothstein SJ, Chardon F, Suzuki A. Biological functions of asparagine synthetase in plants. Plant Sci. 2010;179:141–53.

  25. 25.

    Gaufichon L, Rothstein SJ, Suzuki A. Asparagine metabolic pathways in Arabidopsis. Plant Cell Physiol. 2015;57:675–89.

  26. 26.

    Lea PJ, Sodek L, Parry MA, Shewry PR, Halford NG. Asparagine in plants. Ann Appl Biol. 2007;150:1–26.

  27. 27.

    Lam H-M, Wong P, Chan H-K, Yam K-M, Chen L, Chow C-M, et al. Overexpression of the ASN1 gene enhances nitrogen status in seeds of Arabidopsis. Plant Physiol. 2003;132:926–35.

  28. 28.

    Antal T, Mattila H, Hakala-Yatkin M, Tyystjärvi T, Tyystjärvi E. Acclimation of photosynthesis to nitrogen deficiency in Phaseolus vulgaris. Planta. 2010;232:887–98.

  29. 29.

    Nunes-Nesi A, Fernie AR, Stitt M. Metabolic and signaling aspects underpinning the regulation of plant carbon nitrogen interactions. Mol Plant. 2010;3:973–96.

  30. 30.

    Roycewicz P, Malamy JE. Dissecting the effects of nitrate, sucrose and osmotic potential on Arabidopsis root and shoot system growth in laboratory assays. Philos Trans R Soc B. 2011;367:1489–500.

  31. 31.

    Koch K. Sucrose metabolism: regulatory mechanisms and pivotal roles insugar sensing and plant development. Curr Opin Plant Biol. 2004;7:235–46.

  32. 32.

    Barratt DP, Derbyshire P, Findlay K, Pike M, Wellner N, Lunn J, et al. Normal growth of Arabidopsis requires cytosolic invertase but not sucrose synthase. Proc Natl Acad Sci. 2009;106:13124–9.

  33. 33.

    Jang IC, Oh SJ, Seo JS, Choi WB, Song SI, Kim CH, et al. Expression of a bifunctional fusion of the Escherichia coli genes for trehalose-6-phosphate synthase and trehalose-6-phosphate phosphatase in transgenic rice plants increases trehalose accumulation and abiotic stress tolerance without stunting growth. Plant Physiol. 2003;131:516–24.

  34. 34.

    Pramanik MHR, Imai R. Functional identification of a trehalose 6-phosphate phosphatase gene that is involved in transient induction of trehalose biosynthesis during chilling stress in rice. Plant Mol Biol. 2005;58:751–62.

  35. 35.

    Ge LF, Chao DY, Shi M, Zhu MZ, Gao JP, Lin HX. Overexpression of the trehalose-6-phosphate phosphatase gene OsTPP1 confers stress tolerance in rice and results in the activation of stress responsive genes. Planta. 2008;228:191–201.

  36. 36.

    Fernandes JC, Goulao LF, Amâncio S. Regulation of cell wall remodeling in grapevine (Vitis vinifera L.) callus under individual mineral stress deficiency. J Plant Physiol. 2016;190:95–105.

  37. 37.

    Guerriero G, Hausman JF, Strauss J, Ertan H, Siddiqui KS. Lignocellulosic biomass: biosynthesis, degradation and industrial utilization. Eng Life Sci. 2016;16:1–16.

  38. 38.

    Landi S, Esposito S. Nitrate uptake affects cell wall synthesis and modeling. Front Plant Sci. 2017;8:1376.

  39. 39.

    Kurai T, Wakayama M, Abiko T, Yanagisawa S, Aoki N, Ohsugi R. Introduction of the ZmDof1 gene into rice enhances carbon and nitrogen assimilation under low-nitrogen conditions. Plant Biotechnol J. 2011;9:826–37.

  40. 40.

    Yanagisawa S, Akiyama A, Kisaka H, Uchimiya H, Miwa T. Metabolic engineering with Dof1 transcription factor in plants: improved nitrogen assimilation and growth under low-nitrogen conditions. Proc Natl Acad Sci. 2004;101:7833–8.

  41. 41.

    Pireyre M, Burow M. Regulation of MYB and bHLH transcription factors: a glance at the protein level. Mol Plant. 2015;8:378–88.

  42. 42.

    Xiong H, Li J, Liu P, Duan J, Zhao Y, Guo X, et al. Overexpression of OsMYB48-1, a novel MYB-related transcription factor, enhances drought and salinity tolerance in rice. PLoS One. 2014;9:e92913.

  43. 43.

    Merchante C, Alonso JM, Stepanova AN. Ethylene signaling: simple ligand, complex regulation. Curr Opin Plant Biol. 2013;16:554–60.

  44. 44.

    Bi YM, Wang RL, Zhu T, Rothstein SJ. Global transcription profiling reveals differential responses to chronic nitrogen stress and putative nitrogen regulatory components in Arabidopsis. BMC Genomics. 2007;8:281.

  45. 45.

    Peng MS, Bi YM, Zhu T, Rothstein SJ. Genome-wide analysis of Arabidopsis responsive transcriptome to nitrogen limitation and its regulation by the ubiquitin ligase gene NLA. Plant Mol Biol. 2007;65:775–97.

  46. 46.

    Zheng D, Han X, An Y, Guo H, Xia X, Yin W. The nitrate transporter NRT2.1 functions in the ethylene response to nitrate deficiency in Arabidopsis. Plant Cell Environ. 2013;36:1328–37.

  47. 47.

    Hunter T, Karin M. The regulation of transcription by phosphorylation. Cell. 1992;70:375–87.

  48. 48.

    Santos AA, Lopes KVG, Apfata JAC, Fontes EPB. NSP-interacting kinase, NIK: a transducer of plant defence signalling. J Exp Bot. 2010;61:3839–45.

  49. 49.

    Kiba T, Krapp A. Plant nitrogen acquisition under low availability: regulation of uptake and root architecture. Plant Cell Physiol. 2016;57:707–14.

  50. 50.

    Krapp A, Berthomé R, Orsel M, Mercey-Boutet S, Yu A, Castaings L, et al. Arabidopsis roots and shoots show distinct temporal adaptation patterns toward 5nitrogen starvation. Plant Physiol. 2011;157:1255–82.

  51. 51.

    Sharma SS, Dietz KJ. The significance of amino acids and amino acid-derived molecules in plant responses and adaptation to heavy metal stress. J Exp Bot. 2006;57:711–26.

  52. 52.

    Andrews M, Sprent JI, Raven JA, Eady PE. Relationships between shoot to root ratio, growth and leaf soluble protein concentration of Pisum sativum, Phaseolus vulgaris and Triticum aestivum under different nutrient deficiencies. Plant Cell Environ. 1999;22:949–58.

  53. 53.

    Kaiser WM, Weiner H, Huber SC. Nitrate reductase in higher plants: a case study for transduction of environmental stimuli into control of catalytic activity. Physiol Plant. 1999;105:384–9.

  54. 54.

    Masclaux-Daubresse C, Reisdorf-Cren M, Pageau K, Lelandais M, Grandjean O, Kronenberger J, et al. Glutamine synthetase-glutamate synthase pathway and glutamate dehydrogenase play distinct roles in the sink-source nitrogen cycle in tobacco. Plant Physiol. 2006;140:444–56.

  55. 55.

    Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

  56. 56.

    Robinson MD, Oshlack A. A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol. 2010;11:R25.

  57. 57.

    Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28:511–5.

  58. 58.

    Schmittgen TD, Livak KJ. Analyzing real-time PCR data by the comparative CT method. Nat Protoc. 2008;3:1101–8.

  59. 59.

    Conesa A, Götz S, García-Gómez JM, Terol J, Talón M, Roble M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21:3674–6.

  60. 60.

    Zeng J, He X, Wu D, Zhu B, Cai S, Nadira UA, et al. Comparative transcriptome profiling of two Tibetan wild barley genotypes in responses to low potassium. PLoS One. 2014;9:e100567.

  61. 61.

    Tian T, Yue L, Hengyu Y, Qi Y, Xin Y, Zhou D, et al. agriGO v2.0: a GO analysis toolkit for the agricultural community, 2017 update. Nucleic Acids Res. 2017. https://doi.org/10.1093/nar/gkx382.

  62. 62.

    Bardou P, Mariette J, Escudié F, Djemiel C, Klopp C. jvenn: an interactive Venn diagram viewer. BMC Bioinformatics. 2014;15:293.

Download references

Acknowledgments

We thank Prof. Dongfa Sun (Huazhong Agricultural University, China) for providing Tibetan wild barley accessions.

Funding

We deeply thank Natural Science Foundation of China (31620103912), China Agriculture Research System (CARS-05) and Jiangsu Collaborative Innovation Center for Modern Crop Production (JCIC-MCP) for their financial supports to this research, including experimental implementation, sampling and data analysis.

Availability of data and materials

The datasets used and/or analyzed during the current study are available from the first author on reasonable request.

Author information

XQ and JZ designed and performed research. XQ and GC analyzed data. XQ and GZ wrote the paper. All authors have read, edited and approved the current version of the manuscript.

Correspondence to Guoping Zhang.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Figure S1. Temporal expression of the HvHRT2.1 gene in XZ149 and XZ56 under low N stress. The relative expression was calculated by the expression of LN stress divide by that of control. Different lowercase represents significant differences according to the Duncan’s multiple range, P < 0.05, n = 4. Primers of HvNRT2.1 and GAPDH for real time PCR are listed in Table S1. (DOCX 16 kb)

Additional file 2:

Table S2. Summary of mapping reads of the RNA-Seq (XLSX 11 kb)

Additional file 3:

Table S3. The FPKM value of 520 DEGs in XZ149 and XZ56 (XLSX 113 kb)

Additional file 4:

Table S4. Gene accession numbers and sequences of 520 DEGs (XLSX 478 kb)

Additional file 5:

Table S5. DEGs at 12 d after low N stress in XZ149 and XZ56 (XLSX 160 kb)

Additional file 6:

Figure S2. Hierarchical cluster of DEGs at three time points in XZ149 and XZ56. The samples and treatments are displayed above each column. Genes are displayed by different colors and relative levels of expression are showed by a color gradient from low (green) to high (red). (PNG 142 kb)

Additional file 7:

Table S6. DEGs consistently expressed at three time points under low N stress in each genotype (XLSX 13 kb)

Additional file 8:

Figure S3. GO annotation and enrichment analysis (PNG 65 kb)

Additional file 9:

Table S7. Enzymes and sequences in one KEGG pathway (XLSX 12 kb)

Additional file 10:

Table S8. Traits related with N metabolism in two wild barley genotypes XZ149 and XZ56 under low and normal N levels (DOC 40 kb)

Additional file 11:

Table S1. The primers used in real time PCR (DOC 53 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Keywords

  • Barley
  • Chronic low N stress
  • Low N tolerance
  • RNA-Seq
  • Genotypes
  • Differentially expressed genes