Accepted Article

Global transcriptional profiling of a cold-tolerant rice variety under moderate cold stress reveals different cold stress response mechanisms

Junliang Zhaoa,b,†, Shaohong Zhanga,b,†, Tifeng Yanga,b, Zichong Zenga,b, Zhanghui Huanga,b, Qing Liua,b, Xiaofei Wanga,b, Jan Leachc, Hei Leungd and Bin Liua,b,*

a

Rice Research Institute, Guangdong Academy of Agricultural Sciences, Guangzhou 510640, China

b

Guangdong Provincial Key Laboratory of New Technology in Rice Breeding, Guangdong Academy

of Agricultural Sciences, Guangzhou 510640, China c

Bioagricultural Sciences and Pest Management and Program in Plant Molecular Biology, Colorado

State University, Fort Collins, Colorado 80523–1177 d

Plant Breeding, Genetics, and Biotechnology, International Rice Research Institute, Los Banos, 4031

Laguna, Philippines

*Corresponding author, e-mail: [email protected] †

These authors contributed equally

This article has been accepted for publication and undergone full peer review but has not been through the copyediting, typesetting, pagination and proofreading process, which may lead to differences between this version and the Version of Record. Please cite this article as doi: 10.1111/ppl.12291

This article is protected by copyright. All rights reserved

Accepted Article

Gene expression profiling under severe cold stress (4°C) has been conducted in plants including rice. However, rice seedlings are frequently exposed to milder cold stresses under natural environments. To understand the responses of rice to milder cold stress, a moderately low temperature (8°C) was used for cold treatment prior to genome-wide profiling of gene expression in a cold-tolerant japonica variety, Lijiangxintuanheigu (LTH). A total of 5557 differentially expressed genes (DEGs) were found at four time points during moderate cold stress. Both the DEGs and differentially expressed transcription factor genes clustered into two groups based on their expression, suggesting a two-phase response to cold stress and a determinative role of transcription factors in the regulation of stress response. The induction of OsDREB2A under cold stress was reported for the first time in the present study. Among the anti-oxidant enzyme genes, glutathione peroxidase (GPX) and glutathione S-transferase (GST) were up-regulated, suggesting that the glutathione system may serve as the main reactive oxygen species (ROS) scavenger in LTH. Changes in expression of genes in signal transduction pathways for auxin, abscisic acid (ABA) and salicylic acid (SA) imply their involvement in cold stress responses. The induction of ABA response genes and detection of enriched cis-elements in DEGs suggest that ABA signaling pathway plays a dominant role in the cold stress response. Our results suggest that rice responses to cold stress vary with the specific temperature imposed and the rice genotype.

Abbreviations – ABA, abscisic acid; ABRE, abscisic acid responsive element; APX, ascorbate peroxidase;

ARF,

auxin

response

factors;

DEG,

differentially

expressed

gene;

LTH,

Lijiangxintuanheigu; FDR, false discovery rate; GO, Gene Ontology; GPX, glutathione peroxidase; GST, glutathione S-transferase; qRT-PCR, real-time quantitative PCR; QTL, quantitative trait locus; ROS, reactive oxygen species; SA, salicylic acid; SOD, superoxide dismutase.

Introduction Rice (Oryza sativa) originated from tropical or subtropical areas and is generally sensitive to low

This article is protected by copyright. All rights reserved

Accepted Article

temperatures. Cold damage can occur at different developmental stages in rice. In early spring in temperate or subtropical environments, when rice is at the seedling stage, chilling injury can lead to slow growth, delayed crop maturation, poor establishment, and subsequently, decrease in yield. With global climate changes, rice will be more frequently exposed to extreme weather conditions (low and high temperatures). Therefore, improving the cold tolerance of rice at the seedling stage is an important target in rice breeding programs. The rice response to low temperature stress is complex. Compared with cold sensitive rice seedlings, cold tolerant seedlings exhibit lower reductions in chlorophyll content under cold stress (Dai et al. 1990, Guh 2002, Zhan et al. 2005). Furthermore, the cold tolerant seedlings are more vigorous (Han et al. 2007), and show higher leaf water content (Guh 2002, Kabaki and Tajima 1981), higher anti-oxidative enzyme activity (Guh 2002, Huang and Guo 2005) and higher free polyamine levels (Lee et al. 1997). Genetic analysis suggests that cold tolerance at the seedling stage in rice is a quantitative trait controlled by multiple genes, and numerous quantitative trait loci (QTLs) for cold tolerance at seedling stage in rice have been identified (Andaya and Mackill 2003, Koseki et al. 2010, Zhang et al. 2014). The presence of many QTLs for cold tolerance combined with differences in injury symptoms and physiological and biochemical reactions suggest that diverse mechanisms are involved in cold tolerance at seedling stage in rice. Therefore, a systematic and high throughput genomic approach is needed to identify the genes responsible for cold tolerance and understand the cold response mechanisms before genetic manipulation can be effectively applied to improve cold tolerance in rice. Microarray technology provides a powerful tool for genome-wide identification of the functional genes and dissecting the response mechanisms of rice to cold stress. Microarray-assisted genome-wide gene expression profiling of cold stress has been performed in plants, such as Arabidopsis, rice, wheat, barley, and a number of cold responsive genes and different regulatory pathways in response to cold stress have been identified (Monroy et al. 2007, Seki et al. 2002, Tondelli et al. 2006, Yun et al. 2010). For rice plant, several studies on global gene expression profiling of cold stress at seedling stage have been conducted (Cheng et al. 2007, Rabbani et al. 2003, Yun et al. 2010, Zhang et al. 2012b). The comparative transcriptome profiling of cold stress in two contrasting rice genotypes revealed that the major effect on gene expression was up-regulation in the cold-tolerant genotype and strong repression in cold-sensitive genotype under cold stress (Zhang et al. This article is protected by copyright. All rights reserved

Accepted Article

2012b). Because of the importance of transcription factors in regulation of gene expression, transcription factor genes were a focus in these studies. Among the transcription factor genes, DREB1s, members of the well characterized dehydration responsive element binding proteins that function in abiotic stresses, were induced by cold, while DREB2s were induced by dehydration or salt stress in Arabidopsis (Sakuma et al. 2002). Similar results were observed in rice, i.e., OsDREB1A and OsDREB1B were induced by cold, while OsDREB2A was induced by dehydration and salt (Dubouzet et al. 2003). Genome-wide gene expression profiling of chilling stress in two contrasting rice genotypes also showed that OsDREB1A, OsDREB1B and OsDREB2B were up-regulated during cold stress (Zhang et al. 2012b). Anti-oxidative enzymes such as superoxide dismutase (SOD), catalase and ascorbate peroxidase (APX), which play important roles against cold stress in rice (Guo et al. 2006), were also shown in gene expression studies to be induced under cold stress (Yun et al. 2010). Transcription regulatory pathways involved in cold stress were predicted in various studies by analysis of differentially expressed transcription factors and enrichment of cis-elements in promoters of the differentially expressed genes; these pathways included the ROS-bZIP1 (Cheng et al. 2007), R2R3-MYB-dependent (Yun et al. 2010), and CBF and MYBS3 pathways (Zhang et al. 2012b). These studies facilitate our understanding of the regulation of gene expression in response to cold stress and the molecular mechanisms of cold tolerance in rice. However, many of these studies were conducted under severe cold stress, i.e., 4°C (Rabbani et al. 2003, Zhang et al. 2012a, Zhang et al. 2012b). Previous studies had shown that exposure to different low temperatures resulted in dissimilar injury symptoms, and tolerance to these different temperatures was conferred by separate QTLs (Andaya and Mackill 2003, Zhan et al. 2005, Zhang et al. 2014). Taken together, these studies imply that rice plants have diverse mechanisms to respond to different low temperatures. To uncover the molecular mechanism of cold tolerance at the seedling stage in rice under milder cold stress, we performed a genome-wide gene expression profiling of cold stress using a cold-tolerant japonica variety, Lijiangxintuanheigu (LTH) and high-density Agilent Rice Oligo chips. We used a moderately low temperature (8°C), which frequently occurs during rice growing season, for cold treatment, instead of a severe low temperature (4°C) commonly used in previous studies. We identified a number of novel genes and gene families that responded to 8°C cold stress. The cold responsive genes displayed a two-phase expression pattern based on time. A key difference from the previous studies was the observation that LTH likely uses the glutathione system as a reactive oxygen This article is protected by copyright. All rights reserved

Accepted Article

species (ROS) scavenger. Furthermore, our results suggest that abscisic acid (ABA)-dependent and ABA-independent signaling pathways are activated in rice under cold stress, but the ABA-dependent signaling pathway likely plays a more dominant role in the cold stress response in LTH.

Materials and methods Plant material, cold treatment and sample preparation LTH, a japonica variety from Yunnan, China which exhibits strong cold tolerance at different developmental stages (Ye et al. 2009, Zhang et al. 2012b, Zhang et al. 2014), was used for this study. After incubation at 49°C for 96 h to break dormancy, seeds were soaked in water for 36 h and then incubated at 30°C for 48 h. The germinated seeds were sown in plastic trays (41 x 26.5 x 7.5 cm) filled with fine field soil. Materials were planted in 13 rows in a tray with 11 plants per row. The seedlings were allowed to grow for 14 days at 26°C until the three-leaf stage. For the cold treatment, three-leaf stage seedlings were placed in a Conviron PGV36 growth chamber (Controlled Environments Ltd., Winnipeg) maintained at a constant temperature of 8°C, with a 13-h day length of 126 μmol m–2 s–1 light intensity and a relative humidity of 80 ± 5%. Control seedlings were grown under the same conditions, except the temperature was 26°C. The second and third leaves of seedlings were harvested at 6, 12, 24 and 48 h after cold treatment, quickly frozen in liquid nitrogen and stored at –80°C until use. Three biological replicates were performed for each pair of cold treatment and control at each time point.

Total RNA extraction Total RNA of rice leaves was extracted with Trizol reagent (Invitrogen, Carlsbad, CA, USA) and purified with NucleoSpin RNA Clean-up (MACHEREY-NAGEL, Germany). RNA quality and quantity

were

assessed

by

formaldehyde

denaturing

agarose

gel

electrophoresis

and

spectrophotometry (Nanodrop-1000, Thermo-Fisher, USA), respectively.

Microarray analysis The Agilent Rice Oligo Microarray 4×44K (Agilent-015058, Agilent Technologies, USA) was used for global gene expression analysis. RNA amplification and cDNA labeling were conducted using cRNA Amplification and Labeling Kit (CapitalBio, Beijing, China), according to the manufacturer's This article is protected by copyright. All rights reserved

Accepted Article

protocol. cDNA from control with three biological replicates were equally mixed and labeled with Cy3 dye as common reference in all of the microarray assays. cDNA from all the samples, including cold treatment and control at the four time points (6, 12, 24 and 48 h) were labeled with Cy5 dye, respectively, and paired with the common reference for microarray assay. Microarray slides were hybridized with hybridization buffer containing two paired labeled samples (cold treatment or control sample was paired with common reference) in a final concentration of 5×Denhardt, 3×SSC and 0.2% SDS for 16 h at 65°C. After hybridization, the slides were subjected to four-step stringency washes (2×SSC + 0.2% SDS, twice at 42°C for 4 min; 0.2×SSC, twice at 42°C for 4 min). The slides were scanned using Axon 4000B scanner and analyzed by Genepix Pro 6.0 (Axon Instruments, Watertown, Massachusetts, USA). Three biological replicates were performed for each pair of cold treatment and control.

Data analysis The raw data from microarrays was filtered by Genepix Pro 6.0 and normalized using the Lowess method in Acuity 4.0 (Axon Instruments, USA). The microarray contains 42 990 features representing 21 495 full length cDNA. Each cDNA has a unique Name and A_ID. That is, one A_ID has two features on the microarray, which can be considered as two technical replicates. The expression data for a given A_ID was the mean of the two features. The expression ratio of a given A_ID was expressed as fold change (fold change = expression data of cold treated sample / expression data of control sample). SAM (Significance Analysis of Microarray) software (Stanford University, CA, USA) was used to identify differentially expressed genes (DEGs) between cold treated samples and control samples at each time point (Tusher et al. 2001). DEGs were identified using the combined criteria of three-fold change and a cut-off of q-values (≤ 0.05) in SAM based on three biological replicates. Hierarchical cluster analysis of DEGs was performed using software Gene Cluster 3.0 (Eisen et al. 1998). Gene Ontology (GO) analysis was carried out using the Singular Enrichment Analysis tool offered by agriGO (http://bioinfo.cau.edu.cn/agriGO/index.php) at default settings of Fisher t-test (p < 0.05), false discovery rate (FDR) correction by Hochberg method and five minimum number of This article is protected by copyright. All rights reserved

Accepted Article

mapping entries against species specific pre-computed background reference (Du et al. 2010). Over represented cis-motifs in 2 kb upstream sequence from translational initiation codon of DEGs were analyzed using the Osiris program (http://www.bioinformatics2.wsu.edu/cgibin/Osiris/cgi/home.pl) (Morris et al. 2008). The original microarray data set of this study has been deposited in NCBI’s Gene Expression Omnibus (http://www.ncbi.nih.gov/geo/) under GEO series number GSE58132.

Validation of microarray data with real-time quantitative PCR (qRT-PCR) A total of 23 genes representing three types of expression data in microarray: up regulation, down regulation and no change, were selected for qRT-PCR validation. The same RNA used in microarray assays was used in qRT-PCR validation. Total RNA was reversely transcribed into cDNA using PrimeScript reverse transcriptase (Takara, Japan). The real-time PCR was conducted using Ex-Taq SybrGreen PCR Mix (Takara, Japan) in Roche 480 real-time PCR machine (Roche, Germany). Three biological replicates were performed. The primers used for qRT-PCR were listed in Table S1. The comparison of microarray and qRT-PCR assay data was based on the fold change of genes in certain time point between cold treated samples and control samples. The fold change from microarray and qRT-PCR were log (base 2)-transformed and then compared and plotted to test the similarity.

Results Validation of microarray data by qRT-PCR To validate the expression profile obtained from microarray, 23 genes including OsDREB1A, OsDREB1B, OsDREB2A, OsDREB2B and OsbZIP72 were selected for qRT-PCR analysis. We tested the similarity of the gene expression patterns detected by microarray and those by qRT-PCR (Figure S1). Good correlation between microarray and qRT-PCR data (r = 0.9442) indicated the reliability of the microarray experiments in this study.

Differentially expressed genes under cold stress To identify DEGs in response to moderate cold stress at the seedling stage in rice, a strong cold tolerant japonica variety, LTH was subjected to 8°C, and genome-wide differential gene expression This article is protected by copyright. All rights reserved

Accepted Article

analysis was performed in a time series that included four time points (6, 12, 24 and 48 h). A total of 5557 genes were differentially expressed at least at one time point during cold stress. Among these DEGs, 2293, 3081, 2346 and 2413 genes were differentially expressed at 6, 12, 24 and 48 h after cold treatment, respectively. The number of up-regulated genes was similar to that of down-regulated genes at all time points (Fig. 1A). Based on hierarchical clustering analysis of the DEGs, the transcriptional response to cold stress in LTH occurred in two distinct phases, i.e. an early response phase (6 and 12 h) and a late response phase (24 and 48 h) (Fig. 1B).

Functional analysis of DEGs based on GO classification Enriched GO categories of the DEGs at different time points are listed in Table S2. At a significance of less than 0.05 FDR, 24, 30, 23 and 10 significant enriched GO terms for the up-regulated genes were found at 6, 12, 24 and 48 h after cold stress, respectively. Twenty-two enriched GO terms of the up-regulated genes were common at 6 h and 12 h which was defined as the early phase by clustering analysis. The top three most significant GO terms of the up-regulated genes were “response to stimulus (GO:0050896)”, “transcription (GO:0006350)”, “response to abiotic stimulus (GO:0009628)” for the early phase. Top significant GO terms for late phase were, “regulation of biological process (GO:0050789)”, “biological regulation (GO:0065007)”, “regulation of cellular process (GO:0050794)” for 24 h, and “gene expression (GO:0010467)”, “translation (GO:0006412)”, “response to abiotic stimulus (GO:0009628)” for 48 h. “Transcription (GO: 0006350)” was the common significant enriched GO term for the up-regulated genes at all of four time points. “Photosynthesis (GO: 0015979)” was the most significant enriched GO term of the down-regulated genes at 24 and 48 h. These genes were repressed from 12 h after cold treatment. However, “Translation (GO:0006412)” was the most significantly enriched GO term of the down-regulated genes at 6 and 12 h, but it was the most significant enriched GO term of the up-regulated genes at 48 h.

Differentially expressed transcription factor genes under cold stress To understand the regulatory mechanisms of the transcription factors in regulation of cold stress responses, we investigated the expression patterns of differentially expressed transcription factor genes. Among the 5557 DEGs found in this study, 275 were transcription factor genes that belong to This article is protected by copyright. All rights reserved

Accepted Article

diverse transcription factor families according to the Database of Rice Transcription Factor (http://drtf.cbi.pku.edu.cn/) (Gao et al. 2006). The data for differential expression of transcription factor genes at different time points were listed in Table S3. Based on their expression at different time points, the differentially expressed transcription factor genes were clustered into two groups (groups I: 6 and 12 h; group II: 24 and 48 h) (Fig. 1C). This grouping is consistent with the clustering of DEGs (Figs. 1B and 1C). Most of the differentially expressed transcription factor genes start to be induced at the early time during cold stress. Table 1 listed some important transcription factor families that actively responded to cold stress in this study. Among the 28 differentially expressed AP2-EREBP family genes, 26 genes were induced and only two were suppressed by cold stress. Most AP2-EREBP genes were quickly induced in the early phase during cold stress (Table S3). In the AP2-EREBP family, OsDREB1A (Os09g0522200), OsDREB1B (Os09g0522000), OsDREB2A (Os01g0165000) and OsDREB2B (Os05g0346200) were differentially expressed during cold stress. OsDREB1A, OsDREB1B and OsDREB2B were quickly induced by cold and their expression decreased at 48 h after cold treatment, while OsDREB2A was induced after 24 h. Among the 109 WRKY genes identified in rice genome (Jang and Hwang 2010), 17 genes were differentially expressed. Most of them were up-regulated, but three different patterns were observed: continuous induction, induction only in the early phase, and induction only in the late phase. OsWRKY24 (Os01g0826400), OsWRKY45 (Os05g0322900) and OsWRKY71 (Os02g0181300) were induced continuously. In contrast, OsWRKY12 (Os01g0624700) and OsWRKY30 (Os08g0499300) were down-regulated during cold stress (Table S3). In the NAC gene family, 22 genes showed differential expression. Among these genes, only two were down-regulated; the others were up-regulated, with most induced in the early phase after cold treatment. In addition, the genes in three other transcription factor families, bHLH, bZIP and MyB families, exhibited differential expression during cold stress, but their expression patterns were diverse (Table S3).

Transcription profile of antioxidant enzyme genes under cold stress The expression profile of differentially expressed antioxidant enzyme genes under cold stress was shown in Table S4. One catalase gene (Os02g0115700) was down-regulated. Three superoxide This article is protected by copyright. All rights reserved

Accepted Article

dismutase (SOD) genes (Os06g0143000, Os06g0115400 and Os05g0323900) were down-regulated during the early phase of stress. Among 10 ascorbate peroxidase (APX) genes in the microarray, only one gene (Os09g0538600) was down-regulated during cold stress. In addition, 26 peroxidase genes were differentially expressed and most of them were down-regulated during cold treatment (Table S4 and Fig. 2A). However, it is noticeable that among the 19 differentially expressed glutathione S-transferase (GST) genes, 14 GST genes were up-regulated at least at one time point during cold stress (Fig. 2B). Similarly, a glutathione peroxidase (GPX) gene (Os04g0556300) was also significantly induced during cold stress.

Differential response of phytohormone-related genes to cold stress To understand the regulatory mechanisms of phytohormones in rice response to cold stress, we investigated the expression patterns of auxin, ABA, and salicylic acid (SA) related genes in LTH during cold stress. Nine auxin response factor (ARF) genes and ten AUX/IAA genes were differentially expressed. Among nine ARF genes, eight were up-regulated at 6 h after cold stress, and the other one was induced after 24 h (Fig. 3A). However, among ten AUX/IAA genes, seven were down-regulated at different time points during cold stress (Fig. 3B). In addition, three auxin efflux carrier protein genes were significantly up-regulated at the early phase during cold stress (Table 2). In our study, ABA synthesis related genes showed different expression patterns during cold stress (Table 3). Among the genes related to ABA synthesis, Os03g0810800 was down-regulated at all of four time points, while the two genes, Os03g0184000 and Os07g0204900, which encode two crucial enzymes in ABA synthesis, did not show differential expression during cold stress. In addition, ten genes that responded to both cold stress and ABA in rice reported by Rabbani et al. (2003) were also significantly up-regulated in our study (Fig. 3C). However, most of the 14 genes that only responded to ABA but not cold stress in their study did not show significant up-regulation, or even down-regulation during cold stress in our study (Fig. 3D). Genes related to synthesis and metabolism of SA were also found to be differentially regulated in our study. Three phenylalanine ammonia lyase genes (Os02g0626400, Os04g0518100 and Os04g0518400) and one chorismate mutase gene (Os01g0764400), which are responsible for the phenylalanine pathway in SA synthesis, were induced quickly at 6 h after cold treatment. However, the gene Os09g0361400 which encodes isochorismate synthase, a key enzyme in the pathway for SA This article is protected by copyright. All rights reserved

Accepted Article

synthesis from chorismate, had no change in expression during cold stress. Furthermore, three SA glucosetransferase genes (Os04g0206500, Os04g0206600 and Os09g0518200) and one SA methyltransferase gene (Os01g0701700) were highly up-regulated during cold stress (Table 4).

Differentially expressed protein metabolism-related genes under cold stress To understand the protein metabolism in LTH under cold stress, we investigated the differential expression of the genes related to protein synthesis and degradation, including ribosome, translation initiation factor, and proteasome genes. A total of 109 ribosome genes were differentially expressed and most of them were down-regulated during the early phase; then, their expression levels increased and most of these genes were up-regulated at 48 h (Fig. 4A). Ten translation initiation factor genes were differentially expressed and most of them were up-regulated at the late phase (Fig. 4B). The expression patterns of these translation initiation factor genes were similar to that of the ribosome genes (Figs. 4A and 4B). Nine genes annotated as 26S proteasome subunit genes were differentially expressed. All of them were up-regulated at one or more time points during cold stress. The most significant induction happened at 48 h (Fig. 4C).

Detection of enriched conserved cis-regulatory elements in DEGs Transcription factors play important roles in regulation of plants in response to stresses. Cis-elements determine the sites and binding of transcription factors. To understand the transcription regulatory mechanism of LTH in response to cold stress, over represented cis-motifs in the 2 kb upstream sequence from the translational initiation codon of DEGs were analyzed. Up- or down-regulated DEGs at all four time points were investigated separately. Enriched cis-elements of the up-regulated genes at four time points were detected (Table 5). Our results showed that no enriched transcription factor binding sites were detected in the down-regulated genes. For the up-regulated genes, nine and seven enriched transcription factor binding sites were detected at 6 and 12 h, respectively, but only one enriched binding site was detected at 24 and 48 h. Among the ten enriched binding sites, six are the motifs of ABA responsive element (ABRE). The binding site ACGTABREMOTIFA2OSEM was the common enriched binding site among four time points and the only enriched binding site for 24 and 48 h. Six enriched binding sites were in common at 6 h and 12 h. This article is protected by copyright. All rights reserved

Accepted Article

Discussion LTH exhibits different gene expression patterns under moderate cold stress Previous studies on cold tolerance in rice (Andaya and Mackill 2003, Zhan et al. 2005, Zhang et al. 2014) showed that different sets of QTLs for cold tolerance at seedling stage were identified when different low temperatures were used. This suggests that gene expression of rice seedling under cold stress is dependent upon the temperature regimes imposed. To understand the mechanism of rice seedling in response to moderate cold stress and to identify the genes associated with cold tolerance, we conducted a genome-wide gene expression profiling of a cold tolerant japonica variety, LTH under moderate cold stress (8°C). Our results showed that the number of up-regulated genes was similar to that of down-regulated genes at each time point and the numbers of DEGs across time points were also similar (Fig. 1A); this finding was different from previous studies (Yun et al. 2010, Zhang et al. 2012b). Using the same rice genotype (LTH) but severe cold condition (4°C), Zhang et al. (2012b) observed that cold stress induced a continuous increase in DEGs and the number of up-regulated genes clearly exceeded that of down-regulated genes. Similar result was observed by Zhang et al. (2012a) using two contrast rice genotypes and 4°C for cold treatment. However, using the cold-tolerant rice genotype Nipponbare and a milder cold stress (10°C), Yun et al. (2010) observed that the number of DEGs increased with time, but the number of down-regulated genes was much more than that of up-regulated genes during cold stress. In the present study, DEGs could be clustered into two distinct groups according to their expression at different time points (Fig. 1B), suggesting a two-phase transcriptional pattern of LTH in response to moderate cold stress, i.e. 6 and 12 h as the early phase, 24 and 48 h as the late phase. The GO category enrichment analysis further supports the two-phase pattern of gene expression in LTH during cold stress in the present study. Most of the enriched GO terms of the up-regulated genes at 6 h and 12 h were in common. “Photosynthesis (GO: 0015979)” was the most significant enriched GO term of the down-regulated genes at 24 and 48 h, while “Translation (GO:0006412)” was the most significantly enriched GO terms of the down-regulated genes at 6 and 12 h, but the most significant enriched GO term of the up-regulated genes at 48 h (Table S2). The early response phase at 6 and 12 h featured GO terms involved in stimulus response, characteristic of rice seedlings rapidly responding to cold stress by activating transcription. The late response phase at 24 and 48 h featured terms This article is protected by copyright. All rights reserved

Accepted Article

involved in regulation of metabolism, characteristic of rice downstream response to cold stress. Zhang et al. (2012b) also observed a two-phase response in LTH to lower temperatures (4°C) but the timing of the responses was different. In their study, the early phase was defined at 2 h and the late phase was defined at 8, 24 and 48 h after cold stress. In another similar study, three phases of response patterns were defined, i.e. 2 and 6 h as the early phase, 12 h as the middle phase, 24 and 48 h as the late phase (Zhang et al. 2012a). Thus, the early phase in our study is later and longer. Compared with the other studies, the different expression patterns of DEGs in our study may be attributed to use of a moderately low temperature and/or a different genotype.

Transcription factor genes play key roles in regulation of rice response to cold stress In the present study, the differentially expressed transcription factor genes were clustered into two groups based on their expression at different time points (groups I: 6 and 12 h; group II: 24 and 48 h), consistent with the clustering of DEGs (Figs. 1B and 1C), suggesting the determinative role of transcription factors in regulation of the transcriptional response to cold stress. Furthermore, the functional analysis of DEGs based on GO classification revealed that “Transcription (GO: 0006350)” was the common enriched GO term for the up-regulated genes at all of four time points (Table S2) and most of the differentially expressed transcription factor genes were up-regulated with the induction starting early during cold stress (Fig. 1C). These results indicate the importance of early activation of transcription factor genes in response to cold stress in LTH. The DREB subfamily under AP2-EREBP family has been reported to play important roles against abiotic stresses (Roychoudhury et al. 2013). In this study, we also observed that OsDREB1A, OsDREB1B and OsDREB2B were quickly and significantly induced by cold, but the induction decrease to no change at 48 h after cold treatment, while OsDREB2A was induced 24 h after cold stress (Table S3). However, these results differ somewhat from the previous studies. In Arabidopsis, DREB1s were induced by cold, while DREB2s were induced by dehydration and salt stress (Agarwal et al. 2006). Similar results were observed in rice by Dubouzet et al. (2003). In their study, OsDREB1A and OsDREB1B were induced by cold, while OsDREB2A was induced by dehydration and salt. OsDREB2B was also observed to be induced in rice during cold stress (Matsukura et al. 2010). Using the same rice genotype LTH, but with exposure to severe low temperature (4°C) , Zhang et al. (2012b) also observed the induction of OsDREB1A, OsDREB1B and OsDREB2B, but the This article is protected by copyright. All rights reserved

Accepted Article

induction was continuous during cold stress until 48 h and they did not observe the induction of OsDREB2A. Thus, OsDREB1A, OsDREB1B and OsDREB2B seem to be activated longer under severe low temperature and our observation is the first report of the induction of OsDREB2A by cold stress in rice. Further study is underway to confirm the function of OsDREB2A in cold tolerance. Many genes in WRKY family were differentially expressed under cold stress in the present study. OsWRKY24, OsWRKY45 and OsWRKY71 were induced continuously during cold stress (Table S3). These genes had been previously found to be induced by ABA in aleurone cells (Xie et al. 2005). Induction of OsWRKY45 was also reported in response to various abiotic stresses, including cold stress (Tao et al. 2011). Moreover, over-expression of OsWRKY45 enhanced drought tolerance in Arabidopsis (Qiu and Yu 2009). In contrast, OsWRKY12 and OsWRKY30 were down-regulated during cold stress in our study (Table S3). Other studies indicated that OsWRKY30 was induced by ABA (Shen et al. 2012) and fungal infection (Ryu et al. 2006). The induction of numerous WRKY genes in our study implies the important role of this family in rice cold tolerance. In addition, members of NAC transcription factor family has been reported to function in enhancing rice cold tolerance (Hu et al. 2008).In our study, some members in NAC transcription factor family were also differentially expressed during cold stress. Most of the differentially expressed NAC transcription factor genes were induced in the early phase, suggesting that the early activation of NAC transcription factor genes is important in the regulation of response to cold stress in rice.

The glutathione system may serve as the main ROS scavenger in LTH under cold stress Low temperature often leads to production of ROS which causes oxidative damage to proteins, DNA, and lipids in plants. ROS-scavenging enzymes are major contributors to ROS detoxification, especially under stresses (Suzuki et al. 2012). In rice, Guo et al. (2006) reported that the activity of catalase, APX and SOD was enhanced in cold-tolerant cultivars, but decreased in cold-sensitive cultivars after cold stress. Previous genome-wide expression profiling also revealed that the anti-oxidant enzyme genes were induced under cold stress (Cheng et al. 2007). Over-expression of a rice APX gene could enhance cold tolerance in rice (Sato et al. 2011). However, in our study, many of these anti-oxidant enzyme (APX, SOD, catalase and peroxidase) genes were not differentially expressed, or even down-regulated under cold stress. In contrast, genes encoding GPX and GST were up-regulated during cold stress (Table S4), suggesting that LTH may use the glutathione system as its This article is protected by copyright. All rights reserved

Accepted Article

ROS scavenger.

Auxin, ABA and SA are actively involved in regulation of LTH response to cold stress Plant stress responses are complicated and regulated by phytohormones (Kohli et al. 2013). In the present study, genes related to metabolism or signal transduction of auxin, ABA and SA were differentially expressed, implying the involvement of these phytohormones in the regulation of LTH response to cold stress. For the auxin related genes, most of the differentially expressed ARF genes are induced after cold stress (Fig. 3A), while AUX/IAA genes were mainly suppressed (Fig. 3B). The up-regulation of ARF genes and down-regulation of AUX/IAA genes under cold stress were also observed in Arabidopsis (Hannah et al. 2005). The expression profiling of rice under cold stress also showed that nine AUX/IAA genes were up-regulated and four were down-regulated, while three ARF genes were up-regulated and three were down-regulated (Jain and Khurana 2009). Although our findings were similar to the results observed in these studies, the specific responding genes in the two families were different from these studies, suggesting that the involvement of these two family members in response to cold stress is genotype specific. In addition, we observed significant up-regulation of three auxin efflux carrier protein genes at early phase during cold stress (Table 2). In Arabidopsis, cold stress inhibited the intracellular trafficking of auxin efflux carrier proteins in root (Shibasaki et al. 2009). The up-regulation of auxin efflux carrier protein genes can facilitate auxin transportation and this may account for the cold tolerance of LTH. ABA is a well-known phytohormone that plays an important role in plant responses to environmental stresses. ABA accumulates under abiotic stresses, such as drought and salt (Xiong and Zhu 2003). However, the role of ABA in plant cold tolerance is still obscure. Some studies have shown transient increases of endogenous ABA in response to cold stress (Anderson et al. 1994), and enhanced freezing tolerance in cultured plant cells after exogenous ABA application (Chen and Gusta 1983). But other studies showed that ABA levels were much lower and no accumulation occurred in cold-tolerant rice (Zhao et al. 2008). In our study, Os03g0810800, the gene related to ABA synthesis was down-regulated during cold stress, while the expression of Os03g0184000 and Os07g0204900, which encode two crucial enzymes in ABA synthesis (Fang et al. 2008), did not change during cold stress (Table 3). These results suggest that ABA biosynthesis is not a major event in cold stress, This article is protected by copyright. All rights reserved

Accepted Article

consistent with results reported by Lang et al. (1994). However, many ABA responsive genes were significantly up-regulated in our study (Fig. 3C). Significantly, six out of the ten enriched binding sites of the up-regulated genes, ABADESI1, ABREOSRAB21, ACGTABOX, ACGTOSGLUB1, ABREZMRAB28 and ACGTABREMOTIFA2OSEM are the motifs of ABA responsive element (ABRE) (Table 5), suggesting the active involvement of ABA signal transduction in cold responses. Similarly, none of the known ABA synthesis-related genes was differentially expressed and many ABA-induced genes were induced during cold stress in Arabidopsis (Lee et al. 2005). They proposed that the cold induction of some ABA-inducible genes might be mediated by ABF1 and ABF4 in bZIP family which can bind to the ABA responsive element in ABA-responsive promoters. All the transcription factors reported to date, which bind to ABRE or ABRE-like sequences, belong to the bZIP class of transcription factors (Nijhawan et al. 2008). Our results also showed that 11 out of 15 differentially expressed bZIP transcription factor genes induced during cold stress (Table S3). Among the 11 up-regulated bZIP genes observed in our study, the OsbZIP72 (Os09g0456200) has been proven to be a transcription factor that can bind to ABRE in rice using the yeast hybrid systems (Lu et al. 2009). Thus, although ABA levels may not change greatly under cold stress, it is possible that ABA signal transduction can be involved in cold response. Based on the observations that many ABA responsive genes were significantly up-regulated during cold stress and that the dominated enriched binding sites of the up-regulated genes are the motifs of ABA responsive element (ABRE)，we speculate that ABA signaling pathway may play a pivotal role in the regulatory network of LTH response to cold stress. SA is a plant hormone for plant immunity and is important in plant systemic acquired resistance (Vlot et al. 2009). However, there is growing evidence that suggests SA involvement in responses to abiotic stress (Pál et al. 2013). But controversial results on the functions of SA on plant stress tolerance were reported in different studies (Hara et al. 2012). Our results indicated that the SA synthesis pathway through phenylalanine ammonia lyase was activated under cold stress in LTH (Table 4). In eggplant, exogenous SA treatment improved chilling tolerance by up-regulating the expression of glutathione related genes, such as GPX, GST and glutathione reductase (Chen et al. 2011). Thus, it is possible that the activation of the SA synthesis pathway could lead to induction of GPX and GST gene observed in the present study. Further study should be carried out to elucidate

This article is protected by copyright. All rights reserved

Accepted Article

these issues.

Induction of protein metabolism-related genes in the late phase during cold stress is associated with cold tolerance in LTH In the present study, ribosome genes were down-regulated at the early phase, then their expression levels increased and many of these genes were up-regulated at 48 h (Fig. 4A). The translation initiation factor genes were up-regulated at the late phase, with similar expression patterns of the ribosome genes (Fig. 4B). The 26S proteasome subunit genes were up-regulated at one or more time points during cold stress, and the most significant induction occurred at 48 h (Fig. 4C). Ribosome and translation initiation factor are responsible for protein synthesis, while proteasome is responsible for protein degradation. Thus, the protein metabolism-related genes exhibited similar expression patterns, i.e. induction at late phase, particularly at 48 h during cold stress. This result is consistent with the GO enrichment analysis, in which “Translation” was the most enriched GO term of the up-regulated DEGs at 48 h after cold stress. Ribosomes are speculated to play important roles against abiotic stresses including cold stress, and some ribosome genes are induced under cold stress in plants (Brosché and Strid 1999, Kim et al. 2004, Sáez-Vásquez et al. 2000). It has been reported that ribosome had RNA helicase activity (Takyar et al. 2005). Unwinding RNA helices is an essential step in many crucial physiological processes, including replication, DNA repair, recombination, transcription, pre-mRNA splicing, and translation (Delagoutte and von Hippel 2003). In eukaryotes, the turnover rates of misfolded and damaged proteins as well as numerous regulatory proteins are controlled by the ubiquitin/26S proteasome system (Kurepa and Smalle 2008). Recent reports demonstrated that genes encoding proteasome subunits were induced under cold in Thellungiella halophila (Gao et al. 2009). Fluctuation in proteasomal activity affected plant cell proliferation and expansion rates, suggesting that they play an important role in response to stresses (Kurepa et al. 2009). Generally，the plant protein synthesis is inhibited under cold stress. Thus, the induction of protein metabolism related genes, ribosome, translation initiation factor, and proteasome genes in the late phase of cold stress may be an important feature for cold tolerance in LTH. In conclusion, a genome-wide gene expression profiling of a cold tolerant japonica variety, LTH,

This article is protected by copyright. All rights reserved

Accepted Article

has provided new insights on how gene expression in rice seedling respond to cold stress. In contrast to most studies on rice cold response, we used a moderately low temperature (8°C), which frequently occurs during rice growing season, as a cold stress. Compared with previous studies, we identified novel genes and different regulatory mechanisms in response to cold stress. A comparison of results in this study with those of previous studies, suggests that the mechanisms triggered by cold stress in rice seedlings vary depending on rice genotype and the temperature imposed.

Author contributions JZ and SZ participated in microarray experiments, drafting the manuscript and proposal writing. TY, ZZ, ZH, QL and XW participated in material preparation, RNA extraction and RT-PCR. JL and HL participated in experimental designs and corrected the manuscript. BL conceived of the study, drafted proposal and corrected manuscript. Acknowledgements – This work was supported in part by the Ministry of Science and Technology of China (2006DFB33320), National Natural Science Foundation of China (30671251), Guangdong Natural Science Foundation (10151064001000012 ), Guangdong Scientific and Technological Plan (20110205, 2011A020102008) and China Non-profit Agricultural Special Foundation (201303007).

References Agarwal PK, Agarwal P, Reddy M, Sopory SK (2006) Role of DREB transcription factors in abiotic and biotic stress tolerance in plants. Plant Cell Rep 25: 1263–1274 Andaya V, Mackill D (2003) Mapping of QTLs associated with cold tolerance during the vegetative stage in rice. J Exp Bot 54: 2579–2585 Anderson MD, Prasad TK, Martin BA and Stewart CR (1994) Differential gene expression in chilling-acclimated maize seedlings and evidence for the involvement of abscisic acid in chilling tolerance. Plant Physiol 105: 331–339 Brosché M, Strid Å (1999) The mRNA-binding ribosomal protein S26 as a molecular marker in plants: molecular cloning, sequencing and differential gene expression during environmental stress. Biochim Biophys Acta 1445: 342–344 Chen S, Zimei L, Cui J, Jiangang D, Xia X, Liu D, Yu J (2011) Alleviation of chilling-induced oxidative damage by salicylic acid pretreatment and related gene expression in eggplant seedlings. Plant Growth Regul 65: 101–108 Chen THH, Gusta LV (1983) Abscisic acid-induced freezing resistance in cultured plant cells. Plant This article is protected by copyright. All rights reserved

Accepted Article

Physiol 73: 71–75 Cheng C, Yun KY, Ressom H, Mohanty B, Bajic V, Jia Y, Yun S, de los Reyes B (2007) An early response regulatory cluster induced by low temperature and hydrogen peroxide in seedlings of chilling-tolerant japonica rice. BMC Genomics 8: 175 Cruz RPd, Sperotto RA, Cargnelutti D, Adamski JM, FreitasTerra T, Fett JP (2013) Avoiding damage and achieving cold tolerance in rice plants. Food and Energy Security 2: 96–119 Dai Q, Vergara B, Visperas R (1990) Amelioration of cold injury in rice (Oryza sativa L.): improving root oxidizing activity by plant growth regulators. Philipp J Crop Sci 15: 49–54 Delagoutte E, von Hippel PH (2003) Helicase mechanisms and the coupling of helicases within macromolecular machines Part II: Integration of helicases into cellular processes. Q Rev Biophys 36: 1–69 Du Z, Zhou X, Ling Y, Zhang Z, Su Z (2010) agriGO: a GO analysis toolkit for the agricultural community. Nucl Acids Res 38: W64–W70 Dubouzet JG, Sakuma Y, Ito Y, Kasuga M, Dubouzet EG, Miura S, Seki M, Shinozaki K, Yamaguchi‐ Shinozaki K (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 Eisen MB, Spellman PT, Brown PO, Botstein D (1998) Cluster analysis and display of genome-wide expression patterns. Proc Natl Acad Sci USA 95: 14863–14868 Fang J, Chai C, Qian Q, Li C, Tang J, Sun L, Huang Z, Guo X, Sun C, Liu M (2008) Mutations of genes in synthesis of the carotenoid precursors of ABA lead to pre-harvest sprouting and photo-oxidation in rice. Plant J 54: 177–189 Gao F, Zhou Y, Zhu W, Li X, Fan L, Zhang G (2009) Proteomic analysis of cold stress-responsive proteins in Thellungiella rosette leaves. Planta 230: 1033–1046 Gao G, Zhong Y, Guo A, Zhu Q, Tang W, Zheng W, Gu X, Wei L and Luo J (2006) DRTF: a database of rice transcription factors. Bioinformatics 22: 1286–1287 Guh JO (2002) Mechanisms of chilling tolerance in relation to antioxidative enzymes in rice. Korean J Crop Sci 47: 341–351 Guo Z, Ou W, Lu S, Zhong Q (2006) Differential responses of antioxidative system to chilling and drought in four rice cultivars differing in sensitivity. Plant Physiol Biochem 44: 828–836 Han L, Qiao Y, Zhang S, Zhang Y, Cao G, Kim J, Lee K, Koh H (2007) Identification of quantitative This article is protected by copyright. All rights reserved

Accepted Article

trait loci for cold response of seedling vigor traits in rice. J Genet Genomic 34: 239–246 Hannah MA, Heyer AG, Hincha DK (2005) A global survey of gene regulation during cold acclimation in Arabidopsis thaliana. PLoS Genet 1: e26 Hara M, Furukawa J, Sato A, Mizoguchi T, Miura K (2012) Abiotic stress and role of salicylic acid in plants. In: Ahmad P, Parsad MNV (eds) Abiotic stress responses in plants: metabolism, productivity and sustainability. Springer Science, New York, pp 235–251 Hu H, You J, Fang Y, Zhu X, Qi Z, Xiong L (2008) Characterization of transcription factor gene SNAC2 conferring cold and salt tolerance in rice. Plant Mol Bio 67: 169–181 Huang M, Guo Z (2005) Responses of antioxidative system to chilling stress in two rice cultivars differing in sensitivity. Biologia Plantarum 49: 81–84 Jain M, Khurana JP (2009) Transcript profiling reveals diverse roles of auxin-responsive genes during reproductive development and abiotic stress in rice. FEBS J 276: 3148–3162 Jang JY, Hwang DJ (2010) The WRKY superfamily of rice transcription factors. Plant Pathol J 26: 110–114 Kabaki N, Tajima K (1981) Effect of chilling on the water balance of rice seedlings. Jap J Crop Sci 50: 489–494 Kim KY, Park SW, Chung YS, Chung CH, Kim JI, Lee JH (2004) Molecular cloning of low-temperature-inducible ribosomal proteins from soybean. J Exp Bot 55: 1153–1155 Kohli A, Sreenivasulu N, Lakshmanan P, Kumar PP (2013) The phytohormone crosstalk paradigm takes center stage in understanding how plants respond to abiotic stresses. Plant Cell Rep 32: 945–957 Koseki M, Kitazawa N, Yonebayashi S, Maehara Y, Wang ZX, Minobe Y (2010) Identification and fine mapping of a major quantitative trait locus originating from wild rice, controlling cold tolerance at the seedling stage. Mol Genet Genomics 284: 45–54 Kurepa J, Smalle JA (2008) Structure, function and regulation of plant proteasomes. Biochimie 90: 324–335 Kurepa J, Wang S, Li Y, Smalle J (2009) Proteasome regulation, plant growth and stress tolerance. Plant Signal Behav 4: 924–927 Lang V, Mantyla E, Welin B, Sundberg B, Palva ET (1994) Alterations in water status, endogenous abscisic acid content, and expression of rab18 gene during the development of freezing tolerance This article is protected by copyright. All rights reserved

Accepted Article

in Arabidopsis thaliana. Plant Physiol 104:1341–1349 Lee B, Henderson DA, Zhu JK (2005) The Arabidopsis cold-responsive transcriptome and its regulation by ICE1. Plant Cell 17: 3155–3175 Lee T-M, Lur H-S, Chu C (1997) Role of abscisic acid in chilling tolerance of rice (Oryza sativa L.) seedlings: II. Modulation of free polyamine levels. Plant Sci 126: 1–10 Lu G, Gao C, Zheng X, Han B (2009) Identification of OsbZIP72 as a positive regulator of ABA response and drought tolerance in rice. Planta 229: 605–615 Matsukura S, Mizoi J, Yoshida T, Todaka D, Ito Y, Maruyama K, Shinozaki K, Yamaguchi-Shinozaki K (2010) Comprehensive analysis of rice DREB2-type genes that encode transcription factors involved in the expression of abiotic stress-responsive genes. Mol Genet Genomics 283: 185–196 Monroy AF, Dryanova A, Malette B, Oren DH, Farajalla MR, Liu W, Danyluk J, Ubayasena LW, Kane K, Scoles GJ (2007) Regulatory gene candidates and gene expression analysis of cold acclimation in winter and spring wheat. Plant Mol Bio 64: 409–423 Morris RT, O'Connor TR, Wyrick JJ (2008) Osiris: an integrated promoter database for Oryza sativa L. Bioinformatics 24: 2915–2917 Nijhawan A, Jain M, Tyagi A K, Khurana J P (2008) Genomic survey and gene expression analysis of the basic leucine zipper transcription factor family in rice. Plant Physiol 146: 333–350 Pál M, Szalai G, Kovács V, Gondor O, Janda T (2013) Salicylic Acid-Mediated Abiotic Stress Tolerance. In: Hayat A, Ahmad A, Alyemeni MN (eds) Salicylic Acid. Plant Growth and Development. Springer, New York, pp 1–14 Qiu Y, Yu D (2009) Over-expression of the stress-induced OsWRKY45 enhances disease resistance and drought tolerance in Arabidopsis. Environ Exp Bot 65: 35–47 Rabbani MA, Maruyama K, Abe H, Khan MA, Katsura K, Ito Y, Yoshiwara K, Seki M, Shinozaki K, Yamaguchi-Shinozaki K (2003) 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 133: 1755–1767 Roychoudhury A, Paul S, Basu S (2013) Cross-talk between abscisic acid-dependent and abscisic acid-independent pathways during abiotic stress. Plant Cell Rep 32: 985–1006 Ryu H-S, Han M, Lee S-K, Cho J-I, Ryoo N, Heu S, Lee Y-H, Bhoo SH, Wang G-L, Hahn T-R (2006) A comprehensive expression analysis of the WRKY gene superfamily in rice plants during This article is protected by copyright. All rights reserved

Accepted Article

defense response. Plant Cell Rep 25: 836–847 Sáez-Vásquez J, Gallois P, Delseny M (2000) Accumulation and nuclear targeting of BnC24, a Brassica napus ribosomal protein corresponding to a mRNA accumulating in response to cold treatment. Plant Sci 156: 35–46 Sakuma Y, Liu Q, Dubouzet JG, Abe H, Shinozaki K, Yamaguchi-Shinozaki K (2002) DNA-binding specificity of the ERF/AP2 domain of Arabidopsis DREBs, transcription factors involved in dehydration-and cold-inducible gene expression. Biochem Biophys Res Commun 290: 998–1009 Sato Y, Masuta Y, Saito K, Murayama S, Ozawa K (2011) Enhanced chilling tolerance at the booting stage in rice by transgenic overexpression of the ascorbate peroxidase gene, OsAPXa. Plant Cell Rep 30: 399–406 Seki M, Narusaka M, Ishida J, Nanjo T, Fujita M, Oono Y, Kamiya A, Nakajima M, Enju A, Sakurai T (2002) Monitoring the expression profiles of 7000 Arabidopsis genes under drought, cold and high-salinity stresses using a full-length cDNA microarray. Plant J 31: 279–292 Shen H, Liu C, Zhang Y, Meng X, Zhou X, Chu C, Wang X (2012) OsWRKY30 is activated by MAP kinases to confer drought tolerance in rice. Plant Mol Bio 80: 241–253 Shibasaki K, Uemura M, Tsurumi S, Rahman A (2009) Auxin response in Arabidopsis under cold stress: underlying molecular mechanisms. Plant Cell 21: 3823–3838 Suzuki N, Koussevitzky S, Mittler R, Miller G (2012) ROS and redox signalling in the response of plants to abiotic stress. Plant Cell Environ 35: 259–270 Takyar S, Hickerson RP, Noller HF (2005) mRNA helicase activity of the ribosome. Cell 120: 49–58 Tao Z, Kou Y, Liu H, Li X, Xiao J, Wang S (2011) OsWRKY45 alleles play different roles in abscisic acid signalling and salt stress tolerance but similar roles in drought and cold tolerance in rice. J Exp Bot 62: 4863–4874 Tondelli A, Francia E, Barabaschi D, Aprile A, Skinner J, Stockinger E, Stanca A, Pecchioni N (2006) Mapping regulatory genes as candidates for cold and drought stress tolerance in barley. Theor Appl Genet 112: 445–454 Trzcinska-Danielewicz J, Bilska A, Fronk J, Zielenkiewicz P, Jarochowska E, Roszczyk M, Jończyk M, Axentowicz E, Skoneczny M, Sowiński P (2009) Global analysis of gene expression in maize leaves treated with low temperature: I. Moderate chilling (14 C). Plant Sci 177: 648–658 Tusher VG, Tibshirani R, Chu G (2001) Significance analysis of microarrays applied to the ionizing This article is protected by copyright. All rights reserved

Accepted Article

radiation response. Proc Natl Acad Sci USA 98:5116–5121 Vlot AC, Dempsey DMA, Klessig DF (2009) Salicylic acid, a multifaceted hormone to combat disease. Annu Rev Phytopathol 47: 177–206 Xie Z, Zhang ZL, Zou X, Huang J, Ruas P, Thompson D, Shen QJ (2005) Annotations and functional analyses of the rice WRKY gene superfamily reveal positive and negative regulators of abscisic acid signaling in aleurone cells. Plant Physiol 137: 176–189 Xiong L, Zhu JK (2003) Regulation of abscisic acid biosynthesis. Plant Physiol 133: 29–36 Ye C, Fukai S, Godwin I, Reinke R, Snell P, Schiller J, Basnayake J (2009) Cold tolerance in rice varieties at different growth stages. Crop Pasture Sci 60: 328–338 Yun KY, Park MR, Mohanty B, Herath V, Xu F, Mauleon R, Wijaya E, Bajic VB, Bruskiewich R, de Los Reyes BG (2010) Transcriptional regulatory network triggered by oxidative signals configures the early response mechanisms of japonica rice to chilling stress. BMC Plant Biol 10: 16 Zhan Q, Zhu K, Chen Z, Zeng S (2005) Studies on the QTL for cold tolerance related characters of rice seedling by molecular markers. Hybrid Rice 20: 50–55 Zhang F, Huang L, Wang W, Zhao X, Zhu L, Fu B, Li Z (2012a) 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 Zhang S, Zheng J, Liu B, Peng S, Leung H, Zhao J, Wang X, Yang T, Huang Z (2014) Identification of QTLs for cold tolerance at seedling stage in rice (Oryza sativa L.) using two distinct methods of cold treatment. Euphytica 195: 95–104 Zhang T, Zhao X, Wang W, Pan Y, Huang L, Liu X, Zong Y, Zhu L, Yang D, Fu B (2012b) Comparative transcriptome profiling of chilling stress responsiveness in two contrasting rice genotypes. PLoS One 7: e43274 Zhao X, Dolferus R, Darvey N (2008) Precision breeding of cold tolerant rice. IREC Farmers Newsletter 177: 12–13

Edited by M. Uemura

This article is protected by copyright. All rights reserved

Accepted Article

Figure legends

Fig. 1. Expression patterns of the DEGs in LTH during cold stress. (A) Distribution of the DEGs at different time points. (B) Hierarchical cluster analysis of DEGs. (C) Hierarchical cluster analysis of differentially expressed transcription factor genes. DEGs were identified using the combined criteria of three-fold change and a cut-off of q-values (≤ 0.05) in SAM based on three biological replicates. These include 5557 DEGs and 275 differentially expressed transcription factor genes at four time points (6, 12, 24 and 48 h) during cold stress. The fold change of gene (expression data of cold treatment / expression data of control) was log (base 2)-transformed and subjected to hierarchical clustering analysis using Gene Cluster 3.0. Up-regulated genes are shown in red, down-regulated genes are shown in green and no change in black.

This article is protected by copyright. All rights reserved

Accepted Article

Fig. 2. Heat map view of expression of the differentially expressed peroxidase genes and glutathione S-transferase genes in LTH during cold stress. (A) Peroxidase genes. (B) Glutathione S-transferase genes. 26 peroxidase genes and 19 glutathione S-transferase (GST) genes exhibited differentially expressed at different time points (6, 12, 24 and 48 h) during cold stress were analyzed by hierarchical clustering using Gene Cluster 3.0. The fold change of gene (expression data of cold treatment / expression data of control) was log (base 2)-transformed and subjected to hierarchical clustering analysis using Gene Cluster 3.0. Up-regulated genes are shown in red, down-regulated genes are shown in green and no change in black.

This article is protected by copyright. All rights reserved

Accepted Article

Fig. 3. Heat map view of expression of genes related to phytohormone signaling pathway in LTH during cold stress. (A) Auxin response factor genes. (B) AUX/IAA genes. (C) Genes responding to both ABA and cold stress. (D) Genes responding to ABA but not cold stress. Differentially expressed genes related to phytohormone signal pathway at different time points (6, 12, 24 and 48 h) during cold stress were analyzed by hierarchical clustering using Gene Cluster 3.0. The fold change of gene (expression data of cold treatment / expression data of control) was log (base 2)-transformed and subjected to hierarchical clustering analysis using Gene Cluster 3.0. Up-regulated genes are shown in red, down-regulated genes are shown in green and no change in black.

This article is protected by copyright. All rights reserved

Accepted Article

Fig. 4. Heat map view of expression of genes related to protein metabolism in LTH during cold stress. (A) Ribosome genes. These include 109 differentially expressed genes. (B) Translation initiation factor genes. (C) Proteasome genes. Differentially expressed genes related to protein metabolism at different time points (6, 12, 24 and 48 h) during cold stress were analyzed by hierarchical clustering using Gene Cluster 3.0. The fold change of gene (expression data of cold treatment / expression data of control) was log (base 2)-transformed and subjected to hierarchical clustering analysis using Gene Cluster 3.0. Up-regulated genes are shown in red, down-regulated genes are shown in green and no change in black.

This article is protected by copyright. All rights reserved

Accepted Article

Table 1 Classification of some important cold responded transcription factor genes by family

Transcription

Number of genes changed

Total number of

Percentage of total

factor family

after cold treatment

genes in rice genome

genes within family

ARF

9

48

18.8

bHLH

25

215

11.6

bZIP

15

139

10.8

AP2-EREBP

28

160

17.5

MyB

35

213

16.4

NAC

22

172

12.8

WRKY

17

109

15.6

ZF HD

4

14

28.6

This article is protected by copyright. All rights reserved

Accepted Article

Table 2 Expression patterns of the three auxin efflux carrier genes in LTH during cold stress. Fold change: expression data of cold treatment/expression data of control. q-value: calculated by SAM software based on three biological replicates. ND: Not detected data because of low signal intensity. * Rice Annotation Project Database (RAP-DB): http://rapdb.dna.affrc.go.jp/. RAP ID *

6h Fold

12h q-value

change Os01g0643300

109.4

Fold

24h

q-value

change 0.000

ND

Fold

48h

q-value

change ND

ND

Fold

Annotation

q-value

change ND

ND

ND

Auxin

efflux

carrier

component 3a Os01g0919800

2.8

0.000

6.1

0.000

0.6

0.085

0.6

6.614

Auxin

efflux

carrier

component 6 Os02g0743400

10.5

0.000

12.0

0.000

0.6

0.085

0.5

2.854

Auxin

efflux

component 1

This article is protected by copyright. All rights reserved

carrier

Accepted Article

Table 3 Expression patterns of genes related to ABA synthesis in LTH during cold stress. Fold change: expression data of cold treatment/expression data of control. q-value: calculated by SAM software based on three biological replicates. * Rice Annotation Project Database (RAP-DB): http://rapdb.dna.affrc.go.jp/ RAP ID*

6h Fold

12h q-value

change Os03g0810800

0.3

Fold

24h

q-value

change 0.000

0.2

Fold

48h q-value

change 0.000

0.2

Fold

Annotation

q-value

change 0.000

0.2

0.000

Xanthoxin dehydrogenase

Os03g0184000

1.4

0.2104

0.8

1.625

1.3

10.378

1.5

21.848

Phytoene desaturase

Os07g0204900

1.8

0.000

1.1

18.085

0.6

0.103

0.5

2.854

ζ-carotene desaturase

This article is protected by copyright. All rights reserved

Accepted Article

Table 4 Expression patterns of genes related to salicylic acid synthesis and metabolism in LTH during cold stress. Fold change: expression data of cold treatment/expression data of control. q-value: calculated by SAM software based on three biological replicates. ND: Not detected data because of low signal intensity. * Rice Annotation Project Database (RAP-DB): http://rapdb.dna.affrc.go.jp/ RAP ID*

6h Fold

12h q-value

change

Fold

24h

q-value

change

Fold

48h q-value

change

Fold

Annotation

q-value

change

Phenylalanine pathway Os01g0764400

4.1

0.000

4.8

0.000

0.5

0.019

0.4

0.173

Chorismate mutase

Os04g0518100

2.2

0.000

4.9

0.000

1.9

0.517

1.8

8.171

Phenylalanine ammonia lyase

Os02g0626400

7.8

0.000

14.9

0.000

2.7

0.000

2.1

3.978

Phenylalanine ammonia lyase

Os04g0518400

53.8

0.000

3.1

0.050

0.6

0.322

0.6

5.136

Phenylalanine ammonia lyase

Isochorismate pathway Os09g0361400

0.8

0.210

1.3

1.625

1.0

23.359

1.1

37.567

Isochorismate synthase

Salicylic acid glucosylation Os04g0206600

ND

ND

ND

ND

570.0

0.000

69.6

0.000

Salicylic acid glucosyltransferase

Os09g0518200

15.0

0.000

ND

ND

0.70

0.976

45.0

0.000

Salicylic acid glucosyltransferase

Os04g0206500

67.9

0.000

64.0

0.000

180.7

0.000

99.9

0.000

Salicylic acid glucosyltransferase

Salicylic acid methylation Os01g0701700

ND

ND

137.9

0.000

6.0

0.000

ND

ND

Salicylic acid methyltransferase

This article is protected by copyright. All rights reserved

Accepted Article

Table 5 Enriched transcription factor binding sites of the up-regulated genes at different time points. Enriched transcription factor binding sites of DEGs at different time points were analyzed using Osiris program (http://www.bioinformatics2.wsu.edu/cgibin/Osiris/cgi/home.pl). q-value was calculated by the program. The accession numbers of these sites refer to PLACE (Plant Cis-acting Regulatory DNA Elements,http://www.dna.affrc.go.jp/PLACE/) Time

Transcription Factor

Accession

point

Binding Site

Number in

q-value

Description

1.00E-10

Experimentally determined sequence requirement of ACGT-core of

PLACE 6h ACGTABREMOTIFA2OSEM

S000394

motif A in ABRE of the rice gene ABADESI1

S000007

1.00E-08

Responsive to ABA and desiccation

ABREOSRAB21

S000012

1.00E-08

ABA responsive element (ABRE) of wheat and rice rab21 DE genes

DRECRTCOREAT

S000418

1.00E-06

Core

motif

of

DRE/CRT

(dehydration-responsive

element/C-repeat) DE cis-acting element found in many genes in Arabidopsis and in DE rice ACGTOSGLUB1

S000278

1.00E-05

ACGT motif found in GluB-1 gene in rice required for endosperm-specific expression

ABREZMRAB28

S000133

1.00E-04

ABA and water-stress responses

ACGTABOX

S000130

1.00E-03

"A-box" according to the nomenclature of ACGT elements

MYCATERD1

S000413

1.00E-03

MYC recognition sequence necessary for DE expression of erd1 in dehydrated Arabidopsis

MYCATRD22

s000174

1.00E-03

Binding site for MYC in Arabidopsis dehydration-responsive gene

ACGTABREMOTIFA2OSEM

S000394

1.00E-11

Experimentally determined sequence requirement of ACGT-core of

12 h

This article is protected by copyright. All rights reserved

Accepted Article

motif A in ABRE of the rice gene ABREOSRAB21

S000012

1.00E-07

ABA responsive element (ABRE) of wheat and rice rab21 DE genes

ABADESI1

S000007

1.00E-05

Responsive to ABA and desiccation

ACGTOSGLUB1

S000278

1.00E-05

ACGT motif found in GluB-1 gene in rice required for endosperm-specific expression

ACGTABOX

S000130

1.00E-04

"A-box" according to the nomenclature of ACGT elements

DRECRTCOREAT

S000418

1.00E-03

Core

motif

of

DRE/CRT

(dehydration-responsive

element/C-repeat) DE cis-acting element found in many genes in Arabidopsis and in DE rice POLASIG2

S000081

1.00E-03

poly A signal found in rice alpha-amylase

ACGTABREMOTIFA2OSEM

S000394

1.00E-05

Experimentally determined sequence requirement of ACGT-core of

24 h

motif A in ABRE of the rice gene 48 h ACGTABREMOTIFA2OSEM

S000394

1.00E-05

Experimentally determined sequence requirement of ACGT-core of motif A in ABRE of the rice gene

This article is protected by copyright. All rights reserved