Integration of QTL, Transcriptome and Polymorphism Studies Reveals Candidate Genes for Water Stress Response in Tomato

Water deficit (WD) leads to significant phenotypic changes in crops resulting from complex stress regulation mechanisms involving responses at the physiological, biochemical and molecular levels. Tomato growth and fruit quality have been shown to be significantly affected by WD stress. Understanding the molecular mechanism underlying response to WD is crucial to develop tomato cultivars with relatively high performance under low watering conditions. Transcriptome response to WD was investigated through the RNA sequencing of fruit and leaves in eight accessions grown under two irrigation conditions, in order to get insight into the complex genetic regulation of WD response in tomato. Significant differences in genotype WD response were first observed at the phenotypic level for fruit composition and plant development traits. At the transcriptome level, a total of 14,065 differentially expressed genes (DEGs) in response to WD were detected, among which 7393 (53%) and 11,059 (79%) were genotype- and organ-specific, respectively. Water deficit induced transcriptome variations much stronger in leaves than in fruit. A significant effect of the genetic background on expression variation was observed compared to the WD effect, along with the presence of a set of genes showing a significant genotype × watering regime interaction. Integrating the DEGs with previously identified WD response quantitative trait loci (QTLs) mapped in a multi-parental population derived from the crossing of the eight genotypes narrowed the candidate gene lists to within the confidence intervals surrounding the QTLs. The results present valuable resources for further study to decipher the genetic determinants of tomato response to WD.


Introduction
Drought is among the most common abiotic stress factor affecting plant growth and crop yield, and more frequent episodes of drought are expected to arise with climate change [1]. Extensive research has been dedicated to understanding the mechanisms driving plant adaptation to drought [2]. Water deficit (WD) stress induced by drought can be defined as a period of plant exposure to dry soil, subsequently resulting in reduced growth and yield [3]. A global understanding of the complex interplay between genetic and environmental factors in crop adaptation to WD is therefore a key aim for breeding purposes.
selected from a collection of 360 accessions [31] to uncover the maximum genetic diversity of this panel. The experimental design and plant growth conditions are described in detail in Albert et al. (2018) [6]. Briefly, the eight genotypes were grown in spring-summer 2015 in a greenhouse in Avignon (south of France) under two watering conditions: control and water deficit (WD). The control condition consisted of full irrigation treatment according to evapotranspiration (ETP), while the WD condition was set progressively after the flowering of the second truss of Cervil (earliest genotype): the water supply was reduced by 25% compared to the control for one week, then decreased by 40% until the end of the experiment to apply a mild water deficit. Throughout the experiment, the relative humidity of the peat substrate was controlled with a GRODAN ® moisture probe and monitored in drought pots at around 30%. Three plants per watering regime and accession were grown in the greenhouse, in four-liter (L) plastic pots filled with peat (Klasmann 165) and watered with nutritive solution (2, 4, 6 mmol L −1 of N, P and K, respectively). Temperature, light and air humidity were scored every hour. On average, the temperature was 25 • C during the day and 18 • C during the night, the daily light integral ranged from 5 to 11 MJ m −2 day −1 and the average humidity was 47% during the daytime and 70% at night. Phenotypic measurements were carried out for different traits related to plant development and fruit composition. Plant height (distance from the soil to the 4th truss) and stem diameter (diameter just below the 4th truss) were measured, with one measure recorded for each plant. Flowering time (flw) was measured as the number of days after sowing date. For each genotype, ten fruits were harvested and weighed to measure the average fruit weight (FW) from the 3rd to the 6th truss. Fruits were then pooled in three biological replicates to measure fruit composition traits, from which dry matter weight (DMW) was evaluated by drying fruit pericarps in an oven at 60 • C for four days. Half of each fruit pool was blended to measure pH and soluble solid content (SSC). Fresh pericarps were also sampled from each pool and frozen in liquid nitrogen before being ground into powder for sugar (glucose and fructose, g 100 g 1 FM) and total vitamin C (VitC, mg 100 g 1 FM) measurements, following the protocols used in Albert et al. [6].

Statistical Analyses of Phenotypic Data
To test for a WD effect at the phenotypic level, a two-way analysis of variance (ANOVA) was performed for each trait separately. The level and significance of the G × W interaction was assessed with the following model: y ij = µ + G i + W j + GxW ij + ε ij . In this model, y ij represents the phenotype of genotype i (G i ) in watering condition j (W j ) and GxW ij and ε ij are the genotype × watering regime interaction and residual errors, respectively. For each trait, phenotypic plasticity (PP) was calculated according to the following formula: PP i = (y WD − y control )/y control , where PP i represents the plasticity (WD response) for each genotype (i) and y WD (or y control ) the average phenotype in the WD (or control) condition. The average response of the eight genotypes was further evaluated by computing the relative stress impact: RSI = Y WD − Y control /Y control , where Y represents the average value across all genotypes for a trait, given a condition.

RNA Extraction
For each genotype, total RNA was collected from growing leaves and fruit pericarps (at least five fruits) at the cell expansion stage. Given the differences in their phenological stage, the cherry (SLC) and large fruit accessions (SLL) were sampled for fruit pericarps at 14 and 21 days after anthesis (DAA), respectively. The samples were immediately frozen after collection, then pooled per genotype, organ and condition, with two to three biological replicates. Messenger RNA (mRNA) was extracted using the Spectrum Plant Total RNA Kit and assessed on a Nanodrop 1000. A total of 72 paired-end strand-specific libraries were generated from 1 µg of the total RNA and sequencing was performed on a Hiseq 3000 with the GenoTool platform (INRA Toulouse). The biological replicates in both conditions were disposed in the same lane for each genotype. Detailed information about the RNA extraction protocol and read sequence processing are described in Albert et al. (2018) [6].

Differential Gene Expression Analysis
Differential expression (DE) analysis was performed using a negative binomial generalized linear model with the Bioconductor and R package DESeq2 1.14.1 [32]. The impact of WD on transcriptome variation was evaluated in fruit and leaf samples separately. The analysis was restricted to genes with at least 20 read counts across samples, encompassing 23,552 and 22,864 genes for leaf and fruit samples, respectively. The impact of genotype and watering condition on transcriptome variation was first graphically evaluated through principal component analysis (PCA) on the normalized gene expression. Normalized counts were transformed with the variance stabilizing transformation (VST) before performing the PCA analysis. Thereafter, genes showing differential expression according to the watering condition were identified for each genotype and for fruit and leaf organs separately. A False Discovery Rate (FDR) threshold of 5% [33] was applied to call significantly differentially expressed genes using the Wald test.

Two-Way ANOVA of Transcript-Level Variation
Analysis of variance was performed on the expression level using the expressed genes (same gene set as the DE analysis) to measure, for each single gene, the relative contribution of genotype (G), watering condition (W) and G × W factors to gene expression variation. The ANOVA model used was similar to the one described for phenotypic trait analysis, with the only difference being that y ij this time represents the transcript abundance (VST-transformed values of normalized counts) of genotype i (G i ) in watering condition j (W j ) for each single gene included in the analysis. A FDR threshold of 5% was applied to call significant differences [33]. This analysis was conducted for each organ separately and the proportion of the sum of squares attributed to genotype (G), watering condition (W) and G × W factors were retrieved and used to estimate the relative contribution of each factor.

Gene Ontology Enrichment Analysis
Gene ontology (GO) enrichment analysis was performed using the R-Bioconductor package goseq (version 1.36.0) [34]. Enriched GO terms were investigated for the sets of DEGs in fruit and leaves separately. For each organ, the DEGs were separated into three different categories (up, down and up-down) depending on the pattern of the expression regulation. The up (or down) category involves all DEGs that were consistently up-(or down-) regulated across genotypes or in at least one genotype. The up-down category involved DEGs showing different regulation patterns (i.e., upregulated in one genotype and downregulated in another and vice-versa). Only annotated genes were considered for the GO analysis. Significant GO terms for biological processes (BPs) and molecular function (MF) were selected after multiple testing correction by setting an FDR threshold at a 5% cutoff with the Benjamin and Hochberg method [33]. The SL2.50 version of the reference genome "Heinz" was used and correction for length bias was carried out with the nullp function before GO enrichment testing. The gene space was composed of 27,014 genes, among which 18,837 genes presented at least one piece of GO information.

Co-Localization of the DEGs and Tomato WD QTLs
Diouf et al. (2018) conducted a QTL mapping analysis using a MAGIC population generated from the intercross of the eight genotypes described in the present study under control conditions and WD [18]. The WD stress condition consisted of a 40% reduction in water irrigation. The authors identified 12 WD-responsive QTLs, among which nine were specific to the WD condition and three QTLs were detected using phenotypic plasticity in response to WD for different traits related to plant and fruit growth. The respective genomic regions of these QTLs were compared to the locations of DEGs identified across the eight parental lines. The comparison was made possible using the physical positions of the tomato reference genome (version SL2.50). Besides, the correlations between the expression level and phenotypic allelic effect for each QTL were assessed.

Phenotypic Response to WD
A total of 10 phenotypic traits related to plant development and fruit composition were measured in eight tomato genotypes grown under control and WD conditions. The genotypes significantly differed for all the phenotypic traits measured (p-value < 0.05) and WD treatment significantly affected every trait except flowering time, glucose, pH and SSC. Considering the average response across genotypes, fruit dry matter weight (DMW +21.4%) and fruit fresh weight (FW −20.7%) were the most affected traits (Figure 1). Variability was observed among the genotypes in their response to WD at the phenotype level (Table S1). Ferum, for instance, showed the highest susceptibility to WD compared to the other genotypes for most of the traits ( Figure S1). The G × W interaction was only significant for four traits (height, diameter, FW and glucose).

Phenotypic Response to WD
A total of 10 phenotypic traits related to plant development and fruit composition were measured in eight tomato genotypes grown under control and WD conditions. The genotypes significantly differed for all the phenotypic traits measured (p-value < 0.05) and WD treatment significantly affected every trait except flowering time, glucose, pH and SSC. Considering the average response across genotypes, fruit dry matter weight (DMW +21.4%) and fruit fresh weight (FW −20.7%) were the most affected traits ( Figure 1). Variability was observed among the genotypes in their response to WD at the phenotype level (Table S1). Ferum, for instance, showed the highest susceptibility to WD compared to the other genotypes for most of the traits ( Figure S1). The GxW interaction was only significant for four traits (height, diameter, FW and glucose).

Transcriptome Variability Across the Eight Genotypes
RNA sequencing was carried out using growing leaf and fruit pericarp samples at the cell expansion stage, collected from the eight genotypes under control and WD conditions. Overall, the RNA sequencing process yielded a total of 23,552 genes (67.8% of tomato CDS) and 22,864 genes (65.8% of tomato CDS) with expression levels above background noise, in leaves and fruit, respectively. Principal component analysis (PCA) on the transformed normalized read counts showed a clear clustering of the samples according to genotypes and conditions ( Figure 2). The first two axes of the PCA explained about 52 and 56% of the variation of the gene expression level in fruit and leaf samples, respectively. For both organs, variability in the transcript levels according to genotypes was captured by both the PC1 and PC2 axes. In addition, for leaf samples, the PCA plot separated the conditions following the PC2 axis, Cervil and Levovil appeared as the most discriminant genotypes and, for fruit organs, transcriptomic variation was highly specific for Cervil, as observed in the PCA plot ( Figure 2).

Transcriptome Variability Across the Eight Genotypes
RNA sequencing was carried out using growing leaf and fruit pericarp samples at the cell expansion stage, collected from the eight genotypes under control and WD conditions. Overall, the RNA sequencing process yielded a total of 23,552 genes (67.8% of tomato CDS) and 22,864 genes (65.8% of tomato CDS) with expression levels above background noise, in leaves and fruit, respectively. Principal component analysis (PCA) on the transformed normalized read counts showed a clear clustering of the samples according to genotypes and conditions ( Figure 2). The first two axes of the PCA explained about 52 and 56% of the variation of the gene expression level in fruit and leaf samples, respectively. For both organs, variability in the transcript levels according to genotypes was captured by both the PC1 and PC2 axes. In addition, for leaf samples, the PCA plot separated the conditions following the PC2 axis, Cervil and Levovil appeared as the most discriminant genotypes and, for fruit organs, transcriptomic variation was highly specific for Cervil, as observed in the PCA plot ( Figure 2).

DEGs Under WD Conditions
Considering the eight genotypes together, a total of 4132 and 12,938 DEGs between the two watering conditions were identified in the DE analysis in fruit pericarps and growing leaves, respectively (Table S2 and Table S3). The number of DEGs was variable among genotypes. In fruit, the number of DEGs varied from zero (Criollo) to 2978 (Levovil) and, on average, the SLL accessions showed a higher number of DEGs than the SLC genotypes ( Figure 3A). In particular, Levovil and Ferum, two SLL genotypes, presented the highest number of DEGs between the conditions, although the proportion of up/downregulated genes differed. For leaf samples, the number of DEGs was higher for the different genotypes than in fruit. The total number of DEGs in leaves varied from 2240 (Ferum) to 6177 (LA1420) and the proportion of up/downregulated genes was almost balanced except for Cervil and LA1420. Most of the DEGs identified were organ-specific, highlighting the organ-dependent regulation of gene expression under WD ( Figure 3B). Indeed, depending on the genotypes, 45 to 82% of the DEGs in fruit were not differentially expressed in leaves.

DEGs Under WD Conditions
Considering the eight genotypes together, a total of 4132 and 12,938 DEGs between the two watering conditions were identified in the DE analysis in fruit pericarps and growing leaves, respectively (Tables S2 and S3). The number of DEGs was variable among genotypes. In fruit, the number of DEGs varied from zero (Criollo) to 2978 (Levovil) and, on average, the SLL accessions showed a higher number of DEGs than the SLC genotypes ( Figure 3A). In particular, Levovil and Ferum, two SLL genotypes, presented the highest number of DEGs between the conditions, although the proportion of up/downregulated genes differed. For leaf samples, the number of DEGs was higher for the different genotypes than in fruit. The total number of DEGs in leaves varied from 2240 (Ferum) to 6177 (LA1420) and the proportion of up/downregulated genes was almost balanced except for Cervil and LA1420. Most of the DEGs identified were organ-specific, highlighting the organ-dependent regulation of gene expression under WD ( Figure 3B). Indeed, depending on the genotypes, 45 to 82% of the DEGs in fruit were not differentially expressed in leaves.  A small set of genes were been identified as DEGs in both organs ( Figure 3B), representing 0-13% of the total DEGs according to the genotype. For most of the genotypes, the pattern of gene expression regulation was different between leaves and fruit. For instance, for genotypes with more than 100 consistent DEGs between fruit and leaves-notably Cervil, Levovil and Ferum-the proportion of the genes upregulated in one organ (leaf or fruit) and downregulated in the other was non-negligible, representing 47, 44 and 79%, respectively ( Figure S2). A small set of genes were been identified as DEGs in both organs ( Figure 3B), representing 0-13% of the total DEGs according to the genotype. For most of the genotypes, the pattern of gene expression regulation was different between leaves and fruit. For instance, for genotypes with more than 100 consistent DEGs between fruit and leaves-notably Cervil, Levovil and Ferum-the proportion of the genes upregulated in one organ (leaf or fruit) and downregulated in the other was non-negligible, representing 47, 44 and 79%, respectively ( Figure S2).

Genotype and WD Impact on the Transcript Level
A two-way analysis of variance (ANOVA) was performed using the VST-transformed values of the normalized transcript level to assess the effect of genotype, watering condition and G × W interaction on the regulation of gene expression. This analysis revealed a total of 16,392 and 19,450 genes that were affected by at least one of the abovementioned factors in fruit and leaves, respectively. A high proportion of the genes tested showed a significant genotype effect, highlighting an important effect of the genetic background on gene expression regulation among the eight lines ( Figure S3). A much smaller number of genes were specifically affected by WD or the G × W interaction.

Gene Ontology Enrichment Analysis
GO enrichment analysis was conducted on a set of 3794 (92% of the total fruit DEGs) and 11,804 (91% of the total leaf DEGs) DEGs in fruit and leaves, respectively, yielding a total of 24 significantly enriched GO terms (Table 1). Only genes with a known annotated function were selected for the GO analysis. GO terms associated with "cell redox homeostasis", "metabolic process", "microtubule-based movement" and "protein phosphorylation" were significantly over-represented among the DEGs in leaves. With reference to the molecular function, three GO terms (GO:0016168, GO:0003735 and GO:0008152) related to "chlorophyll binding", "structural constituent of ribosome" and "metabolic process", respectively, were enriched considering DEGs in fruit and leaves separately.

DEG Co-Location with Previously Identified WD-Responsive QTLs
The combination of gene expression and QTL information was used to identify candidate genes under the WD-responsive QTLs detected in the MAGIC population. To illustrate the approach, we focused on three and nine QTLs previously identified in the MAGIC population using the plasticity response under WD or QTLs identified specifically under the WD condition [18]. The selected QTLs were mapped to seven out of the 12 tomato chromosomes ( Figure 4A). The confidence interval regions of these QTLs encompassed 11 to 51 cM and included in their intervals 149 to 988 genes. However, the number of DEGs within the QTL regions varied from 20 to 102 genes, reducing the set of potential candidates by 46-93% according to the QTL ( Figure 4B). The number of DEGs per Mbp was assessed for each of the seven chromosomes carrying the WD-responsive QTLs and these genomic regions covered by the WD-responsive QTLs were enriched with DEGs for some chromosomes ( Figure 4C). For reasons of consistency, QTLs identified on leaf and fruit traits were mapped specifically to the DEGs detected in leaf and fruit organs, respectively.
Genes 2020, 11, x FOR PEER REVIEW 9 of 15 genomic regions covered by the WD-responsive QTLs were enriched with DEGs for some chromosomes ( Figure 4C). For reasons of consistency, QTLs identified on leaf and fruit traits were mapped specifically to the DEGs detected in leaf and fruit organs, respectively.  QTL detection using the parental haplotype probabilities in the MAGIC population allows for the estimation of the allelic effect for each parental line at every QTL position. Correlation analysis was further investigated between the allelic effect and expression level across the eight genotypes. The expression of a total of 46 genes showed significant correlation with the QTL allelic effect, reinforcing their potential implication in regulating fruit phenotype variation under WD (Table S4). The whole process of candidate gene selection is depicted in Figure 5 for the QTL RIP9.1 as an example. Genes within the RIP9.1 QTL interval. Black dots represent non-DEGs, and red dots, DEGs in the present study. Triangles represent the DEGs for which the delta expression level (expression level in WD -expression level in control) was significantly correlated to the allelic effect of the QTL for the eight genotypes. (Bottom) Correlation between the estimated allelic effect at the QTL (x-axis) and the delta log2 expression levels (y-axis) for four candidate genes with their functional annotation.

Discussion
Plant response to drought is a complex mechanism which ultimately results in phenotypic changes that can alter agronomic performance in crops. Such responses strongly depend on the genetic background, leading to the necessity of screening different genotypes/accessions for drought tolerance. Genotype specificity has been depicted in tomato response to WD for different agronomic traits [35][36][37]. The present study highlighted a strong genotype-dependent variation under WD at both the phenotype and transcriptome levels. Genes within the RIP9.1 QTL interval. Black dots represent non-DEGs, and red dots, DEGs in the present study. Triangles represent the DEGs for which the delta expression level (expression level in WD-expression level in control) was significantly correlated to the allelic effect of the QTL for the eight genotypes. (Bottom) Correlation between the estimated allelic effect at the QTL (x-axis) and the delta log2 expression levels (y-axis) for four candidate genes with their functional annotation.

Discussion
Plant response to drought is a complex mechanism which ultimately results in phenotypic changes that can alter agronomic performance in crops. Such responses strongly depend on the genetic background, leading to the necessity of screening different genotypes/accessions for drought tolerance. Genotype specificity has been depicted in tomato response to WD for different agronomic traits [35][36][37]. The present study highlighted a strong genotype-dependent variation under WD at both the phenotype and transcriptome levels.
The WD treatment induced a significant variability in the transcriptome across the eight genotypes. A total of 14,065 DEGs under WD were detected, among which 7393 (53%) were genotype specific. Cervil and Levovil presented the most divergent pattern of gene expression variation at the fruit level (Figure 2), which was consistent with the phenotypic variation since these two genotypes also presented the smallest and largest fruit weight, representing on average 6.1 g and 119 g, respectively. Furthermore, these two genotypes were also identified as highly divergent at the whole genome level when comparing their polymorphism sequences against the reference genome [38]. Moreover, a bi-parental population generated from these two lines yielded the significant discovery of genetic loci controlling trait variation in tomato [6,39,40].
Organ and tissue sampling can significantly alter the gene expression profile in plants [41,42]. Besides, the organ-specific transcript level might be exacerbated by the presence of stress factors, especially at a specific growth stage. The WD treatment in our study was applied from the first inflorescence appearance until the end of the growing season. Fruit sampling was elaborated to limit differences that might arise from growth stages between the SLL (sampled at 21 DAA) and SLC (sampled at 14 DAA) accessions. The low number of fruit DEGs in SLC accessions could have been linked to the sampling strategy we adopted. However, this hypothesis is not supported by the low number of fruit DEGs observed for LA0147 and Stupicke, which are both SLL accessions.
The process of plant adaptation to drought usually starts with cellular sensing and signaling, which activates downstream drought-responsive genes [4]. The activation of the signaling pathway under WD at the early vegetative growth stage is hence expected to lead to better adaptation. This may explain the high number of DEGs in growing leaves.  have shown, in the same eight lines, a higher WD impact on leaves than at the fruit level and outlined a prevalence of osmotic adjustment and photosynthetic adaptation in tomato response to WD [29]. The source-sink relationship is highly altered under the WD stress condition in tomato [43], which could eventually reflect different transcriptome responses.
Comparative analysis of the whole transcriptome variation under WD has been conducted in several species; however, most of the time, only two genotypes are included (usually tolerant vs. sensitive genotypes). These classical designs are very powerful to detect DEGs involved in specific stress responses. Yet, the few studies which assessed the impact of a specific stress and the genetic background simultaneously revealed the significant effect of the genotypes with a frequent GxE interaction at the transcriptome level [6,41]. The eight MAGIC parental lines were selected from a panel of 360 tomato accessions to represent the diversity observed within the cultivated tomato. Our results provide evidence of a strong genotype effect at the transcriptome level, independent of the growing condition. Expression QTL analysis in the MAGIC population should then yield significant results that could help to get more insight into the molecular mechanisms shaping tomato variation.
Integrating QTL information and gene expression variation is a common strategy for candidate gene screening [44]. This strategy yielded significant results in MAGIC populations in maize [45] and cotton [46], where the parental haplotypes were exploited to drastically narrow the potential candidate genes. Here, we propose a strategy combining DEGs and WD-responsive QTLs to identify candidate genes affecting WD stress response in tomato. The whole approach depicted in Figure 5 highlighted 46 candidate genes (CGs), which showed expression variability in the eight parental MAGIC lines under the WD condition for 12 WD-responsive QTLs. Four interesting CGs were indeed identified in the region of the RIP9.1 QTL, among which was Solyc09g010630, a gene coding for a heat-shock protein (Hsp). In the literature, Hsp has been described as playing a role in drought tolerance [47][48][49] and fruit ripening in tomato [50]. Regarding the allelic effect of the RIP9.1 QTL, WD stress accelerated ripening in all parents except Cervil and Criollo ( Figure 5). A total of 82 polymorphisms (75 SNPs and seven INDELs) were identified for Solyc09g010630, among which 54 discriminated the most divergent lines (Cervil and Ferum) regarding the estimated allelic effect of the QTL. The allelic variants of these polymorphisms are in accordance with the pattern of the gene expression and allelic effects across the eight parental lines. Further studies are, however, necessary to test the potential implication of these alleles in tomato fruit response to WD.
Similarly, other interesting candidates were identified. For instance, RIP10.1 is another QTL detected for the plasticity of fruit ripening under the WD condition, which carried a total of 301 genes within its confidence interval region, among which only 37 genes (12%) showed significant differential expression levels in fruit pericarps. Some of these genes presented a significant correlation between their expression level in response to WD and the allelic effect of the QTL across the eight parental lines, encompassing Solyc10g006130, annotated as "ethylene responsive transcription factor (ERF) 3a" and Solyc10g006650, encoding a "flavoprotein wrbA". In tomato, ERFs have been described as being involved in the fruit maturation process and also as drought inducible transcription factors [51,52]. Thus, Solyc10g006130 represents an interesting candidate for studying the interaction between drought and fruit maturation in tomato. Solyc10g006650 also constitutes a good candidate, as the role of flavoprotein in tolerance to osmotic stress has been depicted in Arabidopsis [53].

Conclusions
The integration of several types of omics data (phenotypic, transcriptomic, metabolomic) may help understanding plant adaptation to drought and to apprehend the complex mechanisms involved. The investigation of GxE at the transcriptome level has the potential to target interesting candidate genes for further functional analyses. Collecting and gathering omics data from different organs, tissues, genotypes and conditions is the key step for omics breeding in the coming years in order to develop climate resilient crops. The results presented here are therefore valuable resources for the tomato community for further studies intended to decipher drought response mechanisms.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-4425/11/8/900/s1: Figure S1: Phenotypic variation under WD for each genotype. Figure S2: Expression pattern of the consistent DEGs in fruit and leaves. Figure S3: Venn diagram representing the number of genes whose expression level was significantly affected by the genotype, the watering condition and the genotype × watering condition interaction factors in the ANOVA analysis, in fruits and leaves. Table S1: Results of the two-way interactive ANOVA analysis on the phenotypic traits. Table S2 and Table S3: Differentially expressed genes under WD for each of the eight genotypes in fruit and leaves. Table S4: Candidate genes under the WD-responsive QTLs.
Data Accessibility: The sequences are available at the NCBI SRA repository under project number PRJNA644865. Expression data for both conditions are also available in the TomExpress database http://tomexpress.toulouse.inra.fr. Funding: This work was performed in collaboration with the GeT core facility, Toulouse, France (http://get.genotoul. fr), and was supported by France Génomique National infrastructure, funded as part of the "Investissement d'avenir" program managed by Agence Nationale pour la Recherche (contract ANR-10-INBS-09)." The work was also supported by the ANR project ADAPTOM ANR-13-BSV7-012. I.D. was supported by the WAAPP (West Africa Agricultural Productivity Project) fellowship and hosted as Ph.D. student in the INRA UR1052, GAFL.