3.1 Physiological characteristics of under temperature stress
HTS disrupts the balance of ROS scavenging and lead to oxidative stress (Cuypers et al., 2010; Liu et al., 2021). Therefore, physiological characteristics of Medicago sativa L. under three kind of temperature stresses were determined. The results showed that antioxidant enzyme activities changed and the balance of ROS scavenging was destoryed under all three stresses. CAT (Fig. 1c), POD (Fig. 1d), SOD (Fig. 1e), APX (Fig. 1f) activities significantly increased under T30, T35, T40, and the CAT activity increased higher than other enzymes (Fig. 1c). The results of uperoxide anion content measurement showed that the content under the three stress conditions was higher than the control, and increase at first and then descend with the temperature increasing (Fig. 1a). The accumulation MDA decreased under all three stresses was (Fig. 1b), was negatively correlated with enzyme activity. Endogenous hormone plays a significant role in plant stress response, this study determined six endogenous hormones. The results showed that the endogenous hormone fluctuated under different treatments, the content of CTK (Fig. 1g), IAA (Fig. 1h), ABA (Fig. 1i), GA (Fig. 1k) was higher than the control significantly, especially in T35 and T40 treatment. ETH (Fig. 1j) increased slowly under three treatment, but till lower than CK. With the temperature gradient, the content of BR (Fig. 1l) increased, decreased and then increased, and the changes were significant higher than CK under T40 treatments. The above evidence indicated that the responses to three temperature in alfalfa, the membrane damage caused by T40 stress is more serious, while the activity of enzyme and content of endogenous hormone were highest than other, this can effectively clear up reactive oxygen species, and alfalfa was protected to better adapt to temperature stress.
3.2 Transcriptomic analysis of alfalfa during high temperature stress
In order to understand the gene expression of alfalfa in responses to temperature stress, we performed RNA-seq analysis. High quality reads were obtained. After quality control filtration, the average Q20 was 98.23%, and the average Q30 was 95.11%. The mean GC content was 44.88% (Table S1), indicating that the obtained data were useful for further analysis. The trends of separation obtained from the principal component analysis (PCA) of the samples reflected the degree of differences among them. The PCA results of the RNA-seq samples showed that PC2 separated the control group (CK) from the treatment groups (T30, T35 and T40), and PC1 separated the treatment groups (Fig. 2a). DEGs were identified under T30, T35 and T40 stress to explore the similarities and differences at the transcription level in response to the three stresses in alfalfa. In total, 3418, 5483 and 4341 DEGs were identified in alfalfa during T30, T35 and T40 stress, respectively (Fig. S1). Among them, 1876, 2815 and 2115 genes were up-regulated, while 1542, 2667 and 2226 genes were down-regulated (Fig. 2c). In addition, 1810, 2889 and 2410 genes were specifically differentially expressed under T30, T35 and T40 treatment, respectively, while 544 genes were differentially expressed after all treatments (Fig. 2b). Collectively, these data showed that more genes were up-regulated than down-regulated in response to T30 and T40 stress, more down-regulated genes in response to T35 stress.
To further explore the DEGs involved in the responses to different temperature stress, Gene Ontology (GO) enrichment analysis was performed. All (up and down genes) DEGs in T30_ vs_ CK, T35_ vs _CK, T40_ vs_ CK treatment groups were catagorized into three GO classification (including 228 subclassification) according their functional annotations, specifically, biological process (BP), molecular function (MF) and cellular component (CC) (Fig. S2). In addition, most DEGs appeared in the classification of BP and MF groups, and most DEGs belong to MF group gathered into “binding” and “catalytic activity” subgroups. Under T30, T35 and T40 stress, DEGs significantly enriched in the biosynthetic process, organic substance biosynthetic process, and metabolic process.
KEGG enrichment analysis of EDGs in three comparison groups was performed, there were 59, 102, and 80 KEGG pathways identified in T30 vs CK, T35_ vs _CK, and T40_ vs _CK comparison group. The pathway of oxidative phosphorylation, glycerolipid metabolism, thermogenesis, glycerophospholipid metabolism, tropane, piperidine and pyridine alkaloid biosynthesis were significantly enriched under T30_ vs _CK (Fig. 2d), the pathway of carbon fixation in photosynthetic organisms, glyoxylate and dicarboxylate metabolism was significantly enriched under T35_ vs _CK (Fig. 2e), the pathway of plant hormone signal transduction, carbon fixation in photosynthetic organisms, glyoxylate and dicarboxylate metabolism, photosynthesis-antenna proteins, ascorbate and aldarate metabolism were significantly enriched under T35_ vs CK_ (Fig. 2f).
3.3. Metabolomic analysis of alfalfa during high temperature stress
The positive(POS)and negative(NEG)ion modes were usec for qualitative and quantitative evaluation to improve the coverage of metabolites and detection effect. The PCA results of metabolic samples were consistent with the RNAseq results, which revealed high similarity among the three biological replicates under each treatment. There were differences in the samples analyzed in the four cases: PC2 separated the temperature treatment from the other samples, and PC1 separated the T40 sample samples from the other samples (Fig. 3a). Under comprehensively unbiased non-targeted metabolomic profiling, a total of 2121 metabolites were identified by positive (1223) and negative (898) ion patterns, and grouped into 16 classes (Table. S3). According to the metabolite classification results, there were 7 types of compounds that account for more than 80% of the total metabolites in leaves, including lipids and lipid-like molecules (22.68%), organic acids and derivatives (22.23%), organoheterocyclic compounds (13.20%), benzenoids (12.82%), organic oxygen compounds (7.26%), organic nitrogen compounds (1.84%), nucleosides, nucleotides and analogues (1.51%).
The overall distribution of differential metabolites in each comparison group is shown in Fig. 3c. There were 247 differential metabolites in T30_vs_CK comparing group, of which 173 were significantly up-regulated and 73 were significantly up-regulated, respectively. T35_vs_CK comparison group had 245 differential metabolites, among which 188 and 57 were significantly up- regulated and down-regulated, respectively. There were 286 differential metabolites in the T40_vs_CK comparing group (220 significantly up- regulated and 66 down- regulated) (Fig. 3b).
Among these metabolites, 76(49 up-regulated and 27 down-regulate)metabolites were differentially expressed under three temperature stress, and 97, 71 and 126 DEMs were specifically expressed under T30, T35 and T40 stress, respectively (Fig. 3c). In general, the number of up-regulated metabolites was higher than that of down-regulated metabolites, especially under T40 stress. The DEMs identified during including lipids and lipid-like molecules, benzenoids, organic acids and derivatives, organoheterocyclic compounds, organic oxygen compounds, organic nitrogen compounds, Nucleosides, nucleotides and analogues.
The KEGG enrichment analysis revealed that differential metabolites identified in T30 _vs_CK, T35_vs_CK, and T40_vs_CK treatment groups significantly enriched in11, 23, and 16 pathways related to plant metabolism (Table S3, S4, S5). The significantly enriched KEGG pathways for differential metabolites shared by three treatment groups include ABC transporters; plant hormone signal transduction; alanine, aspartate and glutamate metabolism; glyoxylate and dicarboxylate metabolism; citrate cycle (TCA cycle); glycine, serine and threonine metabolism. Consequently, temperature stress significantly changed some metabolites generation in alfalfa leaves and markedly modified the enrichment route of these metabolites. In addition, we found that another 9 pathways, including protein digestion and absorption; mineral absorption, pantothenate and CoA biosynthesis, carbon metabolism, cysteine and methionine metabolism, galactose metabolism, cAMP signaling pathway, beta-Alanine metabolism, and phenylpropanoid biosynthesis were all noteworthily altered after the leaf was exposed to different temperature. Therefore, temperature treatments influence the biosynthesis of some secondary metabolites, such as sucrose, betaine, biotin, 2-cis-4-trans-abscisic acid, salicylic acid, trans-zeatin, glyceric acid, isocitrate, succinate, glutamic acid etc. In conclusion, these bioactive substances and changes in their enrichment pattern suggested that alfalfa adapts to temperature stress by activating the biosynthesis of some secondary metabolites.
3.4. Integrated analysis of transcriptomic and metabolomic expression levels
Integrated metabolomic and transcriptomic methodology was used to assess the potential associations between DEGs and DEMs to reveal the pathways through which they are jointly involved. There were 30, 33 and 34 transcripts and metabolites common shared under T30, T35, T40 treatments (Fig. 4a-c). The top 10 KEGG pathways with the most participation of gene and metabolite identified in this experiment are illustrated in Fig. 4, and the results showed that the KEGG pathways in which differential metabolites and genes were jointly involved were glycerolipid metabolism (ko00561), glycerophospholipid metabolism, Glyoxylate and dicarboxylate metabolism(ko00630), Plant hormone signal transduction(ko04705), respectively (Fig. 4). Several DEGs and the differential metabolites are involved in all three pathways, although the DEGs did not directly regulate the production of these different metabolites.
Furthermore, the results indicated that the DEGs were significantly enriched in oxidative phosphorylation and glycerolipid metabolism, while the DEMs were significantly enriched in the glyoxylate and dicarboxylate metabolism and Plant hormone signal transduction pathway duringT30 stress (Table. S6). Under T35 stress, DEGs were significantly enriched in the glyoxylate and dicarboxylate metabolism and glycerolipid metabolism pathway, and the DEMs were enriched in the citrate cycle (TCA cycle) and plant hormone signal transduction pathway (Fig. 4e). During T40 stress, both DEGs and DEMs were significantly enriched in the plant hormone signal transduction and Glyoxylate and dicarboxylate metabolism pathway.
3.5 Analysis of shared or specific DEGs and DEMs under different temperature stress
According to the integrated analysis of transcriptomic and metabolomic, the top 10 KEGG were analysed under T30, T35 and T40 stress. A total of 74 DEGs were shared in top 10 KEGG (Fig. 5), among which 32 and 42 DEGs were up-regulated and down-regulated under T30, 34 DEGs up-regulated and 40 down-regulated under T35, 15 DEGs up-regulated 59 down-regulated under T40. The 74 DEG involved 16 KEGG pathway, among which 28 DEGs belongs to plant hormone signal transduction(ko04075),15 DEGs belongs to glyoxylate and dicarboxylate metabolism (ko00630), 10 DEGs belongs to glycerolipid metabolism(ko00561༉, 9 DEGs belongs to carbon fixation in photosynthetic organisms༈ko00710༉. A total of 41 metabolites were differently expression under temperature stress in top 10 KEGG(Fig. 6a), which 10 metabolites shared under the T30, T35 and T40 stress treatments (Fig. 6b).
In addition, KEGG results of the DEMs and DEGs showed that DEMs were associated with glyoxylate and dicarboxylate metabolism (ko00630), plant hormone signal transduction(ko04075), arginine and proline metabolism(ko00330), glycine, serine, and threonine metabolites༈ko00260༉, cysteine and methionine metabolism༈ko00270༉. Interestingly, among the pathways related to DEGs, the plant hormone signal transduction and glyoxylate and dicarboxylate metabolism were also shared with the DEMs, bute the other pathways were different. Specifically, DEGs and DEMs were also related to phenylalanine metabolism, oxidative phosphorylation, and Ascorbate and aldarate metabolism. The results demonstrated that under high temperature stress, DEGs participated in the synthesis of stress-related secondary metabolites, which may reduce stress damage through the antioxidant pathway.
The 18 specific DEMs indicated differences among the responses of alfalfa to different temperature stress (Fig. 6c). A total of 8 metabolites (abscisic acid, salicylic acid, L-hydroxyarginine, phosphocreatine,.alpha.-keto-.gamma.-(methylthio) butyric acid1, betaine, phenylalanine, dexpanthenol) specifically differentially expressed involed 6 KEGG pathway during under T30 (Fig. 6c), including the arginine and proline metabolism(Ko00330), pantothenate and CoA biosynthesis(ko00770), protein digestion and absorption(ko04974). There were 2 DEMs (glycerophosphocholine, Acetylcholine) specifically expressed involved in glycerophospholipid metabolism (ko00564) under T35 stress. There were 8 DEMs(Jasmonic acid, glutamic acid, D-proline, glycine, L-aspartic acid, L-threonine, L-methionine, D-galacturonic acid) were involved in plant signal transduction(ko04075), glyoxylate and dicarboxylate metabolism(ko00630), cysteine and methionine metabolism༈ko00270༉, ascorbate and aldarate metabolism(ko00053), arginine and proline metabolism༈ko00330) during T40 stress. These specific DEGs were involved in the response to high temperature stress by participating in various pathways to induce the synthesis of specific metabolites.
3.5.1 Plant hormone signal transduction
Zeathin biosythesis, caroterniod biosythesis, Phenylpropanoid biosythesis, α-linolenic acid biosynthesis and brassionoteriod biosynthesis were significant changed in the plant hormone signal transduction pathway(Fig. 7a). The metabolite cytoskinine were significantly upregulated in zeanthin biosythesis, and the gene CRE1(MS.gene68824, MS.gene65393) was downregulated. The metabolite jasmoruic acid were significantly upregulated in the α-linolenic acid biosynthesis, with the gene JAR1(MS.gene22809) and COI1 were upregulated, and JAZ(MS.gene89156; MS.gene09703; MS.gene54190; MS.gene21515; MS.gene26683) was downregulated(Fg.7b).The metabolite aslicylic acid were significantly upregulated in the Phenylpropanoid biosynthesis, with gene TGA (MS.gene29028) was downregulated (Fig. 7b). In addition, the gene AUX1 (MS.gene 90982; MS.gene33487) was downregulated in Tryptophan biosythesis, the gene PYR/PYL4 (MS.gene33704; MS.gene33451; MS.gene67442; MS.gene95667) and BZR1/2(MS.gene59331) were upregulated(Fig. 7b).
3.5.2 Glyoxylate and dicarboxylate metabolism
In the metabolism pathway of glyoxylate and dicarboxylate, the citrate, Cis-aconitate(EINECS), I-socitrate, succorate and L-Glutanate were significantly upregulated, while D-glycerate, glutamine synthetase༈GS), ribulose-bisphotosphate carboxylase large chain༈RuBisCO༉, D-glycerate 3-kinase༈GLYK༉ were remarkably downregulated (Fig. 8a). The identified genes involved in encoding EINECS, GS, RubisCO, and GLKY are shown in the Fig. 8. The significant upregulation of EINECS was beneficial to the increment of I-socitrate and succorate, by contrast, the conspicuously down-regulated GS was not conducive to the generation of L-glutamine. (Fig. 8a). Over all, the evident down-regulation of RuBisCO and GLYK was disadvantageous to gettingD-glycerate and 3-phospho-D-glycerate, respectively.
3.6 qRT -PCR verification of RNA–seq data
To validate the RNA-seq results, randomly selected 9 DEGs for qPCR analysis. The results showed that the expression trends of these genes obtained by the two technologies were basically consistent, indicating the validity of the RNA-seq results (Fig. 9).
The blue bar represents FPKM value in RNA-seq data, the orange bar represents the qRT–PCR results and was calculated using the 2−ΔΔCt method.