Comparative proteomic analysis of QTL CTS-12 derived from wild rice (Oryza rufipogon Griff.), in the regulation of cold acclimation and de-acclimation of rice (Oryza sativa L.) in response to severe chilling stress

Background Rice (Oryza sativa L.) is a thermophilic crop vulnerable to chilling stress. However, common wild rice (Oryza rufipogon Griff.) in Guangxi (China) has the ability to tolerate chilling stress. To better understand the molecular mechanisms underlying chilling tolerance in wild rice, iTRAQ-based proteomic analysis was performed to examine CTS-12, a major chilling tolerance QTL derived from common wild rice, mediated chilling and recovery-induced differentially expressed proteins (DEPs) between the chilling-tolerant rice line DC90 and the chilling-sensitive 9311. Results Comparative analysis identified 206 and 155 DEPs in 9311 and DC90, respectively, in response to the whole period of chilling and recovery. These DEPs were clustered into 6 functional groups in 9311 and 4 in DC90. The majority were enriched in the ‘structural constituent of ribosome’, ‘protein-chromophore linkage’, and ‘photosynthesis and light harvesting’ categories. Short Time-series Expression Miner (STEM) analysis revealed distinct dynamic responses of both chloroplast photosynthetic and ribosomal proteins between 9311 and DC90. Conclusion CTS-12 might mediate the dynamic response of chloroplast photosynthetic and ribosomal proteins in DC90 under chilling (cold acclimation) and recovery (de-acclimation) and thereby enhancing the chilling stress tolerance of this rice line. The identified DEPs and the involvement of CTS-12 in mediating the dynamic response of DC90 at the proteomic level illuminate and deepen the understanding of the mechanisms that underlie chilling stress tolerance in wild rice. Electronic supplementary material The online version of this article (10.1186/s12870-018-1381-7) contains supplementary material, which is available to authorized users.


Background
Chilling stress is one of the most common abiotic stresses and can cause severe injury in all stages of rice growth with seedling stage being one of the most vulnerable to chilling stress. The optimal temperature for rice growth is 25-35°C. Low temperature not only causes severe injury to early-season cold-sensitive rice cultivars in the spring but also leads to substantial yield loss of late-season rice in the autumn. With the increase in global climate anomalies, chilling stress occurs frequently, and the range of its influence is widening [1]. Therefore, understanding the mechanisms that underlie chilling tolerance and to breeding chilling-tolerant cultivars for rice production are urgent goals.
Chilling stress triggers a series of changes in physiological and molecular processes [2] and results in the accumulation of reactive oxygen species (ROS) in plant cells. The accumulating ROS cause oxidative stress, which damages the plant cell membrane, decreases enzyme activity, and inhibits the rate of photosynthesis and protein translation [3,4]. Plants have evolved complicated mechanism(s) at multiple levels to adapt to environmental stresses, including chilling stress. On one hand, stress-responsive signaling pathways regulate the expression of several downstream stress-related genes in response to abiotic stresses. For instance, DREB1/CBF belongs to the AP2/ERF superfamily transcription factors. The overexpression of CBF1 and CBF3 strongly induces the expression of COR (cold regulated) genes, thus increasing the tolerance of Arabidopsis to low-temperature conditions [5]. PtrbHLH, a citrus (Poncirus trifoliata) basic helix-loop-helix (bHLH) transcription factor, functions in chilling tolerance by positively regulating POD-mediated ROS scavenging [6]. On the other hand, non-enzymatic and enzymatic antioxidant systems participate in ROS scavenging to protect plant cells from oxidative damage [7]. Enzymatic antioxidants include superoxide dismutase (SOD), guaiacol peroxidase (POD), ascorbate peroxidase (APX), and catalase (CAT) [8]. SODs [9]. However, the mechanisms underlying the rice chilling stress response remain elusive, although a chilling-tolerance gene, COLD1, which encodes a trans-membrane protein that regulates G-protein-dependent low-temperature sensing and is required for the development of cold tolerance, has recently been cloned and functionally characterized in rice [10].
Many previous '-omics' studies, which can reveal highly informative expression patterns of genes/proteins that reflect the differentiation of response among genotypes with contrasting stress tolerance, have been performed to unravel the abiotic stress-responsive mechanisms in plants [11][12][13][14][15][16]. A large proportion of identified DEGs/DEPs were annotated as photosynthetic genes/proteins or ribosomal constituents and were induced in response to environmental stresses. Six DEGs, including RPS10, RPS11, RPL21, RPS23, RPL35a-3, and a 60S acidic ribosomal protein encoding gene were identified and annotated as constituents of the large or small subunits of ribosomes [16]. All of them were repressed in response to early chilling stress. For the DEGs/DEPs associated with photosynthesis, 13 chlorophyll a/b-binding proteins and 10 chloroplast precursor DEGs were found to be down-regulated in rice under chilling stress [17]. Moreover, similar cases have been reported in drought and salt stress studies. For instance, chlorophyll synthetic and binding proteins, including the chlorophyll a/b-binding protein CP24, CHLH, CHLI, and others functioning in photosynthetic pathways, are significantly down-regulated in rice under drought stress [18,19]. In a salt stress study, the levels of the chlorophyll a/b-binding proteins LHCA1, LHCA2, and LHCA4, were also found to be reduced under salt treatment [11]. Therefore, the comparative transcriptomic and proteomic profiling of genes/proteins provides new insights into environmental stress responses in plants.
Proteomic analysis is regarded as an effective strategy for the large-scale screening and identification of proteins. iTRAQ is a robust mass spectrometry technology that allows the quantitative comparison of protein abundance by measuring the peak intensities of reporter ions released from iTRAQ-tagged peptides by fragmentation during MS/MS [20]. iTRAQ-based proteomics has been widely applied in stress tolerance studies in plants [11,[21][22][23][24].
Here, a chilling-tolerant chromosome segment substitution line (CSSL), DC90, was developed by introgression of the genome of DP15 into the recurrent parent 9311 (chilling sensitive) in our previous study [16]. Therefore, the genetic background of DC90 is highly similar to that of 9311. Previous research has demonstrated that a major QTL harbored in the DP15 genomic segment of chromosome 12 accounted for the chilling-tolerance phenotype of DC90. We designated it as CTS-12. In this study, iTRAQ-based quantitative proteomics was performed to elucidate the effects of CTS-12 on the protein levels of DC90 and 9311 under severe chilling stress and to identify chilling tolerance-associated DEPs by comparison with its recurrent parent. By this approach, we aimed to infer the potential mechanisms that are relevant to the chilling stress tolerance of wild rice.

Plant materials, growth, and treatment conditions
The collection of wild rice was done by Rongbai Li and deposited in wild rice resources conservation field of Guangxi University complying with legislation of China.
The CSSL DC90 were then developed by crossing Guangxi common wild rice with 9311 to obtain chilling tolerance indica cultivar by Prof. Li in our previous breeding project. The growth of rice seedlings and chilling-stress phenotyping were performed according to previously described methods [16]. In brief, rice seeds were sowed in a plastic container filled with paddy soil and grown under natural conditions with natural sunlight to the three-leaf stage. Rice seedlings were then exposed to 10/8°C (day/night) for 5 days with a photoperiod of day-13 h/night-11 h by supplementation with artificial light (20,000 Lux and 65% humidity). After chilling treatment, the rice seedlings were allowed to recover at 28/26°C (day/night) for 7 days. Hydroponic culture was performed to grow rice plants for proteomic analysis as previously described [25]. Briefly, DC90 and 9311 plants were grown in 96-well PCR plates in a growth chamber with rice nutrient solution. Since 9311 is chilling sensitive and is unable to survive from 4-day (96 h) chilling treatment (Additional file 1). To study the dynamic changes of the proteins in response to chilling stress and recovery treatment, the following optimized treatment time lengths were used for sampling. The treatment time length, consisted of chilling treatment for 72 h followed by recovery for 60 h. Three replicates of whole-plant samples were harvested at 0 h, 60 h of chilling treatment, and 60 h of recovery. The DC90 samples from the three timepoints were designated CTL1, CTL2, and CTL3, which represented 0 h, 60 h of chilling treatment, and 60 h of recovery, respectively. Meanwhile, CSL1, CSL2, and CSL3 represented the 9311 samples at 0 h, 60 h of chilling treatment, and 60 h of recovery (Fig. 1a). All samples were frozen immediately in liquid nitrogen and stored at − 80°C for total protein extraction.

The determination of chilling stress induced histochemical and physiological changes in leaf tissues
Staining with 3,3′-Diaminobenzidine (DAB) was used to detect H 2 O 2 . The second leaves of three-leaf stage seedlings were cut into pieces and then subjected to vacuum infiltration in DAB solution (1 mg/mL in water, pH 5.8) for 10 min. The infiltrated samples were incubated for 8 h in the dark. Next, the chlorophyll in the samples was removed by boiling in absolute ethanol. The samples were observed and images taken under a stereo-microscope [26]. For the detection of O 2− , leaf pieces were vacuum-infiltrated as described above in 0.1% nitroblue tetrazolium (NBT) with 0.1% Triton X-100 in 10 mM of phosphate buffer (pH 7.2) for 10 min. Chlorophyll was removed, and the samples were observed under a stereo-microscope [27]. Samples (0.5 g) were ground in liquid nitrogen to a fine powder and then resuspended with 10 mL of 50 mM PBS (pH 7.8) containing 5% PVP (w/v), 1 mM DTT, 0.1% Triton X-100 (w/v) and 0.1 mM EDTA at 4°C . After centrifugation at 15000 g for 30 min, the supernatant was used for enzyme activity assays. The activity of CAT, POD, SOD, and APX were measured by the method as described [28]. The quantitation of relative electrolyte leakage and MDA were conducted following the methods as published previously [16,28]. The AsA assay was performed by the method as described [29] with some necessary modifications. Rice seedlings were hydroponically cultured in 96-well PCR plates to the three-leaf stage, AsA (2 mM) was applied to the nutrient solution, and the seedlings were then exposed to 10/8°C (day/night) for 5 days with a photoperiod of day-13 h/ night-11 h.

Protein extraction, trypsin digestion, and peptide desalting
The samples (1 g) were ground in liquid nitrogen into fine power and then resuspended in lysis buffer (8 M urea, 2 mM EDTA, 10 mM DTT and 1% Protease inhibitor cocktail). The homogenates were centrifuged at 13000 g and 4°C for 10 min. The supernatant was precipitated with acetone at − 20°C for 3 h. After centrifugation, the protein pellet was re-dissolved by urea buffer, which contained 8 M urea and 100 mM tetraethylammonium bromide (TEAB). The protein concentration was determined by using the Bradford Protein Assay Kit according to the manufacturer's instructions.
Trypsin digestion and iTRAQ labeling were performed according to the method described in Mu et al. [15] with minor modification. First, 50 μg of protein was reduced with 10 mM DTT at 37°C for 60 min and then alkylated with 25 mM iodoacetamide (IAM) at room temperature for 30 min in the dark. TEAB (100 mM) was added to dilute the protein sample to a urea concentration of less than 2 M. The protein samples were then digested with modified trypsin (sequencing grade) (protein:trypsin = 50:1) at 37°C overnight. A second digestion was performed with protein:trypsin = 100:1 for 4 h.
After trypsin digestion, the peptides were desalted by a Strata X SPE column and vacuum dried. The peptides were then labeled with the iTRAQ 8-plex kit (AB Sciex, Foster City, CA, USA) according to the manufacturer's instructions. The peptides from samples CSL1, CSL2, CSL3, CTL1, CTL2, and CTL3 were labeled with tags 115, 116, 117, 118, 119, and 121, respectively. All labeled samples were multiplexed and vacuum dried.

High-resolution LC-MS/MS analysis
The labeled samples were reconstituted with HPLC solution A (2% ACN, pH 10) and then fractionated by high pH reverse-phase HPLC with a Waters Bridge Peptide BEH C18 (130 Å, 3.5 μm, 4.6 × 250 mm) column. The The analysis of LC-MS/MS was performed by Nano LC 1000 LC-MS/MS with a Proxeon EASY-nLC 1000 coupled to a Thermo Fisher Q Exactive. Trypsin-digested fractions were directly loaded onto a reversed-phase pre-column (Acclaim PepMap®100 C18, 3 μm, 100 Å, 75 μm × 2 cm) at 5 μL/min in 100% solvent A (0.1 M acetic acid in water). Next, the peptide eluents were loaded onto a reversed-phase analytical column (Acclaim PepMap® RSLC C18, 2 μm, 100 Å, 50 μm × 15 cm). The solvent B gradient (0.1% FA in 98% ACN) increased from 15 to 35% for 45 min, from 35 to 98% for 5 min, and finally was maintained at 98% for 5 min at a constant flow rate of 300 nL/min on an EASY-nLC 1000 system. The eluent was sprayed via the NSI source at an electrospray voltage of 1.9 kV and then analyzed by tandem mass spectrometry in Q Exactive. Full-scan MS spectra (from m/z 350 to 1800) were acquired in the Orbitrap with a resolution of 70,000, and ion fragments were detected at a resolution of 17,500. The 20 most intense precursors were selected for subsequent analysis. Tree-based ion trap HCD fragmentation was performed at a collision energy of 30% in the MS survey scan with 45.0 s dynamic exclusion.

Protein identification and quantification
The LC-MS/MS raw data were searched against the Uni-Prot Oryza sativa Japonica Proteome Database with the software Proteome Discoverer (version 1.3, Thermo Scientific). Trypsin was chosen as the enzyme, and two missed cleavages were allowed. Acetylation in N-Term was set as a variable modification and carbamidomethylation (C) as a fixed modification and oxidation. The searches were performed with a peptide mass tolerance of 20 ppm and a product ion tolerance of 0.02 Da with a 1% FDR. The iTRAQ 8-plex was chosen for quantification during the search. Before data exportation, the results were filtered with p < 0.05 (both the significance threshold and the ion score or expected cutoff ).

Criterion for DEP detection
To identify statistically significant DEPs, the relative fold change (RFC) of proteins was determined by the ratio in the treated samples and their corresponding untreated controls according to a previously described method [24]. The RFC of proteins in response to chilling treatment was calculated as the ratio of CSL2/CSL1 (9311) or CTL2/ CTL1 (DC90), while the RFC in response to recovery was determined by the ratio of CSL3/CSL1 (9311) or CTL3/ CTL1 (DC90) (Fig. 1). The proteins with RFC ≥ 1.200 (p < 0.05) were considered up-regulated, and those with RFC ≤ 0.833 (p < 0.05) were considered down-regulated.

Protein hierarchical cluster analysis
Hierarchical cluster analysis was conducted on the quantified proteins from the biological replicates to assess the reproducibility of the MS data using a previously described method [30]. The relative abundance values of the identified proteins were subjected to Log2 normalization. Heatmaps and dendrograms were generated by pHEATMAP (R package) and g-plots [30], respectively. All the statistical analyses were performed in the R environment.

Quantitative real time reverse transcription PCR analysis
Total RNA was isolated from three additional sets of samples (three biological replicates per sample). cDNA synthesis was performed by reverse transcription (RT) with the Thermo Scientific RevertAid First Strand cDNA Synthesis Kit (Cat# K1622) according to the manufacturer's protocol. The sequences of DEP-encoding genes were downloaded from RGAP [31]. Primers were designed for qPCR to detect the relative expression of the selected DEPs at the transcript level, and PCR was performed by a Roche Lightcycler 480 Real-Time PCR System in 10 μL reactions with the SYBR Green PCR Master Mix kit (BIORAD, USA), following the manufacturer's protocol. The relative expression of each gene was calculated according to the 2 -△△CT method [32]. The GAPDH gene (LOC_Os04g40950) was used as an endogenous reference for qPCR.

GO enrichment analysis, functional classification, and expression profile analysis
The GO enrichment (http://geneontology.org) analysis was performed to annotate the DEPs under the categories of biological process (BP), cellular component (CC), and molecular function (MF) based on the GO-slim annotation database with Bonferroni correction by PAN-THER (http://pantherdb.org) [33]. DAVID was used to perform functional clustering of the DEPs with EASE score = 0.05/high-stringency clustering algorithms according to the user manuals of the online tools (https:// david.ncifcrf.gov/home.jsp) [34,35].
Expression profile analysis was performed with the Short Time-series Expression Miner (STEM) clustering algorithm [36] to obtain an overview of the expression changes in DEPs in response to chilling stress. In brief, the expression data of the DEPs was subjected to Log2 transformation, and the Gramene Oryza sativa gene annotation file (http://www.geneontology.org/page/download-annotations) was used as a reference. The STEM clustering method was applied in the analysis. All the other fields were set to the defaults, and the following customized settings were used: the 'Minimum Absolute Expression Change' of the 'Filtering' option was set to 0.02, the 'Correction Method' was 'Bonferroni Correction' , the 'Permutation Test Should Permute Time Point 0' was checked in the 'Model Profiles' option, and the 'Minimum Correlation' was set to 0.6 in the 'Clustering Profiles' option.

Results
Contrasting stress phenotypes of DC90 and 9311 in response to chilling treatment After 5 days chilling stress treatment, no obvious difference between DC90 and 9311 was visible (Fig. 2a, b). During the recovery period, phenotypic differences between DC90 and 9311 were observed to develop gradually (Fig. 2c). By the end of recovery, the 9311 seedlings were completely wilted, whereas the DC90 seedlings were able to survive the chilling stress, as reported previously [16]. Similar phenotypic differences can be obtained under hydroponic culture conditions (Additional file 2).
To understand how CTS-12 mediated the mechanisms underlying chilling stress responses at the proteomic level, hydroponically cultured rice seedlings were subjected to chilling stress and recovery following an optimal treatment scheme (Fig. 1a). The treated samples were subjected to proteomic analysis.

Identification and quantitation of proteins by iTRAQbased LC-MS/MS analysis
By using iTRAQ labeling LC-MS/MS analysis, a total of 4388 proteins were specifically identified from 1,178,076 LC-MS/MS spectra and 19,405 peptides in three independent experiments with a 1% FDR (false discovery rate) (Additional file 3; Additional files 1, 2, 3, 4; Additional file 5). Among these results, 2819 proteins were quantified in all replicates of all samples (Additional files 2, 3, 4).
Hierarchical clustering was conducted to acquire a comprehensive overview of the expression profiles of these proteins by using the pHEATMAP package in R. The results showed that the three replicates of each timepoint of 9311 and DC90 were clustered into the same respective clades (Additional file 6), indicating the reproducibility and reliability of our LC-MS/MS data. Notably, clustering based on mean of protein abundance revealed clear differences among the three timepoints (Fig. 3a), suggesting that the rice plants mobilized numerous proteins and differentially regulated their abundance to cope with chilling stress.

Identification of statistically significant DEPs
Student's t-test was applied to determine if the proteins in the untreated and treated samples were significantly different based on 2819 quantified proteins (p < 0.05). According to the criteria RFC ≥ 1.200 or ≤ 0.833 and p < 0.05, a total of 217 DEPs (144 up−/73 down-regulated) and 434 DEPs (191 up−/243 down-regulated) were identified in seedlings of 9311 and DC90, respectively, under chilling stress treatment (Fig. 3b). More than twice as many DEPs were identified in DC90 as in 9311, providing further evidence that many proteins in DC90 were regulated in response to chilling stress. Under recovery treatment, a total of 224 DEPs (122 up−/102 down-regulated) and 191 (88 up−/103 down-regulated) were identified in 9311 and DC90, respectively (Fig. 3b). Thus, a higher total number of DEPs were identified in DC90 than in 9311 (Fig. 3b, c; Additional files 3, 4, 5, 6). Among the 330 DEPs in 9311, 106 (32.10%) and 113 (34.20%) were specifically identified during the periods of chilling and recovery, respectively. One hundred and eleven (33.60%) were shared by both stages (Fig. 3c). However, in DC90, only 66 out of 559 DEPs (11.80%) were common to both stages, while 65.80% of DEPs (368 out of 559) were specifically identified in the chilling treatment stage (Fig. 3c). These results further supported the notion that different groups of proteins were mobilized in response to the two treatment stages in the chilling-tolerant genotype.
To compare the responsive DEPs between two genotypes during the whole treatment period, the following Fig. 2 Chilling stress phenotypes of DC90 and 9311 at the seedling stage. a, The DC90 and 9311 seedlings before exposure to low temperature; b, The DC90 and 9311 seedlings after 5-day chilling treatment; c, The DC90 and 9311 seedlings after 7-day recovery treatment. Scale bars = 5 cm in a, b, and c strategy was used to select statistically significant DEPs between DC90 and 9311 for further comparative analysis. First, all proteins with p < 0.05 in CSL2/ CSL1, CSL3/CSL1, CTL2/CTL1, and CTL3/CTL1 were selected to obtain 4 protein groups, which included up-regulated (RFC ≥ 1.200), down-regulated (RFC ≤ 0.833), and NA (0.833 < RFC < 1.200) subsets (Additional files 3, 4, 5, 6). Second, the datasets for CSL2/CSL1 and CSL3/CSL1, CTL2/CTL1 and CTL3/ CTL1, including their RFC data, were respectively combined to obtain the protein groups. The proteins with 0.833 < RFC < 1.200 in both stages and the proteins with p < 0.05 in only one treated stage were eliminated from each group (Additional files 4,5,6,7,8). Finally, totals of 206 and 155 were statistically significant DEPs were selected in 9311 and DC90, respectively. Further GO enrichment analysis, functional classification, and expression profile analysis were based on these two lists of proteins.

Validation of gene expression of differentially expressed proteins
Nine common DEPs were randomly selected for qPCR analysis to test the correlation of expression profiles between proteins and their corresponding mRNAs (Additional file 7). We compared the mRNA levels after 60-h chilling and 60-h recovery treatment with the iTRAQ data. The qPCR results showed that the expression profiles of 4 genes, LOC_Os03g15870.1, LOC_Os03g03360.1, LOC_Os08g33370.2, and LOC_Os03g50290.2, were in agreement with the iTRAQ data, whereas 5 genes (LOC_Os03g16050.1, LOC_Os08g37320.1, LOC_Os01g55 830.1, LOC_Os02g03860.1, LOC_Os03g22950.1) displayed wholly different expression patterns from those observed by the iTRAQ analysis (Fig. 4). The inconsistent expression profiles between qPCR and iTRAQ indicated the transcriptional level of genes were not always parallel to their corresponding protein level [37].

GO enrichment analysis and functional classification of DEPs
To further understand the functions of the identified DEPs, GO analysis and functional classification were performed. In the GO analysis, 177 and 142 protein IDs in 9311 and DC90, respectively, were mapped to the PAN-THER database [33] (Additional files 1, 2, 3, 4, 5, 6, 7, 8). The DEPs in both genotypes were significantly enriched in 14 biological processes, 12 cellular components, and 10 molecular function subgroups. Within the biological process subgroup, the cellular process, metabolic process, and primary metabolic process groups were prominent, indicating that the primary metabolic processes are easily affected in response to chilling stress ( Fig. 5a; Additional files 1, 2, 3, 4, 5, 6,7,8). The cell part, intracellular, organelle, and cytoplasm categories were four cellular compartments in which DEPs were highly localized ( Fig. 5b; Additional files 1, 2, 3, 4, 5, 6,7,8). In terms of molecular function, the DEPs from both genotypes were highly enriched in the same or similar sets of terms, such as binding, oxidoreductase activity, structural constituent of ribosome, structural molecule activity, and nucleic acid binding. However, translation elongation factor activity, translation regulator activity, GTPase activity, hydrolase activity, and catalytic activity were specifically represented in DC90 (Fig. 5c; Additional files 1, 2, 3, 4, 5, 6, 7, 8). Strikingly, the majority of over-represented cellular component GO terms were common to both genotypes. Furthermore, the identified DEPs consisted specifically of cytosol, ribosome, and plastid components (Fig. 5b). The functional clustering of DEPs by using DAVID gave results consistent with these findings (Fig. 5d). Among the over-represented functional categories, the structural constituent of Fig. 4 Comparative analysis of the protein and mRNA profiles of 9 representative DEPs. The X-axis represents the timepoints in the chilling and recovery treatments. The Y-axis indicates the normalized relative mRNA and protein levels. The green and blue lines represent the patterns of mRNA expression in 9311 and DC90, and the red and orange lines represent the patterns of protein expression in 9311 and DC90, respectively ribosome, protein-chromophore linkage, photosynthesis, light harvesting, and 14-3-3 protein groups were shared by both 9311 and DC90 (Fig. 5d).
Taken together, our results indicated that chilling stress-responsive proteins were mainly localized in the cytosol, ribosomes, and plastids/chloroplasts and that they might play vital roles in chilling stress-related responses.

Distinctly different dynamic changes of enriched chloroplastic DEPs in two genotypes during the period of chilling and subsequent recovery
To gain deeper insight into the dynamic changes in proteins between DC90 and 9311 during the whole period of chilling and recovery, we clustered the DEPs based on their abundance data using the STEM clustering algorithm [38]. The results showed that all the DEPs of DC90 and 9311 were clustered into 16 profile models (Additional file 9; Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10). The query lists of the two genotypes had 55 DEPs in common, but only 13 (23.64%) had the same profiles (Additional files 1, 2, 3, 4, 5, 6, 7, 8, 9, 10). Among these 55 DEPs, interestingly, 73.36% of which in 9311 either displayed a continuous increase (models 12 and 13) or decrease (models 0, 2, and 3) in abundance over the time period studied or retained the chilling-induced abundance during the subsequent recovery (models 4 and 11).

Different H 2 O 2 scavenging capabilities in DC90 and 9311
Chilling stress induced ROS accumulation in plant [13]. Furthermore, three down-regulated ROS scavenging associated DEPs, prx11 (Q9FYP0), CATB (Q0D9C4), and catalase (Q10S82), were specifically identified in 9311 (Additional files 2, 3,4,5,6,7,8,9,10). The enzymes catalyze the conversion of H 2 O 2 to H 2 O and O 2 [9]. To investigate the possible correlation of the accumulating ROS and dynamic changes in chloroplast photosynthetic and ribosomal proteins, ROS production and scavenging in the samples from the chilling and recovery stages were tested histochemically and physiologically. DAB and NBT staining showed that both DC90 and 9311 began to accumulate ROS at an early timepoint in the chilling treatment (12 h) (Fig. 7a). Beginning with the 36-h timepoint, the DAB-stained spots in DC90 gradually disappeared, whereas H 2 O 2 continued to accumulate until 72 h in 9311, suggesting that different H 2 O 2 scavenging capabilities in DC90 and 9311. This result was confirmed by the analysis of H 2 O 2 and superoxide scavenging-related enzymes, and physiological changes in rice plants (Fig. 7b, c; Additional file 11 C, D, E). In line with proteomic results, the CAT activity was significantly inhibited after 60 h of chilling stress treatment by comparison with that of DC90 and might be a certain reason for the overproduction of H 2 O 2 in 9311, although no difference was observed at the 0-h timepoint or after 60 h of recovery in either genotype. The activity of another H 2 O 2 -scavenging enzyme, POD, was induced at the 60-h chilling and 60-h recovery timepoints in both DC90 and 9311 by comparison with the untreated control, but no significant differences between genotypes were observed at either timepoint. The activity of SOD and APX showed very similar patterns in DC90 and 9311 at each timepoint (Fig. 7a, d; Additional file 11 E).
The results suggested that the overproduction of H 2 O 2 under chilling stress might account for the chilling sensitivity phenotype of 9311. This idea was further supported by the exogenous application of ascorbic acid (AsA) to 9311 under chilling stress (Fig. 7e). AsA is an antioxidant reagent that participates in H 2 O 2 scavenging in plant cells. The exogenous application of AsA rescued the chilling-induced cell death in 9311.
Taken these together, the abundances of both common and specific proteins identified in DC90 and 9311 exhibited completely different chilling-and recovery-induced regulation. The results also suggested that the CTS-12 might mediate the chilling-induced cold acclimation and de-acclimation of the abundance of photosynthetic and ribosomal proteins in the periods of chilling and recovery. Thereby conferred wild rice chilling tolerance under severe chilling stress. Meanwhile, the balance between H 2 O 2 production and elimination in rice plants under chilling stress might play a certain role in biochemical cold acclimation and de-acclimation.

Discussion
In this study, DC90 was used to dissect the mechanisms underlying the cold tolerance conferred by DP15 via comparative proteomic analysis. Two continuous treatment stages were designed to facilitate obtaining a global picture of the dynamic regulation performed during the transition between the chilling and recovery stages.

Chilling-responsive downstream regulation in rice occurs mainly in chloroplasts
The chloroplast is an important organelle for photosynthesis in plant cells. However, recent studies showed that the chloroplast not only plays its universal role in photosynthesis but also participates in the response to environmental stresses [39,40]. Specifically, photosynthetic  genes/proteins were regulated in response to environmental stresses [11,17,18]. When rice plants were subjected to chilling, drought, or salt stresses, photosynthetic activity was greatly inhibited. Genes associated with chlorophyll biosynthesis, chlorophyll a/b-binding proteins, the light-harvesting complex, the PSII and PSI core complexes, and chloroplast precursors were significantly down-regulated [11,17,18]. Ribosomal proteins also perform extra functions in response to environmental stresses [16,[41][42][43]. For example, in soybean, chilling induced the expression of 3 ribosomal genes, RPS6, RPS13, and RPL37 [44]. In rice, 34 RPLs in tissues covering the major stages of rice growth were highly responsive to diverse stresses, indicating that ribosomal genes appear to play roles in stress amelioration in addition to housekeeping [45]. In this study, among 55 common DEPs, two major group proteins were highlighted. One consisted of 23 chloroplast-localized proteins, including core components of photosystems I and II, such as psbA, psbC, psbB, and PSI-N, as well as proteins involved in chlorophyll binding (CAB1R), light harvesting (psaA), and photoprotection (Q943K1). The other consisted of 6 chloroplastic ribosome constituents, Q10NM5, RPL18, RPL5, rps9, Q6KA00, and RPL12-2. STEM analysis showed that both groups of DEPs were significantly regulated in response to chilling and recovery but exhibit entirely different regulatory patterns in 9311 and DC90 ( Fig. 6; Table 2; Additional files 9, 10). Thus, the chloroplast is inferred to be an important target organelle for downstream regulation in response to environmental stresses in wild rice.
Rice adaptation to chilling stress by the proteomic regulation of cold acclimation and de-acclimation of chloroplastic proteins Plant growth and distribution are affected by low-temperature stress. Overwintering plants can adapt to freezing temperature via the orchestrated expression of several sets of genes. This process, called cold acclimation, triggers complex biochemical and physiological changes in plant cells, including alterations in membrane fluidity and lipid composition, the accumulation of compatible solutes and the regulation of gene expression [46][47][48]. De-acclimation is an important reverse regulatory process that enables overwintering plants to percept warm temperatures and resume growth in spring. Although rice is not an overwintering plant, the dynamic changes in the majority of chloroplast photosyntheticand ribosome-associated proteins in DC90 and 9311 during the transition from chilling to recovery appeared to reflect the phenomena of cold acclimation and de-acclimation at the proteomic level (Fig. 6a, b; Table 2; Additional file 10). Intriguingly, the dynamic changes in these proteins during the recovery stage were distinct between the two genotypes. In 9311, both the ribosomeand photosynthesis-associated DEPs were incapable of regaining their normal abundance, indicating that seemed proteomic de-acclimation in 9311 was impaired severely even when the plants were restored to optimal conditions (Fig. 6a). However, strikingly, the majority of DC90 DEPs, in contrast to those of 9311, returned to normal abundance during recovery, indicating that the DC90 plants underwent proteomic cold acclimation and de-acclimation regulatory process during the whole period of the experimental regime (Fig. 6b). Therefore, these processes might facilitate the survival of rice plants under severe chilling stress.
Imbalance of H 2 O 2 production and scavenging may be responsible for impairing the de-acclimation of chloroplastic proteins in 9311 Oxidative burst, which is caused by the generation of large quantities of reactive oxygen species (ROS, e.g., O 2 .-, H 2 O 2 , OH˙, 1 O 2 ), is one of the earliest responses of plant cells to various environmental stresses [49]. ROS accumulation under environmental stresses may cause damage to lipids, DNA, and proteins; ROS also act as signaling molecule to trigger downstream molecular regulation processes, including ROS scavenging and DNA/protein damage repair. Therefore, the balance between ROS production and the activities of antioxidative enzymes determines whether oxidative signaling and/or damage will occur [50]. DAB and NBT staining showed ROS accumulation in both genotypes at an early timepoint in the chilling treatment (12 h) ( Fig. 7a; Additional file 11 A, B). However, our study showed no significant differences in superoxide accumulation and scavenging between DC90 and 9311 at the tested timepoints (Fig. 7a, d), whereas the scavenging of H 2 O 2 was distinctly different between DC90 and 9311. The accumulated H 2 O 2 in DC90 was gradually scavenged beginning at the 36-h timepoint in the chilling treatment, whereas the increasing accumulation of H 2 O 2 occurred from 12 to 72 h in 9311, resulting in cell death at the later timepoints (84-120 h). This result was confirmed by the analysis of the activity of H 2 O 2 scavenging enzyme, catalase. The activity of CAT in 9311 was significantly down-regulated after 60 h of chilling treatment, whereas there was not significantly induced in DC90 (Fig. 7a, b). This in line with three ROS scavenging associated DEPs, prx11 (Q9FYP0), CATB (Q0D9C4), and catalase (Q10S82), which were specifically identified in 9311 in LC-MS/MS. They were all down-regulated under chilling stress. This behavior indicated that DC90 can maintain the balance between H 2 O 2 production and elimination and survive severe chilling stress. The ability to scavenge ROS and reduce their damaging effects may be Fig. 8 Working model proposed to depict the CTS-12-mediated mechanism underlying the chilling tolerance conferred by common wild rice. Yellow ellipses designate the abundance of proteins at the optimal temperature. Green ellipses designate the down-regulated abundance of proteins in the chilling or recovery stage. Red ellipses designate the up-regulated abundance of proteins in the chilling or recovery stage. RPs designates ribosomal proteins; PSPs designates photosystem proteins related to the regulatory mechanism of cold acclimation and de-acclimation at the proteomic level in wild rice.
Cold acclimation and de-acclimation of photosynthetic and ribosomal proteins plays a vital role in the survival of rice under chilling stress In summary, a working model was proposed to depict the proteomic cold acclimation and de-acclimation in rice cells adapting to severe chilling stress and to speculate on the role of CTS-12 in these processes (Fig. 8). When rice plants are exposed to low temperature, the cell membranes percept the chilling stress and trigger a decrease in the rate of protein translation and photosynthetic activity, which manifests as the down-regulation of chloroplast photosynthetic and ribosomal proteins. This down-regulation appeared to be an important process in the cold acclimation of rice plants. When the adverse environment is removed, the capability of de-acclimation of rice plants will survive chilling stress. The inhibition of photosynthetic activity reduces the utilization of absorbed light energy. The excess light energy results in the generation of toxic ROS, which in turn triggers photodamage to photosystem proteins, such as D1, a core component of photosystem II [51]. Damaged D1 is removed by FtsH proteolytic degradation [52][53][54]. Therefore, the turnover of D1 protein is thought to be a key regulatory step in the PSII repair cycle [55]. In this study, ribosomal proteins, which are constituents of the protein translation apparatus, were significantly inhibited by chilling stress, and thereby inhibited the turnover of D1 protein. That is, the balance between ROS scavenging and production must be maintained for the de-acclimation of rice plants when the adverse environment is removed. In 9311, the imbalance of H 2 O 2 production and scavenging may cause the interruption of D1 protein turnover and cell death, and ultimately impair de-acclimation of rice plants when the adverse environmental factors are removed. Thus, we speculate that the major QTL, CTS-12, might directly or indirectly mediate the process of cold acclimation and de-acclimation of photosynthetic and ribosomal proteins in wild rice under severe chilling stress and recovery stage, respectively. Moreover, the differences of H 2 O 2 production and scavenging between 9311 and DC90 arouse us to question what's the possible correlation between H 2 O 2 balance and CTS-12-mediated the cold acclimation and de-acclimation of photosynthetic and ribosomal proteins in rice under chilling stress. This will be addressed in our further investigation.

Conclusion
Using iTRAQ-based proteomic analysis, we compared a chilling tolerance CSSL, DC90 with its chilling sensitive recurrent parent, 9311 under chilling stress and recovery treatment. In our results, 206 and 155 DEPs were identified in 9311 and DC90, respectively. Strikingly, the majority of which were enriched in the 'structural constituent of ribosome' , 'protein-chromophore linkage' , and 'photosynthesis and light harvesting' categories. And STEM analysis based on these DEPs revealed distinct dynamic responses of both chloroplast photosynthetic and ribosomal proteins between 9311 and DC90. Based on the proteomic analysis and physiological investigation of ROS in 9311 and DC90, we proposed that CTS-12, a QTL conferred DC90 chilling tolerance, might mediate the dynamic response of chloroplast photosynthetic and ribosomal proteins under chilling and recovery, and thereby enhancing the chilling stress tolerance of wild rice. Furthermore, the balance of H 2 O 2 production and elimination in rice plants may play a certain role in maintaining the dynamic response of chloroplast photosynthetic and ribosomal proteins in whole period of chilling and recovery stress (cold acclimation / de-acclimation). The identified DEPs and the involvement of CTS-12 in mediating the dynamic response of DC90 at the proteomic level illuminate and deepen the understanding of the mechanisms that underlie chilling stress tolerance in wild rice.