Transcriptomic Analysis Provides Novel Insights into the Heat Stress-Induced Response in Codonopsis tangshen

Codonopsis tangshen Oliv (C. tangshen) is a valuable traditional Chinese medicinal herb with tremendous health benefits. However, the growth and development of C. tangshen are seriously affected by high temperatures. Therefore, understanding the molecular responses of C. tangshen to high-temperature stress is imperative to improve its thermotolerance. Here, RNA-Seq analysis was performed to investigate the genome-wide transcriptional changes in C. tangshen in response to short-term heat stress. Heat stress significantly damages membrane stability and chlorophyll biosynthesis in C. tangshen, as evidenced by pronounced malonaldehyde (MDA), electrolyte leakage (EL), and reduced chlorophyll content. Transcriptome analysis showed that 2691 differentially expressed genes (DEGs) were identified, including 1809 upregulated and 882 downregulated. Functional annotations revealed that the DEGs were mainly related to heat shock proteins (HSPs), ROS-scavenging enzymes, calcium-dependent protein kinases (CDPK), HSP-HSP network, hormone signaling transduction pathway, and transcription factors such as bHLHs, bZIPs, MYBs, WRKYs, and NACs. These heat-responsive candidate genes and TFs could significantly regulate heat stress tolerance in C. tangshen. Overall, this study could provide new insights for understanding the underlying molecular mechanisms of thermotolerance in C. tangshen.


Introduction
Codonopsis tangshen Oliv is a perennial herbaceous plant belonging to the family Campanulaceae. C. tangshen is widely used in traditional Chinese medicine due to its tremendous health benefits. Traditional medicine prepared from the dried roots of C. tangshen is commonly used to strengthen the spleen and lungs, nourish the blood, and maintain body fluids [1]. The main medicinal constituents of C. tangshen include polysaccharides, saponins, alkaloids, sesquiterpenes, glycosides, polyphenolic polyacetylenes, etc. [2,3].
Modern pharmacology and clinical studies also demonstrated that the main bioactive components of C. tangshen have antioxidant, antibiosis, anti-inflammatory, antitumor, and immune enhancement effects [4][5][6][7][8][9]. C. tangshen is mainly cultivated in Chongqing, Hubei, Sichuan, Hunan, Guizhou, and Shanxi provinces [10]. C. tangshen grows best in highaltitude areas with mild and cool climates [11]. However, a temperature above 30 • C could inhibit the growth and development of C. tangshen, thus reducing yield and quality [12], indicating that extremely high temperature is a severe problem for the cultivation of C. tangshen.
It is well documented that short-term or long-term heat stress exposure could significantly decline growth, development, and yield in many plant species. For instance, a 1 • C rise could reduce crop yields by 17% during the growing season [13]. Continuous heat Life 2023, 13, 168 3 of 20 photographed for phenotypic analysis and leaf samples were collected from control and heat-treated plants for physiological and biochemical analysis. All the treatments were repeated in three biological replications.

RNA Extraction and cDNA Library Preparation
Total RNA was extracted from the leaf samples using TRIzol reagent (Invitrogen, Carlsbad, CA, USA). The quality and quantity of RNA samples were monitored using the Nano Photometer spectrophotometer (IMPLEN, Westlake Village, CA, USA) to ensure high-quality samples for transcriptome sequencing. About 1 µg RNA per sample was used to construct the paired-end libraries using NEBNext ® Ultra™ RNA Library Prep Kit for Illumina ® (NEB, Ipswich, MA, USA) according to the manufacturer's instructions. First-strand cDNA was generated using random hexamer primer and M-MuLV Reverse Transcriptase, followed by synthesizing second-strand cDNA using DNA Polymerase I and RNase H. Then, fragmentation was performed using divalent cations under increased temperature. After adenylation of 3' ends of DNA fragments, adaptor sequences with hairpin loop structure were ligated to the fragments. Subsequently, the fragments about 240 bp in length were purified using the AMPure XP system (Beckman Coulter, Beverly, USA). Then, PCR was carried out using Phusion High-Fidelity DNA polymerase, Universal PCR primers, and Index Primer. Subsequently, PCR products were purified, and the quality of the libraries was checked with the Agilent Bioanalyzer 2100 system. Six libraries were prepared using the NEBNext ® Ultra™ RNA Library Prep Kit for Illumina ® (NEB, USA), following the manufacturer's recommendations. Sequencing was performed using an Illumina Hiseq 2000 platform in Novogene Co., Ltd. (Beijing, China).

Measurement of Physiological Properties
Physiological properties, including chlorophyll content, electrolyte leakage (EL), and malonaldehyde (MDA) content were measured and analyzed according to a previous study [16]. For chlorophyll content analysis, 0.1 g leaves were collected and immersed in 95% ethanol. The samples were kept in the dark for 48 h and absorbance of the chlorophyll extract was read at 665 and 649 nm using a spectrophotometer (UV2600, UNIC, Shanghai, China). The chlorophyll content was calculated according to the formula described in a previous study [16].
For EL analysis, 0.5 g fresh leaves were collected and washed with distilled water to remove debris. Then, 15 mL distilled water was added to each tube and kept in a shaker incubator at 200 rpm for 6 h at 28 • C. Subsequently, the initial conductivity (EL1) was measured using a conductivity meter (JENCO-3173, Jenco Instruments, Inc., San Diego, CA, USA). Then, samples were autoclaved at 121 • C for 30 min. Immediately, the tubes were cooled to room temperature and the second conductivity (EL2) was measured. The relative EL was calculated by using the formula: For MDA content determination, 0.3 g frozen leaves were used to prepare the crude enzyme extract. Leaves were ground into powder in liquid nitrogen using a pre-chilled mortar and pestle (4 • C), homogenized with sodium phosphate buffer (PBS), pH 7.4, and centrifuged at 12,000 rpm for 20 min at 4 • C. The supernatant was used to measure the MDA content. A total of 1 mL extract was mixed with a 2 mL reaction solution containing 20% (v/v) trichloroacetic acid and 0.5% (v/v) thiobarbituric acid. The mixture was boiled in a water bath for 30 min, cooled to room temperature, and centrifuged at 12,000× g for 10 min at 20 • C. The absorbance was read at 532 and 600 nm with a spectrophotometer (UV2600, UNIC, Shanghai, China). The MDA content was determined using the following formula: MDA (nmol g −1 FW) = [(A532 − A600) × V × 1000/ε] × W Where A600 and A532 are the absorbances at 600 and 532 nm, respectively, ε is the specific extinction coefficient (155 mM cm −1 ), V is the volume of crushing medium and W is the fresh weight of leaves (g). The physiological traits analysis was performed with three biological replicates, and one-way ANOVA was used to analyze differences between the control and heattreated group.

Processing and Assembling of Illumina Reads
Paired-end reads were generated from the six libraries. Clean reads were obtained by removing low-quality reads and reads containing adapters or poly-N. Subsequently, the Q20, Q30, GC-content, and sequence duplication levels were analyzed. The downstream analyses were carried out using high-quality clean data. Transcriptome assembly was performed using the Trinity version 2.8.5 with all parameters set to default [

Differential Gene Expression Analysis
For each sample, gene expression levels were evaluated by RSEM [32]. Briefly, clean data were mapped back onto the assembled transcriptome, and then a read count for each gene was obtained based on the mapping results. Then, differential expression analysis of two groups was performed with DESeq2 version 3.11 [33] using a model based on the negative binomial distribution. The p values were adjusted to control the false discovery rate with the Benjamini and Hochberg methods. P-adj < 0.05 and |log2foldchange| > 1 were set as the cutoff value to determine significantly differentially expressed genes (DEGs).

GO and KEGG Pathway Enrichment Analysis
The clusterProfiler 4.0 (https://github.com/GuangchuangYu/clusterProfiler, accessed on 11 March 2021) was used for GO and KEGG pathways enrichment analysis of the differentially expressed genes (DEGs) [34]. Adjusted p value < 0.05 (as predicted by clusterProfiler software) was used as the threshold for the GO and KEGG pathway enrichment analysis.

Validation of RNA-Seq Data by Real-Time Quantitative PCR
To validate the RNA-Seq results, real-time quantitative RT-PCR was performed for selected target genes using the previous method [37]. Briefly, total RNA was extracted from frozen leaves using the Trizol reagent (Invitrogen, Carlsbad, CA, USA) according to the manufacturer's instructions. The quality of RNA extract was monitored on 1% agarose gels, and also checked using the NanoPhotometer ® spectrophotometer (IMPLEN, Westlake Village, CA, USA). A total of 2 µg of RNA extract was used to prepare cDNA for each sample using RevertAid First Strand cDNA Synthesis Kit (Invitrogen, Carlsbad, CA, USA). The PCR reaction was performed using a total of 20 µL consisting of 10 µL of SYBR mix, 0.5 µL of each primer, 2 µL of cDNA, and 7µL dd water. The ABI StepOne Plus Real-Life 2023, 13, 168 5 of 20 Time PCR system (Applied Biosystems, Foster City, CA, USA) was used for PCR reaction. Melting curve analysis was used to check the specificity of the amplified PCR product after each PCR reaction. The glyceraldehyde-3-phosphate dehydrogenase (GAPDH) gene was used as a housekeeping control. The primers used in this study are listed in Table S1. The 2 −∆∆CT method was used to determine the relative expression of genes [38]. One-way ANOVA was applied to analyze the statistical differences between different groups.

Effect of Heat Stress on Physiological and Biochemical Traits
Short-term (0.5 h) heat stress treatment did not significantly affect the morphology of C. tangshen. However, long-term heat stress exposure induced morphological and physiological alterations in C. tangshen. After 7 d of heat stress, the leaves of C. tangshen showed wilting and became withered ( Figure 1A,B). Heat stress also significantly decreased chlorophyll content by 20.71% relative to the control ( Figure 1C). Heat stress dramatically increased EL and MDA content by 54.67% and 95.74%, respectively, relative to the control ( Figure 1D,E).

Illumina Sequencing and De Novo Assembly
Six cDNA libraries were constructed and sequenced from the leaves of C. tangshen. After quality control, a total of 44 Table 1). The raw sequencing data have been submitted to the Genome Sequence Archive (GSA) database at the BIG Sub (CRA004261). In addition, de novo assembly was carried out and identified 74,144 transcripts and 32,719 unigenes with an average length of 2346 bp and 2472 bp, respectively. All the unigenes are listed in Table S2.

Functional Annotations of C. tangshen Unigenes
The functional annotations of unigenes were performed by searching against several databases, including KOG (Clusters of Orthologous Groups of Proteins), GO (Gene Ontology), KO (KEGG Ortholog Database), NR (NCBI non-redundant protein sequences), Swiss-Prot (a manually annotated and reviewed protein sequence database) ( Table 2). A total of 19,546 (59.7%) unigenes had hits in the GO database and were classified into three GO categories, including "biological process (BP)," "cellular component (CC)," and "molecular function (MF) (Figure 2).

Analysis of DEGs in Response to Heat Stress in C. tangshen
To obtain insights into the transcriptomic response of C. tangshen to heat stress, we further analyzed the differentially expressed genes (DEGs) list that was identified using the DESeq package using the CK (control samples) and HT (heat-treated samples). A tota of 2691 DEGs were identified, including 1809 upregulated and 882 downregulated. The expression and annotation of DEGs are shown in Table S4 and Table S5, respectively.
Moreover, the DEGs from the three control and heat treatment replicates were clus tered using the hierarchical clustering method ( Figure 5). The DEGs were classified into two groups and showed contrasting expression profiles between CK and heat treatment Genes with similar expression levels were preferentially clustered together. There were no significant differences in the expression of DEGs between replicates.

Analysis of DEGs in Response to Heat Stress in C. tangshen
To obtain insights into the transcriptomic response of C. tangshen to heat stress, we further analyzed the differentially expressed genes (DEGs) list that was identified using the DESeq package using the CK (control samples) and HT (heat-treated samples). A total of 2691 DEGs were identified, including 1809 upregulated and 882 downregulated. The expression and annotation of DEGs are shown in Table S4 and Table S5, respectively.
Moreover, the DEGs from the three control and heat treatment replicates were clustered using the hierarchical clustering method ( Figure 5). The DEGs were classified into two groups and showed contrasting expression profiles between CK and heat treatment. Genes with similar expression levels were preferentially clustered together. There were no significant differences in the expression of DEGs between replicates.

GO and KEGG Enrichment Analyses of DEGs
GO analysis was performed to investigate the functions of DEGs. As shown in Figure 6, Table S6, DEGs were classified into three GO categories, MF, BP, and CC. The highly enriched MF GO terms were unfolded protein binding, chaperone binding, DNA binding transcription factor activity, ATPase activator activity, channel activity, and heat shock protein binding. In the BP category, protein folding, protein refolding, regulation of pH, sodium ion transport, lipid catabolic process, photosynthetic electron transport chain, and "photosynthesis, light reaction" were significantly enriched. However, we did not find significantly enriched GO terms in the CC category.

GO and KEGG Enrichment Analyses of DEGs
GO analysis was performed to investigate the functions of DEGs. As shown in Figure  6, Table S6, DEGs were classified into three GO categories, MF, BP, and CC. The highly enriched MF GO terms were unfolded protein binding, chaperone binding, DNA binding transcription factor activity, ATPase activator activity, channel activity, and heat shock protein binding. In the BP category, protein folding, protein refolding, regulation of pH, sodium ion transport, lipid catabolic process, photosynthetic electron transport chain, and "photosynthesis, light reaction" were significantly enriched. However, we did not find significantly enriched GO terms in the CC category.
HT-03 HT-02 HT-01 CK-03 CK-02 CK-01 Furthermore, to investigate the metabolic pathways regulated in response to heat stress in C. tangshen, DEGs were systematically subjected to the KEGG pathway enrichment analysis ( Figure S2, Table S7). We found 30 enriched metabolic pathways. Interestingly, the upregulated DEGs were mainly enriched into plant-pathogen interaction (67

DEGs Related to Short-Term Heat Stress in C. tangshen
To identify heat stress-responsive candidate genes, we focused on genes related to HSF-dependent pathways, ROS scavenging enzymes, plant hormone signaling transduction, and transcription factors, including bHLH, bZIP, MYB, WRKY, and NAC. When exposed to heat stress, the increased expression of HSPs usually form a complex with the HSF [39]. We found 130 HSPs and 8 HSFs enriched by heat stress-induced upregulated genes in the HSF-dependent pathway. ROS levels increased significantly in plants exposed to high temperatures [40]. We identified six GSTs and one PAO significantly upregulated genes involved in ROS scavenging. Furthermore, genes related to plant hormone signaling transduction were significantly induced by heat treatment, including eleven ERFs, two ABFs, and one PYL. Calcium-dependent protein kinases (CDPKs) genes are critical components in the heat stress-mediated calcium signal transduction pathway [41]. Here, a total of 15 CPKs and 5 calmodulin-like proteins (CMLs) were significantly upregulated by heat stress in C. tangshen (Table 3). In addition, heat-responsive transcription factors (TFs), including bHLH, bZIP, WRKY, and NAC, were identified in C. tangshen. Out of the 10 bZIP genes, 8 were significantly upregulated by heat treatment compared to the control (Table 4). Additionally, 14 DEGs encoding the NAC TFs were identified, of which 11 were significantly upregulated, while the other 3 were remarkably downregulated after heat treatment (Table 5). We also identified nine WRKY genes that were significantly upregulated by heat stress (Table 6). Moreover, three bHLH genes were identified and remarkably downregulated upon heat stress ( Table 7). The clustering of the representative DEGs related to the different gene categories in heat treatment compared with CK is shown in Figure 7.

Validation of the DEGs by qRT-PCR
The RNA-Seq results were validated by analyzing the expression of ten unigenes re-  Figure 7. Clustering of the differentially expressed genes related to the representative gene category in heat treatment compared with CK. Red color denoted genes with high expression levels, while blue color indicated genes with low expression levels. The relative expression levels of DEGs were calculated using the log2 ratio.

Validation of the DEGs by qRT-PCR
The RNA-Seq results were validated by analyzing the expression of ten unigenes related to TFs or stress responses using qRT-PCR. The transcriptional level of three stressresponsive unigenes (HSP, HSF, and GST) was significantly increased in the Ht group compared to the CK group (Figure 8). Six TFs genes were also checked. Five TFs were markedly upregulated by heat stress (ERF, ABF, bZIP, MYB, and NAC), and only one (bHLH) was downregulated by heat stress. Additionally, the expression of a light-harvesting chlorophyll a/b binding protein (LHC) gene was inhibited by heat stress. In summary, the expression profiles of these unigenes from the qRT-PCR were consistent with that from the RNA-Seq, indicating that the RNA-Seq data were reliable and accurate.

Heat Stress Regulates Physiological and Biochemical Traits in C. tangshen
Heat stress is a severe abiotic factor limiting plants' growth and development. In this study, we demonstrated the heat stress-induced biochemical and transcriptomic responses of C. tangshen. Transcriptomic analysis was carried out after short-term heat stress (0.5 h) treatment. Extended heat stress for 7 d caused a noticeable phenotypic change in C. tangshen, including leaf wilting. Electrolyte leakage and lipid peroxidation have been used to study heat stress-induced membrane damage in plants [42,43]. Physiological and biochemical analysis showed that C. tangshen was severely affected by 7 d heat stress treat-

Heat Stress Regulates Physiological and Biochemical Traits in C. tangshen
Heat stress is a severe abiotic factor limiting plants' growth and development. In this study, we demonstrated the heat stress-induced biochemical and transcriptomic responses of C. tangshen. Transcriptomic analysis was carried out after short-term heat stress (0.5 h) treatment. Extended heat stress for 7 d caused a noticeable phenotypic change in C. tangshen, including leaf wilting. Electrolyte leakage and lipid peroxidation have been used to study heat stress-induced membrane damage in plants [42,43]. Physiological and biochemical analysis showed that C. tangshen was severely affected by 7 d heat stress treatment, evidenced by pronounced levels of electrolyte leakage and lipid peroxidation. Additionally, chlorophyll content was significantly reduced by heat stress. Our results agreed with the findings reported in P. ternate [16]. Together, these results indicate that heat stress is a significant limiting factor for the growth and development of C. tangshen. Thus, identifying key heat-responsive genes and understanding the underlying molecular mechanisms could be crucial to improve the thermotolerance of C. tangshen.

Transcriptomic Response of C. tangshen to Heat Stress
Plants reprogram transcriptomic profiles and develop molecular mechanisms to cope with heat stress [44][45][46]. In this study, we conducted RNA-seq to understand the molecular mechanisms of heat stress response in C. tangshen at the transcriptome level. Here, 2691 DEGs were identified, including 1809 upregulated and 882 downregulated. Further analysis demonstrated that most of the DEGs induced by heat stress were related to the HSF-HSP pathway, ROS scavenging enzymes, plant hormone signaling transduction, and TFs, including bZIP, MYB, NAC, WRKY, and bHLH, which have been reported in previous studies [25,27,47,48]. This study is the first report to investigate the transcriptomic profiles of C. tangshen in response to heat stress. Moreover, we identified candidate heat-responsive genes and transcription factors that could play a crucial role in improving the thermotolerance of C. tangshen. Overall, our results provide essential information for the functional characterization of heat-responsive genes in C. tangshen.

ROS-Scavenging Activity and HSF-HSP Network Involved in HSR
Heat stress induces excessive production of reactive oxygen species (ROS) in plant cells, leading to oxidative stress, which further causes DNA damage, protein denaturation, and impairment of cell membranes [49]. ROS-scavenging enzymes, such as GSTs and APXs, play a crucial role in heat stress tolerance in plants [50]. In this study, six GST and two APX genes were significantly induced in response to heat stress, indicating that these genes could be involved in scavenging ROS and protecting C. tangshen from heat stress-induced damages.
The HSF-HSP network is activated to protect plant cells against heat stress-induced damages [51]. In Arabidopsis, HSPs, including Hsp 100, Hsp90, Hsp70, and Hsp60, were significantly induced to provide mechanisms against heat stress [52]. In this study, heat stress induced the expression of HSF and HSP genes in C. tangshen. For instance, 8 HSFs and 130 HSPs genes, including HSP20, HSP70, and HSP90 were, significantly upregulated by heat stress. Previous studies showed that HSP proteins could be crucial in improving heat tolerance in plants. Overexpression of HSP70 from Brassica campestris improved heat tolerance in tobacco [53]. Additionally, overexpression of a peony HSP70 enhanced heat tolerance in transgenic Arabidopsis [54]. In addition, GmHsp90A2 from soybean was identified as a positive player of heat tolerance [51]. Thus, HSP70 and HSP90 could be key genes to improve the thermotolerance of C. tangshen. Meanwhile, further study is required to investigate how HSPs and HSFs function in C. tangshen under heat stress.

DEGs Related to Plant Hormone Signal Transduction
The phytohormones, including ABA, ethylene, and brassinosteroids, play crucial roles in plants' growth development and heat stress response [53][54][55][56][57][58][59]. Ethylene-responsive transcription factors (ERFs) are a transcriptional subfamily downstream of the ethylene signaling pathway [60,61]. ERFs are involved in floral organ development, fruit ripening, abiotic stresses, hormonal signal transduction, and pathogenesis [61,62]. Most ERFs were significantly upregulated in plants treated with heat stress [63]. In this study, most ERF genes were induced in response to heat stress (Figure 9), suggesting that ERFs might function in protecting C. tangshen from heat stress damage. ABFs constitute a significant transcription factor in the ABA signaling pathway that activates ABA/stress-responsive genes in response to abiotic stress [64]. The binding of ABA to PYL proteins inhibits protein phosphatase type 2Cs (PP2Cs) activity, thereby activating SNF1-related kinases (SnRK2s) to phosphorylate and activate ABFs [65]. A study on the heat-tolerant rice cultivar "Annapurna" showed that PP2Cs and ABFs were upregulated at 37 • C temperature stress [22]. In agreement, all the PP2Cs and most ABFs were upregulated by heat stress in C. tangshen, indicating the critical role of the ABA signaling pathway in C. tangshen thermotolerance. These results suggest that different hormone signaling pathways are involved in the thermotolerance of C. tangshen. Several functional studies indicate that ABF genes could play a remarkable role in improving plants' stress tolerance. For instance, overexpression of OsABF1 positively regulated drought tolerance in rice [66]. Additionally, overexpression of IbABF4 from sweet potatoes could improve heat tolerance in Arabidopsis [67]. Thus, the ABF genes identified in this study could play a crucial role in improving the heat tolerance of C. tangshen.

DEGs Related to Transcription Factors
Transcriptional factors, including MYBs, NACs, WRKY, bHLH, and bZIP, play a pivotal role in heat stress response [68][69][70][71][72]. In this study, 17 MYBs, 14 NACs, 9 WRKY, 3 bHLH, and 10 bZIP were significantly regulated in response to heat stress in C. tangshen. Overexpression of crucial stress-responsive TFs could remarkably improve stress tolerance in plants. For instance, MYB305, isolated from Lilium longiflorum, with transactivation ability in yeast and plant cells, was reported to play a positive role in thermotolerance [73]. In Arabidopsis, AtWRKY25, AtWRKY26, and AtWRKY33 were involved in resistance against pathogenic bacteria and abiotic stresses, including heat stress [74]. Likewise, NAC19 and bZIP28 identified from Arabidopsis positively regulated thermotolerance [69,72]. Overexpression of AtMYB68 enhanced drought and heat tolerance in Arabidopsis and Brassica napus [75]. Additionally, overexpression of OsbZIP46 could positively regulate drought tolerance in rice [76]. In a recent study, a maize NAC TF, ZmNAC074, improved heat tolerance in Arabidopsis [77]. Thus, we speculated that the TFs identified in this study could play a significant role in improving the thermotolerance of C. tangshen.
Additionally, a large number of calcium-dependent protein kinases (CDPKs) were upregulated in response to heat stress in C. tangshen. CDPK proteins play crucial roles in plants' abiotic stress response by regulating downstream target proteins [78]. For instance, overexpression of CDPK6 enhanced salt and drought tolerance in transgenic Arabidopsis [79]. In rice, OsCDPK7 positively regulates cold, drought, and salinity stresses [80]. Additionally, ZmCDPK7 was involved in improving heat tolerance in maize [81]. These results indicate that CDPK proteins are a central player in multiple stress responses. Thus, the upregulation of CDPKs could be one of the mechanisms of thermotolerance in C. tangshen. In addition, candidate CDPKs identified in this study can be used for functional analysis in response to heat stress. genes could play a remarkable role in improving plants' stress tolerance. For instance, overexpression of OsABF1 positively regulated drought tolerance in rice [66]. Additionally, overexpression of IbABF4 from sweet potatoes could improve heat tolerance in Arabidopsis [67]. Thus, the ABF genes identified in this study could play a crucial role in improving the heat tolerance of C. tangshen.

Conclusions
This study demonstrated how C. tangshen responded to short-term heat stress at the transcriptome level. Heat stress is a significant factor limiting the growth and development of C. tangshen. Heat stress reprogrammed the transcriptomic profiles of C. tangshen and regulated various genes related to the HSP-HSF network, ROS-scavenging activity, hormone signaling transduction, and TFs. Overall, these findings could provide essential transcriptional information for improving the thermotolerance of C. tangshen. Further functional studies are needed to understand the detailed molecular mechanisms of these candidate genes in response to heat stress in C. tangshen.

Data Availability Statement:
The raw sequencing data have been submitted to the Genome Sequence Archive (GSA) database at the BIG Sub (https://ngdc.cncb.ac.cn/gsub/, accessed on 1 October 2022). The accession number is CRA004261.

Conflicts of Interest:
The authors declare no conflict of interest.