Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Comparative transcriptome sequencing of tolerant rice introgression line and its parents in response to drought stress

BMC Genomics volume 15, Article number: 1026 (2014) Cite this article

Abstract

Background

Rice (Oryza sativa. L) is more sensitive to drought stress than other cereals, and large genotypic variation in drought tolerance (DT) exists within the cultivated rice gene pool and its wild relatives. Selective introgression of DT donor segments into a drought-sensitive (DS) elite recurrent parent by backcrossing is an effective way to improve drought stress tolerance in rice. To dissect the molecular mechanisms underlying DT in rice, deep transcriptome sequencing was used to investigate transcriptome differences among a DT introgression line H471, the DT donor P28, and the drought-sensitive, recurrent parent HHZ under drought stress.

Results

The results revealed constitutively differential gene expression before stress and distinct global transcriptome reprogramming among the three genotypes under a time series of drought stress, consistent with their different genotypes and DT phenotypes. A set of genes with higher basal expression in both H471 and P28 compared with HHZ were functionally enriched in oxidoreductase and lyase activities, implying their positive role in intrinsic DT. Gene Ontology analysis indicated that common up-regulated genes in all three genotypes under mild drought stress were enriched in signaling transduction and transcription regulation. Meanwhile, diverse functional categories were characterized for the commonly drought-induced genes in response to severe drought stress. Further comparative transcriptome analysis between H471 and HHZ under drought stress found that introgression caused wide-range gene expression changes; most of the differentially expressed genes (DEGs) in H471 relative to HHZ under drought were beyond the identified introgressed regions, implying that introgression resulted in novel changes in expression. Co-expression analysis of these DEGs represented a complex regulatory network, including the jasmonic acid and gibberellin pathway, involved in drought stress tolerance in H471.

Conclusions

Comprehensive gene expression profiles revealed that genotype-specific drought induced genes and genes with higher expression in the DT genotype under normal and drought conditions contribute jointly to DT improvement. The molecular genetic pathways of drought stress tolerance uncovered in this study, as well as the DEGs co-localized with DT-related QTLs and introgressed intervals, will serve as useful resources for further functional dissection of the molecular mechanisms of drought stress response in rice.

Background

Rice (Oryza sativa L.) is a staple food for more than half of the world’s population, especially in developing countries. Drought is the most serious environmental stress, limiting crop growth and productivity: drought-induced loss in crop yield probably exceeds losses from all other causes [1]. Drought tolerance (DT), therefore, is a major aim of rice breeding, especially in tropical Asian and African countries [2]. DT is a complex trait, and a number of quantitative trait loci (QTLs) for DT in rice have been identified [3]. However, breeding drought-tolerant rice is hard to achieve by conventional strategies, including marker-assisted selection. Understanding of the molecular mechanisms underlying DT is therefore needed for successful, knowledge-based crop improvement [4].

Determining the mechanisms directly involved in DT remains a challenging task because it involves several metabolic and morphologically adaptive pathways [5, 6]. Abscisic acid (ABA) is an important phytohormone involved in drought stress tolerance in plants, whose mechanism in plant DT is relatively clear. Under drought stress, ABA-mediated stomatal closure is the mechanism used by plants to adapt to water deficiency [7–9]. Reactive oxygen species (ROS), including hydrogen peroxide, which are widely generated under stress, have been proposed to function as second messengers in ABA signaling in guard cells [10–13]. In guard cells, ABA-stimulated ROS accumulation activates plasma membrane calcium channels and triggers stomatal closure [11, 14]. It was reported that jasmonic acid (JA) is also associated with stomatal closure under drought stress, the detailed molecular mechanisms remain elusive [15–17]. Though there is no evidence in relevance of gibberellic acid (GA) and stomatal closure under drought stress, reduction of GA levels and signaling has been shown to contribute to plant growth inhibition under several abiotic stresses, including cold, salt and osmotic stress [18, 19].

Many efforts have been made to identify the genes involved in drought stress tolerance in a number of plant species. Several drought-responsive genes encoding late embryogenesis abundant proteins, dehydration-responsive element binding protein, and protein phosphatase 2C were characterized as key components in the molecular network of DT [20]. In rice, genome-wide gene expression analyses identified many drought stress-responsive genes [21–24]. However, many of these drought stress-responsive genes have unknown functions, and details of their DT molecular mechanisms remain to be determined [20, 25].

Rice is more sensitive to drought stress than other cereals, and large genotypic variation in DT exists within the cultivated rice gene pool and its wild relatives [2]. Selective introgression of DT donor segments into a drought-sensitive (DS) elite recurrent parent by backcrossing is an effective way to improve drought stress tolerance in rice [2, 26]. In a previous study, we developed a DT introgression line (IL), H471, using the DT donor PSBRC28 (P28) and the DS recurrent parent Huang-Hua-Zhan (HHZ) (unpublished data). Compared with completely different genotypes with contrasting performance on the target trait, the ILs can largely reduce the genetic background noise in comparative transcriptomic analysis with the recurrent parent because the selected DT ILs carry a small number of genomic segments from a known DT donor [26]. Taking advantage of the combination of genome DNA re-sequencing and next-generation RNA sequencing (RNA-seq), ILs can effectively identify genes related to DT, thereby increasing our knowledge of the molecular mechanisms of this complex trait in rice [23]. In the present study, DT IL and its parental lines were used to analyze their transcriptome changes under drought stress comparatively, with the aim of extending our understanding of the genetic mechanisms of DT in rice.

Results

Drought stress physiology of three rice genotypes

At the tillering stage, there was no visible difference observed among H471, P28, and HHZ after 1 and 2 days of drought stress; however, obvious leaf rolling in HHZ was observed after 3 days of drought stress. This phenomenon was not observed in H471 and P28 until 4 days of drought stress (Figure 1A). The yield performance of three genotypes were remarkably reduced by drought stress compared with the well-watered control; however, H471 and P28 achieved 26% and 21% higher grain yields compared with HHZ under drought conditions (Figure 1B).To investigate the physiological difference in DT of three genotypes, several indices of drought-induced effects on leaves at the tillering stage were measured. The water loss rate (WLR) from excised leaves was determined for the three genotypes: H471 and P28 showed relatively lower WLR than HHZ over a period of 10 h (Figure 2A). Accordingly, the relative water contents (RWCs) in H471 (76.9%) and P28 (78.8%) were significantly higher than that in HHZ (67.6%) after 3 days of drought stress (Figure 2B).

Figure 1
figure 1

Phenotypes of three genotypes under well-watered and drought stress conditions. (A) Well-watered and 3 days drought stressed phenotypes of P28, HHZ, and H471 at the tillering stage. (B) Grain yield performance of HHZ, P28, and H471 under drought stressed and well-watered conditions. Each column represents mean ± s.d. (nine replicates); **p <0.01 versus HHZ (ANOVA, Dunnett’s multiple comparison test).

Full size image
Figure 2
figure 2

Differential drought physiological traits of H471, P28, and HHZ. (A) Water loss rate of P28, HHZ, and H471. For each replicate, 15 fully expanded leaves at the tillering stage were used in a triplicate experiment. (B) Relative water content of HHZ, P28, and H471 treated with 3 days of natural drought using the fully expanded leaves. (C) The stomatal closure rate in HHZ, P28, and H471 (upper epidermis: 6 randomly selected high power field for HHZ; 5 for P28; 6 for H471; lower epidermis: 8 randomly selected high power field for HHZ; 6 for P28; 6 for H471). (D) Relative electrolyte leakage of HHZ, P28, and H471 after 3 days of drought stress. (E), (F) Activity of ROS-scavenging enzymes (CAT and APX) in three rice genotypes subjected to drought stress. Each column represents mean ± s.d. (three replicates); **p <0.01; *p <0.05 versus HHZ (ANOVA, Dunnett’s multiple comparison test).

Full size image

Stomatal closure is one of the first responses to drought conditions that might control plant dehydration [27]. To investigate the stomatal status of the three rice genotypes at the tillering stage under drought stress, leaf surfaces were examined using scanning electron microscopy (SEM). H471 and P28 had significantly higher stomatal closure rate than HHZ after 1 day of drought stress (p <0.01, Kruskal-Wallis ANOVA, Dunn’s multiple comparison test) (Figure 2C). Relative electrolyte leakage (REL), which is an indicator of cell membrane injury, was also detected. The RELs of H471 and P28 were significantly lower than that of HHZ after 3 days of drought stress (Figure 2D), indicating that H471 and P28 experienced significantly less cell membrane injury than HHZ under drought stress conditions. Additionally, the activities of catalase (CAT) and ascorbate peroxidase (APX) in H471 and P28 were significantly higher than that in HHZ under drought stress conditions (Figure 2E and F), showing their active detoxification by reactive oxygen scavenging regulation in the two genotypes in response to drought [28]. Taken together, these results demonstrated that the DT of IL H471 was significantly improved over that of HHZ, reflecting the introgression of favorable alleles from DT donor P28 into the HHZ background.

Global gene expression profiling of three genotypes under well-watered and drought stress by transcriptome sequencing

At the tillering stage, total RNA from leaves of HHZ, P28, and H471 at 1 and 3 days of drought stress and its corresponding well-watered control were paired-end sequenced using Illumina sequencing technology. A total of 17.8–23.5 million reads of 100 bp in length were generated for each sample, and the number of mapped reads were in the range of 14.4–19.6 million (Table 1). The unique matching ratio was in the range of 77.7–79.7% (Table 1), the unique matching reads were used for further analysis. The high-quality reads from individual libraries were mapped to the rice genome; more than 22,774 mapped genes per library were determined simultaneously.

Table 1 Mapping results of RNA-Seq reads of HHZ, P28 and H471 under 1d and 3d drought stress and control conditions
Full size table

The detected expressed genes in all samples were subjected to cluster analysis. As shown in Additional file 1, the three genotypes under 1 and 3 days of drought stress and control conditions were separated from each other: H471 and HHZ were clustered together in each condition subgroup, which was consistent with their similar genetic background. Based on this result, the transcriptomic response to 1d and 3d drought stress could be classified as mild stress response (MR) and severe stress response (SR), respectively, in consistent with the previous reports on abiotic stress gene profiling [29, 30].

Correlation analysis was used to assess the quality of Illumina sequencing results between two replicates of each sample. The results indicated that the coefficient of correlation between the biological samples was high, supporting the reproducibility of the results (Additional file 2). To validate the Illumina sequencing results, quantitative real-time reverse transcription-PCR (qRT-PCR) was used to assess the expression levels for 27 genes of rice independently. The genes and primer sets used are shown in Additional file 3. RNA samples extracted from three additional replicate sets were used as templates. The high correlation (R2 = 0.93, p <0.01) between RNA-seq and qRT-RCR expression values indicated that there was a good agreement between both approaches (Additional files 4A and 4B).

Intrinsic transcriptome differences of H471 and P28 compared with HHZ prior to drought stress

Phenotypic differences resulting from gene expression variation have been observed in species [30, 31]. To investigate the intrinsic differences in gene expression between the DT genotypes and the DS genotype, the gene expression levels in H471 and P28 were compared with that in HHZ under control conditions. Between P28 and HHZ, 1282 genes were identified as differentially expressed, which is consistent with their large genetic difference. However, there were only 343 genes detected to be differentially expressed between H471 and HHZ under normal growth conditions (Additional file 5). This is consistent with only a few of chromosome fragments being introgressed from P28 to HHZ. The genes with higher basal expression level in H471 compared with HHZ were mainly functionally enriched in oxidoreductase activity, lyase activity, carboxylic acid metabolic process, response to stress, and cofactor binding (Additional file 6–1).

These constitutively differentially expressed genes in P28 and H471 compared with HHZ under control conditions could be classified into two groups based on their expression patterns under drought stress conditions. The first group comprised those that were unresponsive to drought stress, including 503 and 74 genes from P28 and H471, respectively, which are functionally enriched in redox regulation and apoptosis-associated proteins (Additional file 7–1 and 7–2). The second group comprised 779 and 269 genes differentially expressed in P28 and H471 after at least 1 or 3 days of drought stress (Additional file 7–3 and 7–4). Among them, 67 genes (Additional file 7–5) with higher basal expression in both P28 and H471 were functionally enriched in oxidoreductase and lyase activity. Further analysis indicated that 10 of these genes were colocalized in the introgressed regions, implying their positive role in the response to drought stress.

Comparative transcriptome profiling of three genotypes under drought stress

To determine the similarities and differences in drought-induced transcriptomes in the three genotypes, the gene expression alterations in the three genotypes under 1 and 3 days of drought stress compared with their respective controls were analyzed. In H471, P28, and HHZ, 7862, 7717, and 7625 DEGs were identified, respectively, after 1 or 3 days of drought stress.

After 1 day of drought stress, there were 5617, 5849, and 5579 DEGs detected in H471, P28, and HHZ, respectively. Venn diagram analysis indicated that 3945 DEGs (1966 up-regulated and 1979 down-regulated) were commonly regulated in the three genotypes by 1 day of drought stress (Table 2, Figure 3A and B) which accounted for 70.23%, 67.45% and 70.71% of total DEGs in H471, P28 and HHZ, respectively, indicating that the mild drought stress response of the different genotypes is mostly the same. GO analysis showed that the commonly up-regulated genes were functionally enriched in the regulation of the biological process, signaling process, and carbohydrate metabolic process (Additional file 6–2); while the shared down-regulated genes in the three genotypes were mainly associated with cellular protein metabolic process, phosphate metabolic process, and transport (Additional file 6–3).

Table 2 Summary of differentially expressed genes in H471, HHZ and P28 under 1d and 3d drought stress compared with its respective well-watered control
Full size table
Figure 3
figure 3

Venn diagram of differentially expressed genes (DEGs) in H471, HHZ, and P28 under 1 and 3 days of drought stress. DEGs were identified with adjusted p-value <0.001. (A), (B), (C), and (D) show the Venn diagram results for the three genotypes of 1 day up-regulated, 1 day down-regulated, 3 days up-regulated, and 3 days down-regulated genes, respectively, under drought stress.

Full size image

After 3 days of drought stress, 5290, 4957, and 5116 genes were identified as DEGs in H471, P28, and HHZ, respectively (Table 2), showing relatively fewer DEGs could be detected at the SR phase compared with the MR phase. Among them, 3308 DEGs (1484 up-regulated and 1824 down-regulated) were commonly regulated in three genotypes at the SR phase which accounted for 62.53%, 66.73% and 64.66% of total DEGs in H471, P28 and HHZ, respectively (Figure 3C and D). GO analysis of these DEGs indicated that carbohydrate metabolic process, regulation of cellular process, transcription, macromolecule modification, cellular protein metabolic process, and transport were highly enriched in the SR DEGs (Additional files 6–4 and 6–5). Taken together, the results showed that broad functional categories of genes were commonly involved in drought stress response.

Effect of introgression on the transcriptome of H471 in response to drought stress

The backcross introgression strategy is widely used for crop improvement. The introgressions combine the genetic background of the recurrent parent, which could result in novel changes in gene expression. To evaluate the effect of introgression on transcriptome of H471 under drought stress, the genome-wide gene expressions in H471 and HHZ under drought stress were compared. The results indicated that 460 and 380 genes were up- and down-regulated, respectively, in H471 compared with HHZ under drought stress (Additional file 8). Only a small proportion of the DEGs (128 up-regulated, 103 down-regulated) colocalized with the introgressed regions (Additional file 8), indicating introgression could contribute to new DT expression phenotypes in H471 relative to HHZ.

GO enrichment analysis of these DEGs highlighted the following functional categories: signaling transduction, transcription regulation, stress response, hormones (GA, JA) signal transduction, and ROS homeostasis (Additional file 9).

Genes related to signaling transduction: 18 genes encoding receptor kinases (RKs) were differentially expressed in H471 vs. HHZ under drought conditions. Among them, two genes encoding cysteine-rich receptor-like protein kinases (RLKs) (LOC_Os07g43560, LOC_Os07g43570) were up-regulated in H471 vs. HHZ and one was repressed (LOC_Os04g30040). Eight genes encoding leucine-rich repeat transmembrane protein kinases were up-regulated and five were down-regulated in H471 relative to HHZ. Additionally, 26 genes encoding kinase proteins were differentially expressed in H471 vs. HHZ under drought stress. These included four genes encoding calmodulin-dependent protein kinases (CDPKs) (LOC_Os03g20380, LOC_Os03g43440, LOC_Os05g26870, and LOC_Os10g39420). Three genes encoding EF hand calcium-binding proteins (OsCML15, OsCML18, and OsCML31) were up-regulated in H471, and two genes encoding CaLB domain proteins (LOC_Os08g20544, LOC_Os02g27130) were down-regulated in H471 compared with HHZ under drought stress, indicating their different roles in the drought stress response.

Phytohormone related proteins: 16 genes encoding ten JA and six GA related proteins were differentially expressed in H471 vs. HHZ under drought stress (Additional file 9). Kyoto encyclopedia of genes and genomes (KEGG) analysis of these DEGs showed that alpha-linolenic acid metabolism (ko00592), which is associated with JA synthesis, was highlighted in H471 (Additional file 10). Eight JA biosynthesis-related genes, encoding four lipoxygenase (LOX), one alpha-dioxygenase (DOX1), one acyl-CoA oxidase (ACX), one enoyl-CoA hydratase/3-hydroxyacyl-CoA dehydrogenase (MFP2) and one acetyl-CoA acyltransferase (fadA), were differentially expressed in H471 compared with HHZ under drought stress. Meanwhile, two genes (OsJAZ1 and OsJAZ7) involved in the JA signaling pathway were also significantly up-regulated in H471 vs. HHZ under drought stress. Consistently, the JA content in H471 was significantly higher than that in HHZ under drought stress, even though the JA contents in all three genotypes were repressed by drought compared with their respective control (Figure 4, Additional file 11). A GA biosynthesis-related gene, encoding putative gibberellin 20-oxidase 2 (GA20ox2) was significantly up-regulated after 1 day of drought stress, while a GA deactivation-related gene encoding gibberellin 2-oxidase was evidently up-regulated after 3 days of drought stress in all genotypes. Importantly, four putative gibberellin receptor encoded genes were detected to be differentially expressed in H471 compared with HHZ under drought stress. These results indicated that the JA and GA pathway is involved in drought stress tolerance in the DT IL H471.

Figure 4
figure 4

Jasmonic acid (JA) related pathways may play key roles in drought tolerance in H471. (A) JA synthesis and signaling pathways that contain differentially expression patterns between H471 and HHZ under drought conditions. Gene names in red and green color indicate up- and down-regulated in response to drought stress, respectively. (B) Hierarchical clustering of JA-related genes and some co-regulated genes. The color scale indicates the expression value. (C) JA contents in H471, HHZ, and P28 under drought conditions. Each column represents mean ± s.d. (three replicates); **p <0.01; *p <0.05 versus HHZ (ANOVA, Dunnett’s multiple comparison test).

Full size image

Transcription regulation related genes: 36 transcription factor (TF) genes were differentially expressed between H471 and HHZ under drought stress. Among them, three AP2/EREBP genes (LOC_Os01g64790, LOC_Os03g08470, LOC_Os04g57340) were significantly up-regulated in H471, which were previously identified to be highly involved in drought and salt stress [32]; OsWRKY4 (LOC_Os01g53040) and OsWRKY9 (LOC_Os01g18584) were highly expressed in H471 relative to HHZ under drought; these two genes were associated with biotic stress tolerance in a previous study [33, 34]. The others included four MYB TFs, five bHLH TFs, and seven NAC proteins, most of which were up-regulated in H471 compared with HHZ under drought, implying their positive role in the drought stress response.

Genes encoding redox regulation-related proteins: A set of 19 genes related to redox regulation were enriched in H471 compared with HHZ under drought stress. These genes included those encoding ascorbate peroxidase, oxidoreductase, peroxidase precursor, glutathione S-transferase, and glutathione synthetase (Additional file 9), showing that redox regulation is involved in the molecular mechanisms of drought stress tolerance.

Genes related to carbohydrate metabolism and osmotic adjustment: 9 genes encoding UDP-glucoronosyl and UDP-glucosyl transferase (UGT) family proteins and genes encoding three osmotins and one malate synthase were differentially expressed in H471 vs. HHZ under drought. These genes were previously reported to be involved in sugar metabolism and detoxification [35, 36] and osmotic adjustment [37].

Co-regulatory gene networks of H471 in response to drought stress

To explore the genetic networks associated with DT, all DEGs belonged to the above function categories in H471 compared with HHZ under drought stress (Additional file 9) were subjected to co-expression network analysis. Forty of them were found to be co-regulated and formed a complex network (Figure 5). In this network, the genes could be separated into four groups according to their putative functions. Group A was enriched in stress signaling transduction, including 13 genes encoding receptor kinases, protein kinases, and Ca2+ related protein. Four genes in group B were functionally involved in hormone signaling pathways: three genes encoding lipoxygenases (LOXs) and one gene encoding GA20ox2; LOX is involved in JA biosynthesis and signaling [38, 39]; and GA20ox2 is related to the biosynthesis of gibberellin [40]. The genes in group C were enriched in transcription regulation, including two NAC, one MYB, three CCT/B-box zinc finger proteins, and one HLH type TFs. The group D genes were mainly involved in ion transport, including two high affinity K+ transporter 5 s, three ABC transporters, and an ATP synthase subunit C family protein; and redox homeostasis, including peroxidase superfamily proteins, two glutathione S-transferases, and one short-chain dehydrogenase/reductase; and carbohydrate metabolism, including UDP-glycosyltransferase and malate synthase. These co-regulated genes in group D were evidently downstream in the drought stress tolerance mechanism of H471.

Figure 5
figure 5

Co-regulatory networks of genes differentially modulated in H471 under drought stress (false discovery rate 0.05). Co-regulation analysis was based on the calculation of pairwise Pearson correlation coefficient (PCC) of logarithmic expression values (with a cutoff of 0.75) in Rice Oligonucleotide Array Database. Four subgroups, indicated with letters from A–D which also distinguished by shapes, were identified as their putative functions. Different colors indicate genes with unique function class: olive, Ca2+ signal related; gold, receptor kinase; blue, kinase; violet, TF; green, K + transporter; Pink, redox regulation; black, osmotin; red, carbohydrate metabolism. The thickness of the edges is proportional to the PCC. Asterisks indicate the genes located in the introgression segment from P28. Arrows indicate up- and down-regulated in response to drought stress, respectively; and 1, 3 and B with arrows show 1 day, 3 days, both 1 and 3 days drought stress, respectively.

Full size image

Colocalization of DT related DEGs in the introgression fragments and QTL intervals

According to the results of genotyping by resequencing, 26 P28 fragments were introgressed into H471 (Figure 6). We analyzed the previously reported DT-related QTLs based on the Gramene QTL database (http://archive.gramene.org/qtl/) compared with the introgression intervals in H471. The results indicated that 39 known DT QTLs were co-localized in the introgressed regions in H471 (Figure 6).

Figure 6
figure 6

The recombination map of H471 and drought stress-related quantitative trait loci (QTLs) located in or near the introgressed regions. Names of QTLs, differentially expressed genes and introgressed bins are shown on the sides of the chromosome, more details are provided in Additional file 12.

Full size image

Correspondingly, all DEGs, including those differentially expressed in H471 compared with HHZ under control, 1 day or 3 days of drought stress conditions, were mapped onto the rice chromosome: 289 of these DEGs were localized on the introgressed intervals. Furthermore, 205 DEGs were localized onto 15 DT-related QTL intervals that overlap with the introgressed fragments in H471 (Additional file 12). Seven genes on the introgressed region of chromosome 1 were detected to be significantly up-regulated in H471 relative to HHZ after 3 days of drought stress: bHLH TF, lateral root primordium (LRP) protein, beta-1,3-glucanase 3, glycosyl hydrolase, glutathione S-transferase, and bidirectional amino acid transporter. Several drought-related QTLs, including those for leaf rolling, osmotic adjustment, and relative water content [41–43] were colocalized on this region; thus, up-regulation of these seven genes in H471 could be positively related to physiological and metabolic adaptation to drought stress. There was also a cluster of 10 genes located on the introgressed segment of chromosome 4 that had no corresponding DT QTLs. These genes included eight wall-associated kinase (WAK) family proteins and two cysteine-rich RLKs were found to be constitutively down-regulated in H471 compared with HHZ. It was evidenced that WAKs and cysteine-rich RLKs are essential for the normal regulation of cell enlargement and abiotic stress sensing [44, 45], reduced expression of these genes in H471 might negatively affect cell growth and the stress response.

Discussion

The DT IL showed enhanced drought stress tolerance

ILs have been widely used in genetic analysis and molecular breeding [26, 46]. Drought-selected ILs have provided a useful genetic resource to improve our understanding of the genetic and molecular basis of DT in rice [2]. In the present study, the DT IL H471 and donor P28 were more tolerant to drought stress than the recurrent parent HHZ. Under drought stress, H471 and P28 showed delayed leaf rolling, lower RWL, higher RWC, higher stomatal closure rate, relatively lower cell membrane injury, and with significantly higher CAT and APX activities than HHZ, indicating that H471 and P28 underwent complex morphological and physiological changes to resist drought stress [47]. These results showed that H471 and P28 used a typical DT strategy to cope with drought stress by limiting water loss and the enhanced DT of H471 over HHZ was caused by introgression of favorable alleles from donor P28.

Genome-wide transcriptome reprogramming of H471 and its parental lines in response to time-series drought stress

A comparative transcriptome analysis of three genotypes under 1- and 3-day drought stress indicated that the DT and DS genotypes shared a large proportion of drought-induced DEGs, revealing common drought stress-responsive processes. GO analyses showed that commonly up-regulated genes at the MR phase were highly enriched in signaling transduction and transcription regulation, consistent with previous reports that abiotic stress triggered the MR phase of transcriptome alterations of genes related to signaling cascades, including receptor kinases, transcription factors and components of calcium signaling [30, 48]. However, after 3 days of drought stress, the common DEGs were mainly involved in carbohydrate metabolic process, post-translational protein modifications, transport, and redox regulation, implying that downstream genes related to metabolite adaptation, ROS homeostasis, and post-transcriptional regulation were involved in the severe drought stress response, in accordance with previous results [48, 49].

Effect of introgression on the transcriptome in response to drought stress

The DEGs on the introgressed regions in H471 relative to its recurrent parent HHZ could have direct effects on its DT phenotype changes. Identification of those DEGs co-localized with identified DT-related QTLs provides useful data for DT molecular breeding and gene functional dissection. Meanwhile, it was reported that donor introgressions combine recurrent parental alleles that might result in novel changes in expression, and many genes showing non-parental expression patterns were identified outside of the introgressed fragments [50]. In this study, a set of genes in the DT introgression line H471 were exclusively differentially expressed compared with its parental lines, especially under drought stress, showing a unique expression related to the DT phenotype. These expression changes beyond the introgression regions might correspond to the well-known transgressive or nonparental expression in hybrid crop plants [51]. However, it was determined that the genome-wide expression changes in an IL might result from activation of transposon mPing[52] or transgressive siRNA [53]. The molecular mechanisms of this non-parental expression alteration, especially under drought stress, need to be further elucidated.

A complex genetic network including the JA and GA pathway is involved in drought stress tolerance

Drought stress tolerance is a complex trait and involves many genes. Hence, deciphering the molecular mechanisms underlying DT in plants is a challenging task. Genome-wide identification of drought-responsive regulons in contrasting DT genotypes with similar genetic backgrounds could help to dissect novel genetic components involved in DT. We comparatively analyzed the differential gene expression between H471 and HHZ under drought stress: 840 genes were differentially expressed, indicating substantial transcriptome reprogramming in H471 under drought stress compared with HHZ.

A set of genes related to signaling transduction including RKs, CDPKs, and CBPs were up-regulated in H471 compared with HHZ. RLKs play an important role in plant growth and responses to abiotic stresses by activating initial signaling transduction [45]; CDPKs and CBPs are the main components in the stress signaling pathway, which acts by modulating ABA signaling and reducing the accumulation of ROS [54, 55]. Up-regulation of these genes in H471 implied their positive role in DT by enhancing signaling pathways under drought stress.

Plant hormones, including ABA, GA and JA, play key roles in their ability to adapt to changing environments [15, 18, 19, 56]. GO and KEGG analyses revealed that genes related to JA biosynthesis and the signaling pathway were differentially expressed in H471 relative to HHZ under drought stress, implying that JA was involved in drought response. JA plays an important role in plant growth, development and stress response [57]. The JA signaling and biosynthesis genes were found to be significantly regulated by drought in Arabidopsis[58] and rice [59]. In this study, the JA levels in all three genotypes were repressed by drought compared with their control, this result indicated that JA is probably not required at high concentration under drought stress, and an increase in JA content might negatively affect plant growth as suggested by previous reports [58, 60]. However, the JA content was significantly higher in H471 than in HHZ and P28 under drought stress, which correlated with the increased expression of JA biosynthesis genes in H471 relative to HHZ, demonstrating that JA plays an important role in DT. Previous studies indicated that JA affected the transcript levels of genes related to antioxidants under water stress [61]; relative higher JA content in H471 might enhance drought tolerance by modulating antioxidant homeostasis.

It was shown that JA could interact synergistically and antagonistically with other phytohormones [60, 62]. In this study, differentially expressed GA-related genes in H471 compared with HHZ under drought stress provide evidence for the role of GA metabolism and the regulation of the GA signaling pathway on exposure to drought stress which consistent with latest study. Emerging evidence for interaction of the GA-signaling molecule DELLA with components of the signaling pathway for the stress hormone JA [19] suggests that GA signaling might integrate multiple hormone signaling pathways in response to drought stress, however the crossroads of these interactions still remain to be elucidated.

TFs are hub regulators in response to biotic and abiotic stresses. Several TFs, including AP2/EREBPs, WRKYs, bHLHs, and NACs, were highly up-regulated in H471 compared with HHZ under drought stress. All these TFs were previously reported to be involved in transcription regulation of abiotic stresses tolerance [20, 63]. Importantly, AP2/EREBPs, WRKYs, and NACs are key regulators of ABA-mediated stomatal closure and, hence, drought responses [64, 65]. bHLHs play an important role in the JA-mediated regulatory network of the abiotic stress response [66]. Differential expression of these TFs in H471 relative to HHZ implied that both ABA and JA play a central role in regulating drought stress tolerance in H471.Co-expression analysis showed that there were core DT genes forming a complex network in H471 (Figure 5), including genes related to signaling transduction, the JA pathway, and TFs, and those involved in downstream functions, such as carbohydrate metabolism, ion transport, and ROS regulation. In the network, a bHLH TF (LOC_Os01g11910) is connected to a receptor-like kinase (LOC_Os05g48660) and a hub 5-AMP-activated protein kinase (LOC_Os02g04330), interacts with glutathione S-transferase (LOC_Os01g72130), LZ-NBS-LRR protein (LOC_Os11g39020), and peroxidase (LOC_Os07g48040), and is then further connected with two JA biosynthesis lipoxygenases (LOC_Os03g49380 and LOC_Os12g37260). All these genes form a complex interacting network. However, their real interactions need to be confirmed by further experiments, including yeast-two-hybridization screening. These co-expressed gene networks provide useful information for dissecting the molecular mechanisms underlying drought stress tolerance in rice.

Conclusions

In this study, a DT IL and its parental lines, including the DT donor and DS recurrent parent, were used to characterize the differences of leaf transcriptome dynamics under 1 and 3 days of drought stress at the tillering stage using high-throughput RNA sequencing. Drought induced transcriptome reprogramming in a DT IL H471 could result from introgressed chromosome segments from the DT donor P28, and the differentially expressed genes in the H471 relative to the HHZ under drought stress might contribute to the enhanced drought tolerance, finally improved yield performance under stress. Co-expression analysis revealed a complex regulatory network, including genes related to the signaling transduction, JA and GA pathways, transcription regulation, redox control and osmotic adjustment, involved in drought stress tolerance. The data obtained in this study could extend our understanding of the molecular mechanisms of DT in rice.

Methods

Plant materials and experimental treatments

Three rice genotypes were used in this study. HHZ is a widely used indica inbred rice in South China, with high yield and good quality, which is sensitive to drought stress. P28 is an indica rice variety from the Philippines and is well adapted to drought stress. H471 is a BC1F5 drought-tolerant IL with a few chromosomal fragments introgressed from the donor parent P28 into background of the recurrent parent HHZ. Genome-wide single nucleotide polymorphism analysis by re-sequencing showed that H471 differs from HHZ at 26 genomic segments from P28, with sizes in the range of 57–6057 kb (Figure 6). The recombination bins were judged as described in Huang et al. [67], using the genome-resequencing data of HHZ, P28, and H471 (unpublished data).

To evaluate the DT performance of the three genotypes, a pot experiment was arranged in a randomized complete block design with three treatments (well-watered, 1-day drought stressed, and 3-day drought stressed) and six replications or pots in a greenhouse at the Institute of Crop Sciences of Chinese Academy of Agricultural Sciences (Beijing, China). The days were counted after the available water content (AWC) in the soil reached 20% to allow drought measurements at precisely determined intervals. Three healthy seedlings of the three genotypes were transplanted equidistantly into a strip pot (15 cm in height and 50 cm in diameter) filled with 2 kg of sterilized field soil, which contained about 50% AWC of the soil, as measured by soil moisture meters (TZS-W, Zhejiang Top Instrument Co. Ltd). Seedlings of each genotype were planted in six pots giving a total of six plants; all plants were grown with 14 h daylight at 28°C and a 10-h dark period at 25°C under controlled conditions. Withholding water at the tillering stage started the drought stress treatment. The soil water content reached 15%, 10%, and 7.5% after 1 day, 3 days, and 4 days of drought treatment, respectively.

Yield performance was evaluated in experiments under drought stress and well-watered conditions at the experimental farm of International Rice Research Institute (Philippines). Seedlings of three genotypes were transplanted into a three-row plot with 45 plants per plot at a spacing of 15 × 25 cm between rows and plants within each plot; three replications for each genotype. For drought stress treatment, water was drained and irrigation was held at the peak tillering stage until maturity. In the well-watered control, the field was irrigated at weekly intervals until 2 weeks before harvesting. Grain yield was measured on the 10 plants sampled at maturity from the middle row of each plot.

SEM analysis

The second leaves from 1-day drought-stressed plants grown in the pots were used for SEM analyses, as described in You et al. [68] with minor modifications. Fresh leaf samples were pre-fixed for 3 h in 3% glutaraldehyde-sodium phosphate buffer (0.1 M) at room temperature and rinsed three times with 0.1 M sodium phosphate buffer. Post-fixation was performed with 2% OsO4 at 4°C. The samples were dehydrated through an ethanol series and infiltrated with an isoamyl acetate series. The samples were then coated with metal particles for analysis by SEM to observe the guard cells. A Hitachi S750 scanning electron microscope (http://www.hitachi-hitec.com/global/em/) was used to take photographs, and the numbers of guard cells in randomly chosen fields were counted and analyzed statistically.

Physiological traits of the three genotypes under drought stress

Detached leaves are weighed and saturated with water for 24 h, then weighed again and dried for 48 h, and weighed again. RWC was calculated using the following formula: RWC (%) = [(FM - DM)/(TM - DM)] × 100, where FM, DM, and TM are the fresh, dry, and turgid masses of the tissue weighed, respectively. Monitoring the fresh weight loss at the indicated time points (per hour) measured the WLRs of detached leaves [68]. Measuring solute leakage from rice leaf tissue evaluated the REL, according to the method of Arora et al. [69], with minor modifications. Three replicates of 0.5 g fresh leaves were sampled from control and drought-treated plants. After being cut into 1-cm pieces, the 0.5 g leaf samples were immersed in 20 mL distilled water in a test tube for 1 h with the help of a vacuum pump. After standing for 2 h at 25°C, water conductivity was measured. Leaf discs were then killed in the same solution by autoclaving, and total conductivity was measured at room temperature. Percent injury arising from each treatment was calculated from the conductivity data using the equation: % injury = [(% L(t) - % L(c))/(100 - % L(c))] × 100), where % L (t) and % L(c) are percent conductivity for treated and control samples, respectively. Antioxidant enzyme activity, including catalase (CAT) and ascorbate peroxidase (APX), were determined following previously reported methods [70].

Enzyme-linked immuno-sorbent assay (ELISA) was performed to measure JA content. Three replicates of 200 ~ 500 mg fresh leaves were sampled from control and drought-stressed plants and immediately ground with 80% methanol under ice-bath conditions, endogenous JA content was measured according to the manufacture’s protocol (RB Plant-JA, Cat No.DRE-P10695, USA).

RNA extraction, RNA-seq library construction, and sequencing

Three top leaves for each sample (two replicates for each sample) were harvested for each genotype under 1 day and 3 days of drought stress and under well-watered control conditions. All samples were immediately frozen in liquid nitrogen and stored at -80°C. The TRIzol Reagent (Invitrogen, USA) was used to extract total RNA, which was quantified by a Qubit RNA assay kit (Applied Biosystems, CA, USA). An Agilent 2100 Bioanalyzer (Agilent Technologies) was used to check the RNA integrity. The ribosomal RNA (rRNA) was removed from 8 Î¼g of total RNA using a RiboMinusâ„¢ Plant Kit (Invitrogen), followed by a Ribo-Zero Gram-Negative Bacteria kit (Epicentre), according to the manufacturer’s instructions. The TruSeq RNA Sample Preparation kit (Illumina) was used to construct the paired-end fragment library, with minor modifications. Briefly, the rRNA removed RNA was fragmented and the first-strand cDNA was synthesized using random hexamers and SuperScript II Reverse Transcriptase. The RNA template was then removed and a replacement strand was synthesized to generate double-stranded (ds) cDNA. After end repair and 3′ end adenylation, an indexed adapter was ligated to the dsDNA. Fragments of 300–350 bp were excised and enriched by 12 cycles of PCR. The QUBIT and Agilent 2100 Bioanalyzer, assessed the yield and size distribution of the PCR products, respectively. CapitalBio Corporation, Beijing, China subjected the produced libraries to cluster generation on cBot and sequencing on a HiSeq 2000 platform (Illumina) with paired-end 100 base pair reads. The Illumina instrument software performed primary data analysis and base calling. Raw sequence data are available in the NCBI’s Gene Expression Omnibus (GEO) database under the accession number GSE57950.

Transcriptome data analysis

An in-house perl script was used to remove adaptor sequences and low-quality sequences from the raw reads. The retained high-quality pair-end reads of rice for each sample were mapped to the rice genome of RGAP at MSU [71] by TopHat [72] and then assembled using Cufflinks [73] to construct unique transcripts sequences, using the parameter: -g -b -u -o. Cuffcompare [73] was used to compare the assembled transcript fragments of each sample to the reference annotation, constructing a non-redundant transcripts data set among the samples. The number of mapped clean reads for each gene was counted and normalized into the reads per kilo base per million value [74]; Cuffdiff [73] was then used to identify DEGs. Finally, genes with a p-value ≤0.001 were designated as significantly differential expressed between each pair of samples.

Gene function annotations were performed based on the Rice Genome Annotation Project version 7 [71]. AgriGO was used to perform GO enrichment analysis [75]. The Kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment was performed using a hypergeometric test. The analysis incorporated false discovery rate correction using the Benjamini and Hochberg method of multiple hypotheses testing to reduce false negatives [76]. A weighted gene co-expression network analysis [77] was used to construct gene coexpression networks based on a Pearson correlation coefficient >0.75.

Quantitative real-time reverse transcription-PCR (qRT-PCR) analysis

To validate the results of the Illumina sequencing experiment, an independent set of samples (three biological replicates per sample) were collected as described for the RNA sequencing. Total RNA was treated with DNase I (TransGene, Beijing, China) to remove residual genome DNA and cDNA synthesis was performed using EasyScript First-Strand cDNA Synthesis SuperMix (TransGene, Beijing, China) according to the manufacturer’s protocol. qRT-PCR was used to verify a subset of DEGs, using the methods described by Swarbrick et al. [78]. The gene sequences were downloaded from the rice genome of RGAP at MSU [71]; Primer 3 software (http://frodo.wi.mit.edu/) was used to design the primers (Additional file 3). Twenty-seven rice genes with various function categories were selected and tested in 20-μL reactions using the SYBR® Green PCR Master Mix kit (Applied Biosystems, CA, USA), following the manufacturer’s protocol, via an ABI Prism 7900 Sequence Detection System (Applied Biosystems). The relative expression of each gene was calculated according to the method of 2-△△Ct[79]. The Actin 1 gene (LOC_Os03g50890) was used as endogenous references for qRT-PCR, and all analyses were performed with three technical and three biological replicates.

References

  1. Farooq M, Wahid A, Kobayashi N, Fujita D, Basra SMA: Plant drought stress: effects, mechanisms and management. Agron Sustain Dev. 2009, 29: 185-212. 10.1051/agro:2008021.

    Article  Google Scholar 

  2. Lafitte HR, Li ZK, Vijayakumara CHM, Gao YM, Shi Y, Xu JL, Fu BY, Yu SB, Ali AJ, Domingo J, Maghirang R, Torres R, Mackill D: Improvement of rice drought tolerance through backcross breeding: evaluation of donors and selection in drought nurseries. Field Crops Res. 2006, 97: 77-86. 10.1016/j.fcr.2005.08.017.

    Article  Google Scholar 

  3. Hao W, Lin HX: Toward understanding genetic mechanisms of complex traits in rice. J Genet Genomics. 2010, 37: 653-666. 10.1016/S1673-8527(09)60084-9.

    Article  CAS  PubMed  Google Scholar 

  4. Millan T, Clarke H, Siddique K, Buhariwalla H, Gaur P, Kumar J, Gil J, Kahl G, Winter P: Chickpea molecular breeding: new tools and concepts. Euphytica. 2006, 147: 81-103. 10.1007/s10681-006-4261-4.

    Article  Google Scholar 

  5. Price AH, Cairns JE, Horton P, Jones HG, Griffiths H: Linking drought-resistance mechanisms to drought avoidance in upland rice using a QTL approach: progress and new opportunities to integrate stomatal and mesophyll responses. J Exp Bot. 2002, 53: 989-1004. 10.1093/jexbot/53.371.989.

    Article  CAS  PubMed  Google Scholar 

  6. Krasensky J, Jonak C: Drought, salt, and temperature stress-induced metabolic rearrangements and regulatory networks. J Exp Bot. 2012, 63: 1593-1608. 10.1093/jxb/err460.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  7. Schroeder JI, Allen GJ, Hugouvieux V, Kwak JM, Waner D: Guard cell signal transduction. Annu Rev Plant Physiol Plant Mol Biol. 2001, 52: 627-658. 10.1146/annurev.arplant.52.1.627.

    Article  CAS  PubMed  Google Scholar 

  8. Zhu JK: Salt and drought stress signal transduction in plants. Ann Rev Plant Biol. 2002, 53: 247-273. 10.1146/annurev.arplant.53.091401.143329.

    Article  CAS  Google Scholar 

  9. Fan LM, Zhao Z, Assmann SM: Guard cells: a dynamic signaling model. Curr Opin Plant Biol. 2004, 7: 537-546. 10.1016/j.pbi.2004.07.009.

    Article  CAS  PubMed  Google Scholar 

  10. McAinsh MR, Clayton H, Mansfield TA, Hetherington AM: Changes in stomatal behavior and guard cell cytosolic free calcium in response to oxidative stress. Plant Physiol. 1996, 111: 1031-1042.

    CAS  PubMed Central  PubMed  Google Scholar 

  11. Pei ZM, Murata Y, Benning G, Thomine S, Klusener B, Allen GJ, Grill E, Schroeder JI: Calcium channels activated by hydrogen peroxide mediate abscisic acid signalling in guard cells. Nature. 2000, 406: 731-734. 10.1038/35021067.

    Article  CAS  PubMed  Google Scholar 

  12. Kohler B, Hills A, Blatt MR: Control of guard cell ion channels by hydrogen peroxide and abscisic acid indicates their action through alternate signaling pathways. Plant Physiol. 2003, 131: 385-388. 10.1104/pp.016014.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  13. Bright J, Desikan R, Hancock JT, Weir IS, Neill SJ: ABA induced NO generation and stomatal closure in Arabidopsis are dependent on H2O2 synthesis. Plant J. 2006, 45: 113-122. 10.1111/j.1365-313X.2005.02615.x.

    Article  CAS  PubMed  Google Scholar 

  14. Hamilton DW, Hills A, Kohler B, Blatt MR: Ca2+ channels at the plasma membrane of stomatal guard cells are activated by hyperpolarization and abscisic acid. Proc Natl Acad Sci U S A. 2000, 97: 4967-4972. 10.1073/pnas.080068897.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  15. Creelman RA, Mullet JE: Jasmonic acid distribution and action in plants: regulation during development and response to biotic and abiotic stress. Proc Natl Acad Sci U S A. 1995, 92: 4114-4119. 10.1073/pnas.92.10.4114.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  16. Akter N, Sobahan MA, Uraji M, Ye WX, Hossain MA, Mori IC, Nakamura Y, Murata Y: Effects of depletion of glutathione on abscisic acid and methyl jasmonate-induced stomatal closure in Arabidopsis thaliana. Biosci Biotechnol Biochem. 2012, 76: 2032-2037. 10.1271/bbb.120384.

    Article  CAS  PubMed  Google Scholar 

  17. Daszkowska-Golec A, Szarejko I: Open or close the gate – stomata action under the control of phytohormones in drought stress conditions. Front Plant Sci. 2013, doi:10.3389/fpls.2013.00138

    Google Scholar 

  18. Phelix CF, Feltus FA: Plant stress biomarkers from biosimulations: the Transcriptome-To-Metabolomeâ„¢ (TTMâ„¢) technology - effects of drought stress on rice. Plant Biol. 2014, doi:10.1111/plb.12221

    Google Scholar 

  19. Colebrook EH, Thomas SG, Phillips AL, Hedden P: The role of gibberellin signalling in plant responses to abiotic stress. J Exp Biol. 2014, 217: 67-75. 10.1242/jeb.089938.

    Article  CAS  PubMed  Google Scholar 

  20. Shinozaki K, Yamaguchi-Shinozaki K: Gene networks involved in drought stress response and tolerance. J Exp Bot. 2007, 58: 221-227.

    Article  CAS  PubMed  Google Scholar 

  21. Zhou J, Wang X, Jiao Y, Qin Y, Liu X, He K, Chen C, Ma L, Wang J, Xiong L, Zhang Q, Fan L, Deng XW: Global genome expression analysis of rice in response to drought and high-salinity stresses in shoot, flag leaf, and panicle. Plant Mol Biol. 2007, 63: 591-608. 10.1007/s11103-006-9111-1.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  22. Degenkolbe T, Do PT, Zuther E, Repsilber D, Walther D, Hincha DK, Kohl KI: Expression profiling of rice cultivars differing in their tolerance to long-term drought stress. Plant Mol Biol. 2009, 69: 133-153. 10.1007/s11103-008-9412-7.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  23. Wang D, Pan YJ, Zhao XQ, Zhu LH, Fu BY, Li ZK: Genome-wide temporal-spatial gene expression profiling of drought responsiveness in rice. BMC Genomics. 2011, 12: 149-10.1186/1471-2164-12-149.

    Article  PubMed Central  PubMed  Google Scholar 

  24. Moumeni A, Satoh K, Kondoh H, Asano T, Hosaka A, Venuprasad R, Serraj R, Kumar A, Leung H, Kikuchi S: Comparative analysis of root transcriptome profiles of two pairs of drought-tolerant and susceptible rice near-isogenic lines under different drought stress. BMC Plant Biol. 2011, 11: 174-10.1186/1471-2229-11-174.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  25. Gong Y, Rao LQ, Yu DQ: Abiotic Stress in Plants. Chemistry "Agricultural Chemistry". Chapter 6. Edited by: Stoytcheva M, Zlatev R. 2013, doi:10.5772/55865

    Google Scholar 

  26. Li ZC, Mu P, Li CP, Zhang HL, Li ZK, Gao YM, Wang XK: QTL mapping of root traits in a doubled haploid population from a cross between upland and lowland japonica rice in three environments. Theor Appl Genet. 2005, 110: 1244-1252. 10.1007/s00122-005-1958-z.

    Article  CAS  PubMed  Google Scholar 

  27. David TS, Henriques MO, Kurz-Besson C, Nunes J, Valente F, Vaz M, Pereira JS, Siegwolf R, Chaves MM, Gazarini LC, David JS: Water use strategies in two co-occurring Mediterranean evergreen oaks: surviving the summer drought. Tree Physiol. 2007, 27: 793-803. 10.1093/treephys/27.6.793.

    Article  CAS  PubMed  Google Scholar 

  28. Kadkhodaie A, Razmjoo J, Zahedi M: Peroxidase, Ascorbate Peroxidase and Catalase activities in Drought Sensitive, Intermediate and Resistance Sesame (Sesamum indicum L.) Genotypes. Int J Agron Plant Prod. 2013, 4: 3012-3021.

    Google Scholar 

  29. Yun KY, Park MR, Mohanty B, Herath V, Xu F, Mauleon R, Wijaya E, Bajic VB, Bruskiewich R, De los reyes BG: Transcriptional regulatory network triggered by oxidative signals configures the early response mechanisms of japonica rice to chilling stress. BMC Plant Biol. 2010, 10: 16-10.1186/1471-2229-10-16.

    Article  PubMed Central  PubMed  Google Scholar 

  30. Zhang T, Zhao XQ, Wang WS, Pan YJ, Huang LY, Liu XY, Zong Y, Zhu LH, Yang DC, Fu BY: Comparative transcriptome profiling of chilling stress responsiveness in two contrasting rice genotypes. PLoS One. 2012, 7: e43274-10.1371/journal.pone.0043274.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  31. Zhang X, Byrnes JK, Gal TS, Li WH, Borevitz JO: Whole genome transcriptome polymorphisms in Arabidopsis thaliana. Genome Biol. 2008, 9: R165-10.1186/gb-2008-9-11-r165.

    Article  PubMed Central  PubMed  Google Scholar 

  32. Jin P, Huang LY, Wang D, Wu HM, Zhu LH, Fu BY: Expression profiling of rice AP2/EREBP genes responsive to abiotic stresses. Sci Agr Sin. 2009, 42: 3765-3773.

    CAS  Google Scholar 

  33. Kim CY, Lee SH, Park HC, Bae CG, Cheong YH, Choi YJ, Han C, Lee SY, Lim CO, Cho MJ: Identification of rice blast fungal elicitor-responsive genes by differential display analysis. Mol Plant Microbe In. 2000, 13: 470-474. 10.1094/MPMI.2000.13.4.470.

    Article  CAS  Google Scholar 

  34. Zhou YL, Xu MR, Zhao MF, Xie XW, Zhu LH, Fu BY, Li ZK: Genome-wide gene responses in a transgenic rice line carrying the maize resistance gene Rxo1 to the rice bacterial streak pathogen, Xanthomonas oryzae pv. oryzicola. BMC Genomics. 2010, 11: 78-10.1186/1471-2164-11-78.

    Article  PubMed Central  PubMed  Google Scholar 

  35. Liu FX, Xu WY, Wei Q, Zhang ZH, Xing Z, Tan LB, Di C, Yao DX, Wang CC, Tan YJ, Yan H, Ling Y, Sun CQ, Xue YB, Su Z: Gene expression profiles deciphering rice phenotypic variation between Nipponbare (japonica) and 93–11 (indica) during oxidative stress. PLoS One. 2010, 5: e8632-10.1371/journal.pone.0008632.

    Article  PubMed Central  PubMed  Google Scholar 

  36. Ross J, Li Y, Lim E, Bowles DJ: Higher plant glycosyltransferases. Genome Biol. 2001, 2: reviews3004-

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  37. Goel D, Singh AK, Yadav V, Babbar SB, Bansal KC: Overexpression of osmotin gene confers tolerance to salt and drought stresses in transgenic tomato (Solanum lycopersicum L.). Protoplasma. 2010, 245: 133-141. 10.1007/s00709-010-0158-0.

    Article  CAS  PubMed  Google Scholar 

  38. Bell E, Mullet JE: Lipoxygenase gene expression is modulated in plants by water deficit, wounding, and methyl jasmonate. Mol Genet Genomics. 1991, 230: 456-462. 10.1007/BF00280303.

    Article  CAS  Google Scholar 

  39. Gardner HW: 9-Hydroxy-traumatin, a new metabolite of the lipoxygenase pathway. Lipids. 1998, 33: 745-749. 10.1007/s11745-998-0265-z.

    Article  CAS  PubMed  Google Scholar 

  40. Sasaki A, Ashikari M, Ueguchi-Tanaka M, Itoh H, Nishimura A, Swapan D, Ishiyama K, Saito T, Kobayashi M, Khush GS, Kitano H, Matsuoka M: Green revolution: a mutant gibberellin-synthesis gene in rice. Nature. 2002, 416: 701-702. 10.1038/416701a.

    Article  CAS  PubMed  Google Scholar 

  41. Courtois B, McLaren G, Sinha PK, Prasad K, Yadav R, Shen L: Mapping QTLs associated with drought avoidance in upland rice. Mol Breed. 2000, 6: 55-66. 10.1023/A:1009652326121.

    Article  CAS  Google Scholar 

  42. Lilley JM, Ludlow MM, McCouch SR, O’Toole JC: Locating QTL for osmotic adjustment and dehydration tolerance in rice. J Exp Bot. 1996, 47: 1427-1436. 10.1093/jxb/47.9.1427.

    Article  CAS  Google Scholar 

  43. Price AH, Townend J, Jones MP, Audebert A, Courtois B: Mapping QTLs associated with drought avoidance in upland rice grown in the Philippines and West Africa. Plant Mol Biol. 2002, 48: 683-695. 10.1023/A:1014805625790.

    Article  CAS  PubMed  Google Scholar 

  44. Wagner TA, Kohorn BD: Wall-associated kinases are expressed throughout plant development and are required for cell expansion. Plant Cell. 2001, 13: 303-318. 10.1105/tpc.13.2.303.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  45. Osakabe Y, Yamaguchi-Shinozaki K, Shinozaki K, Tran LS: Sensing the environment: key roles of membrane-localized kinases in plant perception and response to abiotic stress. J Exp Bot. 2013, 64: 445-458. 10.1093/jxb/ers354.

    Article  CAS  PubMed  Google Scholar 

  46. Zamir D: Improving plant breeding with exotic genetic libraries. Nat Rev Genet. 2001, 2: 983-989. 10.1038/35103590.

    Article  CAS  PubMed  Google Scholar 

  47. Xu LX, Han LB, Huang BR: Antioxidant enzyme activities and gene expression patterns in leaves of kentucky bluegrass in response to drought and post-drought recovery. J Amer Soc Hort Sci. 2011, 136: 247-255.

    CAS  Google Scholar 

  48. Caldana C, Degenkolbe T, Cuadros-Inostroza A, Klie S, Sulpice R, Leisse A, Steinhauser D, Fernie AR, Willmitzer L, Hannah MA: High-density kinetic analysis of the metabolomic and transcriptomic response of Arabidopsis to eight environmental conditions. Plant J. 2011, 67: 869-884. 10.1111/j.1365-313X.2011.04640.x.

    Article  CAS  PubMed  Google Scholar 

  49. Trzcinska-Danielewicz J, Bilsk A, Fronk J, Zielenkiewicz P, Jarochowsk E, Roszczyk M, Jończyk M, Axentowicz E, Skoneczny M, Sowiński P: Global analysis of gene expression in maize leaves treated with low temperature: I. Moderate chilling (14°C). Plant Sci. 2009, 177: 648-658. 10.1016/j.plantsci.2009.09.001.

    Article  CAS  Google Scholar 

  50. Koenig D, Jiménez-Gómez JM, Kimura S, Fulop D, Chitwood DH, Headland LR, Kumar R, Covington MF, Devisetty UK, Tat AV, Tohge T, Bolger A, Schneeberger K, Ossowski S, Lanz C, Xiong G, Taylor-Teeples M, Brady SM, Pauly M, Weigel D, Usadel B, Fernie AR, Peng J, Sinha NR, Maloof JN: Comparative transcriptomics reveals patterns of selection in domesticated and wild tomato. Proc Natl Acad Sci U S A. 2013, 110: 2655-2662. 10.1073/pnas.1309606110.

    Article  Google Scholar 

  51. Bell GD, Kane NC, Rieseberg LH, Adams KL: RNA-seq analysis of allele-specific expression, hybrid effects, and regulatory divergence in hybrids compared with their parents from natural populations. Genome Biol Evol. 2013, 5: 1309-1323. 10.1093/gbe/evt072.

    Article  PubMed Central  PubMed  Google Scholar 

  52. Yasuda K, Tsukiyama T, Karki S, Okumoto Y, Teraishi M, Saito H, Tanisaka T: Mobilization of the active transposon mPing in interspecific hybrid rice between Oryza sativa and O. glaberrima. Euphytica. 2013, 192: 17-24. 10.1007/s10681-012-0810-1.

    Article  CAS  Google Scholar 

  53. Shivaprasad PV, Dunn RM, Santos BA, Bassett A, Baulcombe DC: Extraordinary transgressive phenotypes of hybrid tomato are influenced by epigenetics and small silencing RNAs. EMBO J. 2012, 31: 257-266. 10.1038/emboj.2011.458.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  54. Reddy AS, Ali GS, Celesnik H, Day IS: Coping with stresses: roles of calcium- and calcium/calmodulin-regulated gene expression. Plant Cell. 2011, 23: 2010-2032. 10.1105/tpc.111.084988.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  55. Asano T, Hayashi N, Kikuchi S, Ohsugi R: CDPK-mediated abiotic stress signaling. Plant Signal Behav. 2012, 7: 817-821. 10.4161/psb.20351.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  56. Peleg Z, Blumwald E: Hormone balance and abiotic stress tolerance in crop plants. Curr Opin Plant Biol. 2011, 14: 290-295. 10.1016/j.pbi.2011.02.001.

    Article  CAS  PubMed  Google Scholar 

  57. Cheong JJ, Choi YD: Methyl jasmonate as a vital substance in plants. Trends Genet. 2003, 19: 409-413. 10.1016/S0168-9525(03)00138-0.

    Article  CAS  PubMed  Google Scholar 

  58. Harb A, Krishnan A, Ambavaram MM, Pereira A: Molecular and physiological analysis of drought stress in Arabidopsis reveals early responses leading to acclimation in plant growth. Plant Physiol. 2010, 154: 1254-1271. 10.1104/pp.110.161752.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  59. Du H, Liu H, Xiong L: Endogenous auxin and jasmonic acid levels are differentially modulated by abiotic stresses in rice. Front Plant Sci. 2013, 4: 397-

    Article  PubMed Central  PubMed  Google Scholar 

  60. Aimar D, Calafat M, Andrade AM, Carassay L, Abdala GI, Molas ML: Drought tolerance and stress hormones: From model organisms to forage crops. Plants and Environment. Edited by: Dr. Hemanth Vasanthaiah. 2011, ISBN: 978-953-307-779-6, InTech, Available from: [http://www.intechopen.com/books/plants-and-environment/drought-tolerance-and-stress-hormones-from-model-organisms-to-forage-crops],

    Google Scholar 

  61. Shan CJ, Liang ZS: Jasmonic acid regulates ascorbate and glutathione metabolism in Agropyron cristatum leaves under water stress. Plant Sci. 2010, 178: 130-139. 10.1016/j.plantsci.2009.11.002.

    Article  CAS  Google Scholar 

  62. Wasternack C, Hause B: Jasmonates: biosynthesis, perception, signal transduction and action in plant stress response, growth and development. Ann Bot-london. 2013, 111: 1021-1058. 10.1093/aob/mct067.

    Article  CAS  Google Scholar 

  63. Lindemose S, O’Shea C, Jensen MK, Skriver K: Structure, function and networks of transcription factors involved in abiotic stress responses. Int J Mol Sci. 2013, 14: 5842-5878. 10.3390/ijms14035842.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  64. Rushton DL, Tripathi P, Rabara RC, Lin J, Ringler P, Boken AK, Langum TJ, Smidt L, Boomsma DD, Emme NJ, Chen X, Finer JJ, Shen QJ, Rushton PJ: WRKY transcription factors: key components in abscisic acid signalling. Plant Biotechnol J. 2012, 10: 2-11. 10.1111/j.1467-7652.2011.00634.x.

    Article  CAS  PubMed  Google Scholar 

  65. Nakashima K, Takasaki H, Mizoi J, Shinozaki K, Yamaguchi-Shinozaki K: NAC transcription factors in plant abiotic stress responses. Biochim Biophys Acta. 1819, 2012: 97-103.

    Google Scholar 

  66. Niu Y, Figueroa P, Browse J: Characterization of JAZ-interacting bHLH transcription factors that regulate jasmonate responses in Arabidopsis. J Exp Bot. 2011, 62: 2143-2154. 10.1093/jxb/erq408.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  67. Huang XH, Feng Q, Qian Q, Zhao Q, Wang L, Wang A, Guan JP, Fan DL, Weng QJ, Huang T, Dong GJ, Sang T, Han B: High-throughput genotyping by whole-genome resequencing. Genome Res. 2009, 19: 1068-1076. 10.1101/gr.089516.108.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  68. You J, Zong W, Li XK, Ning J, Hu HH, Li XH, Xiao JH, Xiong LZ: The SNAC1-targeted gene OsSRO1c modulates stomatal closure and oxidative stress tolerance by regulating hydrogen peroxide in rice. J Exp Bot. 2013, 64: 569-583. 10.1093/jxb/ers349.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  69. Arora R, Pitchay DS, Bearce BC: Water-stress-induced heat tolerance in geranium leaf tissues: a possible linkage through stress proteins?. Physiol Plant. 1998, 103: 24-34. 10.1034/j.1399-3054.1998.1030104.x.

    Article  CAS  Google Scholar 

  70. Bonnecarre’re V, Borsani O, Dı’az P, Capdevielle F, Blanco P, Monza J: Response to photoxidative stress induced by cold in japonica rice is genotype dependent. Plant Sci. 2011, 180: 726-732. 10.1016/j.plantsci.2011.01.023.

    Article  Google Scholar 

  71. Kawahara Y, de la Bastide M, Hamilton JP, Kanamori H, McCombie WR, Ouyang S, Schwartz DC, Tanaka T, Wu J, Zhou S, Childs KL, Davidson RM, Lin H, Quesada-Ocampo L, Vaillancourt B, Sakai H, Lee SS, Kim J, Numa H, Itoh T, Buell CR, Matsumoto T: Improvement of the Oryza sativa Nipponbare reference genome using next generation sequence and optical map data. Rice. 2013, 6: 4-10.1186/1939-8433-6-4.

    Article  PubMed  Google Scholar 

  72. Trapnell C, Pachter L, Salzberg SL: TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009, 25: 1105-1111. 10.1093/bioinformatics/btp120.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  73. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L: Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010, 28: 511-515. 10.1038/nbt.1621.

    Article  CAS  PubMed Central  PubMed  Google Scholar 

  74. Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B: Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat Methods. 2008, 5: 621-628. 10.1038/nmeth.1226.

    Article  CAS  PubMed  Google Scholar 

  75. Du Z, Zhou X, Ling Y, Zhang Z, Su Z: agriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010, 38: 64-70. 10.1093/nar/gkq310.

    Article  Google Scholar 

  76. Benjamini Y, Drai D, Elmer G, Kafkafi N, Golani I: Controlling the false discovery rate in behavior genetics research. Behav Brain Res. 2001, 125: 279-284. 10.1016/S0166-4328(01)00297-2.

    Article  CAS  PubMed  Google Scholar 

  77. Cao P, Jung KH, Choi D, Hwang D, Zhu J, Ronald PC: The rice oligonucleotide array database: an atlas of rice gene expression. Rice. 2012, 5: 17-10.1186/1939-8433-5-17.

    Article  PubMed  Google Scholar 

  78. Swarbrick PJ, Huang K, Liu G, Slate J, Press MC, Scholes JD: Global patterns of gene expression in rice varieties undergoing a susceptible or resistant interaction with the parasitic Strigahermonthica. New Phytol. 2008, 179: 515-529. 10.1111/j.1469-8137.2008.02484.x.

    Article  CAS  PubMed  Google Scholar 

  79. Livak KJ, Schmittgen TD: Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001, 25: 402-408. 10.1006/meth.2001.1262.

    Article  CAS  PubMed  Google Scholar 

Download references

Acknowledgements

This work was supported by the National 863 Project (2014AA106A603), the Bill & Melinda Gates Foundation Project (OPP51587), and the Shenzhen Peacock Plan.

Author information

Author notes

    Authors and Affiliations

    1. Institute of Crop Sciences/National Key Facility for Crop Gene Resources and Genetic Improvement, Chinese Academy of Agricultural Sciences, South Zhong-Guan-Cun Street 12#, Beijing, 100081, China

      Liyu Huang, Fan Zhang, Fan Zhang, Wensheng Wang, Yongli Zhou, Binying Fu & Zhikang Li

    2. Graduate School of Chinese Academy of Agricultural Sciences, South Zhong-Guan-Cun Street 12#, Beijing, 100081, China

      Fan Zhang

    Authors
    1. Liyu Huang
      View author publications

      You can also search for this author in PubMed Google Scholar

    2. Fan Zhang
      View author publications

      You can also search for this author in PubMed Google Scholar

    3. Fan Zhang
      View author publications

      You can also search for this author in PubMed Google Scholar

    4. Wensheng Wang
      View author publications

      You can also search for this author in PubMed Google Scholar

    5. Yongli Zhou
      View author publications

      You can also search for this author in PubMed Google Scholar

    6. Binying Fu
      View author publications

      You can also search for this author in PubMed Google Scholar

    7. Zhikang Li
      View author publications

      You can also search for this author in PubMed Google Scholar

    Corresponding authors

    Correspondence to Yongli Zhou or Binying Fu.

    Additional information

    Competing interests

    The authors declare that they have no competing interests.

    Authors’ contributions

    BF, YZ and ZL designed the experiments and drafted the manuscript. LH, FZ, FZ, and WW conducted the phenotype and RNA-seq experiments. FZ and WW performed the microarray data analysis, and FZ carried out the qRT-PCR confirmation of candidate genes. BF revised the manuscript. All authors read and approved the final manuscript.

    Liyu Huang, Fan Zhang contributed equally to this work.

    Electronic supplementary material

    12864_2014_6721_MOESM1_ESM.ppt

    Additional file 1: Hierarchical cluster analysis of nine sample pools (columns) and all expressed genes, under control and drought stress conditions (rows). A PowerPoint file containing hierarchical cluster analysis of nine sample pools (columns) and all expressed genes, under control and drought stress conditions (rows). The raw data represented here can be obtained from GEO: GSE57950. In the colored panels, each horizontal line represents a single gene and the colored line indicates the expression level (in a log scale) of the gene relative to the median in a specific sample: high expression in red, low expression in green. (PPT 182 KB)

    12864_2014_6721_MOESM2_ESM.ppt

    Additional file 2: Correlation analysis of Illumina sequencing results between two replicates of each sample for H471, P28, and HHZ, under control (ck), 1 day and 3 days of drought stress, respectively. A PowerPoint file containing correlation analysis of Illumina sequencing results between two replicates of each sample for H471, P28, and HHZ, under control (ck), 1 day and 3 days of drought stress, respectively. (PPT 295 KB)

    12864_2014_6721_MOESM3_ESM.xls

    Additional file 3: Information of primers used in qRT-PCR analysis. Excel file containing information of primers used in qRT-PCR analysis. (XLS 22 KB)

    12864_2014_6721_MOESM4_ESM.ppt

    Additional file 4: Comparison of transcription measurements by Illumina sequencing and quantitative real-time reverse transcription-PCR (qRT-PCR) assays. A PowerPoint file containing comparison of transcription measurements by Illumina sequencing and quantitative real-time reverse transcription-PCR (qRT-PCR) assays. (A) The correlation coefficient (R2) between the two datasets is 0.93. (B) Comparative analysis of six candidate genes expression level by qRT-PCR and RNA-seq. qRT-PCR quantification values were compared with HHZ_ck. Error bars indicate the standard deviation. Actin 1 was used as an endogenous control. (PPT 800 KB)

    12864_2014_6721_MOESM5_ESM.xls

    Additional file 5: DEGs in P28 and H471 compared with HHZ under a well-watered condition. Excel file containing a list of DEGs in P28 and H471 compared with HHZ under a well-watered condition. (XLS 338 KB)

    12864_2014_6721_MOESM6_ESM.xls

    Additional file 6: GO enrichment analysis of DEGs in various groups. Excel file containing five lists of GO enrichment of DEGs in various groups: 6–1, GO enrichment analysis of genes with the higher basal expression in H471 compared with HHZ; 6–2, GO enrichment analysis of the commonly up-regulated genes in three genotypes under 1 day of drought stress; 6–3, GO enrichment analysis of the commonly down-regulated genes under 1 day of drought stress; 6–4, GO enrichment analysis of the commonly up-regulated genes in three genotypes under 3 days of drought stress; 6–5, GO enrichment of the commonly down-regulated genes in three genotypes under 3 days of drought stress (XLS 47 kb). (XLS 47 KB)

    12864_2014_6721_MOESM7_ESM.xls

    Additional file 7: Genes list of higher basal expression in P28 and H471 compared with HHZ. Description: Excel file containing five list of genes with higher basal expression in P28 and H471 compared with HHZ. (XLS 337 KB)

    12864_2014_6721_MOESM8_ESM.xls

    Additional file 8: DEGs list in H471 compared with HHZ under drought stress. Excel file containing a list of DEGs in H471 compared with HHZ under drought stress. (XLS 158 KB)

    12864_2014_6721_MOESM9_ESM.xls

    Additional file 9: List of DEGs with various function categories in H471 compared with HHZ under drought conditions. Excel file containing a list of DEGs with various function categories in H471 compared with HHZ under drought conditions. (XLS 57 KB)

    12864_2014_6721_MOESM10_ESM.xls

    Additional file 10: List of alpha-linolenic acid metabolism-related genes in H471 vs. HHZ DEGs. Excel file containing a list of alpha-linolenic acid metabolism-related genes in H471 vs. HHZ DEGs. (XLS 51 KB)

    12864_2014_6721_MOESM11_ESM.ppt

    Additional file 11: Validation of JA related genes by qRT-PCR. A PowerPoint file containing comparison validation of eight JA related genes by qRT-PCR assay. qRT-PCR Quantification values were compared with HHZ_ck. Error bars indicate the standard deviation. Actin 1 was used as an endogenous control. (PPT 924 KB)

    12864_2014_6721_MOESM12_ESM.xls

    Additional file 12: Information on introgressed chromosome segments and the differentially expressed genes in H471 compared with HHZ under control and drought stress conditions and their corresponding DT-related QTLs. Excel file containing information on introgressed chromosome segments and the differentially expressed genes in H471 compared with HHZ under control and drought stress conditions and their corresponding DT-related QTLs. (XLS 61 KB)

    Authors’ original submitted files for images

    Rights and permissions

    Open Access  This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made.

    The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.

    To view a copy of this licence, visit https://creativecommons.org/licenses/by/4.0/.

    The Creative Commons Public Domain Dedication waiver (https://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

    Reprints and permissions

    About this article

    Check for updates. Verify currency and authenticity via CrossMark

    Cite this article

    Huang, L., Zhang, F., Zhang, F. et al. Comparative transcriptome sequencing of tolerant rice introgression line and its parents in response to drought stress. BMC Genomics 15, 1026 (2014). https://doi.org/10.1186/1471-2164-15-1026

    Download citation

    • Received:

    • Accepted:

    • Published:

    • DOI: https://doi.org/10.1186/1471-2164-15-1026

    Keywords

    BMC Genomics

    ISSN: 1471-2164

    Contact us