Comparative transcriptome analysis reveals the function of SlPRE2 in multiple phytohormones biosynthesis, signal transduction and stomatal development in tomato

Transcriptomic, physiological, and qRT-PCR analysis revealed the potential mechanism by which SlPRE2 regulates plant growth and stomatal size via multiple phytohormone pathways in tomato. Paclobutrazol resistance proteins (PREs) are atypical members of the basic/helix-loop-helix (bHLH) transcription factor family that regulate plant morphology, cell size, pigment metabolism and abiotic stress in response to different phytohormones. However, little is known about the network regulatory mechanisms of PREs in plant growth and development in tomato. In this study, the function and mechanism of SlPRE2 in tomato plant growth and development were investigated. The quantitative RT-PCR results showed that the expression of SlPRE2 was regulated by multiple phytohormones and abiotic stresses. It showed light-repressed expression during the photoperiod. The RNA-seq results revealed that SlPRE2 regulated many genes involved in photosynthesis, chlorophyll metabolism, phytohormone metabolism and signaling, and carbohydrate metabolism, suggesting the role of SlPRE2 in gibberellin, brassinosteroid, auxin, cytokinin, abscisic acid and salicylic acid regulated plant development processes. Moreover, SlPRE2 overexpression plants showed widely opened stomata in young leaves, and four genes involved in stomatal development showed altered expression. Overall, the results demonstrated the mechanism by which SlPRE2 regulates phytohormone and stress responses and revealed the function of SlPRE2 in stomatal development in tomato. These findings provide useful clues for understanding the molecular mechanisms of SlPRE2-regulated plant growth and development in tomato.


Introduction
In nature, plant growth and development are precisely manipulated by endogenous phytohormones and environmental stress-induced physiological and biochemical responses. Generally, these responses are regulated by activating or inhibiting the expression of specific genes. Transcription factors (TFs) act as key regulators and interact Communicated by Rachel Wells. with the cis-acting elements of responsive gene promoters to maintain the normal growth and development of plants. Among these transcription factors, the basic/helix-loop-helix (bHLH) family is a large family that controls phytohormone and stress-regulated gene responses in plants.
According to recent studies, bHLH family TFs regulate the expression of many genes involved in phytohormone responses, light responses and multiple plant growth and development processes, such as seed dormancy and germination, flowering time determination, cell fate determination, and stress response processes (Hao et al. 2021;Qian et al. 2021;Sanagi et al. 2021). The bHLH family genes Phytochrome-Interacting Factor 3 (PIF3) and PIF4 regulate light signaling, photomorphogenesis and chlorophyll metabolism in Arabidopsis (Job and Datta 2021;Xu and Zhu 2021). GLABRA3 (GL3) mediates root hair formation as a component of the WER-GL3-TTG1 transcriptional complex (Qiu et al. 2021). The bHLH transcription factor MYC2 acts as a key regulator in the jasmonate response and mediates the plant responses to biotic and abiotic stresses (Du et al. 2017;Verma, et al. 2020). In tomato, SlbHLH95 controls the fruit ripening process and impacts trichome formation Zhang et al. 2020). The bHLH TFs regulated these processes by binding to the cis-elements of target promoters, such as G-box (CAC GTG ) and E-box (CANNTG), and thus manipulated the expression of target genes (López-Vidriero et al. 2021;Yang et al. 2021).
In bHLH family TFs, the basic domain functions in DNA binding to conserved cis-elements of target genes while the HLH domain participates in the formation of a homodimer or heterodimer to regulate the expression of target genes (Toledo-Ortiz et al. 2003). The bHLH TFs include a subfamily of atypical bHLH TFs known as Paclobutrazol Resistances (PREs), since these bHLH genes confer resistance to paclobutrazol-induced plant growth inhibition in Arabidopsis ). These atypical bHLH proteins lack the basic domain of bHLH proteins and primarily function in the competitive inhibition of their interaction partners by forming heterodimers (Hornitschek et al. 2009;Seo et al. 2011). The functions of PRE genes have been reported in several species. In Arabidopsis, initial studies revealed the role of AtPREs in gibberellin (GA), brassinosteroid (BR) and light signaling pathways that regulate plant growth via their interaction with bHLH proteins (Bai et al. 2012;Hao et al. 2012;Hyun and Lee 2006;Lee et al. 2006;Wang et al. 2009). Further studies suggested their function in the auxin response, abiotic stress responses and floral organ development (Shin et al. 2019;Zheng et al. 2017Zheng et al. , 2019. The Increased Lamininar Inclination1 (OsILI1) gene, a homologue of PREs in rice, positively regulates the rice leaf angle in response to BRs (Zhang et al. 2009). Recently, the role of PRE in carotenoid metabolism, receptacle development and abiotic stress responses was demonstrated in Satsuma mandarin (citrus), Fragaria × ananassa Duch. (strawberry) and Malus domestica (apple) (Endo et al. 2016;Li et al. 2022a;Medina-Puche et al. 2021. In tomato, five PRE gene homologues exist in its genome, among which SlPRE2 was first investigated as an important regulator of leaf morphology, internode elongation and fruit development (Zhu et al. 2017(Zhu et al. , 2019. SlPRE5 overexpression plants showed rolling of mature leaves, elongated leaf petioles and reduced chlorophyll contents, similar to the phenotype of SlPRE2 overexpression lines (Li et al. 2022a). However, the SlPRE2 overexpression plants showed reduced water loss rate and stomatal aperture, but their young leaves had elevated water loss rate compared with the wild-type (WT) plants (Zhu et al. 2017), which was not reported in the SlPRE5 overexpression lines. Thus, further investigation of the function and the regulatory mechanisms of SlPREs, especially the SlPRE2, in tomato plant growth and development are still needed. In this study, the expression of SlPRE2 in response to multiple phytohormones and abiotic stresses was investigated. The RNA-seq technique was employed to identify the genes and metabolic and signaling pathways that were regulated by SlPRE2. In addition, the role of SlPRE2 in stomatal development was studied. The obtained data will be a valuable resource for exploring the molecular mechanisms of SlPRE2 regulated tomato plant growth and development.

Expression analysis of SlPRE2 under exposure to phytohormones and abiotic stresses
In recent years, PREs have been reported to be involved in hormone signaling pathways, as first reported in the gibberellin pathway and also found in the BR, abscisic acid (ABA) and auxin pathways Zhang et al. 2009;Zheng et al. 2017Zheng et al. , 2019. SlPRE2 shows a GA-inducible expression pattern, and overexpression of this gene leads to a plant growth and development phenotype that matches the function of gibberellin in tomato (Zhu et al. 2017(Zhu et al. , 2019. In this study, the expression activity of SlPRE2 under GA treatment was reconfirmed (Fig. 1A). SlPRE2 showed increased expression after 3, 6, 9 and 24 h of gibberellic acid (GA 3 ) application to tomato leaves. These results are consistent with the induced expression response of SlPRE2 in GA 3 -treated tomato fruit as previously reported (Zhu et al. 2019). Additionally, indole-3-acetic acid (IAA) and methyl jasmonate (MeJA) treated samples showed elevated SlPRE2 expression level ( Fig. 1B and C). However, the expression of SlPRE2 could be suppressed by ABA and 1-aminocyclopropane-1-carboxylate (ACC) treatment at different time points ( Fig. 1D and E). The responsiveness of SlPRE2 to osmotic abiotic stress was also investigated. SlPRE2 showed reduced expression in both leaves and roots when tomato plants were treated with PEG and sodium chloride ( Fig. 1F-I). These results indicated that SlPRE2 could be regulated by multiple phytohormones and abiotic stresses.

Expression of SlPRE2 during the light cycle
AtPRE6/AtKDR was first reported to show increased expression during the daytime, and SlPRE2 showed high expression under low light conditions and lower expression under high light conditions (Hyun and Lee 2006;Zhu et al. 2017). Since PREs show interaction activity with the light signaling regulators HFR1 and PAR1 (Hao et al. 2012;Hong et al. 2013), it appeared that SlPRE2 was involved in Fig. 1 Expression profiles of SlPRE2 in response to phytohormones and stresses. qRT-PCR was used to investigate the expression of SlPRE2 gene. A-E the expression of SlPRE2 under GA 3 , IAA, MeJA, ABA, and ACC phytohormone treatments, respectively. F-I the expression of SlPRE2 in tomato leaves and roots under PEG 6000 and salt treatment, respectively. For each time point, the relative expression of SlPRE2 in treated samples was compared with the control samples at the same time point. Each value represents mean ± SD from three biological replicates Fig. 2 The expression curve of SlPRE2 in 16 h light and 8 h dark growth condition. The leaf samples with similar growth stage were collected from 1-month-old plants at a series of timepoint for two days. Values are mean ± SD from three biological replicates the light response based on its expression under different light conditions. In this study, the expression of SlPRE2 at a series of time points in 16 h light and 8 h dark growth conditions were investigated. qRT-PCR analysis showed that SlPRE2 expression first decreased firstly and then increased during 16 h of light exposure (Fig. 2). Under the dark growth conditions, the expression of SlPRE2 decreased slightly and then increased until the light period began. A previous study showed that SlPRE2 exhibited a reduced expression level under low-light growth conditions compared with that under high-light conditions (Zhu et al. 2017), consistent with the altered expression of SlPRE2 observed between light and dark growth conditions in this study. SlPRE2, which shows differences expression during the light cycle, might modulate light-regulated plant development via a similar mechanism of AtPRE interaction with other bHLH factors (Hao et al. 2012;Hyun and Lee 2006).

PRE2OE leaves showed opened stomata in young leaves
In a previous study, the overexpression of SlPRE2 led to the rolling of mature leaves in tomato (Zhu et al. 2017). They showed a reduced water loss rate and stomatal closure, which were speculated to be strategies for regulating water loss in the mature leaves of SlPRE2-overexpression lines (PRE2OE). However, the young leaves of PRE2OE plants presented a higher water loss rate than those of the WT, which suggested that PRE2OE may present an opposite stomatal state in young leaves. In this study, the abaxial epidermis of young PRE2OE leaves were peeled and observed using a microscope. In contrast to the relatively closed stomata of mature leaves, the young PRE2OE leaves showed widely opened stomata compared with the stomata of the WT leaves (Fig. 3A). They presented a reduced stomatal pore length and increased stomatal pore width (Fig. 3B). The width/length ratio was higher in young PRE2OE leaves (Fig. 3C). These results indicated that the elevated water loss rate was related to the opened stomata in young PRE2OE leaves.
Transcriptome sequencing, de novo assembly, and annotation of tomato genes in PRE2OE and wild-type plants Total RNA was isolated from the leaves of 2-months-old PRE2OE and WT tomato. RNA sequencing was performed on the BGISEQ-500 platform. For these samples, an average of 6.42 Gb of total reads with an average of 92.47 Q30 bases were acquired (Table S1). Adaptor and low-quality sequences were trimmed. The preliminary assembly showed that 94.71-95.31% of the clean reads were mapped to the Solanum lycopersicum reference genome, and 83.55-84.36% of the clean reads were mapped to Solanum lycopersicum reference genes (Table S2). Then, the expression of each library was quantified as normalized FPKM (fragments per kilobase per million mapped) values. Principal component analysis (PCA) of transcript differences showed significant differences among the WT and PRE2OE samples (Fig. 4A). The high similarity among the biological replicates from the WT and PRE2OE samples demonstrated that the RNAseq results were consistent. A total of 23,622 genes were assembled in these sequenced samples (Fig. 4B, Table S3). To obtain a comprehensive understanding of transcript expression related to the expression of SlPRE2, differentially expressed genes (DEGs) between the PRE2OE lines and C stomatal aperture (ratio of stomatal width:length) in epidermal peel of young leaves from 1-month-old PRE2OE and the wild-type plants. Results are shown as means ± SD from three independent experiments (n > 300). The statistic mean, SD and t test P value were shown in corresponding figures the WT were identified. Hierarchical cluster analysis was used to estimate the variance in DEGs, and 964 significantly DEGs [adjusted P value (Q value) ≤ 0.05 and |log 2 FC|≥ 1] were identified in the PRE2OE plants compared with the WT plants ( Fig. 4C, Table S3). Among the DEGs, 417 genes were upregulated while 547 genes were downregulated in PRE2OE plants (Fig. 4D). The maximally upregulated and downregulated genes showed log2FC values of 7.89 and − 4.17, respectively (Table S3). Therefore, the reliability of the RNA-Seq data was confirmed by the consistency between the RNA-Seq analysis and qRT-PCR results (Fig. 4E).

Transcription factor and promoter analysis of DEGs
Based on gene family analysis, 55 transcription factors with significant differences were identified in the PRE2OE lines. They included ARF, bHLH, bZIP, C2H2, G2-like, GRAS, GRF, MADS-box, MYB, NAC, NF-YB, TCP, WRKY and YABBY family genes (Fig. 5A, Table S4). They may function in plant development, fruit ripening, light responses and stress responses. SlPRE2 was highly expressed in PRE2OE lines (log2FC: 7.89). SlPRE1/STYLE2.1 (Solyc02g087860), a close homologue of SlPRE2, and SlPIF7a (Solyc03g115540), which is involved in the light response, both showed downregulated expression in PRE2OE lines (Rosado et al. 2016). The plant growth and development-related genes Solyc04g081240 (ARF5), Solyc12g010410 (homeobox protein knotted-1 like), Solyc02g062960 (HD-ZIP) and Solyc06g005090 (LOB domain-containing protein) were more highly expressed in PRE2OE than in the WT. Additionally, several genes of the ERF and NAC families exhibited altered expression in PRE2OE plants. For example, SlERF5 was upregulated (Log2FC: 1.17) in the PRE2OE line. SlERF5 has been reported to be regulated by multiple abiotic stresses and positively controls plant tolerance to drought and salt stress (Pan et al. 2012). These results suggested that SlPRE2 might regulate plant development and plant adaptation to environmental stress by controlling the expression of these transcription factors. Furthermore, the promoter sequences of 964 DEGs were analyzed. In the promoters (2000 bp upstream of start code ATG) of these genes, a large number of transcription factor-binding sites (TFBSs) were identified (Fig. 5B, Table S5). In these TFBSs, many bHLH family TFBSs were identified,

GO enrichment and KEGG pathway analysis of DEGs
To gain insight into the functional categories of the DEGs induced by the expression of SlPRE2, a GO enrichment analysis was performed to investigate the distribution of DEGs in the major biological process (BP), cellular component (CC) and molecular function (MF) categories. The most enriched terms in BP, CC and MF were metabolic process (GO:0008152, 329 DEGs), membrane (GO:0016020, 347 DEGs) and catalytic activity (GO:0003824, 438 DEGs), respectively (Fig. 6A). The cellular process (GO:0009987, 287 DEGs), cell (GO:0005623, 334 DEGs) and binding (GO:0005488, 380 DEGs) categories were also highly enriched in PRE2OE lines compared with WT. Detailed results of the GO enrichment analysis are shown in Table S6.
To investigate the pathways of the DEGs regulated by SlPRE2, the DEGs were aligned to Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Table S7). The DEGs identified upon the overexpression of SlPRE2 in tomato were enriched in 109 metabolic pathways including 83 "metabolism" pathways, 17 "genetic information processing" pathways, 4 "environmental information processing" pathways, 3 "cellular processes" pathways and 2 "organismal systems" pathways. Among the 109 pathways, 17 pathways had P value ≤ 0.05. Among these pathways, the most significantly enriched pathways were KEGG pathway enrichment analysis of the DEGs "photosynthesis-antenna proteins", "steroid biosynthesis", "photosynthesis", "carbon fixation in photosynthetic organisms", "arachidonic acid metabolism", and "fatty acid degradation" (Q value ≤ 0.05). These results suggested a role of SlPRE2 in the photosynthesis and fatty acid metabolism pathways. However, "carbon metabolism", "glutathione metabolism", and biosynthesis processes of phytohormones such as brassinosteroid, cytokinin and gibberellin were also enriched in PRE2OE lines (Fig. 6B, Table S7). In addition, the PRE2OE lines showed DEGs that were involved in several pathways primarily focused on amino acid metabolism, the biosynthesis of secondary metabolites, and carbohydrate metabolism. These results indicated that the constitutive expression of SlPRE2 had a major effect on multiple plant growth and development processes.

Expression of genes involved in photosynthesis and chlorophyll metabolism
In previous studies, PRE family genes have been reported to function in phytohormone and light response processes (Castelain et al. 2012;Hyun and Lee 2006;Hong et al. 2013). In tomato, SlPRE2 showed a light intensity-regulated expression pattern, and the overexpression of SlPRE2 resulted in a lower chlorophyll content and curled leaves (Zhu et al. 2017). These results suggested that SlPRE2 may regulate plant photosynthesis or light responses. Consistent with this speculation, the transcriptome analysis showed that several genes involved in photosynthesis had significantly reduced expression in PRE2OE lines compared with the WT (Fig. 7). In PRE2OE lines, the expression of Ferredoxin (FD) (Solyc10g075160), which plays roles in the chloroplast electron transport chain, was significantly decreased compared with that in the WT. Several genes encoding photosystem I reaction center subunits had lower expression level in PRE2OE leaves. Additionally, the expression of genes involved in chlorophyll biosynthesis, including light-dependent NADPH-protochlorophyllide oxidoreductase (LPORs, Solyc07g054210, Solyc10g006900, Solyc12g013710) and geranylgeranyl reductase (GGR , Solyc03g115980) genes, was significantly decreased in PRE2OE leaves compared with the WT, consistent with the reduced chlorophyll content in PRE2OE lines (Fig. 7) (Zhu et al. 2017). Furthermore, the Stay-Green (SGR) protein (Solyc08g080090), which plays a critical role in chlorophyll degradation, presented significantly increased expression in PRE2OE plants, while the chlorophyllase coding genes Solyc12g005300 and Solyc09g065620, which putatively catalyze the hydrolysis of phytol from chlorophyll showed decreased expression. These results indicated the altered metabolism of chlorophyll via the expression of SlPRE2 in tomato.

Genes involved in carbohydrate metabolism
The organic substances produced by photosynthesis are the material basis for the survival of plants, and the catabolism of carbohydrates in respiration releases energy for plant growth and development. Here the overexpression of SlPRE2 reduced the chlorophyll content in tomato leaves (Zhu et al. 2017), and photosynthesis-and chlorophyll metabolism-related genes also showed altered expression in these plants (Fig. 7). This implied that the carbohydrate metabolism process may have changed since the altered carbon assimilation from the photosynthesis process. Based on the KEGG and GO pathway enrichment analysis, many genes associated with carbohydrate metabolism had significant differential expression in PRE2OE lines. In the glycolysis/gluconeogenesis pathway, 12 genes showed significantly changed expression (Fig. 8). Among these genes, the genes encoding phosphoglycerate kinase (Solyc07g066610) and phosphopyruvate hydratase (Solyc10g085550) had lower expression levels while the genes encoding aldehyde dehydrogenase (Solyc05g005700) and pyruvate decarboxylase (Solyc09g005110) had higher expression levels in PRE2OE, which suggested the downregulation of the glycolysis metabolism pathway. Additionally, in the TCA cycle, Solyc12g011000 (citrate synthase), which performs the condensation of acetyl-CoA with oxaloacetate to form citrate, had an increased expression level in the PRE2OE plants. However, Solyc01g106480 (a malate dehydrogenase coding gene), which catalyzes the NAD/NADH-dependent interconversion of the substrates malate and oxaloacetate had a reduced expression level. In addition, the Solyc01g110360 (fructose-bisphosphate aldolase 1, FBA1), Solyc02g062340 (fructose-bisphosphate aldolase 1, FBA2), Solyc02g069620 (ribose-5-phosphate isomerase, RPI), and Solyc10g086720, which are involved in the pentose phosphate pathway, had lower expression in PRE2OE than in the WT. These results indicated that the overexpression of SlPRE2 affected the glycometabolism of tomato plants.

Genes involved in phytohormone biosynthesis and signaling pathways
Phytohormone-manipulated plant growth and development can be regulated by a series of transcription factors, among which the bHLH family genes played an important role. Previously, the functional mechanism of PREs was reported to be involved in the gibberellin, brassinosteroid and light signaling pathways, especially in the cell expansion process (Zhang et al. 2009;Zhu et al. 2017Zhu et al. , 2019. However, recent research has revealed that PREs participate in the auxin, ABA and abiotic stress responses (Zheng et al. 2017(Zheng et al. , 2019. Therefore, the effect of SlPRE2 on the transcription of genes related to phytohormone metabolism and signaling was analyzed via the comparison of genes in PRE2OE lines and the WT. In the dataset, more than 40 genes that were involved in phytohormone metabolism and signaling responses showed significantly altered expression under the overexpression of SlPRE2 (Table 1). These phytohormone pathway-related genes were enriched in different pathways, including the gibberellin, brassinosteroid, and cytokinin biosynthesis and signaling pathways. In the gibberellin biosynthesis process, Solyc03g006880 (encoding GA20ox) and three gibberellin 2-oxidaseencoding genes had reduced expression level in PRE2OE lines. Solyc03g006880 is reported to negatively respond to exogenous GA 3 treatment (Livne et al. 2015). Therefore, the reduced expression of this gene in PRE2OE plants further supported its role in gibberellin biosynthesis (Zhu et al. 2019). In the brassinosteroid biosynthesis pathway, several genes involved in steroid biosynthesis had lower expression, and DET2, a key enzyme in brassinosteroid Fig. 7 The expression profiles of DEGs related to photosynthesis, photosynthesis-antenna protein, and porphyrin and chlorophyll metabolism in the PRE2OE and the wild-type. The scale bar denotes the value of log10(FPKM+1). The blue and red colors in the heatmap represent repressed and induced expression level biosynthesis, had a reduced expression level in PRE2OE lines. However, Solyc02g084740 and Solyc02g085360, which encode CYP90C1 and CYP90B3, showed higher level of expression in PRE2OE plants. Additionally, downstream genes of the BR signal pathway, such as Solyc03g093120, Solyc03g093130 and Solyc03g093110, had reduced expression level, suggesting that the PRE2OE plants may have lower BR biosynthesis and signal responses than the WT plants. In addition, several genes involved in the cytokinin pathway had decreased expression, and Solyc12g008900 (encoding cytokinin oxidase, which catalyzes the degradation of cytokinin) had increased expression in PRE2OE plants, indicating that the constitutive expression of SlPRE2 may inhibit cytokinin biosynthesis. Furthermore, several genes involved in the IAA, ABA and SA signaling pathways also showed altered expression, even though changes in their biosynthesis were not observed. Overall, these results indicated the functions of SlPRE2 in multiple phytohormone responses, especially the gibberellin, brassinosteroid and cytokinin response pathways.

Potential protein interaction network of the DEGs
Atypical bHLH transcription factors lacking the basic domain could be small interfering peptides that affect the transcriptional activity of their interaction partners (Seo et al. 2011;Hong et al. 2013). The PRE proteins are reported to interact with proteins of the bHLH family (such as PAR1, AIFs, and IBH1), MADS-box family, NF-YC family, and HSF family in Arabidopsis and rice (Bai et al. 2012;Hao et al. 2012;Ikeda et al. 2013;Trigg et al. 2017;Zhang et al. 2009). In tomato, SlPRE2 showed interaction activity with SlPAR1 and SlAIFs in a yeast two-hybrid system (Zhu et al. 2019). Among these reported PRE interaction partners, only Solyc01g096370, a bHLH family gene homologue to At4g00870, had a reduced transcription level in PRE2OE lines (log2FC − 0.9965, Q value 0.0002), even though 50 homologues of PRE interaction partners were predicted (Table S8). To analyze the potential protein association networks within the SlPRE2 regulated pathway, the STRING protein-protein interaction network (Version 11.5, https:// cn. string-db. org/) was employed. This approach has been  (Aamir et al. 2017;Cirauqui et al. 2018;Mahadevan et al. 2016). The 964 DEGs identified in PRE2OE vs. WT were used as input for interaction network prediction at the high confidence level of 0.700, employing active interaction sources of text mining, experiments and databases. The results showed that protein interactions involved in carbohydrate metabolism, photosynthesis, and diterpenoid and steroid biosynthesis were overrepresented, consistent with the reported functions of PREs in gibberellin, brassinosteroid and light responses in plant species (Fig. 9). All the protein nodes and their predicted functional partners are provided in Table S9.

SlPRE2-regulated genes involved in stomatal size
According to the transcriptome analysis, four genes involved in stomatal development showed reduced expression in PRE2OE plants (Fig. 10A). The homologues of Solyc03g093940, an MYB transcription factor family gene, including AtMYB60 and AtMYB44, function in stomatal movement in Arabidopsis. The overexpression of these two genes increased stomatal closure (Jung et al. 2008;Oh et al. 2011). Solyc01g098190 (SlNHX4) and Solyc01g096740 (SlTOM), which encode sodium cation antiporter and metal chelator transporter, were predicted to regulate stomatal closure and ion transmembrane transport (Ali et al. 2021;Keisham et al. 2018;Kravchik et al. 2022). Solyc08g066610 (SlEPFL9), which encodes an epidermal  (Lu et al. 2019;Niwa et al. 2013). Among these four stomatal development-related DEGs, Solyc03g093940 and Solyc01g096740 showed high expression in roots, and Solyc08g066610 had high expression in young and mature leaves (Fig. 10B). These results suggested the tissue specificity of these genes in SlPRE2 regulated stomatal development/size and indicated that Solyc08g066610 may play a more crucial role in these processes.
In addition, to validate the transcriptomic data, several DEGs were chosen for qRT-PCR analysis in PRE2OE and WT plants (Fig. 10C-E). The selected genes were involved in phytohormone biosynthesis and signaling, photosynthesis and chlorophyll metabolism, carbohydrate metabolism, and stomatal development. The qRT-PCR results were consistent with the transcriptomic data (Fig. 4E).

Discussion
PREs participate in the regulation of the transcriptional activity of many genes and therefore influence plant development and stress responses in plant species (Bai et al. 2012;Endo et al. 2016;Li et al. 2022b; Medina-Puche et al. Fig. 10 The expression of stomatal development-related genes in SlPRE2 overexpression lines and verification of selected genes by qRT-PCR. A the log2FC of four genes related to stomatal development in PRE2OE compared with the wild-type. B the expression profile of stomatal development-related genes in different tissues and organs. C relative expression of genes related to phytohormone biosynthesis and signaling. D relative expression of genes related to photosynthesis and chlorophyll metabolism, carbohydrate metabolism. E relative expression of genes related to stomatal development. Hypo, hypocotyl; Cotyl, cotyledons; Meri, vegetative meristems; YL, young leaves; ML, mature leaves; Yfb, young flower buds; AF, anthesis flowers; 10DPA, 10 days postanthesis; 20DPA, 20 days post-anthesis; RF, ripening fruit 1 3 2021; Zhang et al. 2009). In this study, the overexpression of SlPRE2 altered the expression of 964 genes (Fig. 4). The promoter analysis of these DEGs suggested that SlPRE2 may activate or inhibit downstream genes through transcription factors of bHLH family (Fig. 5). The SlPRE2regulated downstream genes participated in pathways such as phytohormone biosynthesis and signaling, light signaling, and photosynthesis. In Arabidopsis, the six PRE subfamily genes participate in the gibberellin, brassinosteroid and light signaling pathways by controlling genes involved in these pathways. Recently, AtPREs were shown to regulate AtABI3, AtDREB2A and the flower development-related genes AtARGOS, AtIAA19, AtACS8, and AtMYB24 (Shin et al. 2019;Zheng et al. 2019). Citrus CubHLH1, a homologue of PREs, modulates the expression of carotenoid biosynthetic genes, especially the SlHYb and SlCCD1, in tomato fruit (Endo et al. 2016). High-throughput sequencing revealed the function of FaPRE1 in cell size and fruit ripening in strawberry (Medina-Puche et al. 2019. These results indicated the function of SlPRE2 in plant development through a complicated transcriptional regulatory mechanism. PREs were first reported to regulate plant growth and cell size via the gibberellin, brassinosteroid and light signaling pathways (Hao et al. 2012;Hyun and Lee 2006;Lee et al. 2006;Zhang et al. 2009). In previous studies, SlPRE2 was shown to function in tomato plant development in several ways. SlPRE2 overexpression lines grow fast and show rolling of mature leaves, longer shoots and reduced pigment contents (Zhu et al. 2017). However, the silencing of SlPRE2 altered seed size in tomato fruit, consistent with the relatively high expression of this gene in developing tomato fruit (Zhu et al. 2019). In addition, the SlPRE2 transgenic plants all showed altered sensitivity to exogenous gibberellic acid. Similar to the function of SlPRE2 in tomato, PRE genes function in plant development in several species by modulating multiple phytohormone pathways. AtPRE6 is a transcriptional repressor that negatively regulates the auxin response (Zheng et al. 2017). A previous publication reported that AtPREs were involved in the plant responses to ABA and salinity (Zheng et al. 2019). In rice, ILI1 overexpression leads to leaf inclination and acts downstream of brassinosteroids by directly activating OsBZR1-regulated transcription (Zhang et al. 2009(Zhang et al. , 2014. The CubHLH1 protein showed DNA binding activity against the promoter regions of carotenoid biosynthetic genes and promoted ABA accumulation (Endo et al. 2016). Strawberry FaPRE1 regulates ripening-related genes via auxin and ABAmodulated expression activation, while it does not respond to gibberellin (Medina-Puche et al. 2019). These results extended the understanding of the function of PREs in plant development.
In tomato, transcriptome analysis revealed that SlPRE2 regulated genes related to carbohydrate metabolism, photosynthesis and chlorophyll metabolism. Many genes involved in GA, BR and ABA biosynthesis and signaling were revealed in PRE2OE plants (Table 1, Fig. 10C). Additionally, the expression of SlPRE2 was regulated by exogenous GA 3 , IAA, ABA, MeJA, ACC and also abiotic stresses including drought and salt stresses (Fig. 1). These results suggested the broad participation of SlPRE2 in tomato plant development through these phytohormone and stress response pathways. In addition, overexpression of SlPRE5, a close homologue of SlPRE2 in tomato, resulted in changes in chlorophyll metabolism, photosynthesis, and gibberellin biosynthesis and signaling pathways, which is highly consistent with the function of SlPRE2 (Li et al. 2022a). These resulted further explained the conserved role of SlPREs in tomato development.
In these SlPRE2-regulated phytohormone metabolic and signaling processes, ABA plays a crucial role in the plant response to abiotic stresses such as drought and salinity Li et al. 2020;Vonapartis et al. 2022). Under environmental stress conditions, ABA biosynthesis is enhanced, and the elevated ABA content induces seed dormancy, inhibits plant growth, minimizes water loss from transpiration, and induces stress-related gene expression (Bharath et al. 2021;Nakashima and Yamaguchi-Shinozaki 2013). These response processes facilitate plant adaptation to the stress environment. The PRE2OE plants showed rolling of mature leaves and a lower water loss rate, while the young leaves had a higher water loss rate than the WT leaves (Zhu et al. 2017). It was speculated that the rolling of leaves and closing of stomata might be strategies for retaining water during plant growth. Consistent with this prediction, the young leaves of PRE2OE had widely opened stomata which may be responsible for the higher water loss rate in these young leaves (Fig. 3). PREs are generally thought to be involved in the gibberellin and brassinosteroid response pathways. However, some studies have revealed the relationship between PREs and ABA, suggesting the participation of PREs in the ABA-mediated plant stress response (Medina-Puche et al. 2019;Zheng et al. 2019).
As mentioned above, the expression of SlPRE2 affected the leaf water loss rate and stomatal size, especially in young leaves. Several genes involved in the biosynthesis and signal transduction of GA, BR and ABA, which may function in stomata, were observed to present altered expression in PRE2OE plants (Table 1). Considering that GA, BR and ABA all regulate stomatal development and movement (Kim et al. 2018;Li et al. 2020Li et al. , 2021Nir et al. 2017), SlPRE2 may act as a positive regulator of stomatal size by modulating these phytohormone pathways. On the other hand, consistent with the opened stomata, four genes related to stomatal development exhibited reduced expression in PRE2OE plants (Fig. 10). Among these genes, SlEPFL9 showed specific expression in young leaves and mature leaves, suggesting that SlEPFL9 may be an important response gene in SlPRE2-regulated plant development and that its function in this process should be further investigated.
In summary, the expression of SlPRE2 was affected by multiple phytohormones and abiotic stresses, and SlPRE2 might be involved in these phytohormone-and stressregulated plant growth processes. The overexpression of SlPRE2 changed the expression of genes involved in the biosynthesis and/or signal transduction of several phytohormones. Photosynthesis, chlorophyll metabolism and carbohydrate metabolism were also affected. The stomatal size of young leaves and several genes related to stomatal development were shown to be regulated by the overexpression of SlPRE2. It is likely that the SlPRE2 acts as a regulator of young leaves development and stomatal behavior, and subsequently leads to the rolling mature leaves with reduced water loss rate as the adaptation of plants to environment. Examination of the function of SlPRE2 using an inducible gene expression system and the search for its interacting partners will be helpful for further understanding the role of the regulatory mechanism of SlPRE2 in tomato plant development and abiotic stress response.

Plant materials and growth conditions
Wild-type (WT) Solanum lycopersicum Mill. cv. Ailsa Craig (AC) was used in this study. The SlPRE2 overexpression lines PRE2OE#2 and PRE2OE#8 have been described (Zhu et al. 2017(Zhu et al. , 2019. All plants were grown in greenhouse with 16 h light (27 °C) and 8 h dark (24 °C) cycles and were well watered with fertilizer every three days. All of the experimental samples were collected from plants, immediately frozen in liquid nitrogen and stored at − 80 °C for total RNA isolation.

RNA-seq and data analysis
The total RNA from tomato leaf samples was isolated using RNAiso Plus (Takara) following the manufacturer's instructions and purified using the RNeasy MinElute kit (Qiagen). RNA quality and concentration were measured using 1% agarose electrophoresis, a NanoDrop 2000 microspectrophotometer and an Agilent 2100 RNA Nano 6000 Assay Kit. RNA samples were sent to BGI Genomics for library construction and Illumina sequencing. The samples were sequenced on the BGISEQ-500 platform in 150 bp paired-end runs. High-quality reads were obtained by removing low-quality reads, adapter-contaminated reads, and poly-N reads using SOAPnuke software. All clean reads were obtained and transcriptome assembly was performed using StringTie. Gene functions were annotated using the tomato genome annotation (version 4.0) in the SGN database (https:// solge nomics. net/). PlantTFDB (http:// plant tfdb. gaolab. org/ index. php) was used to annotate transcription factors.
For differentially expressed gene (DEG) analysis, the FPKM method was used to analyze the expression level of each transcript, and a Q value ≤ 0.05 and |log2FC|≥ 1 were used to identify significant differences in gene expression. Gene Ontology (GO) enrichment and KEGG pathway enrichment analysis of the DEGs were performed using KOBAS (http:// kobas. cbi. pku. edu. cn/) and DAVID Bioinformatics Resources (https:// david. ncifc rf. gov/ home. jsp).

Phytohormone and abiotic stress treatments
For phytohormone treatments, 4-weeks-old WT tomato seedlings were treated with 50 μM GA 3 , 10 μM IAA, 50 μM MeJA, 50 μM ABA, and 50 μM ACC. Ethanol (0.05%) was used as a mock control. Leaf samples were collected after treatment for 0 h, 3 h, 6 h, 9 h, 12 h and 24 h. For abiotic stress treatments, 4-weeks-old tomato seedlings were treated with 20% polyethylene glycol 6000 (PEG 6000), and 200 mM sodium chloride (NaCl). Tomato plants with normal growth were used as controls. Leaf and root samples were collected after treatment for 0 h, 3 h, 6 h, 9 h, 12 h and 24 h. Three individual seedlings were collected for one sample, and three biological replicates were collected for each treatment at each time point. Samples were frozen in liquid nitrogen and stored at − 80 °C for total RNA isolation.

RNA isolation and reverse transcription and real-time PCR
For the verification of DEGs in the RNAseq results, the total RNA of frozen samples was isolated using the RNA Easy Fast Plant Tissue Kit (Tiangen, Beijing, China). cDNA biosynthesis was conducted by using M-MLV reverse transcriptase (Promega, Beijing, China). The synthesized cDNA was diluted at a 0.25 dilution ratio and utilized as a working template. The quantitative real-time PCR (qRT-PCR) reaction consisted of 10 μL Fast SYBR qPCR Mix (GenStar, Beijing, China), 1 μL of a 10 mM forward and reverse primer mixture, 2 μL cDNA template, and RNasefree water up to a final volume of 20 μL. qRT-PCR was performed in the CFX96 touch Real-time system (Bio-Rad, USA) with the following program: 95 °C for 3 min, followed by 40 cycles of 95 °C for 15 s, and 60 °C for 15 s, and a melting curve analysis during the 60-95 °C melt. The SlCAC and SlEF1α genes were used as reference genes for relative expression analysis (Expósito-Rodríguez et al. 2008). The