An Integrative Study Showing the Adaptation to Sub-Optimal Growth Conditions of Natural Populations of Arabidopsis thaliana: A Focus on Cell Wall Changes

In the global warming context, plant adaptation occurs, but the underlying molecular mechanisms are poorly described. Studying natural variation of the model plant Arabidopsis thaliana adapted to various environments along an altitudinal gradient should contribute to the identification of new traits related to adaptation to contrasted growth conditions. The study was focused on the cell wall (CW) which plays major roles in the response to environmental changes. Rosettes and floral stems of four newly-described populations collected at different altitudinal levels in the Pyrenees Mountains were studied in laboratory conditions at two growth temperatures (22 vs. 15 °C) and compared to the well-described Col ecotype. Multi-omic analyses combining phenomics, metabolomics, CW proteomics, and transcriptomics were carried out to perform an integrative study to understand the mechanisms of plant adaptation to contrasted growth temperature. Different developmental responses of rosettes and floral stems were observed, especially at the CW level. In addition, specific population responses are shown in relation with their environment and their genetics. Candidate genes or proteins playing roles in the CW dynamics were identified and will deserve functional validation. Using a powerful framework of data integration has led to conclusions that could not have been reached using standard statistical approaches.


Introduction
In the global warming context, temperature fluctuations will be major concerns [1]. Local adaptation is defined as the response of natural populations following their interactions with their environmental conditions. Freezing events, without any preceding chilling period, associated to sudden elevations of temperature, are critical for plant development. Thus, coping with thermal constraints will be a major challenge to maintain agricultural productivity by doing a selection of warm-adapted and cold-resistant species [2]. Looking at natural variations within wild species is crucial to elucidating the molecular bases of the phenotypic adaptation [3].
(L'Hospitalet-près-l'Andorre), which are living between 700 and 1400 m in altitude in the Pyrenees Mountain [38]. We have used one genotype per homogeneous population, and 20 plants of this genotype for each of the three biological replicate for all the analyses. Seeds were sowed in Jiffy-7 ® peat pellets (Jiffy International, Kristiansand, Norway). After 48 h of stratification at 4 • C in darkness, plants were grown at two different temperatures, at 22 or 15 • C, under a light intensity of 90 µmol photons.m −2 s −1 . They were cultivated under a long-day condition (16 h light/8 h dark) with 70% humidity. Rosettes were collected just before bolting (stage 5.10 [41]) from four-or six-week-old plants grown at 22 or 15 • C, respectively. Floral stems were also collected at the same stage of development (first flower, stage 6 [40]): 6 weeks for Col; 7 weeks for Grip and Roch; and 8 weeks for Hosp and Hern. They were collected two weeks later for all the plants grown at 15 • C.

Macrophenotyping
Rosettes of plants grown at 22 and 15 • C were analyzed at the time of sampling. We have measured their diameter and fresh mass. The number of leaves was counted. Before freezing, pictures were taken to measure the rosette areas with the ImageJ software [42]. Regarding the floral stems, the length, the number of cauline leaves, the mass and the diameter at the base of the floral stem were measured before freezing. Data are described in [43].

Histological Staining of CWs
Whole rosettes and base of floral stems were harvested in 50 mL Falcon tubes, and rapidly infiltrated under vacuum with FAA (10% formalin (37% formaldehyde solution, Sigma-Aldrich, Saint-Quentin Fallavier, France); 50% ethyl alcohol; 5% acetic acid; 35% distilled water). They were fixed for 16 h at 4 • C. The dehydration and paraplast infiltration protocol was described in [44]. The whole rosettes were processed and stained as described in [16] to keep the phyllotaxy intact. The stained sections were analysed using a Nanozoomer slide scanner (Hamamatsu, Shizuoka, Japan) to produce whole slide scan at 20× resolution. The scans were analyzed using NDP view (Hamamatsu) and the images were directly extracted from the viewer to assemble the figure.

Extraction of Proteins from Purified CWs
CW purification was performed as described [45]. We have carried out the sequential extraction of proteins from purified CWs as described, successfully using solutions containing CaCl 2 0.2 M and LiCl 2 M and buffered at pH 4.6 [46]. Typically, 0.2 g of lyophilized CWs was used for one extraction and about 500 µg proteins were extracted. The final protein extract was lyophilized. Proteins were quantified with the CooAssay Protein Assay kit (Interchim, Montluçon, France).

Sequential CW Polysaccharides Extraction and Identification
The sequential extraction of CW polysaccharides was performed in four steps and detailed in [16]. In summary, 100 mg of a deproteinized CW fraction were used. Four successive extractions were carried out to obtain extracts enriched in pectins (E1 and E2, respectively using CDTA 50 mM and 50 mM Na 2 CO 3 ) and hemicelluloses (E3 and E4, respectively using 20 mM NaBH 4 and 4 M NaOH). Each extract was hydrolyzed in 2 N TFA for 1 h at 120 • C. After 10× dilution in UHQ water, monosaccharides were analysed by High-Performance Anion-Exchange Chromatography coupled to Pulsed Amperometric Detection (HPAEC-PAD; Dionex, Sunnyvale, California, USA) using a CarboPac PA1 column (Dionex). L-Fuc, L-Rha, L-Ara, D-Gal and GalA (Sigma-Aldrich); D-Glc (Merck, Darmstadt, Germany); D-Xyl (Roche, Mannheim, Germany) were used as standard monosaccharides for identification and quantification. Data are described in [43]. CW polysaccharides reconstruction was performed using formula previously described [16,43] and adapted from [47].

Identification of Proteins by LC-MS/MS
The identification of proteins extracted from CWs was performed by LC-MS/MS at the PAPPSO proteomic platform (pappso.inrae.fr) after tryptic digestion in solution as described [48]. Parameters for MS data processing in the X!Tandem software (JACKHAMMER, 2013.6.15, www.thegpm.org/tandem/) and the X!Tandem Pipeline 3.3.4 [49] are detailed in [48]. Tryptic digestion was declared with no possible miscleavage. Only proteins identified with at least two different specific peptides in the same sample and found in at least two biological replicates were validated. Furthermore, quantification was performed on peptides with standard deviation retention times lower than 20 s and peak width lower or equal to 100 s. Raw data are described in [50].

RNA Sequencing
Protocols of the transcriptomic analysis were as described in [16]. RNA-Seq data and procedure details are available in Supplementary Table S2. Resulting sequences are available at NCBI short read archive (SRA, BioProject PRJNA344545). Raw data are described in [50].

Bioinformatics Annotation of Proteins and Quantification
The prediction of sub-cellular localization and functional domains of proteins was performed with the ProtAnnDB tool (www.polebio.lrsv.ups-tlse.fr/ProtAnnDB/) [51]. A protein was considered as a CWP if two bioinformatics programs predicted it as secreted, no intracellular retention signal was found and no more than one trans-membrane domain was found as described [52]. We have checked the prediction of the cell wall localization of many of the identified CWP families by doing a literature search as well as a mining of the WallProtDB database dedicated to cell wall proteomes (www.polebio.lrsv.ups-tlse.fr/WallProtDB). Quantification was only operated for CWPs using the MassChroQ software [53] as in [54]. The functional annotation of CWPs is given according to WallProtDB. All the LC-MS/MS data have been deposited at PROTICdb (http://proteus.moulon.inra.fr/ w2dpage/proticdb/angular/#/projects/197) and CWP MS data are also available at WallProtDB.

Statistical Analysis and Data Integration
Data analysis and methodology were achieved as described [40]. We have performed multi-block supervised analyses with the sparse version Multi-Block Partial Least Squares Discriminant Analysis (MB-PLS-DA) implemented in the mixOmics R package version 6.1.1. (http://www.bioconductor.org/ packages/release/bioc/html/mixOmics.html) [39]. They have taken into account all phenomics and metabolomics variables associated to a restrictive pool of 40 CWPs and 40 genes having modified transcript levels and encoding CWPs identified in the proteomic studies. To improve the visualization, the clustered image map was done with the high scoring 20 proteins and 20 genes and all the variables of the phenomics and metabolomics blocks.

Results
We have selected four A. thaliana populations (Roch, Grip, Hern and Hosp) among the 30 new populations collected in the Pyrenees Mountains recently characterized [38]. They constitute four independent homogeneous populations originating from various altitudes (Table 1) and they do not belong to the genetic cluster of Col (named III) which was taken as a reference ecotype. Roch and Grip belong to the same genetic cluster (named I), whereas Hern and Hosp belong to a third one (named II) (see [37] for the genetics characterization of the populations). In addition, within each genetic cluster, one population originates from moderate altitude (Roch and Hern) and the other one from high altitude (Grip and Hosp). The altitude clustering is also illustrated by the climate PC1 index that summarizes all the environmental and climatic data associated with each population.  [38]. 3 Climate PC1 (index that describes and explains the climate gradients). 4 Altitude at which Col was originally collected.

The Morphological Phenotypes Are Temperature-and Organ-Dependent
The macro-and micro-phenotypes of rosettes of plants grown at 22 and 15 • C were analyzed. We have collected floral stems at equivalent degree days in order to be at the same developmental stage despite differential stem grow speed. The impact of the 15 • C growth temperature on rosettes was striking. The number of rosette leaves increased for all the Pyrenean population as well as for Col (e.g., a 60% increase for Col, Figure 1B), compared to that of plants grown at 22 • C. Besides, the number of rosette leaves at the time of bolting could be contrasted: Col plants grown at 15 • C had the same number of rosette leaves as Roch grown at 22 • C and Hosp grown at 15 • C. The higher the altitude at which the Pyrenean populations were collected, the higher the number of rosette leaves prior bolting whatever the growth temperature (compare data in Figure 1B and Table 1). Except for Hern, an increase of the number of cauline leaves at the first flower-stage was observed for the 15 • C growth temperature ( Figure 1C). Similarly, Hern had a tendency to show opposite floral stem traits (mass, diameter and length) as compared to the other populations grown at 15 • C (Supplementary Figure S1).
The leaf number, combined with the mass, the diameter, the density and the projected area of rosettes (Supplementary Figure S2) were integrated in the statistical study described thereafter (see § 2.5 and 2.6).
We have completed this macro-phenotyping with a microscopic phenotyping of the petioles of rosettes leaves. The petioles were cut into cross-sections ( Figure 2). Two examples of the whole phyllotaxy present on the tissue arrays are given for Col ( Figure 2U) and Roch ( Figure 2V) grown at 15 • C to show the homogeneity of the staining pattern in a given condition. In a CW context, the astra blue/safranin red staining has a relative specificity for hydrophilic polysaccharides stained in blue and hydrophobic CW compounds stained in red. At both 22 and 15 • C growth temperatures, the petiole of Col displayed a typical staining pattern with astra blue-positive cortical parenchyma CWs vs. safranin red-positive epidermis cuticle and vasculature lignified walls. At 22 • C, and in reference to Col cross sections, the safranin red staining was more intense in the epidermis cuticle of Grip, Hern and Hosp, and somehow in the vasculature of Grip and Hern, probably reflecting the corresponding cuticle reinforcement and vasculature increased lignification, respectively. At 15 • C, the safranin red staining was more intense in all the tissues of Roch, Grip, Hern and Hosp as compared to Col. These included the epidermis, the lignified vasculature as well as the cortical parenchyma. Similarly to Grip, Hern and Hosp grown at 22 • C, the epidermis and vasculature safranin red intense staining probably corresponded to cuticle reinforcement and vasculature increased lignification, respectively. However, the nature of the hydrophobic CW compounds responsible for the cortical parenchyma safranin reactivity was unclear.

The CW Composition of Aerial Organs Varies Depending upon Temperature Growth Conditions
We have analyzed the monosaccharide composition in rosettes and stems grown at 15 and 22 • C. We could reconstruct major CW polysaccharides as described [16,51]. Rosettes had higher ratios of linear to branched pectins than floral stems. Rosettes and floral stems contained lower ratios at 15 • C than at 22 • C ( Figure 3). Furthermore, the ratio variations between the two growth temperatures were higher in floral stems than in rosettes for all the populations with the exception of Hern. We have performed a similar analysis for the other polysaccharides (Supplementary Figures S3 and S4). To summarize, rosettes of plants grown at 15 • C contained more xyloglucans and a higher ratio of branched rhamnogalacturonans I than those of plants grown at 22 • C. Conversely, floral stems of Col, Grip and Hosp contained more xyloglucans when grown at 22 • C than at 15 • C. Interestingly, for rhamnogalacturonans I, homogalacturonans and to a lesser extend xyloglucans, differences between populations were more pronounced in floral stems at 15 • C than at 22 • C. Finally, homogalacturonans (also called linear pectins) were less abundant in rosettes and stems of plants grown at 15 • C than in those of plants grown at 22 • C.

The CW Proteome Varies at Sub-Optimal Growth Temperature
The CW proteomic analysis allowed identifying 364 and 414 different CWPs in rosettes and floral stems, respectively (Supplementary Table S1). Eighteen and 29 new CWPs were identified in this study as compared to the previously published A. thaliana Col and Sha rosette and floral stems proteomes, respectively [48,54]. We have found CWPs specifically accumulating in certain populations and/or at either growth temperature. For instance, AT5G59320 (AtLTP3) and AT5G59310 (AtLTP4) were more abundant at 15 • C, whereas AT4G37800 (AtXTH7) was more abundant at 22 • C in all the Pyrenean populations and Col. AT1G78450 (unknown function) was specific to Roch (see Supplementary Table  S1 for the detailed results).
Even if the principal component analysis (PCA) is an unsupervised method, it clearly highlighted groups of samples and allowed evaluating the good reproducibility of a given experiment [40], when performed on quantified CWPs ( Figure 4). The analysis of the two organs showed a growth temperature-specific response of each population at the CW proteome level. Regarding rosettes, the first component (PC1) separated the samples according to the growth temperature (22 vs. 15 • C) and the second one (PC2) according to the population (Hosp vs. the other populations) ( Figure 4A). For floral stems, PCA clustered together the samples according to the genetic clusters defined above: clusters III (Col) and I (Roch/Grip) vs. cluster II (Hern/Hosp) ( Figure 4B). The growth temperature affected less clusters III and I than cluster II as shown by the better clustering of Col, Roch and Grip compared to Hern and Hosp. Col, Roch, Grip and Hosp grown at 22 • C were located farther right on the PC1 axis as compared to the corresponding populations grown at 15 • C. Conversely, the Hern samples showed a reverse location with regard to PC1 compared to the other populations.

Regulation of Gene Expression Depends on Both Growth Temperatures and Populations
RNA-Seq analyses allowed the identification of 19,763 and 22,570 genes with significant levels of expression in rosettes and floral stems respectively (see raw data and statistical analyses in Supplementary Table S2). As for CW proteomics data, the unsupervised multivariate analysis of the "whole transcriptome" of rosettes ( Figure 5A) and floral stems ( Figure 5B) showed the reproducibility of the data by the clustering of the biological replicates.
The transcript levels of the genes encoding CWPs identified by the CW proteomic analysis were extracted (Supplementary Table S2) in order to correlate proteomics and transcriptomics data. PCA restricted to these transcripts (defined as "CW transcriptome") showed similar separation as that obtained with the "whole transcriptomes" (compare Figure 5A with Figure 5C; Figure 5B with Figure 5D), as observed in a previous study [16].
For rosettes, the analysis brought out the strong effect of the temperature that separated the samples on PC1 ( Figure 5A,B). The PCAs of "whole transcriptomes" and "CW transcriptomes" of rosettes both isolated Hosp grown at 15 • C. A closer analysis indicated that the set of genes explaining the discrimination of Hosp grown at 15 • C was not related to the CW proteomes (Supplementary Table  S2). In turn, the analysis focused on the CWP transcriptomes highlighted a few candidates such as AT2G25510 (unknown function), AT2G14610 (Cys-rich secretory protein) which were overexpressed or AT5G26000 (glycoside hydrolase family 1) which was down-regulated. In floral stems, PCA highlighted the particular behaviors of Hern and Col which were both separated from the other samples ( Figures  4B and 5B,D; Supplementary Table S2). The impact of the growth temperature on gene expression levels was observed in both rosettes and floral stems for the Pyrenean populations and for Col and allowed discriminating genes exhibiting specific expression patterns.

An Integrative Study Highlights the Different Responses of the Rosettes and Floral Stems According to the Growth Temperature
Thereafter, we refer to each omics dataset (phenomics, metabolomics, proteomics and transcriptomics) as a block and we have used a previously described framework for the integration of the omics data [40].
To provide an insight into the results of the integrative process, we have considered the pairwise correlations between each pair of blocks [39]. The higher the values, the stronger the relationships between the blocks ( Figure 6A). In rosettes, these values were higher than 0.83 (value obtained between CW transcriptomics and phenomics). These values were mainly due to a strong temperature effect as indicated in the segregation of two groups in Supplementary Figure S6A. The highest correlation was found between the CW transcriptomics and CW proteomics blocks (0.97). All the blocks showed a clear segregation between the growth temperatures according to the first component that confirmed the strong effect of the temperature ( Figure 6A and Supplementary Figure 5A). Using the sparse multi-block partial least squares discriminant analysis (MB-PLS-DA) analysis, we could identify the variables from each block involved in the discrimination according to the growth temperature (Supplementary Figure S5A and Table S3). The clustered image map also includes the results of a hierarchical clustering performed jointly on both the variables and the samples. Again, the samples were clearly discriminated according to growth temperatures. The list of discriminated transcripts and CWPs is given in Table 2. Among this set of CW transcripts and CWPs, 23% corresponded to proteins acting on CW polysaccharides, 18% to proteases and 18% to proteins possibly related to lipid metabolism. Besides, most changes occurred at the 15 • C growth temperature (30 out of the 40 candidates). The statistical analysis of the data was performed with a sparse MB-PLS-DA and selected the 20 best hit among each block. 1 ID: proteins with interaction domains (with proteins or polysaccharides); LM: proteins possibly related to lipid metabolism; M: miscellaneous proteins; OR: oxido-reductases; P: proteases; PAC: proteins acting on CW polysaccharides; S: proteins possibly involved in signaling; UF: proteins of yet unknown function. 2 "+"corresponds to proteins or transcripts showing an increased level of accumulation at a given growth temperature. The fold changes can be read in Supplementary Tables S1A and S2B. For floral stems, we have observed lower correlations between block pairs (Supplementary Figure S6A). As for rosettes, we have found the highest correlation (0.9) between the proteomics and CW transcriptomics blocks and the lowest (0.63) between the phenomics and metabolomics blocks. Low correlation could be due to the weak segregation between the growth temperatures in the first component obtained in the phenomics and metabolomics blocks. Table 3 provides an overview of the CW transcripts and CWPs whose levels of accumulation were modified upon different growth temperatures. Among these CW transcripts and CWPs, 22% corresponded to oxido-reductases and 22% to proteins possibly related to lipid metabolism. The statistical analysis of the data was performed with a sparse MB-PLS-DA and selected the 20 best hit among each block. 1 AGI codes in bold means that the protein and the corresponding transcripts were both differentially accumulated. 2 See legend to Table 2, for the meaning of the following abbreviations: ID, LM, M, OR, P, PAC, S, UF. 3 "+" corresponds to proteins or transcripts showing an increased level of accumulation at a given growth temperature. The fold changes can be read in Supplementary Tables S1B and S2D.
The correlation observed between block pairs was probably lowered by the specific phenotype of Hern floral stems at 15 • C (Supplementary Figure S6B). The important effects of the growth temperature and the specificity of Hern were also highlighted by the hierarchical clustering of the samples in the clustered image map of the discriminating variables (Supplementary Figure S6C, Table S4). The dendrogram discriminated the samples on the basis on growth temperatures: Roch, Grip, Col and Hosp grown at 15 • C clustered together whereas all the plants grown at 22 • C together with Hern grown at 15 • C formed another cluster.

The Hosp and the Hern Populations Exhibit Particular Phenotypic Responses Regarding Rosette and Floral Stem Development
The integrative analysis has allowed the identification of relevant candidate genes/proteins related to the temperature effect. The next step was to analysis the specific responses of the different populations at the level of the two analyzed organs, rosettes and floral stems.
For rosettes, when we have compared the metabolomics and the phenomics blocks, the correlation varied from 0.87 between the CW proteomics and CW transcriptomics blocks to 0.56 were compared ( Figure 7A). Regarding the individual plots, each block discriminated different populations (Supplementary Figure 5B): (i) in the CW proteomics block, we have observed three clusters, namely Roch/Grip, Col/Hern and Hosp; (ii) we have found a similar clustering in the CW transcriptomics block with a better segregation; (iii) the phenomics block only allowed a separation of Roch at 22 • C; (iv) the metabolomics block did not allow distinguishing the populations and Col. These differences between the Pyrenean populations profiles resulted from the contrasted values of the variables in each block represented in Supplementary Figure S5B.
The clustered image map of the samples for each block involved in the discrimination according to the Pyrenean populations and Col showed a segregation into three groups ( Figure 7B). From top to bottom, the first one contained only Hosp, the second one Roch grown at both temperatures and Grip grown at 15 • C, and the third one was composed of Hern and Col grown at both temperatures and Grip grown at 22 • C. The specific profile of Hosp was supported by a pool of CW transcripts and CWPs sorted in Table 4 (see also Supplementary Table S4B). One fourth of them corresponded to proteins possibly involved in lipid metabolism and 27% to proteins of yet unknown function.  For floral stems, the correlations between components from each block were highly variable from 0.95 between the CW proteomics and CW transcriptomics blocks to 0.25 between the metabolomics and the phenomics blocks (Supplementary Figure S7A). This large range of correlation was due to the good clustering of the samples observed between the CW proteomics and CW transcriptomics blocks and the poor clustering and large variability observed for the metabolomics and the phenomics blocks (Supplementary Figure S7B). This observation was confirmed by the strong segregation of Hern along the first axis and of Hosp along the second one regarding the CW proteomics and CW transcriptomics blocks on the individual plots (Supplementary Figure S7B). On the other hand, only Hern could be separated between the two growth temperature conditions with the phenomics block, but not as clearly as with the CW proteomics and CW transcriptomics blocks. At the end, the Pyrenean populations and Col were well separated by a few discriminating variables (Supplementary Figure S7C).
Altogether, this analysis allowed segregating the samples in three major groups: (i) Hern, (ii) Hosp, and (iii) Col, Grip and Roch. Using the sparse approach, we could sort a restricted pool of transcripts and CWPs differentially accumulated in floral stems of Hosp and Hern vs. Col, Grip and Roch (Table 5,  Supplementary Table S4B). For Hosp, the CW transcripts and proteins mainly corresponded to proteins acting on cell wall polysaccharides (35%), whereas for Hern, proteases and proteins acting on cell wall polysaccharides represented 42% and 35% of the set, respectively.
The statistical analysis of the data was performed with a sparse MB-PLS-DA and selected the best hit among each block. 1 AGI codes in bold mean that the protein and the corresponding transcripts were both differentially accumulated. 2 See legend to Table 2, for the meaning of the following abbreviations: ID, LM, M, OR, P, PAC, S, UF. 3 "+/−"correspond to proteins or transcripts showing an increased (+) or a decreased (-) level of accumulation. The fold changes can be read in Supplementary Tables S1B and S2D.

Discussion
A comparative study of four populations of A. thaliana originated from a Pyrenees Mountains altitudinal gradient and of the well-described Col ecotype was performed to better understand the CW plasticity in response to environmental changes. We have analyzed the plant response to sub-optimal growth conditions, i.e., at 15 • C, with an integrative approach using heterogeneous omics datasets to evaluate the impact on the CW. Altogether, we could show (i) a major effect of the growth temperature, (ii) the organ specificity of the response, and (iii) the particularities of Hosp and Hern, compared to Grip, Roch and Col. Since many CW transcripts and CWPs were highlighted in this work, we have focused on a few of them in the following discussion.

Temperature Has an Important Impact on Plant Development and on the CW
The observed macro-phenotypes have highlighted the ubiquitous impact of the growth temperature on the development of both rosettes and floral stems for all plants. The petioles of the plants grown at 15 • C, especially the Pyrenees Mountains population, presented a more intense staining of hydrophobic compounds in the cuticle epidermis of the cross-sections and around vessels. In the former case, such compounds could be present in the cuticle. This would be consistent with the observations reported in a previous study performed on the Sha ecotype originated from 3400 m in a high valley of Tajikistan [16]. Indeed, when grown at 15 • C, Sha petioles showed a thicker cuticle than Col petioles. In the latter case, these hydrophobic compounds could be lignin and ferulic acid which could contribute to an increase of CW rigidity [18].
Although some trends are common to rosettes and floral stems upon growth at 22 • C or 15 • C, only two CW transcripts/CWPs of the oxido-reductases functional class were common to the short list of those differentially accumulated at a given growth temperature: AT5G15350 which is an early nodulin (AtEN22) homologous to blue copper binding proteins, and AT1G41830 which is a multicopper oxidase (AtSKS6) homologous to the A. thaliana SKU5 (SKEWED ROOT 5). To our knowledge, the role of the former protein has not been studied yet whereas SKU5 was shown to play a role in root directional growth [55] and more recently in the spaceflight response [56].
To summarize, the impact of the growth temperature on the levels of accumulation of CW polysaccharides, CW transcripts or CWPs was observed in both rosettes and floral stems for all plants. However, this response was mostly organ-specific.

The Response to Growth Temperature is Organ-Specific
The CW structure and the composition of rosettes and floral stems were very different in the two growth conditions. For example, at the 15 • C growth temperature, we have found more xyloglucans and more branched rhamnogalacturonans I in the rosettes of all plants whereas we have observed less xyloglucans in the floral stems of Col, Grip and Hos. Modifications in the composition of polysaccharides were already observed in Miscanthus giganteus ecotypes upon cold stress. In particular, the amount of (1,3)(1,4)-β-d-glucans (mixed-linked glucans) was increased [57]. In tobacco pollen tubes, cold was shown to induce changes in the pattern of deposition of homogalacturonans (both methylated and demethylated), callose and cellulose [58]. Such changes are assumed to modify the CW properties. CWPs possibly involved in these modifications could not be identified among the identified candidates maybe because the observed changes have occurred at earlier stages of development. Indeed, the plants have been continuously grown either at 22 • C or at 15 • C. We do not apply a stress as it is usually defined, i.e., a sudden change in the environmental conditions [16]. The observed changes could also occur at the level of the biosynthesis of the CW polysaccharides, which takes place in the Golgi apparatus and which could not be seen in this study.
Considering the temperature as the factor for discriminant analysis, the integrative analysis of rosettes highlighted an enrichment in three CWPs functional classes (proteins acting on CW polysaccharides, proteases and proteins related to lipid metabolism) among the candidate CW transcripts and CWPs (Supplementary Figure S5A, Table 2). Most of the proteins acting on CW polysaccharides were glycoside hydrolases (GHs). Three GH27 transcripts encoding α-galactosidases (AT3G26380, APSE; AT5G08370, AGAL1; AT5G08380, AGAL2) were up-regulated at 15 • C. An ortholog of AT5G08380 in Thlaspi arvense, another Brassicaceae, was shown to be up-regulated during cold acclimation [59,60]. Besides, APSE, and to a lesser extent AGAL1 and 2, have been assumed to be involved in the hydrolysis of the β-l-Arap residues of the glycan moiety of arabinogalactan proteins in muro [61]. Among proteases, one Ser carboxypeptidase (AtSCPL10) was accumulated at the 15 • C growth temperature and transcripts encoding AtSCPL29 and 49 were more abundant. Recently, a protein of the same family (AtSCPL41) was shown to play a role in membrane lipid metabolism [62]. Besides ASPG1 (ASPARTIC PROTEASE IN GUARD CELL 1) was shown to be involved in drought avoidance via ABA-dependent signaling [63]. Several proteins possibly involved in lipid metabolism, three non-specific lipid transfer proteins (LTPs) were found to be more abundant at 15 • C than at 22 • C. AtLTP3 (AT5G59320) was previously shown to be involved in the response to freezing stress [64] and was part of a tandem duplication with AtLTP4 (AT5G59310). The two corresponding proteins have been found in close proximity in the dendrogram (correlation > 0.9, Supplementary Figure S5A). Even if no study has shown the role of these LTPs in muro, their accumulation in rosettes at 15 • C could be related to the high number of hydrophobic compounds found in the cuticle epidermis of the petiole cross-sections mentioned above. However, the atltp3 mutant did not show any cuticular phenotype in normal growth conditions, but higher sensitivity to drought stress, whereas overexpressing plants were more tolerant to freezing and drought [64]. The levels of expression of AtLTP4 and AtLTP8 were increased in the atltp3 mutant, suggesting functional redundancy [64]. The analysis of double or triple mutants could contribute to a better understanding of the impact of low temperature growth on cuticle formation.
In floral stems, the majority of the candidate CW transcripts and CWPs identified at the 15 • C growth temperature belonged to the oxido-reductases and proteins possibly involved in lipid metabolism functional classes. Two oxido-reductases families are well-represented: blue copper binding proteins (AT5G15350, AtEN22; AT4G31840, AtEN13; AT4G12880, AtEN20) and multicopper oxidases (AT5G21100; AT3G13990, AtSKS11; AT1G41830, AtSKS6). As mentioned above, AtEN22 and AtSKS6 are the two candidates common to rosettes and floral stems. Regarding proteins possibly involved in lipid metabolism, AtLTPg4 (AT1G27950) and AT5G62210 encoding a protein of yet unknown function, were highlighted. Interestingly, the AT5G62210 transcript and the encoded protein exhibiting a PLAT/LH2 domain both accumulated at 15 • C. This protein showed 41% homology with that encoded by the pathogen-induced CAPIP2 gene which was described as a response protein against abiotic stresses [65,66]. A previous transcriptomics analysis of epidermal peels of A. thaliana stems has previously highlighted genes encoding LTPs (AtLTP2, 5 and 7 and AtLTPg4, 6, 7 and 21) and thus allowed identifying new proteins possibly playing roles in the biosynthesis of waxes and cutin [67]. The importance of the cuticle of the floral stems is critical as it constitutes a vital hydrophobic barrier, providing mechanical strength and viscoelastic properties as well as protecting the plant from the environmental stresses [25].
Altogether, the growth temperature had a different effect on rosettes and floral stems at the molecular level. Distinct CW transcripts and CWPs of interest could be highlighted. However, in both cases, the proteins possibly related to lipid metabolism were particularly well-represented, suggesting a critical role of the cuticle in the acclimation to low growth temperatures.

The Hern and Hosp Populations Exhibit Particularities in Their Responses to Growth Temperature
Hosp showed different responses at the level of rosettes at the 15 • C growth temperature. Six CW transcripts and CWPs allowing the discrimination of Hosp had unknown function, like AT5G38980 for which both transcripts and proteins were accumulated. AtLTPg6 (AT1G55260) was more abundant and AtLTP6 (AT3G08770) transcripts were up-regulated. Interestingly, the transcripts levels of these two candidates were also up-regulated in the Sha altitudinal ecotype grown at 15 • C [16] and LTPgs have been reported to be implicated in the transport of lipids through the hydrophilic CWs to build up the cuticle [68]. These results are consistent with the higher level of hydrophobic compounds observed in the epidermis cuticle of Hosp petioles compared to Col.
The behavior of Hern and Hosp was different from that of Grip, Roch and Col with regard to floral stems. However, we could not identify any common candidates when the integrative analysis was focused on the temperature parameter. The specificity of Hern floral stems could be noticed by the clustering of all the phenomics and metabolomics variables, except the ratio of rhamnogalacturonan I branching (Supplementary Figure S7). Co-expressed with these variables, six proteins were more abundant in Hern, three of them belonged to the proteins with interaction domains functional class and two to the proteins acting on CW polysaccharides functional class. Among the proteins with interaction domains, the two lectin genes (AT3G15356, AT3G16530) were known to be sensitive to abiotic stress and up-regulated upon hormonal treatments [69,70]. AT3G04720 (PR4) was positively regulated by jasmonate and ethylene [71]. AT4G16260 encoding a β-1,3-glucosidase belonging to the proteins acting on CW polysaccharides functional class was also regulated by ethylene [72]. These results highlight a possible link between the Hern floral stem phenotype and the jasmonate and ethylene hormones. Hormones were not analyzed in this study but could play a major role in the phenotypic plasticity of floral stems. Interestingly, transcripts encoding proteins acting on CW polysaccharides and proteases which were well-represented in this set of candidates were down-regulated in Hern. Proteases play roles in the generation of signals involved in plant development, in the maturation of proteins and in protein turnover [73,74]. These processes which are still poorly understood could be important in the phenotypic differentiation of Hern floral stems at 15 • C. Hern appears as an interesting population to study the contribution of such proteins to CW plasticity in floral stems. Besides, Hern floral stems exhibited lower levels of six transcripts among which two encoded proteases: AT1G28110 (AtSCPL25) and AT4G00230 (AtSBT4.14).
Hosp also showed particularities at the proteomics and transcriptomics level. With less contrasted altered phenotypic traits than Hern, Hosp exhibited a lower number of CW transcripts and CWPs with modified levels of accumulation. Nevertheless, some candidates CW transcripts or CWPs appeared to specifically accumulate in Hosp such as the Ser (AT2G12480, AtSCPL43) and the Asp (AT3G02740) proteases. In addition, the AT1G78820 transcripts encoding a lectin were accumulated twice as much as in Col. Moreover, the transcripts of two candidates already identified in the stem epidermis [66] were specific to Hosp and were found to be more abundant at 15 • C: AtLTPg4 (AT1G27950) and PME41 (PECTIN METHYLESTERASE 41, AT4G02330). Mutants impaired in AtLTPg4 and AtLTPg21 exhibit significant alterations in the cuticular wax composition [75]. Also, the pme41 mutants were shown to be sensitive to chilling stress and this phenotype could be related to the brassinosteroids responses [75]. This PME could play a role in plant chilling tolerance by modifying the mechanical properties of the CW and specifically the esterification rate of homogalacturonans [76]. Beyond, the possible regulation of the PMEs by hormones such as brassinosteroids could provide further insight into their role in the acclimation of plants to low temperature.

Conclusions
The aim of this study was to identify new players of the CW plasticity upon acclimation to sub-optimal growth temperatures by an integrative study using different A. thaliana populations collected in the Pyrenees Mountains. A contrasted CW plasticity has been observed in rosettes and in stems such as the increase of xyloglucan content in the former and their decrease in the latter at the 15 • C growth temperature. The differential CW plasticity responses between organs could also be illustrated by the differences observed between the CW transcriptomes and proteomes depending on both the population and the growth temperature. The global integration of all four blocks of data (phenomics, metabolomics, proteomics and transcriptomics) supervised by a categorical outcome (temperature or ecotype) has allowed finding new candidates for functional studies, such as LTPs which could play roles in the biogenesis of the cuticle. Such candidates would not have been identified using a mono-block analysis.
This work has revealed that low temperature acclimation affects the expression of many genes at the transcriptional or post-transcriptional levels, many of them having yet unknown biological functions. Importantly, this powerful integrative analysis has allowed identifying interesting candidates for exploring the adaptation of natural populations to low temperature acclimation. Further studies using mutants with altered expression of these candidates are now necessary to elucidate their biological significance. In this way, hypothesis-based approaches exploring understudied traits can provide new issues in A. thaliana research. Furthermore, this knowledge could be transferred to plants of agronomical interest.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-4409/9/10/2249/s1, Figure S1: Phenotype of floral stems in contrasted growth conditions, Figure S2: Phenotype of rosettes in contrasted growth conditions, Figure S3: Reconstruction of the main rosette cell wall polysaccharides from monosaccharide analysis, Figure S4: Reconstruction of the main cell wall polysaccharides of floral stems from monosaccharide analysis, Figure S5. Individual plots project the rosette samples for the four blocks used for the sparse MB-PLS-DA analyses, Figure S6: Graphical representation of sparse MB-PLS-DA analyses that discriminate the floral stems of Col, Roch, Grip, Hern and Hosp according to the growth temperature, Figure S7: Graphical representation of sparse MB-PLS-DA analyses that discriminate the floral stem samples according to Col, Roch, Grip, Hern and Hosp, Table S1: Statistical analysis of proteomics data, Table S2: Statistical analysis of rosettes and floral stems transcriptomics data, Table S3: Statistical results of the sparse MB-PLS-DA of rosettes, Table S4