Sulfate starvation has a major impact on leaf and root growth in early stages of tomato development
As previously reported for other tomato cultivars [27–30], we found that the growth of Solanum lycopersicum cv. Moneymaker seedlings was severely impacted by sulfate starvation. As shown in Fig. 1A, lack of external sulfate has a clear negative effect on the growth of the aerial tissue, as compared with the full nutrient condition. This effect is noticeable from the third week after sowing. Multifactorial analysis of leaf fresh weight and total leaf area showed that there is a strong interaction between sulfate availability and time factors (p-value <0.001, two-way ANOVA, Fig. 1B and Figure S1), indicating that the growth response of leaves to sulfate availability depends on the plant developmental stage. Accordingly, no significant differences between control and sulfate-starved plants were observed at 2 weeks after sowing (p-value=0.81), while a strong dependence of external sulfate on leaf fresh weight and area was observed 3 and 4 weeks after sowing (p-value <0.001) (Fig. 1B and Figure S1). The same interaction between sulfate and time was observed for root weight, with differences between control and sulfate-starved plants being significant only after 3 and 4 weeks of sowing (p-value <0.001, Fig. 1C). These results indicate that whole-plant growth is severely affected by the lack of external sulfate from 3 weeks after sowing.
Consistent with leaf weight and leaf area measurements over time, there is a steady accumulation of sulfate in leaves, that is abolished in sulfate-starved plants in weeks 3 and 4 (Fig. 1D). In contrast, the increase in root weight over time of plants grown in control condition is not explained by a concurrent increase in root sulfate levels (Fig. 1E). Moreover, no major differences in sulfate accumulation are observed in the roots of full nutrient or sulfate-starved plants (Fig. 1E). Together, these results indicate that additional factors besides sulfate content might be controlling root growth response to sulfate starvation.
Sulfate starvation promotes major changes in leaf and root transcriptomes in a time-dependent manner
In order to gain insights into the molecular mechanisms underlying leaf and root growth repression in response to sulfate starvation, we performed a transcriptomic analysis at a temporal and organ scale, using RNA-seq. To perform this analysis, total RNA from roots and leaves was extracted at 2, 3 and 4 weeks after sowing, considering three independent biological replicates. The RNA-Seq data was pseudo-aligned to the ITAG3.20 transcriptome annotation using kallisto [31] as detailed in Materials and Methods. On average, each sample had 30 million reads, of which 27 million reads (90%) were pseudo-aligned to the tomato transcriptome (Figure S2).
Principal component analysis (PCA) of the RNA-seq data showed that samples are well separated along the first component, which explains more than 50% of the variability in leaf and root samples (Fig. 2). Interestingly, full nutrient and sulfate-starved samples from leaves are grouped together at 2 weeks after sowing, indicating similar transcript expression profiles at this time point (Fig. 2). In contrast, control and sulfate-starved samples from root tissue at 2 weeks are well separated, suggesting that sulfate starvation alters the plant transcriptome in an organ- and time-dependent fashion, with an earlier response in root tissue.
Accordingly, no genes were found differentially expressed between control and sulfate-starved leaves at two weeks after sowing, while 438 genes are differentially expressed in roots in this time point (log2 Fold Change>1 and q-value<0.05) (Figure S3A and Table S1). In addition, 55% of these 438 genes are subsequently regulated in leaves (241 genes, Figure S3B), which point out an orchestrated whole organism response to sulfate deprivation over developmental time.
A considerable number of genes were differentially expressed by sulfate starvation in roots and leaves: 7024 genes at 3 weeks after sowing and 6163 at 4 weeks after sowing, with most genes being differentially expressed in leaves (Figure S3A and Table S1). This major reprogramming of the plant transcriptome correlates with the observed differences in leaf and root growth of control and sulfate-deprived plants between weeks 3 and 4, suggesting changes in gene expression are partly responsible for the observed phenotypes (Fig. 1A and B).
In order to identify genes that are differentially expressed in response to sulfate starvation, time or a sulfate starvation-time interaction, we performed a multifactorial analysis of RNA-seq data in each tissue using sleuth [32]. In the case of roots, 6052 genes were significantly affected by time, 3755 genes by sulfate starvation and 714 by the interaction of both factors (q-value <0.05, and log2FC >1) (Fig. 3A). In the case of leaves, 6084 genes were affected by sulfate starvation, 7062 genes by time and 4571 genes by the interaction of both factors (q-value <0.05, and log2FC >1) (Fig. 3B). In contrast to root tissue, these results indicate that the response to sulfate starvation in leaves has a strong dependence on the plant age. In fact, about 90% of sulfate-responsive genes in leaves are significantly modulated by time (Fig. 3B).
To validate the results of the RNA-seq analysis with an alternative RNA quantification methodology, we selected 10 sulfate-regulated genes in root and leaves involved in different molecular and biological functions including sulfate metabolism, transport and transcriptional regulation, as well as a group of genes with unknown functions (Table S2). The expression patterns obtained by qPCR analysis showed a significant positive correlation with RNA-seq data (P <0.0001; 0.89 Pearson correlation value; Figure S4), indicating that most of genes analyzed showed a similar expression patterns than the ones identified by RNA-seq analysis.
Comparative analysis of the tomato and Arabidopsis sulfate starvation-responsive transcriptome
In a previous study, we identified 2046 genes regulated by sulfate availability in Arabidopsis thaliana using an integrative meta-analysis approach of transcriptomic data [15]. This meta-analysis included samples of roots and leaves at different stages of development (from seedling to juvenile plants). In the case of tomato plants, our multifactorial analysis identified 7589 sulfate-responsive genes in roots and leaves. In order to compare the sulfate starvation responsive transcriptome of both Arabidopsis and tomato, we first performed an orthology analysis of sulfate responsive genes. The identification of orthologous groups allows cross-referencing of genes from multiple species [33]. We used the PLAZA 4.0 database [34] to identify the ortholog group to which each sulfate-regulated gene belongs. Thus, we identified 4262 orthologous gene families out the 7589 sulfate-responsive genes in tomato and 1338 orthologous gene families out the 2046 sulfate-responsive genes in Arabidopsis (Table S3). We found that 70% of the Arabidopsis orthologous gene families are associated with sulfated-responsive genes were also present in tomato plants (Fig. 4A). This result indicates that the RNA-seq analysis in tomato plants captured most of the sulfate-responsive gene families previously detected in Arabidopsis.
Over-representation analysis of biological functions shared by both species revealed several GO terms that are expected to be conserved in the sulfate response as “cellular oxidant detoxification”, “sulfate transmembrane transport” or “sulfur utilization” (Table S4). Interestingly, we also found a significant enrichment (adjusted p-value <0.05) in GO terms which have not been previously studied in the context of the response to this nutrient such as those related to cell wall or cytokinin transport (Table S4). On the other hand, specific GO terms of the Brassicaceae family such as GSL biosynthesis are among the most over-represented biological processes in the case of sulfate-responsive genes in Arabidopsis (Fig. 6B and Table S4). The intersection between the lists of orthologues genes of both species also detected sulfate-responsive processes that are found exclusively in tomato plants. The most over-represented biological process in this group of genes was "cellular response to phosphate starvation" suggesting that sulfate deficiency also affects internal phosphate levels in tomato plants. Specifically, we found 52 genes in this biological process including SPX genes [35], phosphatases and those encoding enzymes involved in phosphate mobilization by membrane lipid remodeling such as sulfoquinovosyl-diacylglycerol [36, 37] (Figure S5). In order to verify the response of this set of genes to sulfate starvation in tomato plants, we computed the average fold of change between sulfate starved and control samples in both tissues and we then performed a hierarchical clustering analysis. Thus, we found that 60% of the genes involved in the phosphate starvation response (30 out 52 genes) showed lower expression levels in sulfate-starved plants at 3 and 4 weeks after sowing in both tissues (Cluster 3, Figure S5), suggesting that sulfate starvation might alter internal phosphate levels. To address this hypothesis, we determined the internal phosphate levels at 3 weeks after sowing in both tissues. The results shown in Fig. 4D, indicated that sulfate-starved plants exhibited significantly higher phosphate levels that control plants in both tissues, with a greater phosphate accumulation in leaves.
Temporal dynamics of gene expression of tomato roots and leaves are altered by sulfate starvation
In order to get further insights about the temporal patterns of expression of the genes and associated biological processes affected by sulfate starvation, we computed Pearson correlation indexes for each pair of sulfate-regulated genes in roots and leaves. We performed a hierarchical clustering analysis using Dynamic TreeCut [38], and identified two co-expression clusters for genes exclusively regulated by sulfate and six co-expression clusters for genes regulated by sulfate and time in roots and leaves (Fig. 5A, Fig. 5B, Figure S6, Figure S7). We found that the majority of sulfate-responsive genes (60% in roots and 72% in leaves) are contained in Clusters 1 and 2 of genes regulated by both factors (Fig. 5A and Fig. 5B). We thus decided to analyze these clusters in more detail.
As shown in Fig. 5C and Fig. 5D, the expression profiles of Clusters 1 and 2 in roots and leaves show similar trends. Genes in Cluster 1 are repressed over time, with an earlier decrease in mRNA levels under sulfate starvation conditions, between weeks 2 and 3. On the other hand, genes in Cluster 2 slightly increase between weeks 3 and 4 in control conditions, while an important increase in gene expression is shown between weeks 2 and 3 under sulfate starvation conditions (Fig. 5C and Fig. 5D). Although the expression patterns of the genes contained in Clusters 1 and 2 are similar between roots and leaves, the identity of these sulfate-responsive genes differs between these organs (Table S5 and Table S6). As such, we found different enriched biological processes associated with Clusters 1 and 2 in roots and leaves. In roots, Cluster 1 is enriched in genes associated to growth, such as “epidermal cell division” and “plant-type secondary cell wall biogenesis”, as well as processes related to redox activity, such as “hydrogen peroxide catabolic process” and “cellular oxidant detoxification” (Fig. 5C and Table S7). In leaves, we found an enrichment of genes involved in photosynthesis, microtubule-based movement and several GO terms related with cell division such as “mitotic cell cycle checkpoint” or “mitotic spindle organization” (Table S8). GO terms related to photosynthesis are especially abundant in this cluster, so we analyzed in more detail the distribution of these genes across different categories of primary metabolism using Mapman4 annotation framework [39]. In the case of photosynthesis, most genes are involved in light reactions, followed by the Calvin-Benson cycle and sucrose/starch metabolism (Figure S8). Cluster 1 genes are also involved in other aspects of primary metabolism such as nitrogen assimilation and amino acid biosynthesis as well as fatty acid biosynthesis (Figure S8). In the case of Cluster 2, and consistent with the marked induction in gene expression triggered by sulfate starvation shown in this cluster, we found an overrepresentation of biological processes related to the cellular response to S starvation such as 5'-adenylylsulfate reductase genes (Solyc02g032860, Solyc02g080640 and Solyc03g031620), RESPONSE TO LOW SULFUR (LSU) genes (Solyc03g096760, Solyc03g096770, Solyc03g096780) or sulfate transporters (Solyc05g054740, Solyc04g054730, Solyc04g072760, Solyc12g056930) (Fig. 5C and Table S7). Cluster 2 in leaves is enriched in GO terms related with the response to stress and protein phosphorylation (Fig. 5D and Table S8). Interestingly, we did not find an overrepresentation of S-starvation-related processes in Cluster 2 in leaves, as observed in roots. In fact, classical genes involved in the sulfate starvation response such as SULTR, ATPS, APR or LSU belong to Cluster 5 (Figure S7), suggesting an organ-specific expression of these genes in response to external S during development.
In summary, sulfate starvation promotes significant changes in the temporal expression patterns of a large number of genes in roots and leaves, of which many are associated with functions related to cellular growth, photosynthesis and sulfate uptake and metabolism. Moreover, the comparative analysis of major gene co-expression clusters between roots and leaves indicates an organ-specific response of the tomato plants to sulfate starvation.
Identification of transcription factors involved in the sulfate starvation response in tomato
Our results point at major transcriptomic changes occurring in tomato roots and leaves in response to sulfate starvation. As described in Arabidopsis, some of these changes may be attributed to transcriptional regulation by TFs. Although some TFs controlling sulfate responses have been identified in Arabidopsis [40], to date no information about sulfate-responsive TFs and their role in responses to sulfate starvation in tomato. In order to identify key TFs controlling the response to sulfate starvation of tomato, we constructed gene regulatory networks for roots and leaves [41] as explained in Methods.
In the case of roots, we obtained a gene regulatory network consisting of 24 TFs and 172 targets (Figure S9 and Table S9). TFs of the root regulatory network belong to 17 different families (Table S10), with MYB and bZIP TFs being the most abundant families (4 and 3 members, respectively). In the same way, we constructed a TF-target network using the identified sulfate-regulated genes in leaves which was composed by 45 TFs and 566 targets (Figure S10 and Table S11). Specifically, ERF, TCP and WRKY were the most abundant TF families in the leaves network with 7 and 5 members, respectively (Table S10). Notably, the most abundant TF families in the regulatory network of leaves were not shared with the root network (Table S11), suggesting organ-specific regulation of the transcriptional response to sulfate starvation.
We next focused our analysis on TFs whose predicted target genes have overrepresented biological functions (adjusted p-value <0.05), in order to correlate TF regulation to changes in functional processes. As shown in Fig. 6A and Fig. 6B, we found 3 TFs in roots and 10 TFs in leaves that potentially control genes involved in the sulfate starvation response, phytohormone biosynthesis and signaling, senescence, among others biological processes. Most connected factors in terms of outdegree were a MYB TF (Solyc11g071300) in roots and TCP 21 TF (Solyc03g006800) in leaves, showing 40 and 71 targets respectively (Fig. 6A and Fig. 6B). Notably, MYB target genes are involved in S utilization, while TCP21 ones are involved in photosynthesis (Fig. 6A and Fig. 6B, Table S12).
As mentioned above, our results discovered important changes in gene expression in response to sulfate starvation, suggesting the involvement of different sets of TFs acting in both roots and leaves. However, we found one TF shared by both root and leaf gene regulatory networks which is early induced in roots but later in leaves (Fig. 6C). This TF is member of the EIL family and corresponds to “Ethylene Insensitive 3-Like3” (EIL3), Solyc01g006650. The closest Arabidopsis gene in terms of amino acid identity is SLIM1 (Table S13), which is a key regulatory factor of the sulfate starvation response [17]. In addition, GO enrichment analysis showed that the identified tomato EIL3 TF might control the expression of a set of genes related to sulfate starvation, including classical sulfate-responsive genes such as LSUs as well as genes encoding enzymes involved in the reductive steps of sulfate assimilation (Table S14). Furthermore, we found that 60% of the tomato TF target genes are shared with targets predicted for Arabidopsis SLIM1, according to PlantRegMap, including sulfate-response genes such as LSUs, APRs and SDIs (Table S14). Overall, we identified several candidate TFs to mediate the transcriptome reprogramming that occurs during sulfate starvation response in roots and leaves, highlighting the possible key regulatory role of the EIL3 TF in tomato.