Phenotypic of jujube seedlings post HS
We found a putative heat-resistant jujube cultivar (‘JHR17’) in Turpan city in China, and inbred the cultivar in our laboratory. To obtain an overview of heat-tolerance phenotype of jujube cultivar, seedlings with 14 true leaves were subjected to treatment with HS (45 °C). 0 (control) 1, 3, 5 and 7 days, after the treatment, the leaves of all seedlings did not become withered (Fig. 1A), suggesting this jujube cultivar might be of heat-tolerance.
Stoma is an important channel for gas and water exchange between plants and the atmosphere, which can make adaptive adjustment under various stress conditions. The stomatal density and stomatal opening rate of leaves from each group were assessed. It indicated that the stomatal density and stomatal opening rate post heat treatment were significantly increased, and they showed a trend of first increase and then decrease with the extension of heat treatment (Fig. 1B, Table 1). It suggested that the ‘JHR17’ could reduce the damage by passively changing stomatal density and opening rate.
Table 1
Effect of high temperature stress on stomatal density and stomatal opening rate of seedling
| heat treatment time(d) |
0 | 1 | 3 | 5 | 7 |
Stomatal density (number/figure) | 15.44c | 15.93bc | 19.19a | 19.75a | 18.75b |
Stomatal opening rate (%) | 38.49c | 43.50c | 72.31a | 59.81b | 58.00b |
RNA-seq data summary
We prepared the seedling samples (HR0, HR1, HR3, HR5, and HR7) for the HR cultivar on day 0, 1, 3, 5 and 7 post heat treatment at 45 °C, respectively. Totally, 15 cDNA libraries (HR0-a, HR0-b, HR0-c; HR1-a, HR1-b, HR1-c; HR3-a, HR3-b, HR3-c; HR5-a, HR5-b, HR5-c; HR7-a, HR7-b, HR7-c) were constructed for RNA-seq, which composed three biological replicates at each time point.
Through Illumina HiSeq X-ten platform, we generated over 0.402 billion pair-end reads, corresponding to an average of 26.8 million sequence reads per sample. Using HISAT2 [17]with default settings, about 68.9% clean reads were mapped to the jujube reference genome[11].
In order to understand the spatiotemporal expression patterns of all samples, principal component analysis (PCA) was performed. The three samples in the same time point could form independent clusters (Fig. 2A). Moreover, Pearson correlation analysis for all pairs of RNA-seq samples was performed, demonstrating similar results (Fig. 2B), indicating that gene expression in three replications of every sample was homogeneous (Fig. 2A). This suggests that the replicated samples produced data that are acceptable for further analyses.
Exploration of differentially expressed genes (DEGs)
Using software edge R [18], a total of 1642, 4080, 5160 and 2119 differentially expressed genes (DEGs) (FC ≥ 2 or ≤ 0.5, FDR ≤ 0.01) were detected in the comparisons of HR1 vs. HR0, HR3 vs. HR0, HR5 vs. HR0 and HR7 vs. HR0, respectively, with the highest or lowest FC being 214.08 and 2− 14.12 (Fig. 3A, Table S1). It indicated that heat stress could lead to comprehensive transcriptome changes of the cells from the jujube leaf.
In HR1 vs. HR0, 902 up-regulated and 740 down-regulated genes were found. There were 1850 up-regulated and 2230 down-regulated genes in HR3 vs. HR0. In HR5 vs. HR0, 2167 up-regulated and 2993 down-regulated gene were discovered. 1019 up-regulated and 1100 down-regulated genes were identified in HR7 vs. HR0. The numbers of the upregulated and downregulated DEGs were similar in each comparison, indicating that heat stress had effect on promoting and inhibiting the transcription of numerous genes. Moreover, there are 717 common DEGs among four comparisons above (Fig. 3B). In order to verify the validity of these results, five genes from each list, with altered expression levels and FCs, were chosen as representatives to quantify their expression by RT-qPCR. Results shown in Fig. 5A confirm the robustness of global expression data obtained for both upregulated and downregulated genes (upper and lower panel, respectively).
The molecular response to heat stress
To identify the pathways in which the DEGs were mainly involved, Gene Ontology (GO) enrichment analysis was conducted. It showed that 36, 58, 65 and 37 GO terms were identified in HR1 vs. HR0, HR3 vs. HR0, HR5 vs. HR0 and HR7 vs. HR0, respectively (P < 0.01).
Although the leaves of all seedlings did not become withered under high temperature, multiple DEGs were enriched in “response to stress” (GO: 0006950) and “response to heat” (GO: 0009408) terms for all four comparisons (Fig. 4 and Table S2), indicating the ‘JHR17’ could be sensitive to HS normally and might establish a new steady-state balance of biological processes that enable the organism to function, survive at a higher temperature (45 °C) perfectly. We analyzed these DEGs enriched in “response to stress” and “response to heat” terms, it indicated that the expression level of multiple DEGs associated with response to heat stress were obviously up-regulated after HS (Fig. 4A).
Photosynthesis were affected by heat stress
Photosynthesis occurs in chloroplast, and is a sensitive biological process to high temperature in plants. The “chloroplast organization” (GO: 0009658) and “chloroplast RNA processing” (GO: 0031425) term were identified in the HR1 vs. HR0, HR3 vs. HR0, HR5 vs. HR0 and HR7 vs. HR0, implying that the normal physiology of chloroplasts and photosynthesis would be affected by the HS. Truly, the “photosynthesis, light harvesting” (GO: 0009765) and “photosynthesis” (GO: 0015979) terms were found in HR1 vs. HR0, HR3 vs. HR0 and HR5 vs. HR0 (Table S2). Although the “photosynthesis, light harvesting” (GO: 0009765) and “photosynthesis” (GO: 0015979) terms were not enriched in HR7 vs. HR0, some DEGs belong to the two terms were found.
In order to explore how the HS affected the normal physiology of chloroplasts and photosynthesis, the DEGs were enriched in “chloroplast organization” (GO: 0009399), “chloroplast RNA processing” (GO: 0031425), “photosynthesis, light harvesting” (GO: 0009765) and “photosynthesis” (GO: 0015979) terms were checked and analyzed. To surprise, most of the DEGs enriched in the four terms above were up-regulated by the HS. Moreover, this phenomenon was consistent with the results of the qRT-PCR experiments (Fig. 5). It suggested that the HS might not disrupt the physiology of chloroplasts and photosynthesis and promote the photosynthesis of ‘JHR17’.
Metabolism were affected by heat stress
HS always globally affects the metabolism of the plants. The “myo-inositol hexakisphosphate biosynthetic process” (GO: 0010264) was the only common term associated with metabolism identified in HR1 vs. HR0, HR3 vs. HR0, HR5 vs. HR0 and HR7 vs. HR0. However, multiple specific terms associated with metabolism in four comparisons.
In HR1 vs. HR0, “malate metabolic process” (GO: 0006108), “anthocyanin-containing compound biosynthetic process” (GO: 0009718), “plastoquinone biosynthetic process” (GO: 0010236), “negative regulation of nucleotide metabolic process” (GO: 0045980), “vitamin E biosynthetic process” (GO: 0010189) and “monocarboxylic acid biosynthetic process” (GO: 0072330) were found.
In HR3 vs. HR0, “cellular modified amino acid biosynthetic process” (GO: 0042398), “phenylpropanoid metabolic process” (GO: 0009698), “glutamine biosynthetic process” (GO: 0006542), “nucleotide-sugar metabolic process” (nucleotide-sugar metabolic process), “polyamine catabolic process” (GO: 0006598), “anthocyanin-containing compound biosynthetic process” (GO: 0009718), “cutin biosynthetic process” (GO: 0010143), “positive regulation of flavonoid biosynthetic process” (GO:0009963), “glutathione catabolic process” (GO:0006751), “wax biosynthetic process” (GO:0010025), “glutamate biosynthetic process” (GO: 0006537), “glycine betaine biosynthetic process” (GO:0031456), “positive regulation of auxin metabolic process” (GO:0090355), and “positive regulation of tryptophan metabolic process” (GO:0090358) terms about metabolism were found.
In HR5 vs. HR0, “cysteine biosynthetic process” (GO: 0019344), “positive regulation of catalytic activity” (GO: 0043085), “glucosinolate biosynthetic process” (GO: 0019761), “glycogen biosynthetic process” (GO: 0005978), “cellular glucan metabolic process” (GO: 0006073), “hydrogen peroxide catabolic process” (GO: 0042744), “trehalose biosynthetic process” (GO: 0005992), “tryptophan catabolic process” (GO: 0006569), “long-chain fatty acid metabolic process” (GO: 0001676), “indoleacetic acid biosynthetic process” (GO:0009684), “sucrose metabolic process” (GO:0005985), “glutathione catabolic process” (GO:0006751), “phenylpropanoid biosynthetic process” (GO:0009699), “malate metabolic process” (GO:0006108), “glutamine biosynthetic process” (GO:0006542), “starch biosynthetic process” (GO:0019252), and “cellular carbohydrate metabolic process” (GO:0044262) terms associated with metabolism were identified.
In HR7 vs. HR0, “anthocyanin-containing compound biosynthetic process” (GO: 0009718), “glutamine biosynthetic process” (GO: 0006542), “cellular glucan metabolic process” (GO: 0006073), “positive regulation of catalytic activity” (GO: 0043085), “malate metabolic process” (GO: 0006108), “coumarin biosynthetic process”, “plastoquinone biosynthetic process” (GO: 0010236), “cellular modified amino acid biosynthetic process” (GO: 0042398), “polyamine catabolic process” (GO:0006598) and “wax biosynthetic process” (GO:0010025) terms associated with metabolism were identified.
Validation of RNA-seq by qRT-PCR
In order to verify the reliability of the transcriptome data, six DEGs were randomly selected for expression analysis by qRT-PCR experiments (Fig. 5). The expression patterns shown in the qRT-PCR results (Fig. 5B) were consistent with RNA-seq results (Fig. 5A), with PCCs higher than 0.9.