Genome-wide gene expression profiling of introgressed indica rice alleles associated with seedling cold tolerance improvement in a japonica rice background
BMC Genomics volume 13, Article number: 461 (2012)
Rice in tropical and sub-tropical areas is often subjected to cold stress at the seedling stage, resulting in poor growth and yield loss. Although japonica rice is generally more cold tolerant (CT) than indica rice, there are several favorable alleles for CT exist in indica that can be used to enhance CT in rice with a japonica background. Genome-wide gene expression profiling is an efficient way to decipher the molecular genetic mechanisms of CT enhancement and to provide valuable information for CT improvement in rice molecular breeding. In this study, the transcriptome of the CT introgression line (IL) K354 and its recurrent parent C418 under cold stress were comparatively analyzed to explore the possible CT enhancement mechanisms of K354.
A total of 3184 differentially expressed genes (DEGs), including 195 transcription factors, were identified in both lines under cold stress. About half of these DEGs were commonly regulated and involved in major cold responsive pathways associated with OsDREB1 and OsMyb4 regulons. K354-specific cold-induced genes were functionally related to stimulus response, cellular cell wall organization, and microtubule-based movement processes that may contribute to increase CT. A set of genes encoding membrane fluidity and defensive proteins were highly enriched only in K354, suggesting that they contribute to the inherent CT of K354. Candidate gene prediction based on introgressed regions in K354 revealed genotype-dependent CT enhancement mechanisms, associated with Sir2, OsFAD7, OsWAK112d, and programmed cell death (PCD) related genes, present in CT IL K354 but absent in its recurrent parent C418. In K354, a number of DEGs were co-localized onto introgressed segments associated with CT QTLs, providing a basis for gene cloning and elucidation of molecular mechanisms responsible for CT in rice.
Genome-wide gene expression analysis revealed that genotype-specific cold induced genes and genes with higher basal expression in the CT genotype contribute jointly to CT improvement. The molecular genetic pathways of cold stress tolerance uncovered in this study, as well as the DEGs co-localized with CT-related QTLs, will serve as useful resources for further functional dissection of the molecular mechanisms of cold stress response in rice.
Rice (Oryza sativa L.) is an important cereal species, and a staple food for more than half of the world’s population. During early growth stages, low temperature stress retards rice seedling establishment and plant development, directly impacting yield. Temperate japonica cultivars generally have better seedling-stage cold tolerance (CT) and greater seedling vigor than indica cultivars . By using backcross breeding, however, we have recently identified the putative genetic networks underlying some introgression lines (ILs) with significantly improved CT at the seedling stage under a japonica restorer C418 background (Zhang et al., unpublished data), suggesting the diversity hidden in the indica germplasm. Further elucidation of the molecular basis of CT enhancement of japonica rice arising from introgression of favorable alleles from indica rice would provide rice breeders with much useful information.
Cold stress tolerance mechanisms in plants include cold signal perception, activation of transcription factors (TFs) by signal transduction, and expression of cold-responsive genes for mediating stress tolerance . The plant cell membrane represents an ideal location for the primary temperature stress sensor, because one of the immediate effects of temperature stress in plants is a change in membrane fluidity . At low temperatures, greater membrane lipid unsaturation appears to be crucial to optimum membrane function. When the Arabidopsis fatty acid deficient triple mutant fad3-2 fad7-2 fad8, which is devoid of trienoic fatty acids, was exposed to prolonged low temperatures, a decrease in fluorescence parameters, chlorophyll content, and thylakoid membrane content was observed compared with wild-type plants . Similarly, OsFAD8 has a functional role in maintaining levels of trienoic fatty acids and stress tolerance at low temperatures in rice . In one study , linolenic acid levels increased and palmitic acid levels decreased in cold tolerant rice genotypes exposed to chilling; the opposite behavior was observed in cold sensitive genotypes. Transcriptional regulation is the core component of the complex genetic network associated with plant responses to cold. The C-repeat (CRT)/dehydration responsive element (DRE)-binding TF-mediated cold response pathway has been shown to play a predominant role in CT through the process of cold acclimation . Although rice, a plant of tropical and subtropical origin, lacks mechanisms for cold acclimation, it nonetheless possesses components of this CBF cold-response pathway [8, 9]. OsDREB1A and OsDREB1B are induced by cold stress, and constitutive overexpression of these genes in transgenic Arabidopsis and rice leads to induction of stress-responsive genes, increased high cold and salt tolerance, and growth retardation under normal conditions [8, 9]. Other important signaling pathways have also been shown to be involved in CT. For example, MYBS3 represses the DREB1-dependent cold signaling pathway at the transcriptional level. DREB1- and MYBS3-dependent pathways may complement each other and act sequentially to allow adaptation to immediate and persistent cold stress in rice . The rice R2R3-type OsMyb4 TF controls a hierarchical network comprising several regulatory sub-clusters associated with cellular defense and rescue, metabolism, and development. This network is independent of CBF/DREB, and its sub-regulons operate with possible co-regulators, including nuclear factor-Y . In addition, one ROS (ROS-bZIP)-mediated regulon triggered by chilling stress is independent of Abscisic acid (ABA) and CBF/DREB, and its activation promotes rapid response of rice seedlings to chilling stress . Finally, constitutive and non-cold responsive regulons, which have a differential effect on the cold responsive DRE regulon, also play a key role in CT .
DNA microarray analysis is a well developed high-throughput technology that has been used for many genomic application areas, especially whole-genome gene expression profiling. Although high-throughput RNA sequencing has recently become popular as an alternative to microarray analysis, the microarray platform, with its robust sample process and data analysis pipelines, is still the preferred choice for transcriptomic profiling involving a large number of samples in model plants with well-annotated genomes [14, 15]. Many microarray-based studies have been carried out to identify abiotic stress responsive genes in specific rice varieties and transgenic rice [9, 10, 16, 17], and comparative transcriptional profiling of two contrasting rice genotypes under salinity and drought stress have revealed novel genetic regulatory mechanisms involved in stress tolerance [18, 19]. Unfortunately, it is difficult to use the stress-related genes uncovered in those studies to improve modern varieties in rice breeding, because most of them already exist in elite rice varieties.
In this study, a CT IL and its cold sensitive recurrent parent were assessed in terms of their seedling growth and physiological responses to cold stress treatment. An Affymetrix genome array was used to profile global gene expression changes under a cold stress time series. The results of this study should serve as an initial step in a comprehensive understanding of CT enhancement resulting from introgression of favorable alleles from indica rice into a japonica rice background.
Results and Discussion
Phenotypic differences between the introgression line K354 and its recurrent parent C418 under cold stress
The CT IL K354 and its recurrent parent C418, which possesses a cold sensitive phenotype, were used in this study. Under cold stress, a marked difference in survival rate was observed between K354 and C418 (Additional file 1), with K354 showing better CT and recovery ability than C418. Seedling leaves of C418 exhibited obvious wilting, while only a few necroses on leaf tips were observed in K354 seedlings (Additional file 2A). Compared with C418, cold-treated K354 seedlings experienced significantly less extensive cell membrane injury (relative electrolyte leakage) during the 24 h recovery period following the 48 h cold stress treatment (Additional file 2B). Interestingly, although under control conditions K354 exhibited significantly lower superoxide dismutase (SOD) activity and soluble protein concentrations than C418, after 48 h of cold stress these parameters increased dramatically and reached higher values, even with different changing trends, in K354 compared with C418 (Additional file 2C and 2E). In regard to catalase (CAT) activity, K354 showed a larger change than C418 after 48 h of cold stress (Additional file 2D). Based on phenotypic and physiological traits measured in this study, the CT of K354 was indeed improved over that of C418, due to the introgression of favorable alleles from indica donor Bg300 into the japonica background.
Genome-wide gene expression profiling of K354 and C418 under a cold stress time-series
Comparison of genome-wide gene expression between two genotypes with contrasting CT phenotypes in a cold environment may help identify underlying molecular genetic pathways responsible for cold adaptation of rice plants. In this study, we use an Affymetrix whole rice genome array to profile gene expression changes of K354 and C418 under cold stress at five time points (2, 6, 12, 24, and 48 h). The array contained 49,824 rice genome genes (48,564 for japonica and 1260 for indica). Among these genes, 24,888–27,407 (43.4–47.8%) were found to be expressed in at least one sample in both genotypes under stress or control conditions (Additional file 3). A total of 3184 (5.6%) differentially expressed genes (DEGs) were identified, including comprise 2562 cold-responsive DEGs and 710 genotypic-specific DEGs (Table 1). There were 232/83, 286/243, 503/471, 1075/1092, and 1962/1444 DEGs detected in samples of C418/K354 at 2, 6, 12, 24, and 48 h, respectively, revealing that cold treatment induced a continuous increase in DEGs in both genotypes. There was an apparently delayed response to cold in K354 compared with C418, with 232 DEGs in C418 vs 83 in K354 at 2 h cold stress. Similar results were previously observed in rice under salinity stress; in that study, the sensitive genotype IR29 induced a relatively large number of genes compared to tolerant FL478 . In addition, comparison of gene expression levels between K354 and C418 showed that there were 308 up-regulated and 209 down-regulated genes in K354 compared with C418 under control conditions. This suggests that these genes may be responsible for K354’s intrinsic tolerance to cold stress, which may be largely derived from introgression of donor alleles conferring constitutive tolerance to low temperature.
Gene Ontology (GO) enrichment of all detected cold-responsive DEGs with known and putative functions is shown in Figure 1. The predominant DEGs were functionally involved in metabolism (23.5%, p = 5.85E-33), such as macromolecule metabolism, oxidation-reduction, and primary metabolism; cellular processes (18.6%, p = 6.71E-11), including microtubule-based processes; stimulus response (4.9%, p = 1.04E-8); binding activity (31.0%, p = 5.26E-43), including ion, tetrapyrrole, and nucleotide binding; and catalytic activity (27.4%, p = 1.76E-43), such as hydrolase, transferase, and oxidoreductase activity. GO enrichment analysis of the 710 genotype-specific DEGs, however, indicated that the most prevalent DEGs were involved in programmed cell death (PCD) (3.2%, p = 4.68E-7) and electron carrier activity (3.0%, p = 0.0002), revealing genotype-specific transcriptome changes in response to cold environmental stimuli.
On the basis of hierarchical average linkage cluster analysis (Figure 2), a transcript-level time-dependent cold response could be clearly delineated into early response (phase I, 2–6 h), middle response (phase II, 12 h), and late response (phase III, 24–48 h) components. The samples from K354 and C418 formed separate and distinct clusters showing their genotype-specific responses to cold stress treatment. All 3184 DEGs could be classified into six major groups based on their transcription patterns (Figure 2 and Additional file 4 and Additional file 5). Two major DEG groups shared by both genotypes were either induced (cluster V) and repressed (cluster I) by cold. Clusters II and III represented non-cold-responsive DEGs that were expressed at significantly different levels between C418 and K354 under normal growth conditions. Members of Cluster IV were initially down-regulated but later induced in both genotypes during cold stress. Interestingly, Cluster VI comprised a set of genes that were initially induced and later repressed in K354 during stress (Additional file 4 and Additional file 5). Genes in Clusters II, III, and VI may therefore be mainly responsible for the differential cold stress response between two genotypes, as discussed below.
Using the sample cluster results (Figure 2), 15 DEGs were randomly selected for qRT-PCR analysis to confirm the gene expression changes detected by microarray analysis. The high correlation (R2 = 0.75, P < 0.05) between microarray data and qRT-RCR expression values indicate that there was good agreement between both approaches (Additional file 6).
Time-dependent K354 and C418 gene expression changes in response to cold stress
To assess the differential transcriptome alterations of K354 and C418 under cold stress, DEGs at different time points were compared by Venn diagram (Figure 3). Figure 3A represents the results of comparing induced or repressed genes in C418 under cold stress during phases I, II, and III. A total of 229 (14.8%) and 41 (7.1%) genes were continuously induced and repressed, respectively, in C418 during cold stress. Moreover, 32 (2.1%), 9 (0.6%), and 1095 (70.7%) genes were induced, and 15 (2.6%), 5 (0.9%), and 454 (78.1%) genes were repressed by cold in C418 specifically during phases I, II, and III, respectively. Similarly, a small proportion of genes were stably regulated by cold during phases I and II, but a larger number were exclusively regulated by cold during phase III in K354 (Figure 3B). These results indicate that in both genotypes most genes were regulated under later cold stress, implying that cold stress responsiveness increased with increasing treatment duration.
A total of 156 (78.8%), 303 (91.5%), and 1121 (61.4%) genes were commonly induced and 33 (40.0%), 56 (48.3%), and 248 (39.6%) were commonly repressed by cold in both genotypes during phases I, II, and III, respectively (Additional file 7). The large proportion of shared cold-regulated genes found in C418 and K354 suggests that the major genetic pathways associated with cold response are conserved in both genotypes, consistent with their mostly common genetic background, except for the introgressed genetic segments.
Furthermore, a set of 130 and 21 genes were identified as commonly induced and repressed genes, respectively, in both genotypes over the entire course of the cold treatment (Additional file 8). Among these, 28 TF genes from 9 families were induced, including 11 AP2-EREBP genes, 5 MYB genes, 3 WRKY genes, 3 ZIM genes, and 2 C2H2 genes. AP2-EREBP TFs, including DREB1A (Os09g35030), DREB1B (Os09g35010), DREB1C (Os06g03670), DREB1F (Os01g73770), DREB1G (Os02g45450), ERF3 (Os01g0797600), and 5 other AP2 domain-containing proteins, were found to be significantly up-regulated in both genotypes under cold stress. This suggests that the major cold-responsive genetic pathway shared by C418 and K354 is possibly a DREB1-related pathway. In addition, we detected an important subset of commonly induced genes related to signal transduction in both genotypes, including genes encoding five MAPKKKs, four calmodulin-related calcium sensor proteins, and one CAMK. All these results indicate that conserved transcription regulation and signal transduction pathways are involved in cold stress response in both genotypes.
Unique functional categories of genotype-specific DEGs in K354 and C418
To identify genotype-specific DEGs, we used the following criterion: any gene with more than 5-fold change found only in C418 or K354 at one or more cold-responsive stages was regarded as a C418 or K354 specifically-regulated gene set. After removing redundancies, 641 and 399 cold-responsive DEGs were specifically identified in C418 and K354, respectively (Table 2). GO analysis was performed using Zheng’s method  to functionally classify the genotype-specific DEGs regulated by cold and constitutive DEGs. There were 304 genes determined to be exclusively induced by cold stress in K354 (Additional file 9); their dominant GO categories were metabolic, cellular, stimulus response, and cellular component organization processes. The 95 K354-specific genes (Additional file 10) repressed by cold stress were largely involved in oxidation reduction processes and oxidative stress responses.
In C418, 351 and 290 genes were identified as specifically induced and repressed under cold stress, respectively (Additional file 11 and Additional file 12). For the induced gene set, the most prevalent GO categories were oxidation-reduction (GO: 0055114), trehalose metabolic processes (GO: 0005991), and lipid transport (4 genes, GO: 0006869). Genes repressed exclusively by cold were involved in transmembrane transport (GO: 0055085), protein folding (GO: 0006457), oligopeptide transport (GO: 0006857), and L-methionine salvaging from methylthioadenosine (GO: 019509). The opposite transcription pattern observed for genes related to oxidation regulation can be partly explained by genetic differences in CT between the two genotypes. It has been suggested that oxidation regulation is involved in cold response , but whether high levels of antioxidant enzyme activity under cold stress are positively related to CT is presently unknown.
Among the constitutively highly expressed genes in K354 (Additional file 13), 184 genes were found to be non-responsive to cold stress, indicating an active role in CT in K354. All of these 184 genes were included in Cluster II (Figure 2). GO enrichment analysis indicated that these genes were functionally involved in endogenous stimulus (6.5%, p = 0.0006), protein modification (5.4%, p = 0.0013), signal transduction (5.4%, p = 0.0028), cell differentiation (2.2%, p = 0.0031), and amino acid and derivative metabolism (3.3%, p = 0.0035). Strikingly, ten genes with higher basal expression in K354 were related to PCD (GO: 0012501), and all of them were disease-resistance related genes with an NBS domain (Table 3). It has been reported that cell death pathways in plants can be activated not only by successful recognition of a pathogen during hypersensitive responses but also by a range of abiotic stresses . In addition, PCD has been reported to be activated in tobacco cells in response to cold stress . Furthermore, signal-induced plant PCD can be caused by reactive oxygen species generated by biotic and abiotic stresses, and can be perceived by the NBS domain . In addition, in our study the rice omega-3 desaturase gene OsFAD7 (Os03g18070) was 100 times more up-regulated in K354 than in C418 under control conditions. Rice has only one trienoic acid, linolenic acid, which is produced from linoleic acid by omega-3 fatty acid desaturase . Increased linolenic acid levels are essential for the maintenance of membrane fluidity and chloroplast function at low temperatures . Transgenic tobacco plants over-expressing AtFAD7 have been reported to exhibit chilling tolerance . We also observed that up-regulation of a cell wall-associated receptor kinase (OsWAK112d, Os10g10130) was 200 times greater in K354 than in C418 under control and cold stress conditions. Through an extracellular domain, wall-associated kinases physically link the plasma membrane to the cell wall matrix. They also have the potential to directly signal cellular events through their cytoplasmic kinase domain, which is one of the most likely candidates participating in cell wall-cytoplasm signaling in plant defense reactions . All of these results suggest that the higher activity of defensive proteins and membrane fluidity under control conditions might be involved in CT in K354.
The DREB1 regulon plays an important role in cold stress response in both genotypes
Of the 2384 known or annotated TF genes in the rice gnome , 196 (8.2%) were observed to be differentially regulated by cold in this study (Additional file 14). These TF genes include 26 AP2/EREBP, 21 bHLH, 7 HSF, 19 MYB, 12 NAC, and 17 WRKY genes. A total of 16, 42, and 138 TF genes were identified as K354-specific, C418-specific, or commonly regulated by cold, respectively. Interestingly, five commonly regulated TFs (Os01g58420, Os04g49450, Os07g07974, Os09g33550, and Os11g08210 in Additional file 14) in both genotypes were previously reported to be regulated by cold-induced expression of OsMyb4, suggesting that a DREB/CBF-independent OsMyb4 cold-regulated pathway exists in C418 and K354.
The OsDREB1 regulon has been shown to play an important role in conferring cold stress response in rice . OsDREB1s belong to a sub-family of AP2/EREBP proteins unique to plants and sharing a highly conserved AP2 domain. To identify commonly cold-induced genes within the OsDREB1 regulon, we examined cis-regulatory elements in the 1-kb upstream regions of 991 commonly induced genes in both genotypes using the DRE core motif A/GCCGAC . Of these genes, 298 (30.1%) contained at least one DRE core motif in the 1-kb upstream region (Additional file 15). As shown in Figure 4A, expression patterns of the OsDREB1 regulon (298 genes in Additional file 15) were roughly consistent with that of three OsDREB1 genes (OsDREB1A, OsDREB1B, and OsDREB1C). Of these 298 genes, 108 had been previously shown to be co-expressed with OsDREB1 in other rice microarray experiments using a positive correlation coefficient of 0.9 as the cutoff on the TIGR website  (Additional file 15). OsDREB1A and OsDREB1B were induced by cold with a similar expression pattern but OsDREB1C was expressed differentially during phase I (2–6 h) between the two genotypes (Figure 4). Interestingly, analysis of gene expression profiles revealed the existence of an OsDREB1C-specific regulon, which included 22 genes exclusively regulated by OsDREB1C, (Figure 4B and Table 4) accounting for the main difference in the OsDREB1 regulon between two genotypes. These results indicate that the DREB1 regulon plays an important role in cold stress response in both genotypes, even though there are a few genes, such as the OsDREB1C regulon, whose regulation is unique to genotype-specific pathways.
Co-localization of cold-responsive DEGs and introgressed chromosome segments in K354
Based on the genotyping results, there were 17 introgressed chromosome segments in K354 that differed from C418 (Additional file 16). A number of quantitative trait loci (QTLs) related to CT have been mapped by different research groups [31–42], among them, 31 QTLs that were co-localized onto the 15 introgressed chromosome regions in K354. We correspondingly mapped 148 DEGs onto these introgressed regions (Additional file 17).
There were 32 DEGs co-mapped onto an introgressed region of chromosome 12 where two major CT QTLs were previously identified [31, 42] (Additional file 17). Among these DEGs, two genes encoding the transcriptional regulator Sir2 protein (Os12g07950) and protein phosphatase 2C were highly enriched in K354. According to GO analysis, Sir2 is related to epigenetic function (GO: 0040029), and the protein phosphatase 2C is involved in stress signaling . These two genes are candidates for the CT QTL. In addition, most of the DEGs co-localized on the introgressed segment of chromosome 2 were found to be functionally associated with stress response. These genes include genes encoding cytochrome P450s, terpene synthase, NB-ARC domain containing protein, histone deacetylase, and TF OsWRKY42. Importantly, there was only one candidate gene, OsFAD7, found in the introgressed region near RM251 on chromosome 3 (Additional file 17). OsFAD7 has been reported to be highly induced by extreme temperature, and is involved in membrane stability , indicating it plays a crucial role in maintaining cell membrane fluidity in K354 under cold stress. DEGs on the introgressed regions of RM262 on chromosome 2, RM432 on chromosome 7, and RM38 on chromosome 8 were largely related to signal transduction cascades. These genes included genes encoding terpene synthase (Os02g36140, Os08g07100), protein phosphatase 2C (Os07g32380), and receptor-like kinase OsSIK1 (Os08g02996), which are known to be strongly associated with abiotic stress response [43, 45, 46]. Strikingly, a gene cluster containing four PCD-related genes were co-localized on the introgressed region of chromosome 11 (Additional file 17), indicating their distinct role in CT of K354. These DEGs may serve as functional candidate genes for the identification of a QTL for CT. By combining further functional identification and QTL fine mapping, the co-localized DEGs detected in this microarray analysis may provide the basis for gene cloning and elucidation of the molecular mechanisms responsible for CT in rice.
One objective of this study was an attribution of the differences between C418 and K354 transcriptomes to differences in CT between the two genotypes. These differences, whether at the transcription or phenotypic level, can only be caused by introgressed fragments from Bg300, an indica variety. Although japonica varieties are generally known to be more tolerant to cold stress than indica varieties, we have shown in this study that some CT alleles are “hidden” in the indica background. This means that the CT effect of the alleles may not be favorable in indica, but when introgressed into a different background (i.e., japonica), the allele may manifest a different pattern. Because the genomes of C418 and K354 are almost identical, any transcriptomic differences between the two genotypes must be due to the presence of fragments introgressed from Bg300, whether caused by direct allelic replacement or genome-wide transcriptome reconfiguration as a result of the replacement. Comparison of C418 and K354 transcriptomes can therefore reduce background noise from Bg300. Although we did not do so in this study, the question of which process, i.e., direct allelic replacement or genome-wide transcriptome reconfiguration, is responsible for transcriptome changes between C418 and K354 could be resolved by adding Bg300 to the transcriptome comparison.
Many K354-specific expressed genes that were not mapped onto or near introgressed regions may also play a vital role in CT, possibly through regulation by co-localized DEGs. By this we mean that a large number of genes involved in a signaling pathway related to a certain trait may be activated by expression of upstream genes. Due to the low resolution of the limited SSR markers, we can only roughly estimate the length of the introduced fragments. It is worth noting, however, that whole genomic re-sequencing for K354 and C418 is underway. We believe the re-sequencing results will provide more detailed and precise information regarding the introgressed chromosomal segments. Fine-mapping of these segments is necessary to more precisely identify the functional candidate genes.
Using expression profiling with an Affymetrix rice genome array, we found that at the seedling stage the CT IL K354 and its recurrent parent C418 exhibited a diverse whole genome-wide transcriptional response under both control and cold stress conditions. We detected 3184 genes that were differentially regulated under cold stress, accounting for approximately 6% of total genes on the rice microarray chip. A large proportion of genes, including DREB1 and OsMyb4 regulons, were determined to be commonly regulated by cold in both genotypes. Expression of cold-responsive genes specific to the cold-tolerant genotype K354 may be affected by introgressed chromosome segments, and these genes may contribute to the increased CT of K354. Among genes co-localizing with introgressed chromosomal segments from the indica donor parent were genes involved in PCD, which were constitutive but had higher absolute expression levels in the IL than in the japonica recurrent parent, and some differentially induced genes likely to be involved in the regulation of membrane biogenesis. Our results represent a preliminary elucidation of genotype-dependent CT enhancement mechanisms, and may contribute to the improvement of CT in japonica rice using indica germplasm. Functional validation of CT related genes identified in this study, however, is still required to provide more useful information for CT improvement in japonica rice molecular breeding.
Plant materials and cold treatment
A japonica restorer line, C418, which is the male parent of many elite japonica commercial hybrid cultivars widely grown in northern China, and the CT IL K354, a BC2F6 IL derived from a cross between recurrent parent C418 and Bg300 (an indica variety from Sri Lanka), were used in this study. After an extensive study using 100 simple sequence repeat (SSR) markers that were polymorphic between C418 and Bg300, selected from 600 SSR markers across the rice genome , K354 was found to differ from C418 at 17 genomic segments derived from Bg300 (Additional file 16).
Sterilized seeds of C418 and K354 were allowed to germinate in distilled water for 2 d. Well-germinated seeds were sown in soil-filled plates with holes at the bottom. C418 and K354 were planted at the same spacing in the same plates. The plates were placed in a 65 cm × 44 cm × 14 cm plastic case with nutrient solution, following the IRRI standard protocol , and transferred to a growth chamber (Beijing ZNYT, China). Both cold stress and control treatments were carried out with 24 replicates for each line. Approximately 15 vigorous seedlings from each replicate were allowed to grow until the three-leaf stage in the growth chamber. Growing conditions were a constant 25°C day/night temperature, 12 h photoperiod, and 75–80% relative humidity.
At the three-leaf stage, seedlings were transferred to a 4°C chamber. The time point at which seedlings were placed in the 4°C chamber was defined as 0 h. Whole shoot samples were collected at 2, 6, 12, 24, and 48 h and frozen immediately in liquid nitrogen. Control seedlings were not exposed to cold stress, but were instead maintained at a constant temperature of 25°C; whole shoot samples were taken from the controls at 2 h and 24 h. Collections at the two time points were pooled together as a control. Three biological replicates (each from a separate plastic pot) were prepared for microarray analysis. After collection, samples were kept in a −80°C freezer for later total RNA extraction. To detect phenotypic changes in the stressed seedlings, seedlings were moved to control conditions after 48 h of stress. After 24 h recovery, a photograph of each seedling was taken to serve as a phenotypic record.
Physiological traits of the two genotypes under cold stress
Membrane stability was measured using a previously published procedure , with minor modifications for rice leaf tissue. Three replicates of 0.5 g fresh leaves were sampled from control and cold-treated seedlings. 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 aid 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 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 that of SOD and CAT, were determined following previously reported methods . Soluble protein concentrations were measured according to the protocol of Bradford .
RNA extraction and processing for microarray analysis
RNA was prepared following the recommendations in the Affymetrix GeneChip Expression Analysis technical manual. Briefly, total RNA was isolated from samples frozen in liquid nitrogen using TRIZOL reagent, and then processed by CapitalBio Corporation (Beijing, China) according to the following steps. Two micrograms of purified total RNA were used to synthesize double-stranded cDNA. Biotin-tagged cRNA was generated from an in vitro transcription reaction using a MessageAmp II aRNA amplification kit and then fragmented into 35–200 bp long strands following the Affymetrix protocol. The cRNA was then hybridized to an Affymetrix rice genome array containing 48,564 japonica and 1260 indica sequences at 45°C with rotation for 16 h in an Affymetrix GeneChip 640 hybridization oven. The GeneChip arrays were washed and stained on an Affymetrix 450 Fluidics Station, and then scanned using an Affymetrix 3000 GeneChip Scanner.
Array data statistical analysis
GeneChip Operating Software (GCOS1.4) was used to analyze the hybridization data. Following visual inspection of scanned images, satisfactory images were analyzed to generate CEL raw data files using GCOS1.4 default settings. dChip software was used to perform array normalization, following the dChip user’s manual, according to an invariant set approach. For the comparison analysis, SAM (Significant Analysis of Microarray) software was applied using the two-class unpaired method to identify DEGs between cold-stressed samples and control samples for a given genotype, and between two genotypes at each time point. As there is not a fixed standard threshold for significant differential gene expression, we identified DEGs using the empirical criterion of more than 5-fold change and a cut-off of p < 0.05 in SAM based on three independent biological replicates. An average linkage hierarchical cluster analysis was performed on the DEGs using GeneCluster version 3.0 software .
The putative function of each DEG corresponding to a probe on the Genechip was predicted based on its Affymetrix annotation combined with the TIGR definition . GO enrichment analysis was performed using the web-based software toolkit GOEAST . To analyze DEG regulatory elements, promoter regions (−10 to −1000 bp upstream of the start codon) of selected DEGs were scanned for the DRE motif A/GCCGAC using a Perl program from CapitalBio Corporation. Gene sequences, including 1 kb upstream sequences, of all selected DEGs were downloaded from the TIGR website .
Quantitative RT-PCR validation of DEGs
A total of 15 genes were selected for qRT-PCR validation of the expression level changes identified in the microarray analysis. The selected genes included OsDREB1 regulon genes and genes that were significantly differentially expressed between K354 and C418 and near or within regions of donor introgressed fragments in K354. Sequences corresponding to the selected genes were downloaded from the TIGR rice database , and Primer 5.0 software was then used to design qRT-PCR primers based on exonic gene regions. An Actin gene was used as the endogenous control. First-strand cDNA was synthesized by reverse transcription using 2 μg of total RNA in a 100 μL reaction mixture using a Reverse Transcription System (Promega, A3500). Diluted synthesized cDNA (1 μL) was used for qRT-PCR analysis along with 200 nM of each primer mixed with SYBR Green PCR Master Mix (Takara Code: DRR041A). qRT-PCR was performed in triplicate for each data point using an ABI7300 sequence detection system and the same RNA samples used for the microarray analysis. Relative expression values were calculated after normalizing against the maximum expression value. The data were further normalized to facilitate profile matching to data obtained from the microarray experiments.
The entire set of original microarray data has been deposited in NCBI’s Gene Expression Omnibus  under GEO Series number GSE37940.
Differentially expressed genes
Programmed cell death
Quantitative trait loci
Simple sequence repeat
Significant analysis of microarray.
Mackill DJ, Lei X: Genetic variation for traits related to temperate adaptation of rice cultivars. Crop Sci. 1997, 37: 1340-1346. 10.2135/cropsci1997.0011183X003700040051x.
Yadav SK: Cold stress tolerance mechanisms in plants. A review. Agron Sustain Dev. 2010, 30 (3): 515-527. 10.1051/agro/2009050.
Horvath I, Glatz A, Varvasovszki V, Torok Z, Pali T, Balogh G, Kovacs E, Nadasdi L, Benko S, Joo F, Vigh L: Membrane physical state controls the signaling mechanism of the heat shock response in Synechocystis PCC 6803: identification of hsp17 as a “fluidity gene”. Proc Natl Acad Sci USA. 1998, 95 (7): 3513-3518. 10.1073/pnas.95.7.3513.
Routaboul JM, Fischer SF, Browse J: Trienoic fatty acids are required to maintain chloroplast function at low temperatures. Plant Physiol. 2000, 124 (4): 1697-1705. 10.1104/pp.124.4.1697.
Prakash M, Nair G, Kang IS, Moon BY, Lee CH: Effects of low temperature stress on rice (Oryza sativa L.) plastid ω-3 desaturase gene, OsFAD8 and its functional analysis using T-DNA mutants. Plant Cell Tiss Organ Cult. 2009, 98 (1): 87-96. 10.1007/s11240-009-9541-y.
da Cruz RP, Golombieski JI, Bazana MT, Cabreira C, Silveira TF, da Silva LP: Alterations in fatty acid composition due to cold exposure at the vegetative stage in rice. Braz J Plant Physiol. 2010, 22 (3): 199-207. 10.1590/S1677-04202010000300007.
Thomashow MF: So what’s new in the field of plant cold acclimation? Lots!. Plant Physiol. 2001, 125 (1): 89-93. 10.1104/pp.125.1.89.
Dubouzet JG, Sakuma Y, Ito Y, Kasuga M, Dubouzet EG, Miura S, Seki M, Shinozaki K, Yamaguchi-Shinozaki K: OsDREB genes in rice, Oryza sativa L., encode transcription activators that function in drought-, high-salt- and cold-responsive gene expression. Plant J. 2003, 33 (4): 751-763. 10.1046/j.1365-313X.2003.01661.x.
Ito Y, Katsura K, Maruyama K, Taji T, Kobayashi M, Seki M, Shinozaki K, Yamaguchi-Shinozaki K: Functional analysis of rice DREB1/CBF-type transcription factors involved in cold-responsive gene expression in transgenic rice. Plant Cell Physiol. 2006, 47 (1): 141-153.
Su CF, Wang YC, Hsieh TH, Lu CA, Tseng TH, Yu SM: A novel MYBS3-dependent pathway confers cold tolerance in rice. Plant Physiol. 2010, 153 (1): 145-158. 10.1104/pp.110.153015.
Park MR, Yun KY, Mohanty B, Herath V, Xu FY, Wijaya E, Bajic VB, Yun SJ, De Los Reyes BG: Supra-optimal expression of the cold-regulated OsMyb4 transcription factor in transgenic rice changes the complexity of transcriptional network with major effects on stress tolerance and panicle development. Plant Cell Environ. 2010, 33 (12): 2209-2230. 10.1111/j.1365-3040.2010.02221.x.
Cheng C, Yun KY, Ressom HW, Mohanty B, Bajic VB, Jia YL, Yun SJ, De Los Reyes BG: An early response regulatory cluster induced by low temperature and hydrogen peroxide in seedlings of chilling-tolerant japonica rice. BMC Genomics. 2007, 8: 175-10.1186/1471-2164-8-175.
Jan N, ul-Hussain M, Andrabi KI: Cold resistance in plants: A mystery unresolved. Electron J Biotechnol. 2009, 12 (3): 1-15.
Baginsky S, Hennig L, Zimmermann P, Gruissem W: Gene expression analysis, proteomics, and network discovery. Plant Physiol. 2010, 152 (2): 402-410. 10.1104/pp.109.150433.
Malone JH, Oliver B: Microarrays, deep sequencing and the true measure of the transcriptome. BMC Biol. 2011, 9: 34-10.1186/1741-7007-9-34.
Rabbani MA, Maruyama K, Abe H, Khan MA, Katsura K, Ito Y, Yoshiwara K, Seki M, Shinozaki K, Yamaguchi-Shinozaki K: Monitoring expression profiles of rice genes under cold, drought, and high-salinity stresses and abscisic acid application using cDNA microarray and RNA gel-blot analyses. Plant Physiol. 2003, 133 (4): 1755-1767. 10.1104/pp.103.025742.
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.
Walia H, Wilson C, Condamine P, Liu X, Ismail AM, Zeng LH, Wanamaker SI, Mandal J, Xu J, Cui XP, Close TJ: Comparative transcriptional profiling of two contrasting rice genotypes under salinity stress during the vegetative growth stage. Plant Physiol. 2005, 139 (2): 822-835. 10.1104/pp.105.065961.
Lenka SK, Katiyar A, Chinnusamy V, Bansal KC: Comparative analysis of drought-responsive transcriptome in indica rice genotypes with contrasting drought tolerance. Plant Biotechnol J. 2011, 9 (3): 315-327. 10.1111/j.1467-7652.2010.00560.x.
Zheng Q, Wang XJ: GOEAST: a web-based software toolkit for gene ontology enrichment analysis. Nucl Acids Res. 2008, 36 (suppl2): W358-W363.
Neilson KA, Mariani M, Haynes PA: Quantitative proteomic analysis of cold-responsive proteins in rice. Proteomics. 2011, 11 (9): 1696-1706. 10.1002/pmic.201000727.
van Doorn WG, Beers EP, Dangl JL, Franklin-Tong VE, Gallois P, Hara-Nishimura I, Jones AM, Kawai-Yamada M, Lam E, Mundy J, Mur LA, Petersen M, Smertenko A, Taliansky M, Van Breusegem F, Wolpert T, Woltering E, Zhivotovsky B, Bozhkov PV: Morphological classification of plant cell deaths. Cell Death Differ. 2011, 18 (8): 1241-1246. 10.1038/cdd.2011.36.
Koukalova B, Kovarik A, Fajkus J, Siroky J: Chromatin fragmentation associated with apoptotic changes in tobacco cells exposed to cold stress. FEBS Lett. 1997, 414 (2): 289-292. 10.1016/S0014-5793(97)01008-9.
Das A, Kawai-Yamada M, Uchimiya H: Programmed cell death in plants. Abiotic stress adaptation in plants: physiological, molecular and genomic foundation. 2010, The Netherlands: Springer, 371-383. 1
Kodama H, Hamada T, Horiguchi G, Nishimura M, Iba K: Genetic enhancement of cold tolerance by expression of a gene for chloroplast omega-3 fatty acid desaturase in transgenic tobacco. Plant Physiol. 1994, 105 (2): 601-605.
Li H, Zhou SY, Zhao WS, Su SC, Peng YL: A novel wall-associated receptor-like protein kinase gene, OsWAK1, plays important roles in rice blast disease resistance. Plant Mol Biol. 2009, 69 (3): 337-346. 10.1007/s11103-008-9430-5.
Gao G, Zhong Y, Guo A, Zhu Q, Tang W, Zheng W, Gu X, Wei L, Luo J: DRTF: a database of rice transcription factors. Bioinformatics. 2006, 22 (10): 1286-1287. 10.1093/bioinformatics/btl107.
Nakashima K, Ito Y, Yamaguchi-Shinozaki K: Transcriptional regulatory networks in response to abiotic stresses in Arabidopsis and grasses. Plant Physiol. 2009, 149 (1): 88-95. 10.1104/pp.108.129791.
Yamaguchi-Shinozaki K, Shinozaki K: A novel cis-acting element in an Arabidopsis gene is involved in responsiveness to drought, low-temperature or high-salt stress. Plant Cell. 1994, 6 (2): 251-264.
Ouyang S, Zhu W, Hamilton J, Lin H, Campbell M, Childs K, Thibaud-Niseen F, Malek RL, Lee Y, Zheng L, Orvis J, Haas B, Wortman J, Buell CR: The TIGR rice genome annotation resource: improvements and new features. Nucleic Acids Res. 2007, 35: D883-D887. 10.1093/nar/gkl976.
Andaya VC, Mackill DJ: Mapping of QTLs associated with cold tolerance during the vegetative stage in rice. J Exp Bot. 2003, 54 (392): 2579-2585. 10.1093/jxb/erg243.
Baruah AR, Ishigo-Oka N, Adachi M, Oguma Y, Tokizono Y, Onishi K, Sano Y: Cold tolerance at the early growth stage in wild and cultivated rice. Euphytica. 2009, 165 (3): 459-470. 10.1007/s10681-008-9753-y.
Han LZ, Qiao YL, Cao GL, Zhang YY, An YP, Ye JD, Koh HJ: QTLs analysis of cold tolerance during early growth period for rice. Rice Sci. 2004, 11 (5–6): 245-250.
Iwata N, Shinada H, Kiuchi H, Sato T, Fujino K: Mapping of QTLs controlling seedling establishment using a direct seedling method in rice. Breeding Sci. 2010, 60 (4): 353-360. 10.1270/jsbbs.60.353.
Ji ZJ, Zeng YX, Zeng DL, Ma LY, Li XM, Liu BX, Yang CD: Identification of QTLs for rice cold tolerance at plumule and 3-leaf-seedling stages by using QTLNetwork software. Rice Sci. 2010, 17 (4): 282-287. 10.1016/S1672-6308(09)60028-7.
Jiang L, Xun MM, Wang JK, Wan JM: QTL analysis of cold tolerance at seedling stage in rice (Oryza sativa L.) using recombination inbred lines. J Cereal Sci. 2008, 48 (1): 173-179. 10.1016/j.jcs.2007.08.013.
Lou QJ, Chen L, Sun ZX, Xing YZ, Li J, Xu XY, Mei HW, Luo LJ: A major QTL associated with cold tolerance at seedling stage in rice (Oryza sativa L.). Euphytica. 2007, 158 (1–2): 87-94.
Qian Q, Zeng DL, He P, Zheng XW, Chen Y, Zhu LH: QTL analysis of the rice seedling cold tolerance in a double haploid population derived from anther culture of a hybrid between indica and japonica rice. Chin Sci Bull. 2000, 45 (5): 448-453. 10.1007/BF02884949.
Zhang ZH, Su L, Li W, Chen W, Zhu YG: A major QTL conferring cold tolerance at the early seedling stage using recombinant inbred lines of rice (Oryza sativa L.). Plant Sci. 2005, 168 (2): 527-534. 10.1016/j.plantsci.2004.09.021.
Han LZ, Qiao YL, Zhang SY, Zhang YY, Cao GL, Kim J, Lee K, Koh H: Identification of quantitative trait loci for cold response of seedling vigor traits in rice. J Genet Genomics. 2007, 34 (3): 239-246. 10.1016/S1673-8527(07)60025-3.
Andaya VC, Tai TH: Fine mapping of the qCTS4 locus associated with seedling cold tolerance in rice (Oryza sativa L.). Mol Breed. 2007, 20 (4): 349-358. 10.1007/s11032-007-9096-8.
Andaya VC, Tai TH: Fine mapping of the qCTS12 locus, a major QTL for seedling cold tolerance in rice. Theor Appl Genet. 2006, 113 (3): 467-475. 10.1007/s00122-006-0311-5.
Schweighofer A, Hirt H, Meskiene I: Plant PP2C phosphatases: emerging functions in stress signaling. Trends Plant Sci. 2004, 9 (5): 236-243. 10.1016/j.tplants.2004.03.007.
Cao YP, Shi JL, Li Z, Ming F: Isolation of OsFAD2, OsFAD6 and FAD family members response to abiotic stresses in Oryza sativa L. Hereditas. 2010, 32 (8): 839-847.
Mazid M, Khan TA, Mohammad F: Role of secondary metabolites in defense mechanisms of plants. Biol Med. 2011, 3 (2): 232-249.
Ouyang SQ, Liu YF, Liu P, Lei G, He SJ, Ma B, Zhang WK, Zhang JS, Chen SY: Receptor-like kinase OsSIK1 improves drought and salt stress tolerance in rice (Oryza sativa) plants. Plant J. 2010, 62 (2): 316-329. 10.1111/j.1365-313X.2010.04146.x.
Gramene marker database.http://www.gramene.org/markers/,
Yoshida S, Forna DA, Cock JH: Laboratory manual for physiological studies of rice. 1971, Los Baños: International Rice Research Institute
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 (1): 24-34. 10.1034/j.1399-3054.1998.1030104.x.
Bonnecarrè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 (5): 726-732. 10.1016/j.plantsci.2011.01.023.
Bradford MM: A rapid and sensitive method for quantitation of microgram quantities of protein utilizing the principle of protein-dye binding. Anal Biochem. 1976, 72: 248-254. 10.1016/0003-2697(76)90527-3.
Open Source Clustering Software.http://bonsai.hgc.jp/~mdehoon/software/cluster/,
NCBI’s Gene Expression Omnibus.http://www.ncbi.nlm.nih.gov/geo/,
We thank Dr. Judy Lee for editing the manuscript. This study was funded by the “863” project (2012AA101101) of the Ministry of Science and Technology of China, the “948” Project (#2011-G2B) of the Ministry of Agriculture of China, a program of the International S & T Cooperation (S2012ZR0160), and the Bill & Melinda Gates Foundation (OPP51587). F. Zhang was also supported by the Monsanto’s Beachell-Borlaug International Scholars Program.
The authors declare that they have no competing interests.
FZ, BF and ZL designed the experiments and drafted the manuscript. FZ, LH, WW, XZ and LZ conducted the phenotype and microarray experiments. FZ performed the microarray data analysis, and LH carried out the qRT-PCR confirmation of candidate genes. BF revised the manuscript. All authors read and approved the final manuscript.
Fan Zhang, Liyu Huang contributed equally to this work.
Electronic supplementary material
Additional file 1: Phenotypic performance between C418 and K354 in different experiments. Excel file containing the phenotypic performance of C418 and K354 under both control and abiotic stress conditions in different experiments. (XLS 32 KB)
Additional file 2: Phenotype and physiological changes of K354 and C418 under control and cold stress conditions. A PowerPoint file containing photographs of growth status and physiological conditions of K354 and C418 under control conditions (25°C) and cold stress treatment (4°C) for 48 h. (PPT 375 KB)
Additional file 3: Total number of expressed transcripts. Word file containing the total number of transcripts expressed in the two genotypes at different times under control and cold stress conditions. (DOC 30 KB)
Additional file 4: Expression patterns of differentially expressed transcripts in K354 and C418. A PowerPoint file containing expression patterns of 3184 differentially expressed probes based on log-transformed and normalized expression values from control conditions to 48 h under cold stress in (A) K354 and (B) C418. Green, blue, yellow, cyan, red, and black lines indicate clusters I–VI with the same color in Figure 2, respectively. (PPT 938 KB)
Additional file 6: Validation of gene expression by qRT-PCR. A PowerPoint file containing correlation analysis results between microarray and qRT-RCR experiments based on 15 selected genes representing the OsDREB1 regulon and candidate genes of introgressed segments. Gene expression values were transformed to a log2 scale. The microarray data log2-value (X-axis) were plotted against the qRT-PCR log2-value (Y-axis). (PPT 376 KB)
Additional file 7: Comparative diagram of the total number of differentially regulated genes between C418 and K354 under cold stress. A PowerPoint file containing (A) Cold-induced genes; (B) Cold-repressed genes. Blue, red, and green bars indicate K354-specific (blue bar), common (red bar), and C418-specific (green bar) regulated genes by cold stress at phase I (2–6 h), phase II (12 h), and phase III (24–48 h) of cold-response, respectively. (PPT 128 KB)
Additional file 8: Common genes regulated by cold in both genotypes at all cold-responsive stages. Word file containing 130 and 21 genes commonly induced and repressed by cold in both genotypes at all cold-responsive stages. (DOC 222 KB)
Additional file 9: K354-specific cold-induced genes. Excel file containing a list of K354-specific cold-induced genes at all time points. (XLS 122 KB)
Additional file 10: K354-specific cold-repressed genes. Excel file containing a list of K354-specific cold-repressed genes at all time points. (XLS 56 KB)
Additional file 11: C418-specific cold-induced genes. Excel file containing a list of C418-specific cold-induced genes at all time points. (XLS 139 KB)
Additional file 12: C418-specific cold-repressed genes. Excel file containing a list of C418-specific cold-repressed genes at all time points. (XLS 118 KB)
Additional file 13: Constitutively highly-expressed genes in K354. Excel file containing a list of non-cold-responsive genes with constitutively higher transcriptional expression in K354 than in C418. (XLS 77 KB)
Additional file 14: TF genes detected to be differentially regulated by cold. Excel file containing a list of TF genes differentially regulated by cold in K354 and C418. (XLS 92 KB)
Additional file 15: OsDREB1s and commonly cold-induced genes including a DRE core motif. Excel file containing a list of OsDREB1s and 298 genes commonly cold-induced in both C418 and K354 with at least one DRE core motif at 1 kb up-stream. (XLS 162 KB)
Additional file 16: Genomic distribution of 17 introgression fragments in K354. A PowerPoint file containing information on introgressed fragments, polymorphic SSR markers used, and previously reported QTLs related to CT at the seedling stage near the introgressed regions. (PPT 214 KB)
Additional file 17: DEGs in the introgressed regions of K354. Excel file containing a list of differentially regulated genes detected by comparison of gene expression profiles of two genotypes that were within or neighboring the regions of donor introgression. (XLS 44 KB)
About this article
Cite this article
Zhang, F., Huang, L., Wang, W. et al. Genome-wide gene expression profiling of introgressed indica rice alleles associated with seedling cold tolerance improvement in a japonica rice background. BMC Genomics 13, 461 (2012). https://doi.org/10.1186/1471-2164-13-461