Identification, classification, and expression profile analysis of heat shock transcription factor gene family in Salvia miltiorrhiza

Identification of SmHSF members

Genome Warehouse, BIG Data Center, project number PRJCA003150, which is available at https://ngdc.cncb.ac.cn/bioproject/browse/PRJCA003150, provided the assembly and annotation data for Salvia miltiorrhiza Bunge (Lamiaceae) (Song et al., 2020). The S.miltiorrhiza genome data was screened for HSF genes using the hmmsearch function in HMMER 3.0 (Finn, Clements & Eddy, 2011), and HSF members were identified using the hidden Markov model (HMM) corresponding to the HSF domain (PF00447) that was downloaded from the PFAM database (https://pfam.xfam.org/). E value is under 0.005 and default parameters were used. The S.miltiorrhiza genome database was searched using the BLASTP tool to find similar sequences using the Arabidosis HSF amino acid sequences that were obtained from TAIR (http://www.arabidopsis.org). Additionally, the SMART (Letunic & Bork, 2018) and CDD programs (https://www.ncbi.nlm.nih.gov/cdd/) were used to detect DBD domains in all acquired SmHSF proteins. Using TBtools (Chen et al., 2020) and the physical locational data from the S.miltiorrhiza genome, all of the SmHSFs were mapped to the eight chromosomes and three scaffolds of S.miltiorrhiza (Zhang et al., 2021). In addition, the ExPasy (https://web.expasy.org/compute_pi/) program was used to explore the physical characteristics of predicted HSF proteins. Subcellular localization predictions were generated using Plant-mPLoc with default parameters (http://www.csbio.sjtu.edu.cn/bioinf/plant-multi/) (Chou & Shen, 2010).

Structural, motif, and phylogenetic analysis of SmHSF genes

By comparing predicted coding sequences, the exon-intron distribution of each SmHSF gene was drawn using TBtools (Chen et al., 2020). The online tool MEME5.4.1 (https://meme-suite.org/meme/tools/meme) was used to investigate the conserved motifs of the SmHSF protein sequences. The number of motifs being requested is 10, and the range for motif width is 6 to 50 (inclusive). The results were visualized using TBtools. Using ClustalW in MEGAX with default parameters (https://www.megasoftware.net/), multiple sequence alignments of HSF proteins from Arabidopsis thaliana, Oryza sativa (Chauhan et al., 2011), Solanum lycopersicum (Jin et al., 2020), and S. miltiorrhiza were carried out. The sequences used for alignment are listed in Table S1. The alignment results were used to construct a phylogenetic tree using the Maximum Likelihood method with 1000 bootstrap replicates. Additionally, the evolutionary tree was embellished using Evolview (https://www.evolgenius.info/evolview/) (Subramanian et al., 2019).

Cis-acting elements analysis of SmHSFs

The S.miltiorrhiza genome database provided the 2,000-bp sequence upstream and 200-bp sequence downstream from the transcription start site of each SmHSF gene. These sequences were used to find cis-acting regulatory elements with the online tool PlantCARE (https://bioinformatics.psb.ugent.be/webtools/plantcare/html/). A number of significant cis-acting elements were counted for each SmHSF using TBtools (Chen et al., 2020).

Expression profiles of SmHSFs under drought stress based on transcriptome data

The transcriptomic data generated from the leaf and root of S. miltiorrhiza under moderate drought stress have been described previously (Li, 2020). TBtools was used to create a heat map, and after the data Log2(FPKM+1), row normalization was carried out, to make it easier to compare the expression trends of each SmHSF in leaves and roots and after drought stress.

Plant materials and heat treatments

For heat stress treatment trials, tissue culture seedlings with the same or equivalent growth vigor were chosen. The S.miltiorrhiza tissue culture seedlings, which were grown for two months on MS medium, were moved to an artificial climate chamber with a constant temperature of 42 °C and exposed to heat stress continuously for 24 h, and six-time points (0 h, 1 h, 2 h, 6 h, 12 h, and 24 h) were selected for sample collection. All samples had three biological duplicates and were immediately frozen in liquid nitrogen following collection, and then stored at −80 °C for RNA extraction.

RNA isolation, and cDNA synthesis

Total RNA from S. miltiorrhiza samples were isolated using the Plant Total RNA Isolation Kit Plus following the manufacturer’s protocol (Foregene, Ltd., Chengdu, China). cDNA synthesis was done by using 1 ug of the total RNA samples with RT Easy™ II Kit (Foregene, Ltd., Chengdu, China).

Quantitative Real-Time PCR (qRT-PCR)

Reverse transcribed cDNA products were used as templates of quantitative Real-Time PCR (qRT-PCR). The reaction was carried out on a CFX Opus Real-time PCR system using the Real-Time PCR Easy ™ -SYBR Green I (Foregene, Ltd., Chengdu, China) following the manufacturer’s instructions. NCBI-BLAST Primer designed the primers for 35 SmHSF genes. S. miltiorrhiza Actin was used as an endogenous control for the normalization of expression levels of genes (Jiang et al., 2020). The relative expression levels were calculated using the 2^{− ΔΔ Ct} method. Data were analyzed using one-way ANOVA in GraphPad Prism 9 software (*, P < 0.05; **, P < 0.01; ***, P < 0.001; ****, P < 0.0001). To ensure reproducibility and dependability, three biological replications and three technical replications were implemented for each sample. The primers for the SmHSFs used for qRT-PCR analyses are listed in Table S2.

Identification and chromosomal localization of SmHSF gene family

A total of 35 SmHSF members were found in the Salvia miltiorrhiza genome after blastp and hmmsearch. These genes include HSF family domains and have been validated by NCBI-CD search and SMART. They were given the names SmHSF1 to SmHSF35 based on the gene IDs found in the newly sequenced genome. The members of the SmHSF gene family are listed in Table S3 along with their coding DNA and protein IDs and sequences.

The 35 SmHSFs protein sequences ranged in length from 185 bp to 508 bp, with SmHSF11 having the longest and SmHSF24 having the shortest. SmHSF24 had a molecular weight of 21,538.39 Da, whereas SmHSF11 had a molecular weight of 56,360.16 Da. The highest and lowest members were consistent with the protein sequence. And the range of the isoelectric points (pI) was 4.64 (SmHSF25) to 9.93. (SmHSF29). All SmHSFs were determined to be unstable, with the exception of SmHSF35, which was stable, according to the examination of the instability index (the stability of the protein in a test tube). Additionally, the GRAVY ranged from −0.92 to −0.481, while the Aliphatic Index (AI) ranged from 57.61 to 78.99. Plant-mPLoc subcellular localization predictions suggested that all the HSF proteins were located in the nucleus. All of the information above is displayed in following Table 1. Except for SmHSF1, SmHSF2, and SmHSF3, which were found on the scaffolds, the remaining genes were determined to be spread among eight chromosomes according to the mapping of 35 SmHSFs on the Salvia miltiorrhiza chromosome (Fig. 1). Each chromosome contains a different number of S. miltiorrhiza HSF genes, and the position of the genes does not indicate anything about their function.

Phylogenetic, structural, and motif analysis of SmHSF genes

To investigate the phylogenetic relationship of 35 SmHSFs, a phylogenetic tree was constructed by combining SmHSFs with 21 Arabidopsis HSFs (AtHSFs), 25 rice HSFs (OsaHSFs), and 26 tomato HSFs (SlyHSFs) (Fig. 2). SmHSFs were separated into three groups (HSFA/HSFB/HSFC) in accordance with the grouping of HSFs in Arabidopsis. Each group has different member distributions. The distribution of members in these groups is like that in Arabidopsis thaliana, with Group A having the most members (22 genes) and Group C having the fewest members (two genes). Group A of SmHSF consisted of nine subgroups, from A1-A9; group B, with four subgroups (B1–B4); and group C, with only one subgroup (C1).

To further examine the structural diversity of genes, we analyzed the gene structures of 35 SmHSF genes. Analysis revealed that the majority of SmHSFs members belonged to the same subgroup and shared exon and intron structures (Fig. 3B). The number of exons and introns differs among SmHSFs, however. Of the 35 SmHSFs, 32 have two exons, two genes have three exons, and interestingly, SmHSF24 has five exon regions, two of which are incredibly short, and no UTR was discovered. All SmHSF genes have introns, of which 28 genes have one, five have two, one has three, and one has four.

Ten motifs were found in the SmHSF family members using motif analysis with MEME, the distribution of motifs corresponding to the phylogenetic tree of the SmHSFs family is shown in Fig. 3. Comparable motif compositions among the SmHSFs proteins clustered in the same subgroup imply that the members of the subgroup share similar activities. The DBD domain and HR-A/B are both present in the SmHSF family, just like the HSF family in other species. Motif2, motif3, and motif4 could be found in each member. Following sequence alignment, it was discovered that motif 1, motif 2, motif 3, and motif 6 make up the DBD domain of HSF, whereas motifs 4 and motifs 5 represent HR-A and HR-B, respectively. Except for SmHSF13 and SmHSF24, we detected motif1 and motif6 in all members. Only group A members possess motif 5, which has not been discovered in members of groups B or C. Additionally, only a few members of group A possess motifs 7, 8, and 9, where motif 8 is AHA and motif 9 is NLS. The sequence and distribution information of all motifs are in Table S4.

Cis-acting elements analysis of SmHSFs

To further learn about the biological function of SmHSFs, the cis-acting elements in the 2,000 bp upstream and 200 bp downstream sequences from the transcription start sites of SmHSFs were examined using PlantCARE. A variety of cis-acting components, including those related to stress, hormones, and development, constitute the promoter of each SmHSF (Fig. 4). Members of SmHSFs are broadly diversified in terms of development-related components. The ABA-responsive element (ABRE), ABA and drought responsive element (ABRE3a), methyl jasmonate (MeJA)-responsive element (CGTCA-motif), and another MeJA-responsive element (TGACG motif), the ethylene-responsive element (ERE), gibberellin-responsive element (P-box), and cis-acting element involved in salicylic acid responsiveness (TCA-element) were all found in 30, 20, 30, 30, 28, 7, and 23 SmHSFs, respectively. Also, several light-associated cis-acting elements, like ACE and G-Box, as well as components connected to stress. MYB element, MYC element, DRE core (drought, salt, low temperature, and ABA responses), and STRE element (activated by heat shock, osmotic stress, low pH, and nutrient starvation) were identified in the promoters of 33, 33, 10, and 24 SmHSF genes, respectively. Additionally, anaerobic induction elements (AREs), low-temperature responsiveness elements (LTRs), MYB binding site involved in drought-inducibility (MBSs), defense and stress responsiveness (TC-rich repeats), wounding and pathogen responsiveness element (W-box), and wound-responsive element (WUN-motif) were found in 27, 18, 18, 10, 15, and 20 SmHSFs, respectively. However, the heat stress responsiveness element (HSEs) is only found in SmHSF7. These findings show that SmHSFs may be associated with a variety of transcriptional regulations involving hormones, stress responses, and development.

Expression profiles of SmHSFs under drought stress

Tanshinone and salvianolic acid, have been observed to accumulate in the roots and leaves of S. miltiorrhiza, respectively (Sha, 2015; Wang & Wu, 2010). Additionally, drought stress will have some effect on HSF. To further investigate the expression of HSF in these two tissues and the expression changes in response to drought stress, we selected previous transcriptome data for HSF expression analysis. SmHSF23 and SmHSF 26 were not present in the transcriptome data, while other members were expressed to different degrees in roots and leaves as well as during drought stress. There were variations in the expression patterns of SmHSFs in roots and leaves, as shown in Fig. 5. SmHSF21 and SmHSF29 were not expressed in leaves or roots, respectively. There were 11 members, with the expression in roots being lower than in leaves. Of them, SmHSF10, SmHSF14, SmHSF16, and SmHSF34 clearly differed from the others. Notably, members like SmHSF20, SmHSF6, SmHSF8, SmHSF3, and SmHSF1 turned out to express mostly in the roots. Certain members like SmHSF10, SmHSF16, SmHSF30, SmHSF14, and SmHSF28 have varying expression in the roots and leaves and generally have greater expression levels. Following treatment for drought stress, it was found that the expression of 16 genes was up-regulated in the leaves, with SmHSF9 being the most visibly changed. It was discovered that the expression of 18 genes was up-regulated in roots, with SmHSF1, SmHSF2, SmHSF3, SmHSF5, and SmHSF33 standing out. Eleven SmHSF genes (SmHSF1, 2, 8, 11, 12, 13, 15, 16, 22, 27, and 35) are expressed in both roots and leaves under drought stress, and these genes may be activated by drought stress, further studies are required on how the functions.

Expression profiles of SmHSFs under heat stress

In response to heat stress, HSF genes are crucial for plant heat tolerance. To show how HSF genes react to heat stress in our work, the expression patterns of the SmHSFs gene family were identified using qRT-PCR. The results are depicted in Fig. 6. Despite barely being expressed, SmHSF25 rose in response to heat stress. When exposed to heat stress, the expression of SmHSF11, SmHSF14, and SmHSF30 did not change considerably, and the expression of SmHSF33, SmHSF34 were repressed, while the expression of other SmHSF genes was up-regulated to varying degrees. Notably, the expression of SmHSF2, SmHSF7, and SmHSF9 was markedly upregulated in response to heat stress, particularly SmHSF9, which was almost 4000-fold greater than the control at 2 h, showing that SmHSF2, SmHSF7, and SmHSF9 were implicated in the heat stress response pathway. In addition, the expression of SmHSF3, SmHSF4, SmHSF15, SmHSF27, SmHSF31, SmHSF32, and SmHSF35 also changed significantly after heat stress, and these genes are all worthy of further consideration.