Transcriptome sequencing revealed molecular mechanisms underlying tolerance of Suaeda salsa to saline stress

The halophyte Suaeda salsa displayed strong resistance to salinity. Up to date, molecular mechanisms underlying tolerance of S. salsa to salinity have not been well understood. In the present study, S. salsa seedlings were treated with 30‰ salinity and then leaves and roots were subjected to Illumina sequencing. Compared with the control, 68,599 and 77,250 unigenes were significantly differentially expressed in leaves and roots in saline treatment, respectively. KEGG enrichment analyses indicated that photosynthesis process, carbohydrate, lipid and amino acid metabolisms were all downregulated in saline treatment, which should inhibit growth of S. salsa. Expression levels of Na+/H+ exchanger, V-H+ ATPase, choline monooxygenase, potassium and chloride channels were upregulated in saline treatment, which could relieve reduce over-accumulation of Na+ and Cl-. Fe-SOD, glutathione, L-ascorbate and flavonoids function as antioxidants in plants. Genes in relation to them were all upregulated, suggesting that S. salsa initiated various antioxidant mechanisms to tolerate high salinity. Besides, plant hormones, especially auxin, ethylene and jasmonic acid signaling transduction pathways were all upregulated in response to saline treatment, which were important to gene regulations of ion transportation and antioxidation. These changes might comprehensively contribute to tolerance of S. salsa to salinity. Overall, the present study provided new insights to understand the mechanisms underlying tolerance to salinity in halophytes.


Introduction
The halophyte Suaeda salsa (L.) Pall is a leaf succulent annual plant. It is widely distributed in Eurasia and has been cultivated as a seawater vegetable in desert and coastal areas in P. R. China [1]. S. salsa shows strong tolerance to salinity. In comparison to the control, treatments with moderate salinities (varied from 50 mM to 200 mM NaCl) promoted CO 2 assimilation rate, O 2 production, photosynthesis rate and thus plant growth of S. salsa [2][3][4][5][6]. However, high salinities (> 200 mM NaCl) inhibited its growth significantly [6], displaying much higher tolerance to saline stress than other plants. a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 3.75‰, 7.5‰, 15‰ and 30‰ (geometric sequence by two times), by dissolving commercial sea salts (Tangfeng Slats Company Limited, Tangshan, P. R. China) in deionized water. Deionized water was used as the control. Germination assays were conducted in 9 cm Petri dishes. In each dish, three layers of filter paper were placed at the bottom, 10 ml of saline solution was added and then 100 seeds were placed on filter paper. Seeds were germinated in a light humidified incubator at humidity of 80% RH and dark. The temperature was 25˚C at day time (12 h) and 15˚C at night (12 h). Germinated seeds were counted daily for 10 days to calculate germination rate. To avoid disturbance of evaporation on salinity, germinated seedlings were discarded and ungerminated seeds were transferred to new dishes with fresh solution every two days. Each assay was repeated four times.
Saline treatments of S. salsa seedlings S. salsa seeds were germinated in deionized water as described above. When seedlings reached 10 cm approximately, seedlings were transplanted to plastic containers with sands as cultural substrate and Hoagland's solution as media, and then cultured in a greenhouse at 25˚C. The light cycle was 16 hours: 8 hours (light: dark) and light intensity was approximately 14,400 lux. After two weeks, healthy plants were treated with salinity of 30‰ by dissolving commercial sea salts (Tangfeng Slats Company Limited, Tangshan, China) in Hoagland's solution. Hoagland's solution was used as the control. Each treatment included 20-30 individuals and repeated three times independently. To supplement evaporated water, distilled water was added every day to keep total volume of culture media consistent. After 10 days, plants were harvested. Roots and top four leaves were collected and immediately frozen in liquid nitrogen for RNA extraction. For each treatment, samples from 5 individuals were pooled.

Transcriptome sequencing
Total RNA was extracted using Biozol reagent (Bioer, Hangzhou, China) according to the manufacture's protocol. Quality of total RNA was check by 1% agarose gel electrophoresis, NanoPhotometer spectrophotometer (IMPLEN, CA, USA) and Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA). RNA samples with RNA integrity number (RIN) higher than 8.0 were considered qualified. RNA concentration was measured using Qubit RNA assay kit on Qubit 2.0 Flurometer (Life Technologies, CA, USA).
To construct sequencing libraries, mRNA was enriched using NEBNext Poly(A) mRNA Magnetic Isolation Module (NEB, USA). Sequencing libraries were prepared using NEBNext mRNA Library Prep Master Mix Set for Illumina (NEB, USA) and NEBNext Multiplex Oligos for Illumina (NEB, USA). Sequencing library was viewed on 1.8% agarose gel to check insert size and quantified using Library Quantification Kit-Illumina GA Universal Kit (Kapa, USA).
Qualified libraries were clustered on a cBot cluster generation system using HiSeq 4000 PE cluster kit (Illumina) and sequenced on an Illumina Hiseq 2500 platform. Three independent samples were sequenced as three biological replicates.

Bioinformatics analyses
Clean reads were achieved by removing adaptors, reads with N ratio higher than 1% and low quality reads (with > 50% bases having Phred quality score � 15). Clean reads were subjected to assembly of unigene library using Trinity v2.0.6 [25]. HTSwq v0.6.0 was applied to calculate FPKM values (expected number of fragments per kilobase of transcript sequence per millions base pairs sequenced) of each unigene. Relative expression levels of each gene were compared among different groups using DESeq2 R package. Fold change � 2 and q value (adjusted P value) < 0.001 were considered significantly changed.
DEGs were mapped to Gene Ontology (GO) database and Kyoto Encyclopedia of Genes and Genomes database (KEGG) [26] for enrichment of GO categories and KEGG pathways using BLAST software [27].

Real-time quantitative PCR
Results of Illumina sequencing were validated using real-time qPCR. cDNA was prepared using BioRT cDNA first strand synthesis kit (Bioer, Hangzhou, China) with oligo(dT) primer. qPCR was carried out using BioEasy master mix (Bioer, Hangzhou, China) on a Line Gene9600 Plus qPCR machine (Bioer, Hangzhou, China). DEGs and primers used for qPCR are listed in S1 Table. Transcriptional elongation factor-1 alpha was used as the internal control. Relative expression levels of each gene between saline treatment and the control were calculated using the typical 2 -ΔΔCt method [28]. Three biological replicates were included for each treatment and three technical repeats were performed for each sample.

Effects of salinity on seed germination
Previously, most experiments employed NaCl to represent saline stress, which might misestimate the effects of other components in sea salts on seed germination and plant growth [29][30]. In the present study, we used sea salts to prepare saline solutions, which should be more similar to real environments.
Seeds started to germinate in all treatments after two days. Germination rate significantly decreased in response to saline treatments but increased with germination time. After eight days, in comparison to the control, germination rates were significantly lower in treatments with 30‰ and 15‰, but did not significantly change in other saline treatments (Fig 1). These results suggested that germination of S. salsa seeds could tolerate salinity as high as 7.5‰, which was consistent with previous findings on S. salsa and the general consensus that germination of halophyte seeds is optimum under freshwater but decreases with elevating salinity [31,32].

Transcriptome sequencing
The sequencing data were deposited in NCBI with the reference number of PRJNA512222. Illumina sequencing resulted in 65.21 M to 72.28 M of total clean reads and 6.52 G to 7.23 G of total clean bases for each sample ( Table 1). The sequencing depth was roughly calculated by the ratio of sequencing data to the size of assembled transcriptome, approximately ranging from 22.37 × to 24.81 ×. The total clean reads were higher than those in previous studies on S. salsa [33]. Besides, all Q20 and Q30 values were higher than 98.09% and 91.78%, respectively (Table 1). These indices suggested that the sequencing was deep and the as-obtained data should be enough for further analyses.

De novo assembly
Total clean reads of each sample were assembled independently to get the unigene library of each sample. The total number of unigenes ranged from 40,038 to 133,430, with the mean length of unigenes from 798 bp to 1,180 bp (S2 Table). Benchmarking Universal Single-Copy Orthologs (BUSCO) analyses revealed that at least 2.97% BUSCOs were fragmented or missing (S1 Fig). Alternatively, clean reads of all samples were pooled, equal to 1,234.49 M of clean reads, and then subjected to Trinity assembly. Finally, 196,199 unigenes were obtained, with the mean length of 1,673 bp and the N50 length of 2,780 bp ( Table 2). BUSCO analyses revealed that only 1.98% BUSCOs were missing and none BUSCOs were fragmented (S1 Fig). Obviously, the as-obtained unigenes showed longer mean length and N50 length as well as less fragmented and missing BUSCOs than those assembled from each sample (S1 Table), suggesting that unigene library generated from all samples had higher quality. The as-obtained unigene library was used for further analyses.

Functional annotation and analysis
Among the unigene library, more than 65.45% and 45.91%unigenes were longer than 500 bp and 1,000 bp, respectively (S2

Differentially expressed genes and qPCR validation
Compared with the control, 68,599 (45,321 upregulated and 23,278 downregulated) and 77,250 unigenes (48,682 upregulated and 28,568 downregulated) were significantly differentially expressed in treatment with 30‰ in leaves and roots, respectively (Fig 2). To validate the expression levels predicted by transcriptome sequencing, qPCR was performed for 6 unigenes in leaves and 8 unigenes in roots. The results showed similar tendency to the expression levels calculated by FPKM value (Fig 3), suggesting that DEGs prediction was reliable.

Enrichment of GO categories for DEGs
Top 20 GO categories revealed that in leaves and roots, DEGs in comparison between 30‰ and control were mainly assigned to cellular process, metabolic process, cell, cell part, catalytic activity and binding (Fig 4).

Enrichment of KEGG pathways for DEGs
KEGG enrichment of significantly downregulated genes showed 24 and 16 pathways in leaves and roots, respectively, which shared nine KEGG pathways, including valine, leucine and isoleucine degradation (ko00280), linoleic acid metabolism (ko00591), glycerophospholipid metabolism (ko00564), circadian rhythm-plant (ko04712), plant hormone signal transduction (ko04075), pentose phosphate pathway (ko00030), riboflavin metabolism (ko00740), nicotinate and nicotinamide metabolism (ko00760) and ascorbate and aldarate metabolism (ko00053). These enriched pathways represented harmful effects of salinity on S. salsa. KEGG enrichment of significantly upregulated genes revealed 17 and 18 pathways in leaves and roots, respectively. Between them, only one pathway (ko00261 monobactam biosynthesis) was shared. These results might provide new insights for investigation of molecular mechanisms underlying adaptation of S. salsa to salinity stress (Tables 3 and 4).

Changes of energy metabolism in S. salsa upon saline treatment
Photosynthesis process is quite sensitive to saline stress. In response to saline treatments, photosynthesis process was inhibited in maize [34]and mulberry [35], which should negatively affect accumulation of organic matters and plant growth. Similarly, salinity also displayed negative effects on shoot length of S. salsa [36]. In the present study, carbon fixation in photosynthetic organisms (ko00710) was suppressed in saline-stressed leaves, suggesting that photosynthesis process might be inhibited by salinity, which then should reduce metabolism of organic carbons, displaying downregulation of carbohydrate metabolism (ko00500 starch and sucrose metabolism, ko00051 fructose and mannose metabolism, ko00010, glycolysis/gluconeogenesis, ko00052 galactose metabolism, ko00531 glycosaminoglycan degradation, ko00520 amino sugar and nucleotide sugar metabolism, ko00030 pentose phosphate pathway, ko00640 propanoate metabolism and ko00511 other glycan degradation). Furthermore, downregulation of carbohydrate metabolism would depress lipid metabolism (ko00564, glycerophospholipid metabolism, ko00591 linoleic acid metabolism, ko00071 fatty acid degradation, and ko00592 alpha-linolenic acid metabolism), amino acid metabolism (ko00250 alanine, aspartate and glutamate metabolism and ko00350 tyrosine metabolism) and biosynthesis of complicated polymer compounds (Tables 3 and 4). For example, phenylpropanoid biosynthesis (ko00940), which contributes monolignols, the starting compounds for biosynthesis of lignin [37], was downregulated in the present study. Overall, these changes would comprehensively depress S. salsa growth. These results were partially similar to transcriptomic analyses of S. maritima in treatment with 2% NaCl, which decreased the expression levels of genes related to photosynthesis (particularly the light reaction) and Calvin cycle, but increased most of the genes in the anabolic pathways (such as sucrose and starch synthesis) and catabolic pathways (such as the tricarboxylic acid cycle, glycolysis and the oxidative pentose phosphate pathway). The different changes of sucrose and starch synthesis and glycolysis between these two studies might be attributed to the tested salinities. In the present study, 30‰ salinity was used, which was much higher than that prepared by 2% NaCl.
The KEGG pathway riboflavin metabolism (ko00740) and nicotinate and nicotinamide metabolism (ko00760) is important for biosynthesis of vitamin B2 and B3, respectively. Vitamin B2 is an essential coenzyme to oxidordeuctases (such as succinodehydrogenase, xanthine oxidase and NADH dehydrogenase), participating in degradation of sugars. Vitamin B3 is a major component of coenzyme I (NAD) and coenzyme II (NADPH), which participate in the photosynthesis process, degradation of sugars and lipids. Thus, downregulation of these two pathways could be considered as normal regulation in response to saline-inhibited photosynthesis and accumulation of organic matters.

Changes of ion transportation
To resist saline stress, the first strategy is to relieve harmful accumulation of Na + and Cl -. In NaCl-treated S. maritima, plasma membrane H + ATPase (PM-H + ATPase), Na + /H + exchanger, Na + /H + antiporter, vacuolar-type H + ATPase (V-H + ATPase), betaine aldehyde dehydrogenase (BADH) and choline monooxygenase (CMO) were upregulated, which were involved in ionic adjustment [13]. Similarly, Na + influx transporter and the tonoplast Na + /H + antiporter were upregulated in Suaeda fruticose under saline condition [38]. In the present study, partially similar results were revealed (Table 5). Na + /H + exchanger, V-H + ATPase and CMO were significantly upregulated in leaves and roots. These results suggested that S. salsa might share similar mechanisms with other Suaeda species underlying ionic balance (Na + and H + ) in response to saline stress. The H + -ATPase generates an electrochemical membrane potential in the plasma membrane and Na + /H + exchanger may transport Na + out of cells via this H + gradient [39,40]. However, PM-H + ATPase and Na + /H + antiporter did not change in response to saline treatment in the present study and also in S. glauca [12], indicating that different Suaeda species might differ in regulation of proteins responsible for ionic balance. Besides, potassium channel, two-pore potassium channel, chloride channel protein, S-type anion channel were significantly upregulated in saline-treated leaves in the present study (Table 5). Similar results have been reported in Zostera marina in which various K + channels and transporters showed higher uptake capacity of K + in response to saline treatment than those in the control [41]. Maintenance of K + supply is a crucial feature of salt tolerance, since transportation of K + could effectively adjust osmotic pressure among subcellular compartments and cytosol caused by increased vacuolar Na + concentration [42]. Activity of anion channel is corresponded to chloride channel protein [43]. Upregulation of chloride channel proteins could effectively reduce Claccumulation and increase tolerance to salinity in plants [44].

Activation of antioxidant mechanisms in response to saline treatment
Generally, saline treatments induce oxidative stress to plants [45,46]. The halophyte S. salsa should have some mechanisms to avoid detrimental effects of oxidative stress. In the present study, we examined the expression level of superoxide dismutase (SOD), which is the first antioxidant enzyme in response to environmental stress. The results revealed different patterns between two types of SOD. In comparison to the control, expression level of Fe-SOD increased significantly for 3.4 and 4.8 times in saline-treated leaves and roots, respectively, but expression level of Cu/Zn-SOD did not change significantly (Table 5). Consistently, saline treatment drove greater increase of expression level of Fe-SOD than that of Cu-Zn SOD in rice varieties [47]. Although Fe-SOD and Cu/Zn-SOD have similar catalyzing functions, but their amino acid sequences appear to be unrelated [48]. Thus, they may respond to oxidative stress triggered by different stress factors. More investigations are required to distinguish their biological functions.
Comparison between saline treatment and the control suggested that the KEGG pathway glutathione (GSH) metabolism (ko00480) was upregulated in S. salsa leaves. Although GSH functions in nutrient metabolism and regulation of cellular events (including gene expression, DNA and protein synthesis, cell proliferation and apoptosis, signal transduction, cytokine production and immune response, and protein glutathionylation), the major role of GSH occurs in the antioxidant defense [49]. In response to saline treatments, GSH metabolism was activated in reed [50] and Arthrospira platensis [51]. Thus, activation of GSH metabolism might also contribute antioxidant capacity to protect S. salsa from saline-induced oxidation.
GSH is biosynthesized from glutamate, cysteine, and glycine [49]. In roots, saline treatment upregulated the KEGG pathway glycine, serine and threonine metabolism (ko00260) but downregulated alanine, aspartate and glutamate metabolism (ko00250). How the changes of these two pathways affected GSH metabolism still required more investigations.
The phenylpropanoid pathway was also upregulated in saline treatment (Table 3), which is the upstream pathway to flavonoid biosynthesis [60]. Thus, phenylpropanoids and flavonoids might function together to protect S. salsa from abiotic stresses [37,61]. Similar results were observed in various plants, such as rice [62], wheat [63] and S. glauca [12].

Regulation of plant hormones
The KEGG pathway plant hormone signal transduction (ko04075) is responsible to transduce signals of hormones (including auxin, cytokinine, gibberellin, abscisic acid, ethylene, brassinosteroid, jasmonic acid and salicylic acid) to downstream performance (including cell enlargement, elongation and division, germination, shoot imitation, plant growth and stress response). KEGG enrichment analyses of upregulated or downregulated genes independently both revealed significantly involvement of plant hormone signal transduction pathway in saline-treated leaves, suggesting that this pathway might play important and complicated roles in saline resistance. The FPKM values of different unigenes mapped to the same protein were summed and statistically re-analyzed by Students' T-tests (Fig 7 and S5 Table). Only one gene in the signaling transduction of cytokinine, gibberellin and salicylic acid was upregulated. Based on the changing tendency of only one gene, it was unprecise to conclude that these three hormone signals contributed to saline resistance in S. salsa. In the brassinosteroid signaling transduction, three genes significantly differentially expressed, but their changing tendency were contradictory. Compared with the control, BAK1 was upregulated for 1.7 times, BZR1/2 was upregulated for 1.5 times, but BRI was downregulated for 1.4 times in saline treatment. These changes were quite weak and could be ignored.
Importantly, the signaling transduction of auxin, abscisic acid, ethylene and jasmonic acid involved several genes which mostly showed similar tendency in the present study, suggesting that these signaling pathways might respond to saline treatment in S. salsa and will be discussed in details.

Potential roles of auxin signaling pathway in saline tolerance
The chemical nature of auxin is indole-3-acetic acid. In the present study, the KEGG pathway indole alkaloid biosynthesis (ko00901) was activated in saline-treated leaves. Within this pathway, the unigenes encoding methylesterase (EC 3.1.1.78) and GDSL esterase/lipase (EC 3.1.1.80) were both upregulated in saline treatment in comparison to the control ( Table 5), suggesting that concentrations of indole alkaloids might be elevated. Similar results were reported previously, that indole alkaloids were induced by abiotic stresses, including UV-B irradiation [64], PEG-induced drought [65] and salinity [66]. In the auxin pathway which transduces signals of indole-3-acetic acid, auxin response factor (ARF), auxin responsive GH3 (GH3) and small auxin upregulated RNA family protein (SAUR) were also upregulated for1.7, 2.8 and 2.2 times, respectively, in saline-treated leaves (Fig 7 and S5 Table). Taken together, these results suggested that auxin might promote tolerance of S. salsa to saline stress through auxin signaling pathway.

Changes of abscisic acid signaling in response to salinity
In the abscisic acid (ABA) signaling pathways, abscisic acid receptor (PYR/PYL), serine/threonine-protein kinase (SNRK2) and ABA responsive element binding factor (ABF) were downregulated in saline-treated leaves, compared with the control (Fig 7 and S5 Table), suggesting that ABA signaling pathway was inhibited by saline treatment. These results were different from those in non-halophytes but consistent with previous findings in a related halophyte species S. maritima. As previously reported, saline treatment increased ABA concentrations in non-halophytes tobacco [67], tomato [68], grapevine [69]. Exogenous ABA induced saline tolerant in wheat by elevating expression of peroxidase, which might reduce the active oxygen triggered by salinity [70]. Besides, ABA induced stomatal closure to minimize water loss [71]. However, in the halophyte S. maritima, saline treatments significantly decreased ABA concentration in comparison to the control [72]. Slight decrease of ABA concentration was also observed in NaCl-treated Prosopis strombulifera roots and leaves [73]. These results suggested that inhibition of ABA signaling might be a common strategy for halophytes to resist salinity. Jin et al. [12] also reported changes of expression levels of ABA related genes in saline treatment in S. glauca, but the authors did not explain the detailed mechanisms. The underlying mechanisms required further investigations, but one thing was clear. Decreased level of ABA signaling should ensure stomatal opening and water absorption, which are important for S. salsa growth.

Changes of ethylene signaling in response to salinity
In S. salsa leaves, four genes responsible for ethylene signaling were significantly unregulated in saline treatment compared with the control, including mitogen-activated protein kinase (MPK6), ethylene-insensitive protein 2 (EIN2), ethylene-insensitive protein 3 (EIN3) and ethylene responsive factor 1/2 (ERF1/2) (Fig 7 and S5 Table). These results were consistent with the upregulation of ERFs in saline-treated S. glauca [12] and also in accordance with the general understanding that ethylene signaling is indispensable for tolerance to saline stress in Arabidopsis and many other terrestrial plants [74].
Downstream effectors of ethylene signaling could include reactive oxygen species (ROS) scavengers (such as SOD) and ion transporters [74]. In the present study, expression levels of Fe-SOD and ion transporters (including Na + /H + exchanger, V-H + ATPase, CMO, potassium channel, two-pore potassium channel and chloride channel protein) were all upregulated in response to saline treatment, which might be mediated by the activation of ethylene signaling pathway.
Moreover, ubiquitin mediated proteolysis (ko04120) is a downstream effector pathway of EIN3. It has been reported that this pathway is involved in flooding, heat, ultravioletradiation, oxidative stresses and diseases in plants [75]. In the present study, ubiquitin mediated proteolysis was upregulated in saline-stressed leaves by KEGG enrichment analysis (Table 3), suggesting that ubiquitin mediated proteolysis process might be activated to degrade salinitydenatured or oxidation-denatured proteins [76]and/or inactivate the repressors in plant hormone regulation pathway [77]. This process might be initiated by ethylene accumulation in S. salsa.

Changes of jasmonic acid signaling in response to salinity
Jasmonic acid (JA) could enhance the tolerance of peanut [78], tomato [79]and barley [80]to saline stress by increasing antioxidant enzyme activities [78], lowering Na + concentration and improving stomatal conductance [80]. In the present study, three key genes in this pathway, including jasmonic acid-amido synthetase (JAR1), jasmonate ZIM domain-containing protein (JAZ) and transcription factor MYC2 (MYC2), were significantly upregulated in saline treatment compared with the control, suggesting that activation of JA signaling might improve saline tolerance in S. salsa.

Changes of cell wall-associated genes
In response to environmental stresses, modification of cell wall is a common defense mechanism. In both S glauca and S. maritima, genes related to cell wall dynamics, including laccase, expansins, leucine-rich repeat extensins (LRX), wall-associated receptor kinase proteins and/ or O-acyltransferase WSD1 were upregulated in saline treatment compared with the control. In the present study, expression levels of expansins, LRX and O-acyltransferase WSD1 did not change significantly between saline treatment and the control, but laccase, and another two wall associated proteins cellulose synthase and cellulose synthase interactive [81] were significantly upregulated in the saline treatment ( Table 5), suggesting that cell wall remodeling also took place in S. salsa, which might contribute to saline tolerance in S. salsa as reported in S. glauca and S. maritima [12,13].

Conclusions
In response to treatment with 30‰ salinity, S. salsa displayed suppressed photosynthesis process, carbohydrate, lipid and amino acid metabolisms, which might inhibit growth of S. salsa. To resist saline stress, Na + /H + exchanger, V-H + ATPase, choline monooxygenase, potassium and chloride channels were upregulated to ensure ionic balance. Expression levels of Fe-SOD and genes in glutathione metabolism, biosynthesis of L-ascorbate, flavonoids and phenylpropanoids were all upregulated to produce more antioxidants, thus relieving harmful effects of saline-induced oxidation. Plant hormones played essential roles in saline resistance in S. salsa. Importantly, auxin, ethylene and jasmonic acid signaling transduction pathways were upregulated and abscisic acid signaling transduction was inhibited by saline treatment, which might activate antioxidant mechanisms and ionic adjustment.  Table. FPKM values of key genes involved in biosynthesis of L-ascorbate. Data represent mean ± standard deviation (n = 3). � significantly different from the control (P < 0.05). (DOCX) S4 Table. FPKM values of genes involved in flavonoid biosynthesis. Data represent mean ± standard deviation (n = 3). � significantly different from the control (P < 0.05). (DOCX) S5 Table. FPKM values of genes involved in plant hormone signal transduction pathway. Data represent mean ± standard deviation (n = 3). � significantly different from the control (P < 0.05). (DOCX)