Journal of Genetics & Genomic Sciences Category: Genetics Type: Research Article

Four Cold Responsive Genes That Show Consistent Differential Expression between the Japonica and Indica subspecies of rice

Donghai Mao1*
1 Institute Of Subtropical Agriculture, Chinese Academy Of Sciences, Beijing, China

*Corresponding Author(s):
Donghai Mao
Institute Of Subtropical Agriculture, Chinese Academy Of Sciences, Beijing, China
Tel:+ 86 73184619748,

Received Date: Apr 04, 2019
Accepted Date: Apr 08, 2019
Published Date: Apr 15, 2019


Rice (Oryza sativa L.) is highly sensitive to cold stress during the seedling stage and of the two cultivated subspecies, japonica rice varieties tend to be more tolerant to low temperatures than are indica rice varieties. In this study, we performed comparative transcriptome analysis and found 171 genes that respond to cold under four different conditions in the two rice subspecies. qRT-PCR assays and natural variation analysis were further carried out to identify the differentially-expressed Cold Responsive genes (CORs) involved in the difference in cold tolerance between the japonica and indica rice subspecies. We identified four genes, OsDREB1A, OsABCB5, OsPHIL2 and OsREM4. 1 that showed consistently different expression patterns in the two rice subspecies in response to low temperature. Based on the results of nature variation analysis of these four genes, we predicted they respond to cold through the CBF/DREB1or plant hormone pathway, respectively. The four differentially expressed genes may be key genes lead to different cold tolerance between two rice subspecies. These results not only contribute to a comprehensive understanding of the differences in transcriptome levels of low temperature response genes between rice subspecies, but also provide candidate genes for molecular breeding of cold-tolerant rice cultivars.


Cold responsive genes; Rice; Subspecies; Transcriptome


Rice (Oryza sativa L.) is the most important staple crop that feeds nearly half of the world’s population. Because originated in tropical or subtropical areas, rice is sensitive to low temperatures throughout its entire life cycle. Cold stress can inhibit photosynthesis by reducing the chlorophyll concentration, inducing the accumulation of Reactive Oxygen Species (ROS) and by causing severe damage to various cellular components including membrane lipids, structural proteins and enzymes [1-3]. Low temperatures have the potential to negatively affect the growth and development of rice plants during any developmental stage, from germination to grain filling, resulting in large losses in rice production [4].

Plants have developed several strategies to respond to low temperatures that are known collectively as ‘cold acclimation’ [5]. Numerous studies of cold acclimation have been reported in recent years, and they showed that many pathways participate in the cold-response process [6,7]. The Dehydration-Responsive Element Binding protein (DREB1)/C-repeat binding factor (CBF)-dependent cold-responsive pathway is the most well-known among all the pathways [8]. Another pathway with a central role in the rice response to cold stress is the MYBS3-dependent pathway. Interestingly, DREB1A responds early and transiently, while MYBS3 responds relatively slowly to cold stress in rice, suggesting that theMYBS3-dependent pathway is important for long-term adaptation to persistent cold stress [9]. In addition to these two pathways, the plant hormone systems play an indispensable role in cold stress tolerance [10]. For example, Jasmonate (JA) positively regulates cold tolerance in plants because it can positively modulate the CBF-dependent or –independent pathway, leading to accumulation of cryoprotective compounds such as polyamine, glutathione and anthocyanins [11]. DELLA proteins, the key components of the Gibberellic Acid (GA)-signaling pathway, are involved in the CBF1-mediated cold stress response [12]. The cross-talk between Cytokinin (CK) and Abscisic Acid (ABA) are also important in cold stress signaling. Components of CK signaling, Histidine Kinases (AHKs) and the cold-inducible A-type ARRs, play negative regulatory roles in cold stress response via inhibition of the ABA response [13]. Components of Brassinosteroid (BR) and ethylene signaling were also found to modulate freezing tolerance in Arabidopsis plants [14,15].

There are two major subspecies of rice, Oryza sativa ssp. Indica and Oryza sativa ssp. japonica. Generally, japonica rice varieties are more tolerant to low temperatures than varieties of indica rice [16]. Therefore, cold tolerance in indica rice can be improved by introducing genes from japonica rice through molecular breeding. Recent advances in next generation sequencing technology provide a powerful strategy for discovering new genes associated with cold tolerance in rice [17-21]. In this study, we combined transcriptome sequencing with genome re-sequencing to identify the Cold-Responsive-genes (CORs) that are differentially-expressed between the two rice subspecies. Using this approach, we found four COR genes that are differentially-expressed between indica and japonica varieties. Moreover, these genes show important natural variation between the two subspecies which predicted to result in changes of cold responsive cis-elements. Therefore, these CORs may be useful for improving cold tolerance in rice through molecular breeding.


Global mRNA expression profiles in rice seedlings in response to cold stress

In general, japonica rice varieties display a higher potential for cold tolerance than do indica varieties [16] (Figure S1). In order to fully understand the genetic basis of the difference in the cold response between the two subspecies, we subjected seedlings of two cultivars, Teqing and 02428, to low temperature treatments at 4°C and 10°C. As shown in figure S1, the seedlings showed contrasting cold tolerance in response to the two treatments, suggesting that these two cultivars, Teqing and 02428 are representative of the cold sensitivity of indica rice and the cold tolerance of Japonica rice.

 Figure S1: Phenotypes of the seedlings of the indica and japonica rice cultivars before (A) and after exposure at 4°C for 4 days (B) and 10oC for 5 days(C).Teqing and 02428 are indicated with red arrows). Graphical representation of the relative survival rates of 7 indica and 6 japonica seedlings following exposure to cold stress at 4°C for 4 days (D) or10°C for 5 days (E). Mean values with one or two asterisks were found to be significantly different by Student's t test (*P≤0.05; **P≤0.01; n≥6).

A total of 10 sequencing libraries were constructed from RNA extracted from the shoot tissues of three-leaf stage rice seedlings of the indica rice cultivar Teqing and the japonica cultivar 02428exposed to low temperature conditions(4°C and 10°C) for 0hr, 3hr and 24hr. Because there was no difference between the two 0-hr cold treatment samples, two common controls, I0 and J0, were used for the low temperature treatments (4°C and 10°C) in the indica and japonica rice cultivars, respectively.
The Solexa Genome Analyzer was used to perform high-throughput Tag-seq analysis on these libraries. As shown in table 1, we obtained5.8 to 6.2 million raw sequence tags per library, of which 251, 288 to 301,590 tag sequences were distinct. After filtering out the 3' adaptor sequence, empty reads, low-quality tags, unexpected-length tags and single-copy tags, approximately 5.7 to 6.0 million clean tags remained in each sample. Finally, we obtained 111,174 to 136,023 unique tags for each of the ten libraries. Saturation analysis was applied to estimate whether or not new unique tags can be detected with increases in the total number of tags. The numbers of unique tags increased with the total number of tags and reached a plateau shortly after 1 million tags; no new unique tag was identified as the total number of tags approached 2 million. Therefore, the ten libraries are full representations of transcripts under the different treatments.


Raw Tags

Distinct Tags

Clean Tags

Distinct Clean Tags

Map to Gene Tags

Unique Map to Gene Tags

Map to Genome Tags

Unknown Tags







2821188 (49.31%)

473352 (8.27%)

281027 (4.91%)







2747337 (48.29%)

530238 (9.32%)

262248 (4.61%)







2780519 (48.40%)

704593 (12.26%)

210563 (3.67%)







2788766 (46.10)

576051 (9.52%)

478784 (7.92%)







2552274 (44.70%)

496393 (8.68%)

487862 (8.54%)







2747870 (47.04%)

548967 (9.40%)

427267 (7.31%)







2574611 (45.50%)

566161 (10.01%)

453345 (8.01%)







2810809 (46.87%)

755255 (12.59%)

352231 (5.87%)







2781934 (46.82%)

607244 (10.22%)

418005 (7.04%)







2603310 (45.57%)

5226813 (9.22%)

475573 (8.32%)

Table 1: Numbers of sequence tags in the 10 transcriptome libraries from rice cultivars Teqing and 02428 exposed to low temperatures.

We then annotated the sequence tags based on the Os-Nipponbare-Reference-IRGSP-1.0 genome assembly ( and only clean and unambiguous tags that matched perfectly or had only a single mismatch were analyzed further [22]. Based on these criteria, 45.50%-49.31% of the clean tags in all the ten libraries mapped to unique positions on the reference genes. 8.27%-12.59% of the clean tags from these same libraries were matched to the reference genome database. However, 3.67%-8.54% of the clean tags database from the 10 libraries did not map to the reference genome.

Identification of COR genes in the indica and japonica rice cultivars

Based on hierarchical clustering of the 10 transcriptomes, the 10 experimental samples were classified into two clades (Figure 1). The first clade contained the samples exposed to the low temperature treatments (4 or 10°C) for a short time (3hr) and the 0-hr controls without cold treatment, while the second clade contained the samples exposed to the cold treatments (4 or 10°C) that lasted for 24hr. Moreover, in the second clade, there are two sub-clades that correspond to the two temperatures, and each sub-clade contained both indica and japonica rice cultivars. These results suggest that the differences between the duration of cold treatment were more significant than the differences between either the two low temperature treatments or between the two cultivars, and that the cultivar difference was the smallest.
In order to find CORs in the two cultivars, we compared gene expression between cold treatment samples and the control samples for the same cultivar. The transcripts with at least two-fold differences (|log2Ratio| ≥1 and FDR ≤ 0.001) were referred as CORs. As shown in figure 1, there were a total of 1610 CORs in Teqing and 1307 CORs in 02428 under the 4°C/3hr treatment. After 3-hour exposure to 10°C, 1654 CORs were expressed in Teqing, and 2325 CORs were expressed in 02428. Following 24-hour cold treatment, more CORs were identified in these two cultivars: for the 4°C/24hr treatment, 4870 and 5613 CORs were expressed in Teqing and 02428, respectively, while there were 4712 and 5125 CORs expressed in Teqing and 02428 under the 10°C/24hr treatment, respectively.
Figure 1: (A) Hierarchical cluster analysis of gene expression based on the log ratio fold-change data. I0 and J0 indicate the control samples of Teqing (I) and 02428 (J) before cold treatments, respectively. The number before the letter (I or J) indicates the temperature, the number after the letter indicates the hours of exposure. For example, 10I24 means the Teqing sample exposed in 10oC for 24 hours. (B) Number of up- and down-regulated Cold-Responsive genes (CORs) from Teqing and 02428 rice seedlings exposed to cold (4°C and 10°C) for 3 and 24hr.

Gene ontology and pathway enrichment analysis of the rice CORs

In order to reveal functions of the CORs identified above, we performed Gene Ontology (GO) enrichment analysis. As shown in figure 2, the CORs in the two cultivars were enriched in the “molecular function” terms “transferase activity”, “protein binding”, “hydrolase activity”, “DNA binding”, “protein binding” and “hydrolase activity”. The CORs in Teqing and 02428 were both enriched in the “biological process” terms “metabolic process” and “cellular process”. Genes involved in “intracellular organelle”, “membranes” and “cytoplasm” were also significantly enriched among the “cellular component” GO terms. However, there were no GO terms enriched only in Teqing or 02428.

 Figure 2: Gene Ontology (GO) analysis of COR genes in the “molecular function” (A), “biological process” (B) and “cell component” (C) GO classes in Teqing and 02428 at 4°C and 10°C for 3hr and 24hr. The y-axis and x-axis indicate the percentage of genes in a GO category and the names of the clusters, respectively. The numbers before the letter (I or J) indicate the temperatures, the numbers after the letter (I or J) indicate the hours of exposure, J and I mean Teqing and 02428, respectively. For example, ‘CORs in 4I3’ means differentially expressed genes in Teqing between before (0hr) and after (3hr) cold treatment (4I3 vs I0, |log2Ratio| ?1 and FDR ≤ 0.001).

According to the KEGG (Kyoto Encyclopedia of Genes and Genomes) analysis, there were 15 significantly enriched KEGG pathways in both cultivars (Figure 3), included “plant hormone signal transduction”, “phenylpropanoid biosynthesis”, “spliceosome”, “ABC transporters”, “glycerophospholipid metabolism”, “Glycosylphosphatidylinositol (GPI)-anchor biosynthesis”, “phenylalanine metabolism”, “limonene and pinene degradation”, “isoflavonoid biosynthesis”, “glutathione metabolism”, “glycerolipid metabolism”, “arginine and proline metabolism”, “aminoacyl-tRNA biosynthesis”, “fatty acid elongation in mitochondria” and “vitamin B6 metabolism”. Among these pathways, the pathway with the most significant representation was “plant hormone signal transduction”. These results suggest that plant hormones could play an important role in regulating cold tolerance in rice.
Figure 3: KEGG pathway analysis of CORs in Teqing and 02428.
BG refers to all the genes in the background. ‘CORs in 4I3’ means differentially expressed genes in Teqing between before (0hr) and after (3hr) cold treatment (4I3 vs I0, |log2Ratio| ≥1 and FDR ≤ 0.001). By analogy, the numbers before the letter (I or J) indicate the temperatures, the numbers after the letter (I or J) indicate the hours of exposure; J and I mean Teqing and 02428, respectively.* indicates significantly enriched KEGG pathways with correct P value ≤0.05.

Identification of consistent differentially-expressed COR genes in multiple rice cultivars

CORs that are present in many rice varieties and that respond to different low temperature conditions may be key genes in cold response pathways. To screen such genes, we used Venn diagrams to classify CORs in Teqing and 02428. Based on the results shown in figure 4, we conducted a further Venn diagram analysis and found a total of 171 CORs that were shared by both cultivars and consistently responded to the four different cold treatments. Thus, they may represent crucial cold responsive genes that act in specific pathways.

 Figure 4: Venn diagrams showing the shared COR genes in Teqing and 02428 after exposure to cold treatments. (A) 563 COR genes respond to cold in both cultivars at 4°C for 3 hours; (B) 3361 COR genes respond to cold in both cultivars at 4°C for 24 hours; (C) 716 COR genes respond to cold in both cultivars at 10°C for 3 hours; (D) 3219 COR genes respond to cold in both cultivars at 10°C for 24 hours (CORs only in Teqing are shown in blue circle with the number written in white, CORs only in 02428 are shown in yellow circle with the number written in black, CORs shared by both cultivars are shown in center part of the Venn diagrams with the number written in white). (E) among all the common CORs we get from picture A, B, C and D, 171 COR genes respond to cold in both cultivars in all these cold treatments (Words above the pictures indicate the source of CORs. For example, common CORs in 4°C/3hr means CORs shared by Teqing and 02428, when they are treated in 4°C for 3 hours).

To find the key genes representing the specific pathways, we screened out 13 of the 171 genes we found that were included in predicted cold-response pathways (Supplemental Data 1). To confirm the cold-induced expression of these genes, gene-specific primers were designed for all13 genes and qRT-PCR assays were performed (Figure 5). The expression patterns of these 13 genes agreed with the RNA-seq results, confirming the high quality of the datasets. Among these genes,OsABCB5(LOC_Os01g50100, encoding a member of ABCB subgroup 5 in ABC transporter super family), OsPHIL2 (LOC_Os02g52010, encoding a member of Phi-1-like phosphate induced protein) and OsDREB1A (LOC_Os09g35030, encoding a member of the CBF/DREB1 protein family) had higher expression levels in Teqing than in 02428, while OsCDPK7 (LOC_Os04g49510, encoding a member of calcium-dependent protein kinase family) and OsREM4.1 (LOC_Os07g38170, encoding a member of the remorin family) had significantly lower expression levels in Teqing than in 02428 in response to the cold treatments.

 Figure 5: qRT-PCR analysis of 13 selected COR genes identified by RNA-seq in leaves of Teqing and 02428 plants before (0hr) and after 3 h or 24 h of cold treatment at 4°C or 10°C. The data of expression levels were separately compared with the Teqing sample before cold treatment (0hr). OsActin was used as an internal standard (Student’s t test, *P≤0.05, **P≤0.01, n=3).

To determine whether the observed differences in expression for the five genes are associated with cold tolerance in other cultivars of the two subspecies, we included an additional five cold tolerant japonica rice cultivars (Nipponbare, Ta Hung Ku, Taibei309, Baber, Zhonghua 11) and five cold sensitive indica rice cultivars (9311, Lal-Aman, Pao-Tou-Hung, Ai-Chiao-Hong, Guan-Yin-Tsan) to the study. Details of the cold responses of these 10cultivars are shown in figure S1. All of the japonica rice cultivars showed significantly stronger cold tolerance than the indica cultivars, as expected. We then we assayed the expression levels of the five genes (OsABCB5, OsPHIL2, OsDREB1A, OsCDPK7 and OsREM4.1) in the 10 cultivars. Among these five genes, the expression profiles of OsABCB5, OsPHIL2, OsDREB1A and OsREM4. 1still showed clear distinctions between the indica and japonica subspecies, while OsCDPK7 did not show the indica-japonica difference in the multi-cultivar comparison (Figure 6). In the initial comparison between Teqing and 02428, OsABCB5, OsPHIL2 and OsDREB1Awere expressed at higher levels in Teqing than in 02428; similarly, in the comparison of the five cultivar pairs, the expression of the three genes was also higher in indica than in japonica. OsREM4.1 showed higher expression levels in the japonica rice cultivars, confirming the results from the comparison of Teqing and 02428. After cold treatment, expression of OsABCB5, OsPHIL2, OsDREB1A and OsREM4. 1 was significantly higher than in the control, and the expression of OsPHIL2 was up-regulated >400-fold after 4°C treatment in the indica cultivars.


Figure 6: Relative expression levels determined by qRT-qPCR of the OsABCB5 (A), OsREM4.1 (B), OsPHIL2 (C), OsDREB1A (D) and OsCDPK7 (E) genes in leaves of 6 indica and 6 japonica plants before (0hr) and after 3 h or 24 h of cold (4°C or 10°C) treatment. Values are the averages of three samples ± SE. Error bars indicate the standard deviation. Mean values with one asterisk were found to be significantly different by Student’s t test (*P≤0.05; NS, not significant; n=3).

Analysis of natural variation in these COR genes in multiple rice cultivars

To investigate the natural variations in the sequences of OsABCB5, OsPHIL2, OsDREB1A and OsREM4.1, we extracted the sequences of their 2.0-kb promoter regions and coding regions from the re-sequencing dataset of a diversity panel of 575 rice cultivars that included 294 indica and 239 japonica lines. The original sequencing datasets have been deposited in the Genome Sequence Archive of the Bejing Institute of Genomics, Chinese Academy of Sciences ( under Accession numbers CRA000778, CRA000779 and CRA000995 (Figure 7).

 Figure 7: Schematic diagrams showing the position of SNP/Indels in the2.0-kb promoter regions and Coding Regions (CDS) of the OsABCB5 (A), OsPHIL2 (B), OsDREB1A (C) and OsREM4.1(D) genes. The phylogenetic tree was generated by MEGA5.0 software using neighbor joining method, with the segment length of interior branches indicating bootstrap values (1000 replications). The scale shows nucleotide substitutions per site. In gene schematic diagrams, gray bars indicate exons, and black lines indicate introns, the red vertical lines indicate the SNPs in coding regions while the black vertical lines indicate the SNPs in promoter regions. In the left table, Aus, IndicaJaponica and MIX represent the ecotypes of cultivated rice, such as aus, indicajaponica and mixture of other ecotypes, and the numbers in each block show the number of varieties. In the right table, the position of TSS is 0, negative and positive numbers indicate the positions of polymorphic sites in the promoter and coding region, respectively.

 Based on nucleotide polymorphisms, we could divide the sequences of the 575 lines into four haplotypes for OsABCB5; among these, HapII was the most represented in japonica lines, and Hap I was the most represented in indica lines. There were 12 SNP/Indel variations identified in the OsABCB5 promoter region between indica and japonica lines. Of these, four were predicted to be located within low temperature responsive cis-elements. For example, Hap II and Hap I differed in the nucleotides located at positions -642 G/A, -1633 A/G, -1638 G/A, and -1808 G/A upstream of the start codon; these four SNPs were related to theMYBST1, WRKY71OS/ARR1AT, WBBOXPCWRKY1 and LTRECOREATCOR15motifs, respectively, that are predicted to function in the cold response.

Based on the 79 polymorphic nucleotide sites, a total of five haplotypes of OsPHIL2were identified; among them, Hap I and Hap III were unique to japonica and indica rice lines, respectively. In the promoter region of OsPHIL2, we found 8SNP variations between Hap I and Hap III, and of these, no SNP was related to the cold response directly, but both -639A/G and -1136C/T were predicted to be associated with cytokinin.
By analyzing the candidate regulatory region of OsDREB1A in the 575 accessions, we found four haplotypes; Hap I and Hap II were unique to japonica and indica lines, respectively. In the promoter region, we identified 11SNPs/In Dels, among them, seven SNPs are located at -180, -741, -771, -720, -879, -1009 and -1070bp upstream of the TSS (transcription start site), and correspond to G-A, T-C, T-C, T-C, T-C, T-C, C-T and G-A changes, respectively. The single In Delis 11bpin length (CACAAGCAAAG) and is located 1346-1356bp upstream of the TSS. The SNPs located at positions -180 and -1070 were predicted to be related to auxin, SNPs located at -771 and -1009 were predicted to be related to MYB, SNPs at -741, -879, and the In Del were related to WRKY71OS, ARR1AT,and MYB1AT, respectively. Among these motifs, WRKY71OS and MYB1AT were related to MYB while ARR1AT was related to cytokine in.
We analyzed the promoter and CDS region of OsREM4.1by DNA sequencing and identified 48 variations. Analysis of these variations indicated that two major haplotypes were present in the 575 accessions; Hap I was present in 217 indica rice lines, and Hap II was present in 119 japonica rice lines. We found four In Dels and five SNPs in the candidate regulatory region. There is a SNP located at -483upstream of the TSS, corresponding to a G-C change and related to the ABRELATERD1motifwhichis related to ABA.


The japonica subspecies of Oryza sativa has adapted to a temperate climate, while the indica subspecies has adapted to tropical and subtropical environments. As a result, japonica rice cultivars are generally more tolerant to cold stress than indica rice cultivars (Figure S1). In our study, we used comparative transcriptome analysis on two rice cultivars to identify a group of COR genes in these two cultivars. From the CORs, we further identified four key genes, OsABCB5, OsPHIL2, OsDREB1A and OsREM4.1, that showed ubiquitous and distinctive activity in indica and japonica cultivars in an integrative study using KEGG and GO annotations, Venn diagram analysis, qRT-PCR, and natural variation analysis.

The DREB1/CBF pathway is a well-known cold-responsive pathway in plants [23]. The pathway includes the CBF1, CBF2 and CBF3 genes, which are also known as DREB1B, DREB1C and DREB1A, respectively, in Arabidopsis [24]. However, not all three transcription factors act as positive regulators in the CBF/DREB1 pathway. In Arabidopsis, freezing tolerance of the cbf2 mutant is enhanced by up-regulating the expression of CBF1/DREB1B and CBF3/DREB1A, which means that CBF2/DREB1C plays a negative role in the cold response [25]. OsDREB1A has been reported previously to be induced by cold stress [26,27]. Consistent with these results, we found that the expression of OsDREB1Aincreased after cold treatment. Lourenco et al., have claimed that in OsHOS1 knockdown plants, the reduced expression of OsHOS1 promoted the accumulation of OsICE1 protein and resulted in the up-regulation of OsDREB1A. However, the transgenic plants did not show increased cold tolerance [28]. In our study, the expression level of OsDREB1A was higher in indica rice compared to japonica rice. Therefore, we speculate that OsDREB1Amay play a negative role in cold signaling. In our natural variation analysis, we discovered a SNP related to MYB in the promoter region of OsDREB1Ain HapI (japonica). Interestingly, a MYB-like transcription factor has been reported to be involved in the cold regulation of CBF genes, and MYB15 over expression reduced the expression of CBF genes under cold treatment [29]. Thus, the MYB negative regulation of OsDREB1A may result in the stronger cold tolerance observed in japonica rice cultivars.

Among the four genes we screened, OsREM4.1 is related to the plant hormone ABA, which is known as a cold stress-responsive hormone. Two SNP variations within ABA responsive or MYB elements were located within the promoter region of OsREM4.1. Interestingly, OsREM4.1 was previously characterized and found to be up-regulated by ABA signaling; the OsREM4.1 protein can bind to OsSERK1 and inhibit the formation and activation of the OsBRI1-OsSERK1 receptor complex, which is crucial for BR signaling [30]. In Arabidopsis thaliana, application of exogenous ABA can enhance cold tolerance and induce the expression of several COR genes, including RAB18, RD29A and KIN2 [31,32]. Several studies showed that plants treated with BR grew better at low temperature compared to optimal conditions [33]. In addition, Li et al., described two key components of BR signaling, BIN2 and BZR1, regulate plant freezing tolerance; BIN2 and its homologs play negative roles, while BZR1 plays a positive role in the cold response [14]. In our study, OsREM4.1hada higher level of expression in japonica rice, so we predicted that OsREM4.1 may play a positive role in regulating cold tolerance via crosstalk between the ABA and BR pathways.

Another gene we found that is associated with the BR pathway is OsPHIL2, which encodes a Phosphate-Induced protein (PHI). Interestingly, two BR-response genes in Arabidopsis, EXORDIUM (EXO) and EXORDIUM-LIKE1 (EXL1), were previously reported to have the same PHI conserved region [34]. Expression of the EXO gene is up-regulated by BR and EXL1 has a similar expression pattern to EXO in Arabidopsis in that it can be induced by BR [35,36]. We therefore predicted that OsPHIL2 may be similarly connected to BR signaling pathway. Considering that OsPHIL2 is expressed at higher levels in indica cultivars when exposed to low temperature, OsPHIL2 may negatively regulate cold stress via BR signaling.

ABC transporters constitute one of the largest protein families, and are present in organisms ranging from bacteria to humans [37]. In plants, ABC proteins were originally identified as transporters involved in the final detoxification process [38]. OsABCB5 is predicted to encode an MDR-like ABC transporter. There are several studies that have reported that MDR-like ABC transporters are involved in the auxin signaling pathway. Noh et al., showed that two MDR-like genes, At MDR1 and AtPGP1, are required for auxin transport and auxin-mediated development in Arabidopsis. Also, mutation of two related MDR-like genes (MDR1 and PGP1) resulted in reduced polar auxin transport [39]. Another study found that the Arabidopsis MDR-like ABC transporter4 (AtPGP4) is involved in auxin-mediated root development [40]. Although the role of auxin in cold stress remains unclear, there are reports that potentially link cold stress to auxin by showing that cold stress inhibits the inflorescence gravity response in Arabidopsis thaliana [41,42]. An early study demonstrated that temperature affects the speed of exogenous auxin transport in a variety of plant species [43]. Considering that auxin is at the center of hormonal crosstalk, auxin may mediate the response to cold by interacting with other hormones [44]. These results indicate that OsABCB5may play a crucial role in the cold response in rice seedlings by transporting auxin.

Based on the expression profiles and predicted pathways of the four key genes involved in the cold response in rice, we can conduct further functional studies to determine the biological and biochemical functions of the proteins encoded by these genes.


Plant materials and Evaluation of Phenotypes

In this study, we used Teqing and 02428 as the cold-sensitive indica and cold-tolerant japonica cultivars, respectively. We also chose an additional 5 japonica and 6 indica cultivars for cold temperature treatment. Seeds were allowed to germinate for 3 days in distilled water, after which they were transferred to moist filter paper. The seedlings were grown in 0.5-fold Kimura B solution at pH 6.0 under a 13-h light (28°C)/11-h dark (25°C) photoperiod and 70-80% relative humidity [45]. Three-leaf rice seedlings were transferred to a climatic chamber (Percival Scientific, Perry, IA) and maintained at either 4°C or 10°C or for 4 days and 5 days, respectively. After the plants had recovered at 28°C /25°C day/night cycles for 7 days, the survival rates were determined based on the percentage of surviving seedlings [46]. All treatments were repeated at least three times, and the average of the replicates was used to characterize cold tolerance. The significant differences between subspecies were analyzed using ANOVA or Student’s t-test.

RNA isolation

Total RNA was isolated from plant tissue using TRIzol® reagent (Invitrogen, Carlsbad, CA, USA). The RNA quality was assessed visually by electrophoresis on agarose gels, and the concentration was determined with a NanoDrop spectrophotometer (ND-1000, NanodropTechnologies, Wilmington, DE, USA).

Library construction and Illumina cDNA sequencing

Whole transcriptome analysis was performed using RNA isolated from three-leaf rice seedlings that were cold-treated at 10°C or 4°C for 0, 3 and 24 hours under continuous illumination. We extracted total RNA, used oligo (dT) magnetic beads to purify the mRNA, and then used Oligo (dT) as a primer to synthesize cDNA. The 5' ends of the cDNA tags can be generated with two different restriction end nucleases, NlaIII or DpnII, Usually, the bead-bound cDNA fragments are digested with NlaIII, which recognizes and cuts atthe CATG sites. The 3' cDNA fragments connected to the oligo (dT) beads were washed away and the Illumina adaptor 1 was ligated to the 5'sticky ends of the digested bead-bound cDNA fragments. The junction of the Illumina adaptor 1 and the CATG site is an MmeI recognition site, and it cuts 17bp downstream of the CATG site, producing tags containing adaptor 1. After removing the 3' fragments with magnetic beads precipitation, Illumina adaptor 2 was ligated to the 3' ends of the tags, resulting in tags with different adaptors on each end to form a tag library. Following linear PCR amplification of the sample libraries, the double-stranded cDNA fragments are purified by polyacrylamide gel electrophoresis prior to high-throughput cDNA sequencing using the Illumina HiSeqTM 2000 instrument.

Identification of COR genes

The genes with differential expression between 0 hours and 3/24 hours of cold treatments are Called Cold Responsive genes (CORs). Before comparing the differential expression of genes in response to cold treatment, normalized gene expression levels were obtained by normalizing the number of raw clean tags in each library to the number of Transcripts Per Million clean tags (TPM). A rigorous algorithm method was performed to detect differential expression of genes across samples. A combination of FDR ≤ 0.001 and the absolute value of log 2 Ratio ≥1 were used as the threshold to determine the significance of differentially-expressed genes. GO and pathway enrichment analyses were based on the agriGO ( and KEGG pathway databases ( [47-49]. Cluster analysis was performed with CLUSTER3.0 and viewed using the TREEVIEW software program (

Quantitative reverse transcription polymerase chain reaction (qRT-PCR) analysis of gene expression

An independent experiment was performed to examine the expression of putative COR genes using qRT-PCR. Total RNA was extracted from rice leaves before and after exposure to 4°C or 10°C for 3h and 24h using the TRIzol Reagent (Invitrogen).cDNA was synthesized and qRT-PCR was performed according to the previous report [50]. All reactions were performed with three biological and three technical replicates. OsActin was used as an internal control. The sequences of primers used for the qRT-PCR assays are listed in table S1. The gene expression levels were calculated using the 2-??Ct method. Student’s t-test was used to determine whether there were significant differences between samples.








Phosphate-induced protein 1 conserved region domain containing protein





ABC transporter, ATP-binding protein










dehydration-responsive element-binding protein





ethylene-responsive transcription factor





OsSAUR19 - Auxin-responsive SAUR gene family member





Gibberellin receptor GID1L2





BBTI13 - Bowman-Birk type bran trypsin inhibitor precursor





CAMK includes calcium/calmodulin depedent protein kinases





cytochrome P450





OsFBX391 - F-box domain containing protein





Plastocyanin-like domain containing protein





TE_MEK_ste7_MAP2K.5 - STE kinases include homologs to sterile 7, sterile 11 and sterile 20 from yeast








Table S1: Sequences of the oligonucleotide primers used in the qRT-PCR assays.

Haplotype analysis and SNP discovery of the four key COR genes

Haplotype analysis was performed on the genomic sequences that included the OsABCB5, OsPHIL2, OsDREB1A and OsREM4. 1 promoter and coding regions from the 575 rice accessions. A neighbor-joining tree based on the nucleotide sequences of the four key genes from the rice varieties was constructed using MEGA 5.0. For a genome-wide association study, we used our own SNP database to performed SNP analysis. Additionally, we used the Genome Information Database System ( to predict the functions of the SNPs that we found.


DM conceived and designed the experiments. XH and DM performed the experiments. XH and YT analyzed the data, and XH wrote the manuscript.


We thank Dr. Caiyan Chen for his help in conceiving the project. This work was supported by grants from the National Natural Science Foundation of China (31101211, 31371603).


  1. Aghaee A, Moradi F, Zare-Maivan H, Zarinkamar F, Irandoost HP, et al. (2011) Physiological Responses of Two Rice (Oryza Sativa L.) Genotypes to Chilling Stress at Seedling Stage. African J Biotechnol 10: 7617-7621.
  2. Song SY, Chen Y, Chen J, Dai XY, Zhang WH (2011) Physiological Mechanisms Underlying OsNAC5-Dependent Tolerance of Rice Plants to Abiotic Stress. Planta 234: 331-345.
  3. Uemura M, Joseph RA, Steponkus PL (1995) Cold Acclimation of Arabidopsis Thaliana (Effect on Plasma Membrane Lipid Composition and Freeze-Induced Lesions). Plant Physiol 109: 15-30.
  4. Ye H, Du H, Tang N, Li X, Xiong L (2009) Identification and Expression Profiling Analysis of TIFY Family Genes Involved in Stress and Phytohormone Responses in Rice. Plant Mol Biol 71: 291-305.
  5. Thomashow MF (1999) Plant Cold Acclimation: Freezing Tolerance Genes and Regulatory Mechanisms. Annu Rev Plant Physiol Plant Mol Biol 50: 571-599.
  6. Rihan HZ, Al-issawi M, Fuller MP (2017) Advances in Physiological and Molecular Aspects of Plant Cold Tolerance. J PLANT Interact 12: 143-157.
  7. Shi Y, Ding Y, Yang S (2015) Cold Signal Transduction and Its Interplay with Phytohormones during Cold Acclimation. Plant Cell Physiol 56: 7-15.
  8. Ito Y, Kurata N (2006) Identification and Characterization of Cytokinin-Signalling Gene Families in Rice. Gene 382: 57-65.
  9. Su CF, Wang YC, Hsieh TH, Lu CA, Tseng TH, et al. (2010) A Novel MYBS3-Dependent Pathway Confers Cold Tolerance in Rice. Plant Physiol 153: 145-158.
  10. Kurepin LV, Dahal KP, Savitch LV, Singh J, Bode R, et al. (2013) Role of CBFs as Integrators of Chloroplast Redox, Phytochrome and Plant Hormone Signaling during Cold Acclimation. Int J Mol Sci 14: 12729-12763.
  11. Hu Y, Jiang L, Wang F, Yu D (2013) Jasmonate Regulates The Inducer Of Cbf Expression-C-Repeat Binding Factor/Dre Binding Factor1 Cascade and Freezing Tolerance in Arabidopsis. Plant Cell 25: 2907-2924.
  12. Achard P, Gong F, Cheminant S, Alioua M, Hedden P, et al. (2008) The cold-inducible CBF1 factor-dependent signaling pathway modulates the accumulation of the growth-repressing DELLA proteins via its effect on gibberellin metabolism. Plant Cell 20: 2117-2129.
  13. Jeon J, Kim NY, Kim S, Kang NY, Novák O, et al. (2010) A Subset of Cytokinin Two-Component Signaling System Plays a Role in Cold Temperature Stress Response in Arabidopsis. J Biol Chem 285: 23371-23386.
  14. Li H, Ye K, Shi Y, Cheng J, Zhang X,et al. (2017) BZR1 Positively Regulates Freezing Tolerance via CBF-Dependent and CBF-Independent Pathways in Arabidopsis. Mol Plant 10: 545-559.
  15. Zhang Z, Huang R (2010) Enhanced Tolerance to Freezing in Tobacco and Tomato Overexpressing Transcription Factor TERF2/LeERF2 Is Modulated by Ethylene Biosynthesis. Plant Mol Biol 73: 241-249.
  16. Lv Y, Guo Z, Li X, Ye H, Li X, et al. (2016) New Insights into the Genetic Basis of Natural Chilling and Cold Shock Tolerance in Rice by Genome-Wide Association Analysis. Plant Cell Environ 39: 556-570.
  17. Bai B, Wu J, Sheng WT, Zhou B, Zhou LJ, et al. (2015) Comparative Analysis of Anther Transcriptome Profiles of Two Different Rice Male Sterile Lines Genotypes under Cold Stress. Int J Mol Sci 16: 11398-11416.
  18. Shen C, Li D, He R, Fang Z, Xia Y, et al. (2014) Comparative Transcriptome Analysis of RNA-Seq Data for Cold-Tolerant and Cold-Sensitive Rice Genotypes under Cold Stress. J Plant Biol 57: 337-348.
  19. Dametto A, Sperotto RA, Adamski JM, Blasi ÉA, Cargnelutti D, et al. (2015) Cold Tolerance in Rice Germinating Seeds Revealed by Deep RNAseq Analysis of Contrasting Indica Genotypes. Plant Sci 238: 1-12.
  20. Byun MY, Cui LH, Lee J, Park H, Lee A, et al. (2018) Identification of Rice Genes Associated With Enhanced Cold Tolerance by Comparative Transcriptome Analysis With Two Transgenic Rice Plants Overexpressing DaCBF4 or DaCBF7, Isolated From Antarctic Flowering Plant Deschampsia antarctica. Front Plant Sci 9: 601.
  21. Sperotto RA, de Araújo Junior AT, Adamski JM, Cargnelutti D, Ricachenevsky FK, et al. (2018) Deep RNAseq Indicates Protective Mechanisms of Cold-Tolerant Indica Rice Plants during Early Vegetative Stage. Plant Cell Rep 37: 347-375.
  22. Kawahara Y, de la Bastide M, Hamilton JP, Kanamori H, McCombie WR, et al. (2013) Improvement of the Oryza Sativa Nipponbare Reference Genome Using next Generation Sequence and Optical Map Data. Rice 6: 4.
  23. Ito Y, Katsura K, Maruyama K, Taji T, Kobayashi M, et al. (2006) Functional Analysis of Rice DREB1/CBF-Type Transcription Factors Involved in Cold-Responsive Gene Expression in Transgenic Rice. Plant Cell Physiol 47: 141-153.
  24. Liu Q, Kasuga M, Sakuma Y, Abe H, Miura S, et al. (1998) Two Transcription Factors, DREB1 and DREB2, with an EREBP/AP2 DNA Binding Domain Separate Two Cellular Signal Transduction Pathways in Drought- and Low-Temperature-Responsive Gene Expression, Respectively, in Arabidopsis. Plant Cell 10: 1391-1406.
  25. Novillo F, Alonso JM, Ecker JR, Salinas J (2004) CBF2/DREB1C is a negative regulator of CBF1/DREB1B and CBF3/DREB1A expression and plays a central role in stress tolerance in Arabidopsis. Proc Natl Acad Sci USA 101: 3985-3990.
  26. Dubouzet JG, Sakuma Y, Ito Y, Kasuga M, Dubouzet EG, et al. (2003) OsDREB genes in rice, Oryza sativa L., encode transcription activators that function in drought-, high-salt- and cold-responsive gene expression. Plant J 33: 751-763.
  27. Challam C, Ghosh T, Rai M, Tyagi W (2015) Allele Mining across DREB1A and DREB1B in Diverse Rice Genotypes Suggest a Highly Conserved Pathway Inducible by Low Temperature. J Genet 94: 231-238.
  28. Lourenço T, Sapeta H, Figueiredo DD, Rodrigues M, Cordeiro A, et al. (2013) Isolation and characterization of rice (Oryza sativa L.) E3-ubiquitin ligase OsHOS1 gene in the modulation of cold stress response. Plant Mol Biol 83: 4-5.
  29. Agarwal M, Hao Y, Kapoor A, Dong CH, Fujii H, et al. (2006) A R2R3 Type MYB Transcription Factor is Involved in the Cold Regulation of CBF Genes and in Acquired Freezing Tolerance. J Biol Chem 281: 37636-37645.
  30. Gui J, Zheng S, Liu C, Shen J, Li J, et al. (2016) OsREM4.1 Interacts with OsSERK1 to Coordinate the Interlinking between Abscisic Acid and Brassinosteroid Signaling in Rice. Dev Cell 38: 201-213.
  31. Kurkela S, Borg-Franck M (1992) Structure and expression of kin2, one of two cold- and ABA-induced genes of Arabidopsis thaliana. Plant Mol Biol 19: 689-692.
  32. Lång V, Palva ET (1992) The Expression of a Rab-Related Gene, Rab18, Is Induced by Abscisic Acid during the Cold Acclimation Process of Arabidopsis Thaliana (L.) Heynh. Plant Mol Biol 20: 951-962.
  33. Kamuro Y, Takatsuto S (1991) Capability for and Problems of Practical Uses of Brassinosteroids. Acs Symp Ser 292-297.
  34. Schröder F, Lisso J, Müssig C (2012) Expression Pattern and Putative Function of EXL1 and Homologous Genes in Arabidopsis. Plant Signal Behav 7: 22-27.
  35. Coll-Garcia D, Mazuch J, Altmann T, Müssig C (2004) EXORDIUM Regulates Brassinosteroid-Responsive Genes. FEBS Lett 563: 82-86.
  36. Schröder F, Lisso J, Lange P, Müssig C (2009) The Extracellular EXO Protein Mediates Cell Expansion in Arabidopsis Leaves. BMC Plant Biol 12: 1-12.
  37. Henikoff S, Greene EA, Pietrokovski S, Bork P, Attwood TK, et al. (1997) Gene Families: The Taxonomy of Protein Paralogs and Chimeras. Science 278: 609-614.
  38. Martinoia E, Grill E, Tommasini R, Kreuz K, Amrhein N, (1993) ATP-Dependent Glutathione S-Conjugate “export” Pump in the Vacuolar Membrane of Plants. Nature 364: 247-249.
  39. Noh B, Bandyopadhyay A, Peer WA, Spalding EP, Murphy AS (2003) Enhanced Gravi- and Phototropism in Plant Mdr Mutants Mislocalizing the Auxin Efflux Protein PIN1. Nature 423: 999-1002.
  40. Santelia D, Vincenzetti V, Azzarello E, Bovet L, Fukao Y, et al. (2005) MDR-like ABC Transporter AtPGP4 Is Involved in Auxin-Mediated Lateral Root and Root Hair Development. FEBS Lett 579: 5399-5406.
  41. Fukaki H, Fujisawa H, Tasaka M (1996) Gravitropic Response of Inflorescence Stems In Arabidopsis thaliana. Plant Physiol 110: 933-943.
  42. Wyatt SE, Rashotte AM, Shipp MJ, Robertson D, Muday GK (2002) Mutations in the Gravity Persistence Signal Loci in Arabidopsis Disrupt the Perception and/or Signal Transduction of Gravitropic Stimuli. Plant Physiol 130: 1426-1435.
  43. Morris DA (1979) The Effect of Temperature on the Velocity of Exogenous Auxin Transport in Intact Chilling-Sensitive and Chilling-Resistant Plants. Planta 146: 603-605.
  44. Rahman A (2013) Auxin: A Regulator of Cold Stress Response. Physiol Plant 147: 28-35.
  45. Yamaji N, Xia J, Mitani-Ueno N, Yokosho K, Feng Ma J (2013) Preferential delivery of zinc to developing tissues in rice is mediated by P-Type heavy metal ATPase OsHMA2. Plant Physiol 162: 927-939.
  46. Mao D, Yu L, Chen D, Li L, Zhu Y, et al. (2015) Multiple Cold Resistance Loci Confer the High Cold Tolerance Adaptation of Dongxiang Wild Rice (Oryza Rufipogon) to Its High-Latitude Habitat. Theor Appl Genet 128: 1359-1371.
  47. Du Z, Zhou X, Ling Y, Zhang Z, Su Z (2010) AgriGO: A GO Analysis Toolkit for the Agricultural Community. Nucleic Acids Res 38: 64-70.
  48. Ogata H, Goto S, Sato K, Fujibuchi W, Bono H, et al. (1999) KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Research 28: 27-30.
  49. Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, et al. (2008) KEGG for Linking Genomes to Life and the Environment. Nucleic Acids Res 36: 480-484.
  50. Tan L, Zhu Y, Fan T, Peng C, Wang J, et al. (2019) OsZIP7 functions in xylem loading in roots and inter-vascular transfer in nodes to deliver Zn/Cd to grain in rice. Biochem Biophys Res Commun 512: 112-118.

Citation: Mao D (2019) Four Cold Responsive Genes That Show Consistent Differential Expression between the Japonica and Indica subspecies of rice. J Genet Genomic Sci 4: 011.

Copyright: © 2019  Donghai Mao, et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

© 2022, Copyrights Herald Scholarly Open Access. All Rights Reserved!