Integrated metabolomic and transcriptomic analysis reveal the effect of mechanical stress on sugar metabolism in tea leaves (Camellia sinensis) post-harvest

Plant materials and mechanical stress treatment

We selected the tea cultivar ‘Tieguanyin’ (Camellia sinensis var. Tieguanyin) and harvested it according to the standards of collecting one bud and three leaves. The tea leaves were collected in the garden of Fuqing Nanhushan Tea Co., Ltd. (25°71′N, 119°25′E, Fuzhou, China). The detached tea shoots were divided into two groups. The samples from the first group were treated with continuous turnover for 15 min to simulate mechanical stress treatment (V15); the second group was withered without vibrating stress (temperature 22.0 ± 0.4 °C, humidity 83.8 ± 2%) for 15 min (C15). The fresh leaves were designed as C0 used to control of C15 and V15 (Fig. 1). Three independent biological replicates were performed for each sampling and the second leaves were immediately collected and frozen in liquid nitrogen and then stored in a −80 °C refrigerator for further analsysis.

The mechanical stress test was performed on a mechanical vibrating Yaoqing bench. The technical parameters of the bench were set as follows: table size 1,100 × 1,200 mm, fixed amplitude 25 mm, power 2.2 KW, vibration frequency 330 r/min, maximum load 100 kg, timing control (0–60 min), maximum current 5 A, and voltage 220 V.

Widely-targeted metabolomics analysis

The metabolite extraction and analysis were conducted by Metware Biotechnology Co., Ltd. (Wuhan, China) (Chen et al., 2013). The tea samples were freeze-dried (−50 °C, 10 Pa, 12 h) in a Bisafer-10D vacuum freeze dryer (Saifei Biotechnology Co., Ltd., Nanjing, China) and ground to a powder using a mixing mill (MM 400; Retsch, Haan, Germany) with a zirconia bead for 1.5 min at 30 Hz. A total of 100 mg of powder from each sample was dissolved in 0.6 ml of 70% aqueous methanol solution containing 0.1 mg/l lidocaine as an internal standard and extracted overnight at 4 °C. The mixtures were centrifuged at 10,000 g for 10 min to obtain extracts that were absorbed and filtered for UPLC-MS/MS analysis. The quality control (QC) samples prepared by mixing all samples to check the stability of the instrument.

UPLC and ESI-triple quadrupole-linear ion trap (Q TRAP)-MS/MS conditions

The aqueous solution was analyzed using an LC-ESI-MS/MS system (HPLC: Shim-pack UFLC Shimadzu CBM30A; MS: Applied Biosystems 4500 Q-TRAP; Thermo Fisher Scientific, Inc., Waltham, MA, USA). The system was equipped with a Waters UPLC HSS T3 C18 column (2.1 × 100 mm, 1.8 µm; Waters Company, Milford, MA, USA) with 4 μL at 40 °C. Acetonitrile (0.04% acetic acid) was used as mobile phase B and acidified water (0.04% acetic acid) was used as mobile phase A. The following linear gradient were used for the eluent B: 5% (0 min), 0–95% (0–11 min), 95% (11–12 min), 95–5% (12–12.1 min), and 5% (12.1–15 min). Further, an ESI-triple quadrupole-linear ion trap (QTRAP)-MS detector (Thermo Fisher Scientific, Inc., Waltham, MA, USA) was used to connect the effluent alternatively. Analyst 1.6.1 software (AB SCIEX Pte. Ltd, Framingham, MA, USA) was used for data acquisition, system, and processing. The ESI operating parameters were: turbo spray; ion source; source temperature 550 °C; ion spray voltage (IS) 5,500 V; ion source gas I, 50 psi; ion source gas II, 60 psi; curtain gas, 25 psi; and the collision activated dissociation (CAD) was set at “high”. Quality calibration and instrument tuning were carried out using 10 and 100 μmol/L polypropylene glycol solutions in QQQ and LIT modes, respectively. Decluster potential (DP) and collision energy (CE) were optimized for individual MRM jumps. Monitoring of specific MRM transition sets by metabolites was eluted at each time period (Chen et al., 2013).

The metabolites were identified by comparing their retention times, the mass-to-charge (m/z) values, and the fragmentation patterns with those of the authentic standards or by searching against the public databases (MassBank, KNApSAcK, HMDB, MoTo DB, and METLIN). The characteristic ions of each substance were screened by the triple quadrupole, and the intensity of each characteristic ion was detected and collected. Then the chromatographic peak integration and calibration were calculated using MultiQuant software to obtain the relative contents of each substance in the samples (Fraga et al., 2010).

RNA isolation, cDNA library preparation, and RNA-seq

Ten leaves from each sample were randomly selected for total RNA extraction. Total RNA was extracted using the Trizol kit (TaKaRa, Dalian, China) according to the manufacturer’s protocol. The sequencing library was prepared using the method of Huang et al. (2020a). After a quality check, the construction of cDNA libraries was carried out according to the protocol for the NEBNext® UltraTM RNA Library Prep Kit for Illumina® (NEB, Ipswich, MA, USA). The sequencing library preparations were created on the Illumina HiSeq 2000 platform (Illumina, San Diego, CA, US) and paired-end reads were determined (2 × 150 bp). The poly-N, adapter and low-quality reads (Qphred ≤ 20) were removed to obtain the clean data and calculate the Q20, Q30, and GC content of the clean reads (Data S1) (Zheng et al., 2022).

Gene annotation and differential expression

The filtered sequences were subjected to the tea plant (Camellia sinensis cv. Suchazao) reference genome (Wei et al., 2018) localization analysis using Hisat2 v2.0.4 software (Anders, Pyl & Huber, 2015). The FPKM cacluation for each gene used HTSeq v0.9.1 as a basis for gene expression identification. Transcriptome sequence data are deposited in the NGDC database (GSA accession: CRA007810 and PRJCA011050).

Differential expression analysis was performed on biological replicate samples using the DESeq2 package (Love, Huber & Anders, 2014). Readcount was normalized using DESeq and adjusted for p-values using the Benjamini and Hochberg method to control for false discovery rates. All the genes with an adjusted P ≤ 0.05 as determined by DESeq2 were considered to be differentially expressed. Finally, we predicted the pathways and functions of the differentially expressed genes and made heatmaps using TBtools (Chen et al., 2020).

Quantitative RT-PCR analysis

In order to validate the transcriptome data, 13 DEGs were randomly selected and their expression levels were detected using quantitative RT-PCR (qRT-PCR). The primers were designed using NCBI Primer Blast and listed in Table S1. GAPDH (glyceraldehyde-3-phosphate dehydrogenase) was used as the housekeeping gene (Zhou et al., 2019). The total RNA extraction, cDNA synthesis, and quantitative RT-PCR analysis were performed according to our previous study (Yue et al., 2019). The qRT-PCR was performed using SYBR Premix Ex Taq™ II (TaKaRa, Dalian, China) and determined on a CFX96 Touch Real-Time PCR System (Bio-Rad, Hercules, CA, USA). The relative expression of the genes was calculated using the 2^−ΔΔCt method (Cao et al., 2020). Detection was performed in triplicated for each DEG.

Variations of sugar metabolites under mechanical stress

In this study, widely targeted metabolomics analysis was performed on C0, C15, and V15. As a result, a total of 19 sugar compounds were identified (Fig. 2 and Table S2). The majority of the sugar metabolites accumulated in C15 and V15.

Compared with C0, sugars including N-Acetyl-D-glucosamine (GlcNAc), glucose-1-phosphate (Glu-1P), D-fructose-6-phosphate (D-Fru-6P), D-galactose, D-glucose-6-phosphate (D-Glu-6P), panose, D-sedoheptulose-7-phosphate, D-sucrose, galactinol, D-melezitose O-rhamnoside, D-melezitose, maltotetraose, and melibiose were significantly increased in C15 (P < 0.05). The sugars of GlcNAc, Glu-1P, fructose 6-phosphate-disodium salt, D-Fru-6P, D-galactose, D-Glu-6P, panose, D-sedoheptuiose 7-phosphate, DL-arabinose, and D-mannose were significantly increased in the V15 treatment. Additionally, compared with C15, D-sucrose, galactinol, D-melezitose-O-rhamnoside, D-melezitose, maltotetraose, melibiose, and sucralose were decreased in V15. However, GlcNAc, Glu-1P, fructose 6-phosphate-disodium salt, D-Fru-6P, D-galactose, D-Glu-6P, panose, D-sedoheptuiose 7-phosphate, glucosamine, D-glucose, DL-arabinose, and D-mannose were increased, indicating that the sugar metabolites were regulated by mechanical stress treatment during post-harvest processing.

Differential expression analysis of genes related to the sugar metabolism pathway under mechanical stress

To further explore the transformation mechanism of these sugar metabolites, we performed transcriptome sequencing of the samples. There was a high correlation between the gene expression levels of each sample (R² > 0.95) (Fig. 3). A total of 13 DEGs were randomly selected for qRT-PCR detection showing that their exression patterns were similar to the RNA-seq results (Fig. S1), indicating good reproducibility and reliability of the RNA-seq data of the tested samples. As showed in Fig. 4, a total of 69 differentially expressed genes related to sugar metabolism pathways were obtained from transcriptome sequencing (FC ≧ 1.2). These genes are associated with the starch and sucrose, galactose, and glycolysis pathways of the sugar metabolism pathway.

Figure 5 shows that the expression of genes related to the starch and sucrose metabolism pathways in tea leaves were down-regulated during post-harvest treatments, particularly under mechanical stress. The down-regulated genes in V15 included invertase (INV), trehalose 6-phosphate synthase (TPS), trehalose-6-phosphate phosphatase (TPP), α-trehalase (TRE), phosphoglucomutase (PGM), glucose-6-phosphate isomerase (GPI), starch synthase (SS), 1, 4-α-glucan branching enzyme (GBE), BAM, and AMY when compared to C15. Correspondingly, the sugar metabolites of UDP-glucose, Glu-1P, and D-Glu-6P were increased in tea leaves in response to mechanical stress. In contrast, the up-regulation profiles of hexokinase (HXK), UDPglucose 6-dehydrogenase (UGPase), and 4-α-glucanotransferase (DPE) were associated with the increase of D-Fru-6P, UDP-glucose, and D-glucose in V15.

In the glycolysis metabolism pathway (Fig. 6), the expression levels of HXK, 6-phosphofructokinase 1 (PFK-1), phosphoenolpyruvate carboxykinase (PCK), and pyruvate kinase (PK) genes were up-regulated in V15 compared with C15. The up-regulated expression of these genes reflected the activation of glycolysis in V15. Moreover, mechanical stress decreased phosphoenolpyruvate (PEP) levels in leaves. PEP is an intermediate during glycolysis and its level is inversely proportional to the activity of PFK-1. Additionally, decreasing of PEP could enhance the activity of PFK (Turner & Plaxton, 2003).

The key DEGs and metabolites involved in the galacoste metabolism pathway were identified. As showed in Fig. 7, the genes, including aldose 1-epimerase (AEP), α-galactosidase (AGAL), UDP-glucose 4-epimerase (GALE), UDP-glucose-hexose-1-phosphate uridylyltransferase (GALT), and INV were down-regulated in V15 when compared to C15, and the levels of D-galactose, D-mannose, and D-sorbitol were increased in V15. However, the expression levels of UGPase, galactinol synthase (Gols), and raffinose synthase (RS) were up-regulated in V15, which was related to the increase of D-Glu-6P, UDP-glucose, and D-galactinol levels.

Adaptive strategies of tea leaf metabolome to mechanical stress

The carbohydrates in plants act as energy and substrate sources for a series of physiological activities. These are mainly derived from the degradation of starch under the catalysis of amylase when facing environmental stressors (Shen et al., 2015; Saddhe, Manuka & Penna, 2020). The abiotic stress-induced sucrose accumulation in Oolong tea leaves was regarded as a response to adverse environmental conditions. Under these conditions, sucrose is further hydrolyzed into hexose, including glucose and fructose derivatives, via the catalysis of INV, sucrose synthase (SUS), and HXK, demonstrating a progressive trend in the tea manufacturing process (Wu et al., 2022). In this study, the expression of BAM genes were dramatically up-regulated after indoor-withering and mechanical stress-enriched glucose levels in V15. Our results are similar to the reports that BAM were involved in the response to abiotic stressors during plant development (Ma et al., 2022). The expression of HXK, UGPase, and DPE2 genes were higher in V15, and Glu-1P, D-Glu-6P, and D-Fru-6P amount in V15 were more abundant than those of other groups (Fig. 2). This suggests that the monosaccharides accumulated in tea leaves to deal with the stressors, such as dehydration or the mechanical stress from Oolong tea manufacturing. Furthermore, GlcNAc, fructose 6-phosphate-disodium salt, D-galactose, panose, DL-arabinose, and D-mannose were most abundant in V15 because the compounds are elements of the biosynthesis of cellulose and hemicellulose, and are spontaneously involved in the response to mechanical wounding (Scheller & Ulvskov, 2010; Kadokawa et al., 2015). In addition, our data showed that most of the sugars in Oolong tea increased under mechanical stress, which improved the tea aroma and taste. This is consistent with previous studies that turning-over could improve the aroma and taste of Oolong tea (Zeng et al., 2020).

Transcriptome analysis associated with sugar metabolism pathways