https://orcid.org/0000-0002-2655-5085
Abstract
Water availability influences all aspects of plant growth and development; however, most studies of plant responses to drought have focused on vegetative organs, notably roots and leaves. Far less is known about the molecular bases of drought acclimation responses in fruits, which are complex organs with distinct tissue types. To obtain a more comprehensive picture of the molecular mechanisms governing fruit development under drought, we profiled the transcriptomes of a spectrum of fruit tissues from tomato (Solanum lycopersicum), spanning early growth through ripening and collected from plants grown under varying intensities of water stress. In addition, we compared transcriptional changes in fruit with those in leaves to highlight different and conserved transcriptome signatures in vegetative and reproductive organs. We observed extensive and diverse genetic reprogramming in different fruit tissues and leaves, each associated with a unique response to drought acclimation. These included major transcriptional shifts in the placenta of growing fruit and in the seeds of ripe fruit related to cell growth and epigenetic regulation, respectively. Changes in metabolic and hormonal pathways, such as those related to starch, carotenoids, jasmonic acid, and ethylene metabolism, were associated with distinct fruit tissues and developmental stages. Gene coexpression network analysis provided further insights into the tissue-specific regulation of distinct responses to water stress. Our data highlight the spatiotemporal specificity of drought responses in tomato fruit and indicate known and unrevealed molecular regulatory mechanisms involved in drought acclimation, during both vegetative and reproductive stages of development.
Introduction
Environmental conditions have fundamental impacts on plant development and crop productivity, which will be amplified by climate change (Razzaq et al., 2021). Increased drought duration and an overall decrease in water availability are major contributors to plant abiotic stress under changing environmental conditions (Pokhrel et al., 2021). Thus, there is a pressing need to better understand the molecular mechanisms employed by plants to adapt to water deficit conditions and in all tissues and organs. Over the last decade, many of the physiological and molecular responses of plant vegetative organs to water stress have been characterized in model and crop species (Kaur and Asthir, 2017; Tardieu et al., 2018; Bechtold and Field, 2018; Sun et al., 2020). Generally, plant responses to water stress involve reduced shoot growth and photosynthesis, abscisic acid (ABA) accumulation, stomatal closure, the production of reactive oxygen species (ROS), and changes in metabolites involved in osmotic adjustment and antioxidant defense systems (Fàbregas and Fernie 2019; Sun et al., 2020). However, strategies employed by plants to respond to immediate acute water stress differ from those applied for moderate and prolonged water stress (Harb et al., 2010; Skirycz et al., 2011), so response mechanisms operating at different water deficit severities and time scales must be considered for a comprehensive understanding of performance under water stress (Rymaszewski et al., 2017; Chen et al., 2021). Most of the information related to plant responses to water deficit has been derived from studies of short-term responses to acute stress, which are relatively easy to investigate in controlled lab environments. In contrast, studies of plant acclimation to long-term water stress are labor intensive and require specialized equipment for frequent monitoring and adjusting soil water content (Granier et al., 2006; Skirycz et al., 2011). Accordingly, the physiological and molecular mechanisms involved in plant responses to prolonged drought have been less investigated and have primarily used the model system Arabidopsis (Arabidopsis thaliana). It has been shown that acclimation to long-term water stress includes transcriptional, metabolic, and physiological adjustments leading to an adjusted cellular homeostasis (Harb et al., 2010; Skirycz and Inzé, 2010), but it is not clear whether such information is translatable from model to crop species (Luo et al., 2019; Krukowski et al., 2020).
An example of an important crop for stress studies is cultivated tomato (Solanum lycopersicum), which is sensitive to drought, as manifested by perturbed growth and a reduction of both fruit yield and quality (Diouf et al., 2018). Several studies have described the physiological responses of tomato to water deficit, primarily by analyzing specific morphological, biochemical, and metabolic changes, or changes in expression of individual genes or gene families (Muñoz-Espinoza et al., 2015; Alaguero-Cordovilla et al., 2018; Galdon-Armero et al., 2018; Hao et al., 2019; Gao et al., 2020). There have been few comprehensive analyses of global molecular responses of tomato plants to water stress. Transcriptomic studies have shown that acute water stress causes downregulation of genes related to photosynthesis and cell division, as well as upregulation of ABA responses and genes in pathways that decrease energy dissipation and reduce oxidative damage (Gong et al., 2010; Iovieno et al., 2016; Liu et al., 2017). Investigations into the long-term responses of tomato to moderate water deficit have also revealed substantial transcriptome variations (Albert et al., 2016, 2018; Diouf et al., 2020). However, in tomato, as in other plant species, studies of water stress responses have focused almost exclusively on vegetative tissues (Ma et al., 2020).
Reproductive development, especially during its early stages, is particularly sensitive to water stress and this can be a major factor limiting crop yield (Barnabas et al., 2008; Su et al., 2013; Turc and Tardieu 2018). A number of studies have reported the effects of drought stress on growth, metabolism and quality of fleshy fruit, including tomato (Ripoll et al., 2014; Quinet et al., 2019; Hou et al., 2020; Medyouni et al., 2021). Recent comparative transcriptome analysis investigated responses of fruit to water stress in different tomato genotypes and uncovered quantitative trait loci for fruit quality traits (Albert et al., 2018; Diouf et al., 2020). However, these studies only considered the fruit pericarp during early fruit growth. Fruits are complex organs whose development involves coordinated gene expression among different tissues throughout development (Kang et al., 2013; Shinozaki et al., 2018). Tomato fruit develop from the fertilized ovaries and are composed of distinct interconnected tissues with unique and specialized roles. The inner fruit tissues are enclosed by the pericarp, which surrounds two or more locular cavities divided by the septum. The locules contain the seeds that are attached to the placenta, which acquires a jelly-like consistency during the later stages of fruit development. Vascular bundles that connect the fruit to the stem vasculature converge in the central part of the fruit or columella (Gillaspy et al., 1993; Lemaire-Chamley et al., 2005). Although all plant tissues are potentially sensitive to stress, the molecular responses implemented to cope with adverse conditions differ among organs, tissues, and developmental stages (Martignago et al., 2020). Accordingly, individual cell types and tissues present specific transcriptomic and metabolic responses to abiotic stress (Dinneny et al., 2008; Geng et al., 2013; Miao et al., 2017; Sarabia et al., 2018). Thus, in order to obtain a comprehensive picture of the regulatory mechanisms underlying the acclimation of the tomato fruit to drought, it is important to consider the specific transcriptional responses to water stress in the different constituent tissues.
Here we present a spatially resolved tomato fruit transcriptome (https://tea.sgn.cornell.edu/), that includes different tissue types from growing and ripe fruits subjected to long-term water stress, under different water deficit intensities. We characterized the distribution and timing of gene expression shifts in response to water stress acclimation and revealed tissue-specific changes and coexpression networks associated with the regulation of the fruit response to water deficit. We also collected expression data from tomato leaves and examined transcriptome signatures in response to water stress that are common to, or differ between, leaves and fruit. This study expands our understanding of the molecular mechanisms involved in acclimation to water deficit of vegetative and reproductive tomato organs and can offer insights into the responses to water stress in other fleshy fruit.
Results
Effect of long-term water stress on tomato vegetative growth and fruit development
To study the effect of long-term water stress on tomato plants and fruit, we used an experimental system in which plants were grown under well-watered conditions (control) at 40%–45% volumetric water content (VWC), or subjected to prolonged drought under one of three different levels of water deficit: S1 (35% VWC), S2 (30% VWC), and S3 (25% VWC). Water level deficit at 20% VWC or below was not applied as it caused wilting, considerable flower abortion and very poor fruit set. Leaves and fruit were collected over the course of 4–12 weeks after the start of the water limitation treatment (Figure 1A). We observed that the water stress conditions affected vegetative and reproductive traits, with the S2 and S3 treatments having the strongest impact (Figure 1, D and E; Supplemental Figure S1). A reduction in the size of plants and leaves was observed under all water deficit conditions. Additionally, there was an increase in flower abortion, resulting in lower fruit yield, while fruit set in flowers that reached anthesis was only slightly affected (Supplemental Figure S1A). Fruit size and weight were substantially reduced by lower water availability, with the S2 and S3 conditions causing >50% decrease in fruit weight (Supplemental Figure S1B). Notably, the growth inhibition had a differential effect on fleshy fruit tissues at the mildest stress, with a more pronounced size reduction in the placenta than in pericarp and septum (Supplemental Figure S1C).
Experimental design and phenotype of tomato plants and fruit exposed to prolonged water stress. A, Schematic illustration of the drought treatments and sample collection. Plants were grown under well-watered conditions (Control, C) or under water limiting conditions at three intensity levels (S1, S2, S3). C: 40%–45% VWC; S: 35% VWC; S2: 30% VWC; S3: 25% VWC. B, Circular dripping system and soil moisture sensor used to monitor and control soil water deficit. C, Equatorial fruit section showing the harvested tissues: pericarp, septum, columella, placenta, jelly, and seeds. D, Representative pictures comparing fruit development under control (C) and water limiting conditions (S1, S2, S3). From left to right: 20 dpa and ripe fruit. From top to bottom: whole fruit, longitudinal section and equatorial section. Scale bar = 1 cm. E, Representative pictures comparing the vegetative phenotype of tomato plants grown under control (C) and water limiting conditions (S1, S2, S3). Left: whole plants. Right top: leaves. Right bottom: leaflets. Scale bar = 2 cm.
Changes in gene expression in response to different levels of water stress intensity
Transcriptome analysis of fruit and leaves harvested from plants subjected to different watering regimes was performed using RNA-seq. The gene expression data were incorporated in the Tomato Expression Atlas (TEA; http://tea.solgenomics.net/; Fernandez-Pozo et al., 2017). Visualization of gene expression changes in this database is illustrated in the Supplemental Note and Supplemental Figure S2. A total of 52 samples were generated, each with three biological replicates, comprising leaves and different tissues (pericarp, septum, columella, placenta, locular jelly, and seeds) dissected from fruit at 20 d post anthesis (dpa), which corresponded to the exponential growth phase, and from fully ripe fruit (Figure 1C). We note that while seeds are organs, for the purpose of this study they are included as a category of fruit tissue. Cleaned RNA-Seq reads were aligned to the tomato reference genome version SL3.0 with ITAG 3.2 gene models (https://solgenomics.net/organism/Solanum_lycopersicum/genome), with an average of 11 million mapped reads per library (Supplemental Table S1). Pairwise Pearson correlation coefficient values between biological replicates ranged from 0.92 to 1.00 (Supplemental Table S2). In each tissue, the number of genes with expression levels above the background (set at RPKM ≥1) varied between 15,036 and 20,616, which represented 42% and 58% of the 35,769 genes in the tomato genome, respectively (Supplemental Table S3).
Principal component analysis (PCA) and hierarchical clustering revealed a separation of samples mostly by fruit developmental stage (PC1 51%). Samples also segregated according to tissue type and watering condition (Figure 2). Hierarchical clustering revealed that some tissues, such as columella and placenta at 20 dpa and pericarp and septum at the ripe stage, were grouped according to water deficit intensity, rather than by tissue type, demonstrating the predominant effect of the stress response in transcriptome remodeling (Figure 2B). A distinct clustering by water stress intensity was revealed for all tissues and developmental stages, with a clear separation between “control/S1” samples and “S2/S3” samples (Figure 2B).
Global gene expression changes in tomato leaf and fruit tissues in response to water stress. A, PCA and (B), hierarchical clustering analysis using log2-normalized read values of all genes in the different tissue types, stages, and watering conditions. Leaves and fruit at 20 dpa and at the ripe stage were harvested from control (C) or water stressed plants (S1, S2, S3). Lf: leaf; Pe: pericarp; Se: septum; Co: columella; Pl: placenta; Je: jelly; Sd: seeds. C, Number of DEGs in leaves and in fruit at different water deficit conditions. D, Venn diagram showing the percentage of DEGs shared between the three water deficit conditions. DEGs from all fruit tissues were compiled to obtain the total number of DEGs for each fruit developmental stage. S1: control vs. S1; S2: control vs. S2; S3: control vs. S3.
Differential gene expression analysis of 39 pairwise comparisons comprising all “control vs. stressed” tissue combinations (Supplemental Figure S3A), identified 13,464 unique differentially expressed genes (DEGs), corresponding to 38% of the tomato genes. Among these, 5,391 (40%) and 4,129 (30%) genes were down- or upregulated, respectively, in at least one of the 39 comparisons. Genes regulated in opposite directions (up or down) depending on tissue and/or water stress intensity represented 30% (3,944 genes) of the total DEGs (Supplemental Figure S3B). Gene ontology (GO) term enrichment analysis of the DEGs showed that exclusively upregulated genes were associated with biological processes related to “stress” and “proteolysis,” whereas exclusively downregulated genes were predominantly related to “sulfate transport,” “glutamate receptor,” and “xylan biosynthesis” (Supplemental Figure S3C). We selected genes consistently upregulated or downregulated in more than half of the pairwise comparisons, including leaf and fruit tissue samples, as candidate genes with a putative conserved role in water stress acclimation (Supplemental Table S4). Among these genes were two genes related to sulfate homeostasis: a gene encoding a putative high affinity sulfate transporter (Solyc04g054730) and a gene encoding a tetratricopeptide repeat protein (Solyc06g007970) that has 79% amino acid similarity to Arabidopsis SULPHUR DEFICIENCY-INDUCED 1(AtSDI1), which is involved in the utilization of stored sulfate under sulfur deficiency (Howarth et al, 2009).
The number of DEGs increased with stress intensity in both leaves and fruit, indicating a dose response to water stress (Figure 2C; Supplemental Table S5). The proportion of DEGs shared between the three water deficit conditions was 10% in leaves, and 14% and 20% in 20 dpa and ripe fruit, respectively. The S2 and S3 stress conditions shared a larger proportion of DEGs, with 32% in leaves and 42% and 50% in 20 dpa and ripe fruit, respectively (Figure 2D). GO term enrichment analysis of the DEGs in leaves and fruit at the different water deficit levels, indicated that the biological pathways affected were also dependent on the water stress level. While some processes, such as “secondary metabolite synthesis” in leaves, or “cell wall organization” in 20 dpa and ripe fruit were common to all stress conditions, other pathways were significantly altered only at a specific stress level. Overall, the S2 and S3 stress conditions had more similar affected pathways (Supplemental Figure S4).
Taken together, these results indicated that our experimental analysis discriminated between the different stress intensities both at the physiological and molecular levels. Additionally, they suggest that a distinct transcriptional reprograming, involving substantially more genes, occurs under the stronger (S2 and S3) compared to the mildest (S1) water stress.
Transcriptional tissue and fruit stage specificity in the response to water stress
Leaves shared only 13% of the drought-responsive DEGs with 20 dpa fruit and 9% with ripe fruit, under the S1 stress condition. The percentage of DEGs common to leaves and fruit increased with water stress intensity, reaching 34% in 20 dpa and 22% in ripe fruit at the strongest water deficit (S3). The proportion of shared DEGs between the two fruit stages also increased with the intensity of the water deficit and reached 20% under S3 water stress (Figure 3A). To highlight similarities and differences in water stress responses between leaves and fruit, and between fruit stages, a GO term enrichment analysis was performed for DEGs that were specific or common to leaves and fruit, or to each fruit stage (Figure 3B). The shared DEGs were enriched with GO terms related to “proteolysis,” “lipid storage,” and “lipid transport,” “cell wall,” “ion transport,” “defense,” and “abiotic stress responses.” In contrast, other biological processes were associated with a specific organ and/or stage, such as “cytokinin response” in leaves, “flavonoid metabolism” in 20 dpa fruit, “DNA repair” in ripe fruit, or “photosynthesis” in both 20 dpa and ripe fruit.
Specificity of transcriptional responses to water stress in tomato leaves and fruit. A, Venn diagrams showing the overlap between DEGs from leaves, 20 dpa and ripe fruit at different water stress intensity levels (S1, S2, S3). DEGs from every fruit tissue were compiled for each fruit developmental stage. B, Representative GO terms enriched in DEGs specific or common to leaves, 20 dpa and ripe fruit in response to prolonged water stress (S1, S2, S3). Point size represents FDR corrected P-value for upregulated (red) and downregulated (blue) DEGs.
We next examined gene expression changes associated with distinct fruit tissues (Figure 4). More than half of the total DEGs were specific to only one tissue at all the water deficit levels, and <1% of DEGs were shared between all fruit tissues. The tissue specificity of the transcriptome response was more prominent at the S1 stress level, with DEGs shared between tissues increasing with the stress intensity in both growing and ripe fruit (Figure 4A). The number of DEGs in each fruit tissue varied depending on water stress intensity and developmental stage. At the lowest water stress intensity (S1) the greatest expression changes were observed in the placenta (∼50% of the total fruit DEGs) in 20 dpa fruit, while in ripe fruit, seeds showed the greatest expression changes accounting for >60% of the total fruit DEGs (Figure 4B). A GO enrichment analysis of the DEGs in each tissue highlighted downregulated and upregulated biological processes associated with specific fruit tissues and/or water stress intensity levels (Figure 4, C and D; Supplemental Table S6). In 20 dpa fruit, downregulated DEGs were particularly enriched for those involved in “cell cycle,” “cell wall,” “vesicle transport,” and “microtubule movement” in the placenta, the tissue with the most pronounced gene expression changes at the S1 stress level. Other GO terms, such as “auxin response,” were enriched in the septum, columella, and placenta at the S2 and S3 stress conditions. In ripe fruit, downregulated DEGs were enriched for those related to “photosynthesis,” especially in the pericarp at the S1 stress level (Figure 4C). Analysis of upregulated DEGs in 20 dpa fruit revealed that pathways such as “flavonoid biosynthesis” in the pericarp, or “response to water” in the septum, columella, and placenta, were affected at the S1 stress level. GO terms related to jasmonic acid (JA) signaling were associated with all stress levels in the columella, or with the S2 and S3 stress in the jelly and seeds. In ripe fruit, enriched GO terms in upregulated DEGs included “xyloglucan metabolism” in the columella under all stress conditions. Strikingly, GO terms related to chromatin remodeling and gene silencing were specifically enriched in upregulated genes in the seeds (Figure 4D).
Tissue-specific gene expression changes in the tomato fruit in response to water stress. A, Percentage of total DEGs shared between fruit tissues at different water stress intensity levels (S1, S2, S3). B, Percentage of total DEGs in each fruit tissue. DEGs for each water stress condition (S1, S2 or S3) and developmental stage (20 dpa or ripe fruit) were compiled to obtain the total number of DEGs. C and D, GO enrichment analysis of downregulated and upregulated DEGs in each fruit tissue for each water stress condition. The different water stress intensity levels (S1, S2 and S3) are shown from left to right for each tissue. Pe: pericarp; Se: septum; Co: columella; Pl: placenta; Je: jelly; Sd: seeds. Only representative GO-terms, associated with the highest number of genes, are shown. All GO-terms significantly enriched (FDR < 0.05) in each tissue are listed in Supplemental Table S6.
To further analyze the tissue specificity of the response to drought, for each tissue we compared the percentages of drought-responsive genes among all expressed genes and among DEGs in control conditions (Supplemental Figure S5, A and B). No major differences were found indicating that tissue variations in control conditions did not affect the drought response. We also evaluated how the water stress response differed between tissue-specific genes and genes with homogeneous expression across different tissue types (Supplemental Figure S5, C and D). The percentage of genes that were affected by water stress was overall higher among tissue specific genes than in genes with similar expression in all tissues. This may reflect that processes underlying developmental variation are more likely to respond to stress than the conserved cellular machinery. We also observed that under S1 conditions, the percentage of drought-related genes among genes with similar expression in all tissues was much higher in placenta at 20 dpa and in seeds of ripe fruit than in the rest of tissues. This underscores the stronger sensitivity of these two tissues to a milder water stress.
Gene coexpression networks in tomato tissues exposed to water stress
To further examine the regulation of tissue and stage-specific transcriptional responses to water stress, we performed weighted gene coexpression network analysis (WGCNA). This resulted in 44 distinct coexpressed modules (labelled M1-M44), capturing 16,100 genes, containing between 31 (M44) and 2,628 (M1) coexpressed genes per module (Figure 5A; Supplemental Table S7). Several modules showed leaf, fruit tissue- and/or fruit stage-specific expression patterns in response to water stress. Two gene modules, M1 and M38, were specific to leaves and showed downregulation in response to water stress. Only two modules, M39 and M18, showed tissue specificity irrespective of fruit stage, while several of the modules were associated with a specific fruit stage. Some of the modules showed upregulated or downregulated expression patterns depending on the stress level (e.g. M10) or the tissue (e.g. M4). Noticeably, several modules exhibited both high fruit tissue- and stage-specificity. These included M25 in jelly and M29 in the columella, at the 20 dpa stage. Figure 5B shows the expression pattern of LIPID TRANSFER PROTEIN 2 (LTP2; Solyc10g075110), a representative gene from module M25 encoding a nonspecific lipid transfer protein (nsLTP), which was specifically upregulated in the jelly of 20 dpa fruit at the strongest water stress (S3), illustrating both stress intensity and tissue specific responses. LTP2, shows 67% sequence similarity with Arabidopsis LTP2, which is also induced by water stress (Chae et al., 2010). It has been suggested that nsLTPs play important roles in membrane stabilization, cell wall organization, or signal transduction during stress responses (Liu et al., 2015a; Missaoui et al., 2022).
WGCNA of the response of fruit and leaves to water stress. A, Heat map showing cluster-tissue/stage associations of the 44 WGCNA generated coexpression modules. Each row corresponds to a module and each column corresponds to a tissue and watering condition. From left to right for each tissue and developmental stage: C (control), S1, S2, S3. Color denotes ME values. Examples of stage and tissue-specific modules are shown in red (upregulated), blue (downregulated), and black (mixed pattern). B, Tissue- and stage-specific expression profile of a gene (LTP2, Solyc10g075110) representative of module M25. Images displaying expression levels (RPKM) were extracted from the TEA database (https://tea.solgenomics.net/). C, Network showing connections between core genes in the columella-associated module M29. Each colored circle (node) represents one gene. Diamond shapes represent TFs. Orange node color indicates greater connectivity within the network. The 500 strongest gene–gene interactions are shown. Only edge correlation values > 0.8 are shown and values > 0.9 are indicated in red. The expression profile of underlined genes is shown in (D). D, Expression profiles of hub genes in the columella-associated module M29. Expression data (RPKM) are means (±SD) of three biological replicates. Pe: pericarp; Se: septum; Co: columella; Pl: placenta; Je: jelly; Sd: seeds.
We created a module-GO network that highlighted the tissue-specific regulation of distinct pathways in response to water stress (Supplemental Figure S6, Supplemental Table S8). Modules that were downregulated by water stress and specific to leaves (M1) or to the growing fruit stage (M12), showed an enrichment in GO terms related to “photosynthesis” and “cell growth,” respectively. In contrast, an enrichment in terms related to “starch biosynthesis” was found in M4, which contained genes that were upregulated in fruit internal tissues (Figure 5A). These included genes encoding enzymes involved in starch biosynthesis, which generally showed patterns of upregulation in response to water stress in tissues that accumulate starch during fruit development (Supplemental Figure S7A). Consistent with this observation, starch quantification in pericarp and placenta at 20 dpa confirmed that water stress increased starch levels (Supplemental Figure S7B). In ripe fruit, an enrichment in “cell wall biogenesis and organization” was found in modules associated with the columella (M29, M31).
Notably, several of the tissue-specific modules indicated the modification of different hormone biosynthesis and signaling pathways in response to water stress (Supplemental Figure S6, Supplemental Table S8). Overall, the GO term enrichment analysis suggested a downregulation of the cytokinin pathway (M38) in leaves, an upregulation of the JA pathway in seeds and jelly (M2, M25), alterations in auxin transport in columella and placenta (M10) in 20 dpa fruit, and upregulation of ethylene biosynthesis in the septum and pericarp of ripe fruit (M33).
The gene connectivity in each module was analyzed, with a particular emphasis on identifying transcription factors (TFs) as potential hub genes in tissue- and/or stage-specific modules (Supplemental Tables S9 and S10). This analysis led to the identification of 70 TFs belonging to different TF families, several of which had previously reported roles in abiotic stress responses. To refine the search for hub TFs and better visualize gene interactions in coexpression modules, we modeled a network view using Cytoscape (Shannon et al., 2003). Three of the most striking examples of hub TFs were found in modules M29 (columella), M38 (leaf), and M25 (jelly and leaf; Figure 5, C and D; Supplemental Figure S8). In M29, we identified the tomato homolog of C-REPEAT BINDING FACTOR 3 (CBF3; Solyc03g026270), a TF belonging to the APETALA2/Ethylene-Responsive Factor (AP2/ERF) family, as a central hub connected to genes encoding cell wall remodeling enzymes, such as xyloglucan endotransglucosylase/hydrolases (XTHs; Figure 5C). CBF3 expression showed specific upregulation by water stress in the columella at the ripe stage, with an identical pattern to XTH genes (Figure 5D). M38, enriched in cytokinin-related genes, included a MYB family TF, MYB6 (Li et al., 2016), as the hub gene (Supplemental Figure S8A). MYB6 belongs to the R2R3-MYB subgroup containing transcriptional repressors of the phenylpropanoid pathway (Ma and Constabel, 2019) and of secondary cell wall deposition (Wang et al., 2019a). MYB6 expression was downregulated in response to water stress in parallel with LONELY GUY (LOG), a cytokinin biosynthesis gene, and genes coding for histidine phosphotransfer proteins (AHPs), involved in cytokinin signaling (Hutchison et al., 2006; Supplemental Figure S8B). M25 was enriched in genes related to JA biosynthesis and signaling, and contained AP2-ERF, MYB, and YABBY genes acting as potential hub TFs (Supplemental Figure S8, C and D, Supplemental Table S10).
We used the Analysis of Motif Enrichment (AME) tool (McLeay and Bailey, 2010) to detect enriched cis-regulatory elements 2 kb upstream of gene sequences from the above selected coexpression modules (Supplemental Figure S9, Supplemental Table S11). The results from this analysis supported a role for CBF3 as a central hub TF in M29. Indeed, promoters of the M29 genes were enriched with the DRE core motif (Maruyama et al., 2004; Zou et al., 2011) bound by CBFs, as well as binding motifs for calmodulin-binding transcription activator (CAMTA) proteins, which are CBF regulators (Shi et al., 2018). Motif enrichment in M38 further suggested that the cytokinin pathway is involved in the acclimation response to drought in leaves; specifically, promoters of the M38 genes showed a strong enrichment with AP2-ERF-binding motifs, including the cytokinin-related Cytokinin Response Factors (CRFs) motif (Rashotte and Goertzen, 2010; Raines et al., 2016; Xie et al., 2019). Promoters of genes in M25 were enriched in binding motifs for MYB and AP2-ERF, which is consistent with members of these TF families acting as hub genes in this module (Supplemental Figure S8C, Supplemental Table S10). Additionally, JA signaling-related MYC-binding sites (Howe et al., 2018) were also enriched in M25, suggesting that JA may modulate responses to severe water stress in the jelly and leaves.
Epigenetic-related gene expression changes in tomato seeds during water stress
Analysis of DEGs in seeds from ripe fruit revealed a specific enrichment with GO terms linked to chromatin organization and histone modification (Figure 4D). To gain more insight into transcriptomic changes associated with epigenetic regulation, we examined the expression profiles of 167 genes involved in histone marks, DNA methylation changes, and small RNA biogenesis, including components of the RNA-directed DNA methylation pathway (Bai et al., 2012; Aiese Cigliano et al., 2013; Gallusci et al., 2016). Most of the epigenetic-related gene expression changes occurred in mature seeds, where a strong upregulation of 25 genes was observed, including those encoding histone modifying enzymes, DNA demethylases and methyltransferases, and proteins involved in sRNA biogenesis (Figure 6A; Supplemental Table S12). We hypothesized that these genes play a role in chromatin remodeling and epigenetic regulation in mature seeds that would affect specific loci involved in water stress responses, thus improving resistance to water deficit in the next generation. To test this, we compared the progeny of control plants and the progeny of plants exposed to the S3 water deficit and observed that following a period of water deficit, seedlings grown from seeds generated under drought conditions were less affected by wilting and had a better recovery after rewatering than control plants (Figure 6B).
Epigenetic-related gene expression changes in ripe fruit and transgenerational effect of water stress on seedling performance. A, Heat-map showing gene expression in each fruit tissue under control or water stress conditions. The different watering conditions are shown from left to right for each tissue (Control, S1, S2, S3). Normalized RNA-seq expression was transformed into Z-score. Color key: blue, low level of expression; red, high level of expression. Gene identifiers and expression data (RPKM) are listed in Supplemental Table S12. B, Differences in recovery after water stress between seedlings grown from control and from water stressed (S3) seeds. Water stress was applied to 25-d-old seedlings by soaking pots in water to reach full pot capacity, and subsequently withholding water until all control plants showed a strong wilting phenotype. All plants were then rewatered at full pot capacity and pictures were taken 6 h after rewatering. Eight randomly placed control and S3 seedlings were used for each experiment, and the experiment was repeated 3 times with similar results.
Jasmonate- and ethylene-associated responses in tomato fruit following long-term water stress
Tissue-specific transcriptional changes, as well as gene coexpression networks, revealed variations in JA- and ethylene-related gene expression during the response to long-term water stress (Figure 4D; Supplemental Figure S6, Supplemental Tables S6 and S8). To further explore the spatiotemporal dynamics of JA and ethylene responses, we examined the expression profiles of 100 genes involved in their biosynthesis and signaling pathways (Figure 7; Supplemental Table S13), and observed significant changes in the columella and jelly. Several genes encoding key enzymes in the JA biosynthetic pathway, such as 13-lipoxygenase (13-LOX), allene oxide synthase (AOS), allene oxide cyclase (AOC) and 12-oxophytodienoate reductase (OPR), as well as Jasmonate ZIM-Domain (JAZ) repressor proteins involved in JA signaling (Chini et al., 2017), were strongly upregulated in these two tissues during water stress in 20 dpa fruit (Figure 7, A and C). Accordingly, we observed that water stress caused a strong increase in the levels of JA in both these tissues (Figure 7D). A similar increase in gene expression at 20 dpa was found for 1-aminocyclopropane-1-carboxylate synthase (ACS) and 1-aminocyclopropane-1-carboxylate oxidase (ACO) genes, encoding enzymes involved in ethylene biosynthesis.
Jasmonate- and ethylene-related gene expression changes and hormone levels in tomato fruit under water stress. Heat-maps show gene expression in each fruit tissue at control or water stress conditions for 20 dpa fruit (A) and ripe fruit (B). From left to right for each tissue: C, S1, S2 and S3. Pe: pericarp; Se: septum; Co: columella; Pl: placenta; Je: jelly; Sd: seeds. Normalized RNA-seq expression was transformed into Z-score. Color key: blue, low level of expression; red, high level of expression. Jasmonate (JA)-related genes are in green and ethylene (ET)-related genes are in orange. Gene identifiers and expression data (RPKM) are listed in Supplemental Table S13. C, Pathway for JA biosynthesis showing in red genes involved in the enzymatic steps and upregulated by water stress in columella and/or jelly tissues of 20 dpa fruit. D, JA was measured in columella and jelly fruit tissues at the 20 dpa stage from control and water stressed (S3) plants. Data are means (±SD) of six biological replicates. Statistically significant differences were determined using Student’s t test: **P<0.01). E, Ethylene was measured in fruit at the ripe stage from control (C) and water stressed (S1, S2 and S3) plants. Data are means (±SD) of eight biological replicates. Statistically significant differences were determined using Student’s t test: *P<0.05; **P<0.01.
In contrast, expression changes of JA- and ethylene-related genes in response to water stress showed a less consistent pattern during ripening (Figure 7B). Although gene expression related to JA biosynthesis was less pronounced at this stage, several JAZ genes were still strongly upregulated, mainly in the columella. In the ripe fruit, ACS and ACO genes involved in ripening-associated ethylene production (Liu et al., 2015b) displayed both tissue-specific and opposite regulation in response to water stress: ACS2 and ACO1 were upregulated, whereas ACS4 and ACO4 were downregulated. To investigate this phenomenon further, we quantified ethylene production and observed that water stress increased ethylene biosynthesis in ripe fruit, in a dose-dependent manner (Figure 7E). The expression profiles of ACS2 and ACO1 in response to water deficit suggest that they are involved in the induction of ethylene production by drought stress.
Effect of water stress on carotenoid biosynthesis
Given the key role of ethylene in regulating the ripening program (Li et al., 2021), the observed changes in ethylene production in response to water stress suggested a potential corresponding change in ripening-related traits, such as carotenoid content. The dramatic increase in ethylene production that occurs at the onset of tomato fruit ripening is known to control the rapid accumulation of the carotenoid compounds lycopene and β-carotene, which contribute to the coloration of the ripe fruit (Alba et al., 2005). We observed that fruit exposed to water stress appeared more red-colored than the control fruit, which showed a more orange coloration, suggesting an increase of lycopene (red pigment) and/or a decrease in β-carotene (orange pigment). Chroma Meter measurements confirmed the coloration changes (Figure 8A), showing a decrease of the lightness factor (L*), or darkening, in stressed ripe fruit compared to the control. No significant differences were found for the Chroma (C*) values between stressed and control fruit except for the S3 condition. Analysis of the Hue angle (°h) values confirmed that stressed fruit displayed a more intense red color while control fruit showed a more orange coloration. Quantification of phytoene, phytofluene, lycopene and β-carotene in pericarp tissues using high pressure liquid chromatography (HPLC) further showed that, lycopene, as well as its precursors phytoene and phytofluene, accumulated to a significantly higher level in stressed fruit, whereas β-carotene levels were lower, confirming that water limitation induced a change in carotenoid profiles (Figure 8B). Accordingly, genes encoding phytoene synthase (PSY; Solyc01g015285) and zeta carotene desaturase (ZDS; Solyc01g099640), which are involved in lycopene biosynthesis, showed increased expression in response to water stress. In contrast, genes encoding lycopene β-cyclases (β-LCY3, Solyc06g074240; and β-LCY1, Solyc04g040190) involved in the biosynthesis of β-carotene from lycopene, were downregulated by water stress in ripe fruit (Figure 8C; Supplemental Figure S10).
Effect of water stress on fruit color and carotenoid biosynthesis in the pericarp of ripe fruit. A, Effect of water stress on fruit color. Box plots representing differences in fruit color between plants exposed to different water stress intensity. Fruit color was determined by measuring the lightness factor (L*), the color intensity factor (Chroma, C*), and the actual color appearance factor (Hue angle, °h). Each dataset represents at least 20 biological replicates. Statistically significant differences between control (C) and each of stress condition (S1, S2, S3) were determined using the Student’s t test. The bold line in the center of the boxplots represents the median, the box edges represent the 25th (lower) and 75th (upper) percentiles, and the whiskers extend to the most extreme data points that are no more than 1.5× the length of the interquartile range. B, Effect of water stress on carotenoid content. Data are the mean (±SD) of four biological replicates. Statistically significant differences between control (C) and each of stress condition (S1, S2, S3) were determined using the Student’s t test: *P<0.05; **P<0.01; ***P<0.0001. C, Effect of water stress on carotenoid-related gene expression. Expression data (RPKM) are the mean (±SD) of three biological replicates. PSY: phytoene synthase; ZDS: ζ-carotene desaturase; b-LCY: lycopene ß-cyclase. Statistically significant differences between control (C) and each of stress condition (S1, S2, S3) were determined using the Student’s t test: *P<0.05; **P<0.01.
Changes associated with cuticle synthesis
The hydrophobic cuticle that coats the outer epidermis layer of land plants is critical for limiting water loss from aerial organs, including fruit (Yeats and Rose, 2013; Martin and Rose, 2014). Notably, water deficiency has been reported to induce cuticle biosynthesis in fully expanded immature and ripe tomato fruit (Romero and Rose, 2019), while expanding fruit were not included in this study. To investigate this phenomenon further, we analyzed the changes in the predominant cuticular constituents, cutin and waxes, in expanding fruit at 20 dpa and in ripe fruit. Our results were consistent with those of Romero and Rose (2019) in that levels of waxes and cutin monomers increased in response to water stress in mature fruit (Supplemental Figure S11). However, expanding fruit at 20 dpa from control and water deficient plants, did not show significant differences in cuticular compounds. One explanation for this differential effect on cuticle deposition at different developmental stages is that 20 dpa corresponds to the stage when cuticle deposition is at its highest and it may be that at this point the biosynthetic machinery is already at maximum capacity (Yeats et al., 2010), so water stress has no additional effect. Consistent with this idea, we observed that the expression of cuticle-related genes was not significantly affected, or in some cases was slightly downregulated, by water stress at 20 dpa (Supplemental Table S8, Supplemental Figure S12).
Discussion
Although transcriptomic and genomic approaches have generated large amounts of information regarding plant responses to stress (Zenda et al., 2021; Zhang et al., 2022), there have been few transcriptomic analyses aimed at understanding the molecular mechanisms of such responses in reproductive organs. Here, we present a spatially resolved transcriptome study of tomato fruit grown under long-term water stress. We created gene expression datasets from leaves and six fruit tissues at three water deficit intensity levels, which were then used to characterize tissue-specific transcriptional responses to water stress, and to contrast the responses in leaves and fruit. We also uncovered examples of tissue and stage-specific regulatory mechanisms associated with acclimation to long-term water stress in the fruit.
Water stress responses in tomato leaves and fruit comprise both common and organ-specific changes varying by stress severity
We observed physiological and transcriptional responses that were dependent on the intensity of water stress in both fruit and leaves, and a major difference between the mildest (S1) and more severe water deficits (S2 and S3; Figures 1 and 2). One of the most noticeable differences was reduction in organ size with more severe stress (Figure 1), which is in accordance with previous studies showing that water stress can inhibit the cell cycle process, thereby reducing plant cell division and organ size (Setter and Flannigan, 2001; Su et al., 2013).
Our results suggest that, in leaves, controlled long-term water stress affects biological processes related to cell wall modification and JA-related gene expression changes (Supplemental Figure S4). Downregulation of genes related to photosynthesis was also associated to M1, a coexpression gene module specific to leaves (Figure 5A; Supplemental Figure S6). These are similar to reported effects of acute stress by water withholding, which include downregulation of genes related to photosynthesis and cell wall modification (Iovieno et al., 2016), and changes in JA content (Muñoz-Espinoza et al., 2015). Additionally, acute stress has been shown to cause a progressive increase in ABA and upregulation of ABA signaling-related gene expression (Iovieno et al., 2016). However, in our conditions we did not observe major gene expression changes related to ABA. In contrast, GO terms related to the cytokinin- and auxin-responses were enriched in downregulated DEGs in leaves (Figure 3B; Supplemental Figure S4). This suggests that various levels of stress severity and duration might trigger distinct hormonal responses in the leaves.
A GO enrichment analysis of DEGs common to fruit and leaves highlighted biological processes such as proteolysis, ion transport, lipid metabolism and transport, and cell wall organization, affected by water stress in both organs (Figure 3B). Genes whose expression was consistently altered in response to water stress irrespective of organ, tissue and/or stage, included those involved in the regulation of sulfur homeostasis (Supplemental Table S4). Sulfate, and its derivatives such as glutathione, play important roles in the intrinsic responses of plants to abiotic stress (Gallardo et al., 2014; Chan et al., 2019). Therefore, gene expression changes related to sulfur homeostasis may be needed to re-equilibrate sulfate flux between and within plant tissues, as an important mechanism to confer tolerance to water stress. In addition to a gene encoding a sulfate transporter, we identified a gene with similarity to AtSDI1. In Arabidopsis, sdi1 mutants retain higher root and leaf sulfate concentrations (Howarth et al., 2009), suggesting the participation of the AtSDI1 tomato homolog in the utilization of stored sulfate pools under water stress conditions.
In addition to pathways that were commonly altered in multiple tissues/organs, our analysis highlighted responses that were specific to leaves or fruit. These organ-specific responses to water stress acclimation included changes related to different hormonal pathways, photosynthesis, and secondary metabolism among others, with the level of water deficit intensity influencing which processes were more significantly affected in each organ and stage (Figure 3; Supplemental Figure S4). For example, downregulation of genes related to cytokinin responses was specific to leaves under S2 and S3 water deficit levels. On the other hand, biological processes related to photosynthesis were, surprisingly, more broadly affected by water stress in fruit than in leaves (Figure 3; Supplemental Table S6). The apparent upregulation of photosynthesis-related genes in the 20 dpa fruit by water stress is intriguing. Although it appears that photosynthetic carbon assimilation does take place in green immature tomato fruit, the importance of fruit photosynthesis is still a matter of debate. Tomato fruit lack stomata, which suggests that fruit photosynthesis must depend almost exclusively on re‐assimilation of respiratory CO2 (Quinet et al., 2019). Under optimal conditions, any reduction in the rate of photosynthesis in the fruit can be balanced by the upregulation of leaf photosynthesis and an increase in photoassimilate import (Simkin et al., 2020). Our data support the hypothesis that under water stress, when leaf photosynthesis is compromised (Chaves et al., 2009; Raja et al., 2020; Muhammad et al., 2021), an increased rate of photosynthesis in the fruit may act as a mechanism to maintain fruit growth (Lytovchenko et al., 2011; Bloemen et al., 2013; Simkin et al., 2020). The overall downregulation of photosynthesis-related genes in ripe fruit under water deficit conditions suggests that the distinctive decrease in photosynthesis-associated transcripts that accompanies tomato ripening and chloroplast to chromoplast transition (Simkin et al., 2020), is accentuated by water stress. This could be a side effect of the increased ethylene production under water stress as this hormone is linked to chromoplast development due to its key regulatory role in carotenoid biosynthesis (Klee and Giovannoni, 2011). Chromoplast carotenoids accumulate in green tissues under stress conditions, and it has been suggested that ROS mediate the increase in carotenoid biosynthesis and the transformation of chloroplasts into chromoplasts (Bouvier et al., 1998; Choi et al., 2021). Therefore, we could also hypothesize that the transition from chloroplast to chromoplast is enhanced under water stress as a response mechanism involving carotenoid accumulation to protect against oxidative stress.
Water stress causes tissue and stage-specific transcriptional responses in tomato fruit affecting cell growth and metabolism, and epigenetic pathways
Our analysis uncovered a highly tissue-specific transcriptome response in the fruit: more than half of the total DEGs identified in fruit showed differential expression in only one tissue (Figure 4A). The degree of the response to the water stress intensity level was also tissue dependent. In the growing fruit, the largest expression shift was observed in the placenta under the mildest water deficit (Figure 4B). Interestingly, among different fleshy fruit tissues, the placenta was also the most sensitive to growth inhibition under this condition (Supplemental Figure S1C). The fact that this differential growth reduction in the placenta was not observed in the ripe stage, or at higher levels of water stress deficit, could be due to growth changes progressively affecting other tissues as the fruit develops, or to greater growth alterations caused by more severe stress and that affect equally all tissues. This suggests that the placenta responds to water deficit earlier than other tissues by repressing growth and activating an energy-saving program. Accordingly, GO enrichment analysis indicated that cell-cycle and cell growth-related processes were associated with genes downregulated in the placenta under the mildest stress, further indicating that substantial transcriptional reprogramming of pathways directly or indirectly related to cell growth (e.g. vesicle transport, endocytosis, and microtubule organization) is selectively triggered by this stress condition in this tissue (Figure 4C; Supplemental Table S6).
We identified a specific transcriptome signature in the seeds that accounted for the largest gene expression changes under the mildest stress in ripe fruit (Figure 4B), including changes related to epigenetic regulation and genome stability (Figure 4D; Figure 6A; Supplemental Tables S6 and S12). Epigenetic regulatory mechanisms, such as DNA methylation, post-translational histone modifications and noncoding RNAs play an important role in the regulation of abiotic stress responses (Kim, 2021). The progeny from “stressed seeds” showed improved seedling recovery from water stress (Figure 6B), and we propose that epigenetic modifications in the mature tomato seed may provide a mechanism to impart stress memory, preparing offspring for potential future stress. Such transgenerational effects of drought on seedling performance have been observed in other species, and stress-induced changes in the epigenome of the progeny are suggested as their potential cause (Lämke and Bäurle, 2017; Hatzig et al., 2018; Liu et al., 2021). Our study highlights members of the DEMETER-LIKE DNA demethylases (DMLs) and the SET domain group (SDG) histone methyltransferase gene families, as well as members of the small interfering RNA (siRNA) pathway (Supplemental Table S12), regulated by water deficit in seeds, as candidate genes involved in transgenerational inheritance of epigenomic modifications after water stress in tomato.
The gene expression data also indicated upregulation of flavonoid biosynthesis and glycosylation in the pericarp at the lowest stress, and upregulation of starch biosynthesis in several tissues under all stress conditions in the 20 dpa fruit (Figure 4D; Supplemental Figures S6 and S7A). This suggests an enhancement of these biochemical pathways as a component of fruit responses to water stress. Increased flavonoid biosynthesis and glycosylation may enhance the capacity to scavenge ROS and oxidative tolerance in the growing fruit during drought (Le Roy et al., 2016; Li et al., 2018). The upregulation of starch biosynthesis, further confirmed at the metabolic level (Supplemental Figure S7B), is consistent with reports of increased starch accumulation in tomato under drought and salt stress (Yin et al., 2010; Ripoll et al., 2016) and may contribute to maintaining soluble sugar levels in ripe fruit under stress conditions.
Tissue-specific regulatory networks in response to water stress
The WGCNA further elucidated the tissue-associated pathways and provided insights into their regulation during the response of fruit and leaves to water stress (Figure 5; Supplemental Figure S6, Supplemental Tables S7 and S8).
A module (M29) associated with the columella, a tissue with abundant vascular strands, included strongly upregulated genes involved in cell wall modification, such as XTH genes, during water stress in mature fruit (Figure 5, C and D). This suggests that specific changes in cell wall structure of the vasculature occur in the columella during drought acclimation. An increase in XTH gene expression is often associated to plant abiotic stress responses, suggesting a role of these enzymes in cell wall strengthening to protect plant tissues from water loss (Tenhaken, 2015). Furthermore, a strong connectivity was found between the stress-related TF CBF3 (Rai et al., 2013) and several cell wall–related genes, including those coding for XTHs, suggesting that CBF3 upregulates the expression of these genes. CBFs participate in drought and cold stress responses (Zhao and Zhu, 2016); however, their putative regulatory role in the remodeling of cell walls as a response to abiotic stress, as well as their potential target genes, remain elusive. Our results suggest the participation of CBF3 in a regulatory network upregulating the expression of XTH genes in the columella during water stress, although the direct binding of CBF3 to the promoters of these genes has yet to be verified experimentally.
Other coexpression modules highlighted transcriptional changes in particular components of different hormonal pathways associated with specific tissues. These may constitute a mechanism to tailor hormone-mediated responses to water stress according to tissue type functions. For example, module M38, a downregulated, leaf-specific gene coexpression module, included cytokinin biosynthesis genes and cytokinin signaling genes (Supplemental Figure S6, Supplemental Figure S8, A and B, Supplemental Tables S7 and S8). Moreover, promoters of genes in this module were significantly enriched in CRFs-binding motifs (Supplemental Figure S9, Supplemental Table S11). These results support the repression of the cytokinin response as one of the strategies tomato leaves use to cope with drought stress (Nguyen et al., 2016; Cortleven et al., 2019).
Our data also showed that JA metabolism and signaling were strongly stimulated in a tissue specific manner, with water stress increasing the levels of JA in columella and jelly (Figure 7, A–D). Module M25, containing several JA biosynthesis and JAZ encoding genes, suggested a role for JA in the response of jelly tissue to the more severe stress (Supplemental Figure S8, C and D, Supplemental Tables S7 and S8). JAZ proteins interact directly with various TFs, such as members of the bHLH (MYC), MYB, and WRKY families, to regulate diverse JA responses (Howe et al., 2018). Notably, promoters of genes in M25 were enriched with binding motifs for those TFs, including MYC- and MYB-binding sites (Supplemental Figure S9, Supplemental Table S11). Furthermore, M25 contains MYB24, a hub TF (Supplemental Table S10) with homologs in Arabidopsis and apple (Malus domestica) that have been shown to interact with JAZ proteins and to control JA-mediated processes, such as floral organ development and flavonoid synthesis (Song et al., 2011; Wang et al., 2019b; Zhang et al., 2021). JA signaling and specific JAZ-TF complexes have been reported to contribute to abiotic stress responses; however, their specific role in regulating drought responses is relatively unexplored (Fu et al., 2017). Our results suggest that MYB24 and JAZ proteins may interact as part of a regulatory network dictating tissue-specific responses (i.e. in leaves and jelly) to water stress acclimation in tomato and this will be the target of further experiments. Moreover, the presence of several ERF TFs as well as ACO1 in M25 suggest crosstalk and/or synergistic action between JA and ethylene in regulating these responses (Zhu, 2014; Guo et al., 2018; Yang et al., 2019). The combined action of JA and ethylene in defense responses against pathogens is well known (Zhu, 2014) and it has been suggested that they may also interact during the responses to abiotic stress (Kazan, 2015). Upregulation of the expression of both JA- and ethylene-related genes occurred in response to water stress in the internal fruit tissues, mainly in the columella, at the growing fruit stage (Figure 7A). This further suggests that these two hormones may act in concert to coordinate responses to water stress in tomato fruit tissues. Given the role of JA in stress-induced, long-distance signaling (Howe et al., 2018), the stimulation of JA synthesis and signaling in the columella (Figure 7D) may be part of a mechanism to trigger and coordinate the responses to water stress in the rest of the fruit.
Taken together, these results emphasize that during the response and acclimation to drought, individual fruit tissues exhibit tailored growth, epigenetic and hormonal responses that are regulated to different degrees by the water stress intensity.
Water stress causes transcriptional and metabolic changes in ripening-related pathways
Ethylene is a key regulator of tomato fruit ripening and stress responses (Li et al., 2021), but studies on its role in the response of fruit tissues to water deficit are relatively rare (Adato and Gazit, 1974; Nakano et al., 2003; Hou et al., 2020). While we observed an increase in ethylene production in ripe fruit exposed to prolonged water stress (Figure 7E), the changes in expression of ethylene biosynthesis genes varied substantially between fruit tissues and exhibited contrasting patterns (Figure 7B). This underscores the complex role of ethylene in regulating and potentially coordinating the responses to water stress in various specialized tissues of complex organs. The expression patterns of ACS2 and ACO1, which are responsible for ripening-related ethylene production (Cara and Giovannoni, 2008; Liu et al., 2015b), suggest that they are involved in ethylene synthesis associated with water stress. Based on the expression data and ethylene production, we propose that the upregulation of ACO1 in the columella plays a major role in the enhanced synthesis of ethylene in response to water stress (Supplemental Table S13).
Tomato fruit ripening involves ethylene-dependent transcriptional changes leading to accumulation of lycopene (Liu et al., 2015c). We observed an increase in lycopene levels concomitant with a decrease in β-carotene abundance in the pericarp of ripe fruit that paralleled the ethylene increase in response to water stress (Figure 8B). Similarly, genes encoding enzymes involved in lycopene and β-carotene biosynthesis were up- and downregulated, respectively, by water stress (Figure 8C; Supplemental Figure S10). These changes are consistent with drought-stimulated ethylene production since PSY1 and ZDS are both induced by ethylene, while the LCY genes are repressed by ethylene (Alba et al., 2005; Su et al., 2015; McQuinn et al., 2020). Therefore, we hypothesize that the regulation of genes involved in the carotenoid pathway is mediated by the enhanced ethylene production in stressed fruit.
Summary
Our study shows that water deficit causes distinctive transcriptional remodeling in fruit and leaves and that fruit responses involve tissue-specific changes affecting diverse biological pathways. We found examples of transcriptional changes associated with internal fruit tissues, underscoring the importance of tissue-specific profiling to elucidate mechanisms underlying the coordinated responses of complex organs, such as fruit, to water stress. We also obtained insights into potential regulatory mechanisms during water stress acclimation in fruit and leaves, highlighting TFs and their putative targets and interactors, that may act coordinating tissue-specific cell wall remodeling and hormone responses to drought stress. These data, which are publicly accessible via the TEA database (https://tea.sgn.cornell.edu/), constitute a rich spatially resolved transcriptome resource for tomato, a model system for fleshy fruit, and can be used to guide functional studies to elucidate mechanisms of acclimation to water stress.
Materials and methods
Plant material and water stress treatments
Tomato (S. lycopersicum “M82”) plants were grown in a greenhouse at Cornell University (Ithaca, NY, USA) with a light/dark photoperiod of 16 h/8 h. Soil moisture sensors (MAS-1, Decagon Devices, Pullman, USA) were used in every pot to monitor pot VWC and to regulate a circular dripping irrigation system operating in each pot (Figure 1B). Plants were first grown under well-watered conditions (40%–45% VWC) and water stress treatments started progressively 8 weeks after sowing, when plants started to flower. Soil moisture content was gradually decreased by 5% every week until the desired soil moisture content was reached for each water stress condition: S1 (35% VWC), S2 (30% VWC), and S3 (25% VWC). Flowers were tagged and pollinated manually at anthesis, starting 12 weeks after sowing. Leaves and fruit at 20 dpa and at the red ripe stage (Shinozaki et al., 2018), were collected starting 15 weeks after sowing (Figure 1A). Leaves and fruit tissues were harvested from 14 (control and S1 conditions) to 28 (S2 and S3 conditions) randomly selected individual plants at the same time each day. Leaves consisted of terminal leaflets from the third leaf from the top of the main shoot.
To test the recovery from water stress, pots with 25-d-old seedlings were soaked in water to reach full pot capacity, and then water was withheld until all control plants showed a strong wilting phenotype. Plants were then rewatered to full pot capacity and pictures were taken 6 h after rewatering. Eight randomly placed control and S3 seedlings were used for each experiment, and the experiment was repeated 3 times.
For flower number, fruit traits, ethylene, JA, starch, carotenoids, and cuticle measurements and for gene expression analysis, the Student’s t test (*P < 0.05; **P < 0.01; ***P < 0.0001) was used to assess statistical significance between control and each of the water stress treatments (S1, S2, S3).
RNA isolation, sequencing, and analysis
Plant samples were ground in liquid nitrogen and total RNA was isolated using an RNeasy Mini Kit (Qiagen). Strand-specific libraries were constructed according to Zhong et al. (2011) and sequenced on an Illumina NextSeq 500 from a single end for 75 bp. Raw RNA-Seq reads were processed using Trimmomatic (Bolger et al., 2014) to remove adapter and low-quality sequences, and then aligned to the ribosomal RNA database (Quast et al., 2013) using Bowtie (Langmead et al., 2009) and the mappable reads discarded. The resulting high-quality cleaned reads were aligned to the tomato Heinz genome reference (Tomato Genome Consortium, 2012) SL 3.0 with gene models iTAG3.2 (https://solgenomics.net/organism/Solanum_lycopersicum/genome) using HISAT (Kim et al., 2015). Following alignments, raw counts for each tomato gene were obtained using an in-house script and normalized to reads per kilobase of transcript per million mapped reads (RPKM). Pairwise Pearson correlation coefficient values between biological replicates were calculated with log2-transformed RPKM values (i.e. log2 [RPKM + 1]) using the cor function in the R program (https://www.r-project.org). Genes with averaged RPKM among replicates ≥1 were considered expressed. PCA was performed to compare the log2-transformed RPKM values of the expressed gene profiles among tissue types, stages, and water stress treatments, using the prcomp function in R v3.5.1. The log2-transformed (RPKM + 1) values of all the genes were used for hierarchical clustering using dist (method = “Euclidean”) function and hclust (method = “average”) function in R. Raw read count data were imported to edgeR (Robinson et al., 2010) to identify DEGs between each pair of water stressed and control samples with a cutoff of absolute fold change ≥2 and false discovery rate (FDR) <0.05. Gene annotations, iTAG3.2 GO terms, were obtained from the PANTHER classification system version 10 (Mi et al., 2016). GO enrichment analysis was performed for GO Biological Process Complete using the Bioconductor R library clusterProfiler (Yu et al., 2012), applying a hypergeometric test with FDR correction (adjusted P < 0.05).
Weighted gene coexpression network analysis
Coexpression network modules were identified using the WGCNA package (v1.66) in R (Langfelder and Horvath, 2008). Expressed genes with a high coefficient of variation of averaged RPKM (CV > 0.5) among all samples (16,112 genes) were used for the analysis. A Pearson correlation matrix of pairwise gene expression was constructed using the relative RPKM values, and the adjacency (type = “signed”) function with the soft threshold (power) of 22 (corresponding to an R2-value of 0.903). A topological overlap (TO) matrix was calculated using the correlation matrix and the TOMsimilarity function. Gene modules, with a minimum of 30 genes, were obtained by hierarchical clustering based on the TO similarity using a dynamic hybrid tree-cut algorithm (R package dynamicTreeCut, v1.63-1). Modules with dissimilarity <0.1 were merged to avoid over-splitting modules. Module eigengene (ME), which summarizes the expression profile of each module as the first principal component, was calculated for individual modules. ME-based gene connectivity (kME) was calculated using the signedKME function as a Pearson correlation between the expression of individual genes and the MEs. The hub genes were defined as those with kME > 0.9 within the assigned module. Cytoscape v. 3.7.1 (https://www.cytoscape.org) was used to visualize networks showing connections between core genes within WGNCA modules.
Fruit color, ethylene, and carotenoids measurements
The L, a, and b color space values of red ripe fruit were measured with a Konica Minolta chromameter (CR‐400; Ramsey, NJ, USA) and used to determined Chroma (C*) and Hue angle (h*) coordinates (Arias et al., 2000). The chromameter was calibrated against a CR-A43 white calibration plate. At least 20 fruit per watering regime were measured and for each fruit, three independent measurements were taken in the equatorial part and averaged.
Ethylene was quantified in ripe fruit (Breaker + 8 d) placed in 250 mL open jars overnight at room temperature to reduce harvest stress. Jars were then sealed, incubated at room temperature for 4 h, and 1 mL of headspace was injected into an Agilent 6850 Series II Network GC System gas chromatograph equipped with a flame ionization detector (Agilent, Santa Clara, CA, USA) and an activated alumina column. Ethylene concentrations were calculated using an ethylene standard of known concentration, and normalized by fruit mass. At least eight fruits were used per watering regime.
Carotenoid extraction and quantification by HPLC were performed as previously described (Vrebalov et al., 2009).
Starch quantification
Starch measurement was performed according to Vargas-Ortiz et al. (2013). Briefly, soluble sugars were extracted three times in 50 mM HEPES KOH (pH 7.4), 5 mM MgCl2, 80% ethanol (v/v), at 80°C. The insoluble pellets were dissolved in 0.25 mL 10 mM KOH at 95°C for 1 h. Starch was then hydrolyzed in 50 mM sodium acetate buffer (pH 5.5) at 37°C overnight, by the addition of 10 units of α-amylase and 10 units of amyloglucosidase. Soluble extracts were enzymatically assayed for glucose in a microplate reader (Synergy 2, BioTeck Instruments).
JA quantification
JA was analyzed as described in Patton et al. (2020) with slight modifications. Fresh tissue (150 mg) was extracted in 1 mL of extraction buffer (2:1:0.005 of iso-propanol, HPLC grade H2O, and hydrochloric acid), spiked with a deuterated standard of JA (10 ng/uL). Samples were then centrifuged at 14,000 RPM for 20 min at 4°C and then the supernatant transferred to new tubes containing 1 mL of dichloromethane. After vortexing and centrifuging, the dichloromethane layer containing the compounds of interest and the internal standards was removed, dried in a SpeedVac Concentrator (Savant Instruments Inc., Farmingdale, NY), and resuspended in 200 µL of methanol. Next, 5 µL of the extracted sample was injected into a Dionex UHPLC system (Thermo Scientific, Waltham, MA, USA) with a Kinetix C18 column of particle size 1.7 µm, length 150 × 2.1 mm, 100 Å (Phenomenex, Torrance, CA, USA). JA was detected in an Orbitrap-Q Exactive mass spectrometer (Thermo Scientific, USA) using signature ions and retention times. Data were analyzed using the Xcalibur 3.0 program (Thermo Fisher Scientific Inc. Waltham, MA, USA). Relative amounts of JA were quantified by comparing endogenous concentrations with the 10 ng of JA internal standard.
Cuticle measurements
Cuticular waxes were extracted as previously described (Fich et al., 2020) with 50 mg of tetracosane added as internal standard to 20 dpa samples and 100 mg to ripe fruit samples. Cuticle discs were excised from dewaxed fruit and the cutin component isolated with an enzyme solution as previously described (Fich et al., 2020). Cutin samples were then depolymerized by methanolysis (Philippe et al., 2016) and then waxes and cutin monomers were derivatized and injected into a gas chromatography system as described in Romero and Rose (2019). Four biological replicates were used for each stress condition.
Data availability
Raw RNA-Seq reads have been deposited into the NCBI BioProject database under accession PRJNA853266. Gene expression profiles can be accessed at the Tomato Expression Atlas (TEA; http://tea.solgenomics.net/). All other data supporting the findings are available in the article and the Supplemental Information files.
Accession numbers
The accession numbers of the genes mentioned in this manuscript can be found in Supplemental Tables S9, S12, and S13. Sequence data from this article can be found in the Sol Genomics Network (SGN) database (https://solgenomics.net/organism/Solanum_lycopersicum/genome; ITAG Release 3.2).
Supplemental data
The following materials are available in the online version of this article.
Supplemental Figure S1. Effect of different levels of water deficit on flower number and fruit traits.
Supplemental Figure S2. High-resolution transcript profiling using the TEA database.
Supplemental Figure S3. Differential gene expression analysis of the response of tomato fruit and leaves to prolonged water stress.
Supplemental Figure S4. GO enrichment analysis of the DEGs from fruit and leaves at different water deficit conditions (S1, S2, S3).
Supplemental Figure S5. Effect of the tissue specificity of gene expression in control conditions on the drought-related expression changes at different water deficit conditions (S1, S2, S3).
Supplemental Figure S6. Module-GO network.
Supplemental Figure S7. Effect of water stress on starch biosynthesis in tomato fruit.
Supplemental Figure S8. Connection networks and expression profiles of core genes in the leaf-associated module M38 and jelly/leaf-associated module M25.
Supplemental Figure S9. TF-binding core motifs significantly enriched in the promoters of genes belonging to modules M29, M38, M25.
Supplemental Figure S10. Effect of water stress on the expression of carotenoid-related genes in septum, columella, placenta, jelly, and seeds at the ripe stage.
Supplemental Figure S11. Effect of water stress on the fruit cuticle.
Supplemental Figure S12. Effect of water stress on the expression profiles of genes involved in cutin biosynthesis.
Supplemental Table S1. Total mapped reads in each sample.
Supplemental Table S2. Pearson’s correlation coefficient between biological replicates.
Supplemental Table S3. Number of expressed genes in each sample.
Supplemental Table S4. DEGs showing a consistent pattern in more than half (20) of the 39 pairwise comparisons (“control vs. stressed”).
Supplemental Table S5. Number of DEGs in each tissue at different water stress intensity levels (S1, S2, S3).
Supplemental Table S6. GO term enrichment analysis of the DEGs in tomato leaves and fruit tissues in response to prolonged water stress.
Supplemental Table S7. Module assignment and ME-based gene connectivity in the WGCNA.
Supplemental Table S8. GO term enrichment analysis in WGCNA modules.
Supplemental Table S9. Hub genes in WGCNA modules.
Supplemental Table S10. Hub TFs in WGCNA showing tissue specific responses to water stress.
Supplemental Table S11. TF-binding motifs significantly enriched in the promoters of genes in Modules 29, 38, and 25.
Supplemental Table S12. Genes related to epigenetic regulation showing significant differential expression (Fold change > 2 or Fold change < 0.5, P < 0.05) in at least one fruit tissue and one water stress condition, in ripe fruit.
Supplemental Table S13. JA- and ethylene-related genes showing significant differential expression (Fold change > 2 or Fold change < 0.5, P < 0.05) in at least one fruit tissue and one water stress condition, in fruit at 20 dpa.
Supplemental Note and References. Data visualization and references cited in Supplemental data.
Acknowledgments
We thank Ethan Eddy for help with JA quantification.
Funding
This project was supported by grants to J.K.C.R., J.J.G., C.C., and Z.F. from the Plant Genome Research Program of the US National Science Foundation (IOS-1339287), as well as to J.K.C.R. and J.J.G. by an award (59-8062-9-003P) from the Agricultural Research Service of the United States Department of Agriculture (USDA), and by grants to J.K.C.R. (2020-03667) and to C.C. (2019-67013-29240) from the Agriculture and Food Research Initiative of the USDA.
Conflict of interest statement. None declared.
P.N., J.K.C.R., S.I.S., and C.C. designed the water stress and tissue collection experiments. P.N. performed the experiments and analyzed data. S.I.S. contributed to setting up controlled water deficit conditions. P.N. and Y.X. performed cDNA library construction. P.N., K.B., Y.S., Y.Z., A.P., and Z.F. performed data bioinformatic analysis. A.P. and L.A.M. contributed to database development. G.P. and J.K.C.R. performed cuticle analysis. L.C. contributed to carotenoid quantification. J.V. contributed to ethylene measurements. C.L.C. performed JA quantification. P.N. and C.C. wrote the manuscript. J.K.C.R., J.J.G., Y.S., and Z.F. revised and edited the manuscript. All authors reviewed and approved the final manuscript.
The author responsible for distribution of materials integral to the findings presented in this article in accordance with the policy described in the Instructions for Authors (https://academic.oup.com/plphys/pages/general-instructions) is: Carmen Catalá ([email protected]).
