Genome-Wide Characterisation of the Heat Shock Transcription Factor Gene Family in Betula platyphylla Reveals Promising Candidates for Heat Tolerance

1. Introduction

Plants frequently face numerous abiotic stresses in their surroundings, including conditions like heat, drought and heavy metals stress throughout in their life cycle [1,2]. These environmental factors can negatively influence various aspects of plant health, encompassing growth rates, organ development, yield, and overall quality [3]. Over the course of their extensive evolutionary history, plants have developed a sophisticated transcriptional regulatory system to tackle external challenges, which includes a network governed by Transcription Factors (TFs). These factors are essential regulatory components within cells, overseeing gene transcription levels, engaging in cell differentiation, and mediating responses to environmental stimuli [4]. Transcription factors play a pivotal role in managing both abiotic and biotic stressors, including families such as NAC [5], TCP [6], Dof [7], GATA [8], bHLH [9], GRAS [10], and ERF [11], among others. Notably, HSF has been recognized as a key regulator in stress responses. Among various abiotic stressors, high temperature stress (HS) stands out due to its considerable adverse impact on plants [12,13,14]. Heat stress transcription factors represent a distinct class that becomes activated when plants endure increased temperatures. These activated factors can attach to specific DNA sequences known as Heat Shock Response Elements (HSE), which triggers the transcription of heat shock protein genes, aiding in the protection of plant cells against environmental stressors [15]. Investigating HSF is integral to developing new plant varieties capable of withstanding abiotic stress. In light of the growing global greenhouse effect and the rapid escalation of surface temperatures on Earth, understanding the mechanisms of HSF tolerance to abiotic stress is vital for the selection and breeding of heat-resistant species.

Similar to many transcription factors, the HSF family exhibits a modular architecture. While there is considerable variation in their sequences and sizes, the patterns of promoter recognition are highly conserved across eukaryotes. As illustrated in Figure 1A, the five key structural domains of HSFs include: the DNA-binding domain (DBD), the oligomerization domain (OD), the nuclear export signal (NES), the nuclear localization signal (NLS), and the C-terminal transcriptional activation domain (CTAD) [16]. The DBD, located near the N-terminus, is the most conserved region among HSFs and consists of approximately 100 amino acids, which form a hydrophobic core through three α-helices and four β-strands. This domain contains a conserved Helix-Turn-Helix (HTH) motif, a hallmark of catabolic activating proteins, which enables HSFs to bind specifically to Heat Shock Elements (HSEs: 5'-AGAAAnnTTCT-3'), thus modulating the expression of heat stress-responsive genes [17]. The OD domain, present in all HSFs, plays a critical role in oligomerization. It contains two hydrophobic heptapeptide repeat regions, HR-A and HR-B, each forming a characteristic leucine zipper structure. Under heat stress, HSFs oligomerize into homotrimers via this helical structure, facilitating their binding to HSEs and regulating heat stress responses [18]. The OD structure is also pivotal in classifying HSF family members [19]. Based on sequence homology and structural features of the OD, plant HSFs are categorized into three evolutionarily conserved classes: A, B, and C. Classes A and C feature a compact HR-A/B region, whereas class B HSFs exhibit a similar region but lack the amino acid insertions between HR-A and HR-B present in classes A and C—21 and 7 amino acids, respectively [20]. Most HSFs also contain both a nuclear localization signal (NLS) and a nuclear export signal (NES) at the C-terminus. The NLS, adjacent to the C-terminus of the OD, is enriched in basic amino acids such as lysine and arginine, facilitating nuclear entry. Some HSFs additionally possess a leucine-rich NES, which mediates protein export from the nucleus. Together, the NLS and NES work synergistically to regulate the subcellular distribution of HSFs between the nucleus and cytoplasm [21]. The CTAD, located at the C-terminus, is the least conserved region of HSFs and contains a short AHA motif composed of aromatic, hydrophobic, and acidic amino acids [16]. Although less conserved, the CTAD plays a crucial role in transcriptional activation. The AHA motif specifically enables class A HSFs to function as transcriptional activators, while class B and C HSFs, which lack this motif, generally do not exhibit transcriptional activation activity [22].

The Heat Shock Factor (HSF) family of transcription factors consists of highly conserved genes found across both eukaryotes and prokaryotes. The first HSF gene was identified and cloned in Saccharomyces cerevisiae [23], followed by the discovery of the first plant HSF in tomato [24]. With the advent of high-throughput sequencing technologies, an increasing number of HSF genes from diverse species have been characterized, including Arabidopsis thaliana [25], rice [26]and poplar [27]. These studies have revealed that HSF genes play critical roles in regulating plant growth, development, and responses to environmental stress [28]. For instance, the quadruple mutant of Arabidopsis thaliana (hsfa1a/b/d/e) exhibited severely impaired development, underscoring the pivotal role of HSFA1 in regulating plant growth and development [29]. Similarly, AtHSFA4c has been implicated in root development, with mutations in this gene leading to a reduction in lateral root formation in Arabidopsis thaliana [30]. AtHSFA2 modulates plant development in response to AtABI3 and is essential for seed development [31]. Over-expression of the sunflower HSFA9 gene in tobacco resulted in a substantial accumulation of heat shock proteins (HSPs), enhancing seed longevity and improving the water retention capacity of nutrient organs [32]. The HSFB2a and HsfB2b genes play crucial roles in the development of Arabidopsis thaliana hypocotyls [33,34], while AIHSFB1 accelerates senescence by repressing genes involved in growth and photosynthesis; conversely, its repression promotes growth and delays senescence [35]. In contrast, PpHSF5 in peach has been identified as a potential inhibitor of both root and aerial organ development [36]. These findings collectively suggest that HSFs are integral to the regulation of plant growth and development. High temperature is one of the most common abiotic stresses and a major focus of HSF-related research in stress physiology. The HSFA subfamily plays a central role in the heat stress response [37]. Tomato was the first plant to undergo detailed research on heat stress transcription factors, and overexpression of SlHsfA1 resulted in significantly improved heat tolerance compared to SlHsfA1-silenced plants. Similarly, SlHsfB1 has been shown to enhance heat resistance in tomato [38]. In Arabidopsis thaliana, HSFA6b not only contributes to the heat stress response but also functions as a positive regulator of ABA signaling, improving salt and drought tolerance [39]. Overexpression of HsfA2 in Arabidopsis confers enhanced tolerance to multiple environmental stresses, while knockdown of the gene reduces both basal and acquired thermotolerance as well as tolerance to oxidative stress [40]. In Capsicum annuum, the overexpression of CaHsfA2 [41], as well as the GmHsf-34 gene in soybean [42], significantly increased tolerance to drought and heat stress in Arabidopsis thaliana. Furthermore, heterologous expression of CsHsfA2 in tea plants enhanced the heat tolerance of transgenic yeast [43]. Recent studies have shown that transient expression of 35S: Gooseberry LcHsfA2a in seedlings resulted in significantly lower levels of hydrogen peroxide (H₂O₂) and enhanced heat tolerance following heat stress [44]. These studies highlight the diverse and essential roles of HSFs in plant stress responses. It is important to note that even within the same species, HSFs can exert markedly different, and sometimes opposing, effects on plant physiology. For example, in tomato, HSFB1 acts as a co-factor with HSFA in the heat stress response, while in Arabidopsis thaliana, AtHSFB1 inhibits the activities of AtHSFA2 and AtHSFA7 [45]. These findings uncover the diversity of HSF gene functions across different plant species, suggesting that the precise roles of individual HSFs need to be thoroughly investigated.

B. platyphylla, a deciduous tree belonging to the genus B. platyphylla in the Betulaceae family, is renowned for its environmental adaptability. This species is particularly well-suited for forest restoration in areas impacted by logging or wildfires due to its resilience. Its wood is highly valued in the production of plywood, furniture, and premium-quality paper [46]. While B. platyphylla is generally considered a hardy tree, it faces significant challenges posed by high summer temperatures. Since 1980, the growth of B. platyphylla has been increasingly affected by rising global temperatures. If current trends in global warming persist or intensify, the species' growth and development may face even greater constraints [47]. Given these challenges, exploring the heat tolerance mechanisms of B. platyphylla and enhancing its heat resilience through molecular breeding techniques has become a critical area of research. The recent release of the B. platyphylla genome sequence offers a valuable resource for studying the HSF (Heat Shock Factor) gene family [48]. Using bioinformatics approaches, 21 members of the HSF gene family were identified in B. platyphylla. Their fundamental characteristics were systematically analyzed, including protein multiple sequence alignment, gene structure, conserved domains, phylogenetic relationships, chromosomal distribution, intraspecific and interspecific collinearity, promoter cis-acting elements, and subcellular localization. Expression profiles of BpHSF genes were examined across different tissues (roots, stems, and leaves) and under high-temperature stress. Furthermore, qRT-PCR analysis revealed the expression patterns of genes with significant upregulation under heat stress. The results presented here provide a robust theoretical foundation for advancing the functional exploration of BpHSF genes in B. platyphylla.

2. Results

2.1. Identification and Physicochemical Characterization of the BpHSF Gene Family

Two complementary methods, namely two-way BLAST comparison and Hidden Markov Model (HMM) analysis, were utilized to characterize the HSF gene family in Betula platyphylla. The BLAST analysis identified 27 candidate proteins, while the HMM analysis identified 25. Functional domain verification using SMART software was then performed, excluding proteins lacking or exhibiting incomplete conserved domains. By intersecting the results from both methods, a final set of 21 unique HSF family proteins was identified. These proteins were designated as BpHSF, with the prefix "Bp" representing B. platyphylla, and were further named BpHSFA1a to BpHSFC1c based on their phylogenetic relationships with known Arabidopsis thaliana HSFs. Detailed information, including the physicochemical properties of the 21 BpHSF proteins, is presented in Table 1. The lengths of the BpHSF proteins ranged from 117 amino acids (BpHSFC1c) to 577 amino acids (BpHSFA3), with molecular weights varying from 18.69 kDa (BpHSFC1c) to 64.05 kDa (BpHSFA3). The theoretical isoelectric points (pI) spanned from 4.34 (BpHSFC1c) to 8.78 (BpHSFB1a). Notably, except for BpHSFB1a, BpHSFB4a, BpHSFB4b, and BpHSFB2c, all BpHSF proteins were classified as acidic (pI< 7). The average hydrophilicity values of the BpHSF proteins were negative, indicating a predominantly hydrophilic nature. Additionally, the average instability index of these proteins was calculated to be 55 (above the threshold of 40), suggesting they are likely to be unstable. Subcellular localization predictions revealed that all BpHSF proteins are localized within the nucleus.

2.2. BpHSF Protein Multiple Sequence Comparison

The HSF gene family protein sequences exhibited five characteristic conserved domains arranged sequentially from the N-terminus to the C-terminus: DBD, OD, NLS, AHA, and NES. Evidently, the DNA-binding domain (DBD) remains the most conserved, with all DBD regions sharing a secondary structure comprising three α-helices and four β-sheets. Variations in insertions or deletions are observed across these sequences. For instance, an 8-amino acid insertion is present in BpHSFB1a between α1 and β1(Figure 1B). Additionally, substantial amino acid deletions were identified in BpHSFC1b and BpHSFC1c. Consequently, these deletions may render certain biological functions inactive in B. platyphylla species.

2.3. Chromosomal Localization of BpHSF Genes

Chromosomal localization analysis revealed that the 21 BpHSF genes exhibit distinct distribution patterns across 14 chromosomes, with no BpHSF genes detected on chromosomes 4, 7, 9, and 10 (Figure 2). The gene counts vary among chromosomes, with chromosomes 3 and 8 containing the highest number of BpHSF genes, each hosting four genes. Specifically, BpHSFA9, BpHSFA1b, BpHSFA8a, and BpHSFC1 were located on chromosome 3, while BpHSFA8b, BpHSFA2b, BpHSFB2a, and BpHSFC1c were mapped to chromosome 8. In contrast, chromosomes 1 and 12 harbored only a single BpHSF gene, representing the lowest gene count observed.

Figure 2. Chromosomal localization of the BpHSF gene. B. platyphylla has 14 chromosomes. Scale bars are in Mb, and chromosome numbers are shown above the corresponding chromosomes. Chr: Chromosome.

Figure 2. Chromosomal localization of the BpHSF gene. B. platyphylla has 14 chromosomes. Scale bars are in Mb, and chromosome numbers are shown above the corresponding chromosomes. Chr: Chromosome.

2.4. BpHSF Evolutionary Tree Analysis

To investigate the evolutionary lineage of the HSF gene family in plants, HSF genes from B. platyphylla, Oryza sativa, and Arabidopsis thaliana were compared and analyzed phylogenetically using neighbor-joining phylogenetic trees. This analysis included 67 protein sequences: 21 from B. platyphylla, 25 from Oryza sativa, and 21 from Arabidopsis thaliana. The results classified BpHSF genes into three primary groups, corresponding to three subfamilies (A, B, and C) based on the classification systems of model plants Arabidopsis thaliana and tomato HSF protein families (Table S1). Subfamily A was further divided into nine subgroups (A1~A9), subfamily B into four (B1~B4), while subfamily C included C1 and C2. In B. platyphylla, subfamily A contained eleven BpHSF members, subfamily B contained seven, and subfamily C contained three. The distribution of genes across each subgroup varied, with no members of subfamily A7 detected in B. platyphylla. Additionally, the C2 subgroup appears unique to the OsHSF gene family in Oryza sativa.

Figure 3. Phylogenetic analysis of BpHSF genes was conducted using HSF sequences from B. platyphylla (BpHSF), Arabidopsis thaliana (At), and Oryza sativa (Os). As a result different HSF classes are distinguished by pentagons in various colors, each representing a distinct class.

Figure 3. Phylogenetic analysis of BpHSF genes was conducted using HSF sequences from B. platyphylla (BpHSF), Arabidopsis thaliana (At), and Oryza sativa (Os). As a result different HSF classes are distinguished by pentagons in various colors, each representing a distinct class.

2.5. Protein Motifs and Gene Structure of BpHSF Genes

The MEME tool was employed to identify potential motifs within the BpHSF family, revealing the presence of 10 conserved motifs, designated as motifs 1 through 10 (Figure 4B). Sequence details for each conserved motif were provided in Table S2, with the specific amino acid sequences represented by character strings at each position (Figure S1). Most motifs were approximately 50 amino acids in length, with BpHSFA3 displaying the longest sequence and BpHSFC1c the shortest. To further explore the structural diversity within the BpHSF gene family, a phylogenetic tree was constructed solely for BpHSF proteins, and their gene structures were analyzed (Figure 4A and 4C). Phylogenetic analysis revealed that BpHSF members within the same subfamily share similar motifs, whereas those from different subfamilies show distinct motif patterns, suggesting conservation in protein structure and potential function. In subfamily A, motif patterns were generally consistent, with the exception of BpHSFA8a, which lacks motif 8. Within subfamily B, both BpHSFA5 and BpHSFA4 contain motif 8, while other members exhibit a uniform motif structure. Strikingly, all BpHSF gene family members contain motif 4 and motif 1 universally, and all, except for BpHSFC1b and BpHSFC1c, also contain motif 2 in their amino acid sequences. Some motifs were unique to specific subfamilies: for instance, motifs 10 and 6 are exclusive to subfamily C, while motif 9, with an identical amino acid sequence, was found only in subfamily B.

Exon-intron structure analysis of the 21 BpHSF genes revealed that these genes contain between 1-6 exons and 1-5 introns, with the exception of BpHSFC1c in subfamily C, which lacks intron. The untranslated region (UTR) structure of the BpHSF genes was incomplete: six BpHSF genes possess both 5' UTR and 3' UTR sequences, while four genes (BpHSFA9, BpHSFB1a, BpHSFB4a, and BpHSFC1a) contained only a 5' UTR. While, 11 genes lack UTR sequences altogether.

2.6. Covarianceanalysis of BpHSF

To investigate the co-evolution of BpHSF family genes, we conducted a gene duplication analysis using MCScanX combined with TBtools-II. The analysis identified one tandem duplication pair (BpHSFA8a and BpHSFC1a) on chromosome 3 in the B. platyphylla genome (Figure 5). Additionally, four segmental duplication events were observed, involving the gene pairs BpHSFA1a/BpHSFA1b, BpHSFA6/BpHSFA2a, BpHSFA2b/BpHSFA9, and BpHSFB2a/BpHSFB2b. These findings suggest that gene duplication events have been crucial to the evolution and expansion of the BpHSF family, with segmental duplications particularly contributing to diversity and functional adaptation within the family.

To elucidate the evolutionary trajectory of the BpHSF genome, we conducted a comparative genomic analysis involving six representative species: four dicotyledons (Arabidopsis thaliana, potato, tomato, and poplar) and two monocotyledons (sorghum and Oryza sativa) (Figure 6). The results revealed syntenic relationships between the B. platyphylla genome and all six species, comprising 14 syntenic pairs with A. thaliana, 32 with poplar, 22 with tomato, 19 with potato, 11 with pineapple, and 8 with sorghum. Notably, B. platyphylla exhibited the highest number of collinear gene pairs with poplar (32 pairs), suggesting extensive gene exchange and a closer homologous relationship during their evolutionary history. Although synteny was observed across all representative species, the analyses demonstrated a pronounced trend of higher synteny with dicotyledons compared to monocotyledons, underscoring closer evolutionary affiliations among the dicot specie.

2.7. Analysis of BpHSF Promoter Cis-Acting Elements

Transcription factors are essential in mediating plant responses to biotic and abiotic stresses by regulating various cis-regulatory elements within gene promoter sequences. To explore the potential biological functions and regulatory network of BpHSF, we analyzed cis-regulatory elements within the promoter regions of BpHSF genes (spanning 2000 bp upstream of the start codon). This analysis identified a total of 533 cis-regulatory elements (Table S3), categorized into four groups: light response, stress response, hormone response, and growth and development response (Figure 8A). Light-responsive elements were the most prevalent, comprising 237 elements and accounting for 44.5% of the total, followed by hormone-responsive elements (208, or 39%). Growth and development-responsive elements were the least frequent, potentially indicating a more limited direct role for the BpHSF gene family in plant growth and development.

Evidently, a range of cis-regulatory elements were identified across the BpHSF genes, including the light-responsive G-box, the abscisic acid-responsive element (ABRE), the anaerobic-inducible element (ARE), the drought-inducible MYB-binding site (MBS), and low-temperature-responsive elements (LTRs). These findings suggest a pivotal role for the BpHSF gene family in various environmental stress responses (Figure 8B). However, no cis-elements directly associated with high-temperature responses were detected. Among the 21 BpHSF genes, four (BpHSFA1a, BpHSFA3, BpHSFB1a, and BpHSFB4b) contained more than 30 cis-regulatory elements each. Additionally, stress-responsive elements were particularly abundant in BpHSFB2a, BpHSFB4b, BpHSFA3, BpHSFB1a, and BpHSFA2, highlighting their potential roles in mediating plant responses to environmental challenges such as heat, salinity, and drought stress.

2.8. Secondary and Tertiary Structure Prediction and Subcellular Localisation of the BpHSF Protein

To further investigate the protein structures of BpHSF genes. The results showed that all of these proteins include α-helix, β-turn, irregular coil and extended chain, but the proportion of the secondary structure elements is variable, specifically, the proportion of irregular coil is the highest, varied from 31.12% to 59.1%, followed by α-helix (13.45% to 56.40%), and β-turn is the lowest ranging from 1.86% to 10.73%. α-helix accounted for a high percentage, which can predict that the protein has a high stability (Figure S2). Conserved structural domains of BpHSF proteins were visualized, six 3D structural models of BpHSF proteins were predicted (Table S4). The 3D structural models clearly showed the α-helical and β-folded structures of the DBD structural domains (Figure 9B). Further analysis revealed that BpHSF proteins of the same subfamily have similar three-dimensional structures. The tobacco transient expression analysis revealed that BpHSFA2a, BpHSFB1a and BpHSFC1a were localized within the nucleus (Figure 8B).

2.9. Expression of BpHSF Gene in Different Tissues and Its Response to High Temperature Stressf

The expression profiles of the 21 BpHSF genes were analyzed using the FPKM method based on transcriptome data from three tissues (roots, stems, and leaves) of B. platyphylla (Table S5, Figure 9A). The results revealed significant heterogeneity in the expression levels of BpHSF genes across different tissues, indicating that their expression is regulated in a tissue-specific manner. Conspicuously, three BpHSF genes (BpHSFA8a, BpHSFC1b, and BpHSFC1c) exhibited consistently low expression in all three tissues, suggesting their limited or nonessential roles in these tissues physiological processes. In contrast, other BpHSF genes displayed distinct tissue-specific expression profiles. For instance, BpHSFB1b and BpHSFB2a showed high expression in leaves, while BpHSFA2a, BpHSFA4, BpHSFA5, and BpHSFB1a were highly expressed in roots. Additionally, BpHSFA9, BpHSFA2b, and BpHSFA6 demonstrated elevated expression in stems. These differential tissue-specific expression patterns suggest that BpHSF genes may have important roles in the growth, development, and adaptive responses of B. platyphylla.

To further elucidate the potential roles of BpHSF genes in the heat stress response, we analyzed their dynamic expression patterns using previously acquired transcriptome data (Table S5, Figure 9B). Following 24 hours of heat stress, a subset of BpHSF genes, including BpHSFB2a, BpHSFB2b, BpHSFA2a, BpHSFA5, BpHSFB1b, BpHSFA9, BpHSFA2b, and BpHSFC1a, exhibited significant upregulation, strongly implicating these genes in the early response to high-temperature stress in Betula platyphylla. In contrast, the transcript levels of BpHSFC1c were undetectable, suggesting that this gene is either unresponsive to heat stress or not expressed in B. platyphylla leaves. As the heat stress persisted, several BpHSF genes initially upregulated, such as BpHSFA2a, BpHSFA2b, BpHSFB1b, and BpHSFC1a, underwent significant downregulation. Notably, the pronounced downregulation of BpHSFA2a may reflect an adaptive mechanism by which the plant conserves energy after the initial stress response. Meanwhile, a distinct group of BpHSF genes, including BpHSFA8a, BpHSFA1a, BpHSFA8b, and BpHSFA3, demonstrated marked increases in expression after seven days of continuous heat stress, suggesting their pivotal roles in regulating the late-stage heat stress response.

2.10. BpHSF Expression Validation Through qRT-PCR

Eight heat stress transcription factor genes (BpHSFB2a, BpHSFB2b, BpHSFA2a, BpHSFA5, BpHSFB1b, BpHSFA9, BpHSFA2b, and BpHSFC1a), identified as highly expressed under high-temperature conditions through transcriptome analysis, were further validated using qRT-PCR. The results (Figure 10) revealed that all eight genes reached peak expression levels after 24 hours of heat treatment, followed by varying degrees of decline as stress exposure elevated. Particularly, BpHSFA2a and BpHSFA2b showed significant up-regulation, at 24 hours. BpHSFA2a, in particular, exhibited an exceptional increase, being unregulated 80-100-fold compared to the control, but its expression rapidly returned to baseline levels after 48 hours, indicating a crucial role in the early stages of heat stress response. Further analysis indicated that members of the A subfamily (BpHSFA2a, BpHSFA2b, BpHSFA5, and BpHSFA9) had higher expression levels at 24 hours compared to those of the B and C subfamilies, suggesting an important role for A subfamily genes in the heat stress response. Conversely, BpHSFC1a, a member of subfamily C, exhibited the lowest expression level at 24 hours, being only 1-1.5-fold that of the control, implying a limited role in high-temperature response. These findings are consistent with transcriptomic data and provide a robust experimental foundation for understanding the regulatory dynamics of different HSF subfamilies under heat stress.

3. Discussion

HSFs are a class of transcription factors that play crucial roles in mediating plant responses to heat and various environmental stresses. The HSF gene family has been identified and analyzed in numerous terrestrial plant species. However, the number and composition of HSF gene families vary among species, likely reflecting differences in ecological adaptation, evolutionary history, and environmental stress response mechanisms. HSF gene families have been studied in different plants e.g, eggplant contains 24 SmHSFs [49], walnut has 29 JrHSFs [50], kale-type oilseed rape has 64 HSFs [51], while bread wheat has 78 HSFs [52], and wheat has 82 TaHSFs [53]. Despite these studies, the identification and characterization of the HSF gene family in B. platyphylla have not yet been reported. With the release of the B. platyphylla genome, we have now identified 21 BpHSF genes, which are unevenly distributed across 14 chromosomes. Interestingly, chromosome length did not necessarily correlate with the number of HSF genes. In angiosperms, the HSF gene family is divided into three subfamilies: A, B, and C. Subfamily A includes A1-A9, subfamily B includes B1-B5, and subfamily C contains C1. Based on gene structure and phylogenetic relationships, we classified the 21 BpHSF genes into three subfamilies: A, B, and C. Specifically, subfamily A contains 11 genes, subfamily B contains 7 genes, and subfamily C contains 3 genes. Although the total number of HSF genes in B. platyphylla is similar to that of Arabidopsis, differences were observed in the distribution of specific HSF subclasses. Crucially, the number of HSFA1 subclass members is higher in B. platyphylla, while no members of the HSFA7 subclass were identified. Additionally, there is only one member of the HSFA9 subclass, compared to four in chili pepper [41]. Subfamily C includes three members, of which two belong to class C2, consistent with findings in Oryza sativa (OsHSF) [54] and Minnan (PbHSF) [55], both monocotyledonous plants. Interestingly, the HSFB5 subclass was absent in all of the above studies, further supporting the accuracy of our identification of BpHSF genes in B. platyphylla. Analysis of the physiochemical properties of the proteins encoded by the BpHSF genes indicated that all of the proteins are hydrophilic and predominantly localized in the nucleus. This observation is consistent with findings in carnations [56,57] and turnips [58], suggesting a similar pattern of subcellular localization for HSF proteins.

The highly conserved DNA-binding domain (DBD) consists of approximately 100 amino acid residues across various plant species. Strikingly, the DBDs of BpHSFB1b, BpHSFC1a, and BpHSFC1b are shorter compared to other BpHSFs, possibly due to the genomic or genetic variation in B. platyphylla. Members of the same subfamily within the phylogenetic tree share consistent motif alignments, with motifs 1, 2, and 4 being present in nearly all BpHSFs. The similarity in motif composition among members of the same subfamily suggests that these proteins have conserved structural and functional roles, potentially sharing a common ancestor and being genetically more closely related. In plants, introns are crucial in regulating gene expression. While introns themselves do not encode proteins, they contribute to regulating gene expression, generating alternative mRNA isoforms, and influencing gene stability. Therefore, studying the structural features of BpHSF genes can provide deeper insights into their functions. Our analysis revealed that nearly all BpHSF genes contain between one and five introns, whereas BpHSFC1c lacks introns, which might imply a limited role in gene expression regulation. Additionally, BpHSF genes exhibit significant variation in the location and length of introns, suggesting functional divergence between subfamilies.In the prediction of three-dimensional structures, BpHSF proteins within the same subfamily manifested highly similar structural characteristics. These finding are consistent with studies of HSF genes in papaya [59] and poplar [60], which also showed similar motif compositions and tertiary structures. Previous studies have reported that most HSFs are localized in the nucleus, such as Rhododendron (RsHsf15, RsHsf16, and RsHsf19)[61], peanut (AhHsf20)[62], Salvia divinorum (SmHsf1 and SmHsf7) [63], wheat (TaHsfA2-10) [63], and buckwheat (FtHsf18 and FtHsf19) [65]. Consistent with these findings, our study demonstrated that BpHSFA1a, BpHSFB1a, and BpHSFC1 are all localized in the nucleus.

Studying gene duplication events and synteny relationships provides critical insights into the expansion, functional diversification, and evolutionary history of the BpHSF gene family. In the BpHSF gene family, only one tandem duplication event was identified, whereas five gene pairs exhibited segmental duplication events. These findings suggest that segmental duplication may have been the primary driver behind the expansion of the BpHSF gene family in B. platyphylla. A similar pattern has been observed in passion fruit [66] and hawthorn [67]. Covariance analysis among representative species revealed that B. platyphylla and Populus shared the highest number of covariant gene pairs (32), indicating frequent gene exchanges or shared genomic structures between these species during evolution. This evolutionary analysis of BpHSF genes provides an essential foundation for further functional and evolutionary studies. Cis-regulatory elements (cis-elements) in the promoter regions of gene family members are specific DNA sequences that regulate gene transcription by binding to transcription factors and other regulatory proteins, thereby precisely controlling the level and pattern of gene expression. Promoter analysis revealed that the promoter regions of BpHSF genes harbor different numbers and types of cis-elements, suggesting varied transcriptional regulation among these genes. The identified cis-elements include those responsive to stress, light, hormones, and growth and development, indicating the diverse regulatory roles of BpHSF genes. Interestingly, no heat shock elements (HSE) were detected in the BpHSF promoter regions, which aligns with findings in cannabis [68] and wild grape [69]. However, the presence of HSEs in the barley HSF gene family [70] suggests that the mechanism of HSF gene response to heat stress may vary among plant species. It remains unclear whether the expression of the BpHSF genes in B. platyphylla is directly regulated by heat stress, given the absence of typical HSE elements. Possibly the BpHSF genes responded to heat stress indirectly through other cis-elements. Therefore, further experimental studies are needed to verify the regulation of BpHSF gene expression under heat stress conditions and to elucidate potential alternative regulatory mechanisms involved in environmental stress responses.

Through RNA-seq analysis and qPCR validation, we observed that the BpHSF genes of B. platyphylla reached peak expression after 24 hours of high-temperature treatment, with subfamily A genes (e.g., BpHSFA2a and BpHSFA2b) showing particularly significant responses. Among them, BpHSFA2a exhibited the most pronounced response. These findings are consistent with previous studies in Camellia sinensis (tea tree) [43], Triticum aestivum (wheat) [71], and Capsicum annuum (chili pepper) [72], supporting the idea that HsfA subfamily genes play a crucial role in adopting plants to high-temperature stress by rapidly activating heat-responsive genes [73]. In addition, a recent study found that the transcription factor LpHSFA3, acting as a positive regulator of LpHSFA2a, enhances the heat tolerance of perennial ryegrass [74]. Similarly, DcaHsfA2a and DcaHsfA2b in carnations have been identified as potential candidates for heat tolerance breeding [75]. Under high-temperature stress, Zea mays (maize) genes ZmHsf-01 (A2), ZmHsf-03 (B), ZmHsf-04 (A2), ZmHsf-23 (A6), ZmHsf-24 (A9), and ZmHsf-25 (B) also showed significant upregulation [76]. In soybean, overexpression of the GmHsf-34 (A2) gene increased heat tolerance in Arabidopsis thaliana [77]. While some studies have suggested HSFC1 as a key gene for enhancing plant tolerance to heat and fungicidal stress [78,79]. Our findings indicated that the response of BpHSFC1a was not significant under high-temperature conditions, which diverges from previous reports. For example, PeHSF-C1a in passion fruit has been reported to be involved in the heat stress response [80]. Studies in Populus and Arabidopsis thaliana [81,82] have shown that HSF genes are often activated within 15 to 30 minutes of high-temperature treatment; similarly, genes such as LcHsfA2a, LcHsfA2b, and LcHsfA3 are highly expressed in loblolly pine after just one hour of heat exposure [44]. In contrast, BpHSF genes in B. platyphylla reached their peak expression after 24 hours, suggesting that white birch may employ a unique adaptive mechanism that favors a sustained response to cope with prolonged high-temperature stress. In this study, we found that birch HSF genes exhibited diverse expression responses to high-temperature stress, similar to the signaling and biological processes observed in most plants under heat stress. Thus, the HSF gene family plays an important role in coordinating plant responses to high-temperature stress, acting as a key component of stress signaling (Figure 11).

When ambient temperatures rise to approximately 40°C, intra-cellular calcium ion (Ca²⁺) levels increase, allowing Ca²⁺ to enter the cell and bind to calmodulin (CaM), forming the Ca²⁺/CaM complex. This complex activates various downstream responses, including reactive oxygen species (ROS) signaling and the activation of transcription factors (TFs). Under heat stress, protein denaturation leads to the dissociation of HSF from heat shock proteins (HSPs), resulting in HSF activation. Activated HSF binds to heat shock elements (HSEs) in gene promoters, initiating the expression of HSPs. HSPs assist cells in resisting protein denaturation, maintaining protein stability and function, thereby protecting cells from heat-induced damage. The accumulation of HSPs, along with other stress responses (e.g., ABA-induced stomatal closure), collectively enhances plant tolerance to high-temperature stress, promoting structural and functional regulation and improving overall stress resilience

4. Materials and Methods

4.1. BpHSF Gene Family Identification and Characterization

The complete genome and gene annotation (gff3) file of B. platyphylla were downloaded from Phytozome and analyzed using TBtools-II software [86]. These sequences were compared using TBtools-II, and a bidirectional BLAST search was performed against the full gene sequence database of B. platyphylla, applying a filtering threshold of E ≤ 10⁻⁵. The resulting sequences were considered as potential HSF family members in B. platyphylla.Simultaneously, the Hidden Markov Model (HMM) file for HSF(Pfam no.PF00447) was obtained from the InterPro(https://www.ebi.ac.uk/interpro) database, and HSF-related sequences were screened using HMMER 3.0 software to identify those containing HSF structural domains within the B. platyphylla genome. Protein sequences from B. platyphylla were similarly screened for HSF domains using the Hidden Markov Model. Candidate protein sequences from both methods were refined by removing duplicates and incomplete sequences. The remaining sequences were intersected using Venn plots in TBtools-II. Finally, functional structural domains were verified using the SMART online tool (http://smart.embl-heidelberg.de/), and sequences with missing or incomplete domains were excluded. As a result, 21 B. platyphylla HSF transcription factor family members were identified.The structural and physicochemical properties of the B. platyphylla HSF transcription factor proteins were analyzed using the Expasy ProtParam(https://web.expasy.org/protparam/) tool, which provided information on the number of amino acids(AA), molecular formula(MF), molecular weight(MW), theoretical isoelectric point (pI), , instability index(II), aliphatic index(AI), and Grand average of hydropathicity(GRAVY). Subcellular localization (SL)prediction of the B. platyphylla HSF transcription factor gene family using the online tool WoLF PSORT (https://wolfpsort.hgc.jp/?utm_source=chatgpt.com).

4.2. Multiple Sequence Alignment, Chromosomal Localization and Phylogeny

CLUSTALW(https://www.genome.jp/tools-bin/clustalw) Perform multiple sequence alignment on the BpHSF protein and visualize the structural domains of the BpHSF protein sequence using Jalview software. In order to generate a chromosome map of the HSF gene family in birch, TBtools II software was used to integrate the GFF3 file of birch with the BpHSF gene ID obtained through comparison. ClustalW using default parameters(https://www.genome.jp/tools-bin/clustalw) Multiple sequence alignment was performed on the full-length HSF amino acid sequences of birch, rice, and Arabidopsis Then, MEGA11 software was used for phylogenetic analysis, and the neighbor connection (NJ) method was applied for 1000 guided replicates. A phylogenetic tree generated using iTOL visualization(https://itol.embl.de/).

4.3. BpHSF Conserved Motif Prediction and Gene Structure Analysis

The conserved motifs of gene families were predicted using MEME software (http://memesuite.org/tools/meme) with the motif parameter set to 10, and the MAST XML output file was downloaded for visualization using TBtools-II. The B. platyphylla HSF gene sequences and their corresponding mRNA sequences were extracted, and the exon-intron structures were analyzed and visualized using the GSDS tool (http://gsds.cbi.pku.edu.cn/) [87].

4.4. Analysis of BpHSF Covarianceand Cis-Acting Elements

To investigate the evolutionary expansion and selective pressure acting on the B. platyphylla HSF gene family, gene covariance analysis was conducted. Tandem duplication events among B. platyphylla HSF genes were identified using TBtools-II and MCScanX tools. Subsequently, segmental duplication events and inter-species gene covariance were analyzed using TBtools-II, MCScanX, and BLASTP.

Promoter regions, defined as 2 kb upstream of the transcription start site of all BpHSF genes, were extracted using TBtools-II and submitted to the PlantCARE database (https://bioinformatics.psb.ugent.be/webtools/plantcare/html/) for cis-acting element analysis. All cis-elements in the promoter regions were screened, analyzed for redundancy, and those with significant functional importance were selected and visualized using TBtools-II.

4.5. BpHSF Protein Structure Prediction and Subcellular Localization

SOPMA (https://npsa-prabi.ibcp.fr/cgi-bin/npsa_automat.pl?page=npsasopma.html) is used for protein secondary structure prediction, using default parameters. The SWISS-MODEL database (https://swiss model.expasy.org/) was used to predict the tertiary structure of proteins by homology modeling.

To analyze subcellular localization, transient expression assays were performed in tobacco. Tobacco plants were grown for 4-5 weeks under conditions of 14 h light/10 h dark, 25°C, and 70% relative humidity. Full-length CDSs of BpHSFA2a, BpHSFB1a, and BpHSFC1a were amplified from B. platyphylla leaf cDNA, using primers detailed in Table S6, and cloned into the pCAMBIA2300-eGFP vector. The resulting expression plasmids were transformed into Agrobacterium tumefaciens strain GV3101, with appropriate antibiotics added during culture. The constructed Agrobacterium was then injected into leaves of Nicotiana benthamiana, with the pCAMBIA2300-eGFP empty vector serving as a control. Following injection, tobacco plants were kept in darkness for 12 hours before being transferred to normal light conditions. After 40 hours, fluorescence in the injected areas was observed using a laser confocal microscope (Olympus FV1200), and images were captured. The imaging data were subsequently processed using Adobe Photoshop 2024.

4.6. Plant Materials and Treatments

The experiment was conducted at the field base of Fujian Agriculture and Forestry University. Healthy B. platyphylla seeds were first selected,then soaked, and germinated before being sown in 12-hole seedling trays filled with specialized seedling substrate. The seedlings were kept moist by daily moderate watering to maintain optimal substrate moisture for approximately two weeks. The B. platyphylla seedlings were subsequently grown in a greenhouse at 25°C, with a relative humidity of 60-70%, a photoperiod of 16 h light/8 h dark, and regular moderate watering. After 24 weeks, uniformly grown plants were transplanted from seedling trays and acclimated in Hoagland nutrient solution for one week.

High-temperature stress treatment was conducted using a multi-chamber, multi-temperature gradient thermostat incubator (model: MTI-201B/202B, five chambers, 31 L per chamber). Samples of B. platyphylla roots, stems, and leaves were collected at room temperature (25°C). For the high-temperature treatment, the temperature was set to 40°C, with a humidity of 65% light intensity of 100-110 µmol/(m²·s), and sampling time points at 0, 12, 24, and 48 hours. B. platyphylla leaves collected at 0 h served as the control group. All samples were immediately frozen in liquid nitrogen and stored at -80°C for further analysis. Each treatment was performed in three biological replicates, with 12 B. platyphylla seedlings randomly selected per replicate.

4.7. Expression Profiles Based on RNA-Seq Data

Publicly available Betula platyphylla RNA-seq raw sequencing data were obtained from the NCBI Sequence Read Archive (SRA), including roots (accession no. SRR8953232), stems (accession no. SRR8953231), leaves (accession no. SRR28840994), and samples subjected to high-temperature conditions (accession no. SRP361895) at time points of 6 h, 24 h, 2 days, 7 days, and 14 days, with 0 h (25°C) serving as the control. The raw sequencing data were initially filtered using Fastp v0.23.1.[88] with default parameters to obtain clean data for subsequent analysis. The filtered data were then aligned to the B. platyphylla genome using Hisat2 v2.1.0 [89]. Transcript quantification for each sample was performed using featureCounts v2.0.3 [90]. Gene expression was calculated as the FPKM value for each BpHSF member, log2 transformed, and then visualized as a heatmap using the pheatmap package in R software.

4.8. RNA Extraction and Quantitative Real-Time PCR Analysis

Total RNA was extracted from each tissue sample using the RNA prep PureTotal RNA Extraction Kit for Polysaccharide-Polyphenol Plants (DP441). The RNA was then reverse-transcribed into cDNA using the HiScript III 1st Strand cDNA Synthesis Kit (+gDNA wiper). Primers for RT-qPCR were designed using Primer3Plus (https://www.primer3plus.com) to target eight highly expressed BpHSF genes (BpHSFB2a, BpHSFB2b, BpHSFA2a, BpHSFA5, BpHSFB1b, BpHSFA9, BpHSFA2b, and BpHSFC1a), with all the primer sequences listed in Table S6. 18S rRNA and α-Tubulin were used as internal reference genes [91]. The RT-qPCR reaction mixture consisted of 7 µL ddH₂O, 1 µL cDNA, 2 µL specific primers, and 10 µL SYBR Premix Ex Taq™ II. The RT-qPCR cycling conditions were as follows: pre-denaturation at 95°C for 30 s, followed by denaturation at 95°C for 5 s, annealing at 60°C for 30 s, amplification at 95°C for 5 s, annealing at 60°C for 60 s, and a melting curve analysis (95°C for 5 s, 60°C for 60 s, and 50°C for 30 s). The relative expression levels of BpHSF genes were calculated using the 2⁻^ΔΔCT method [92]. Statistical analyses were performed using one-way ANOVA followed by Duncan's multiple comparison test in SPSS 22.0 software. Graphical representation of data was carried out using GraphPad Prism 10.0.

5. Conclusions

In this study, 21 BpHSF genes were identified in Betula platyphylla, with an uneven distribution across the chromosomes. Comprehensive analyses were performed on the BpHSF genes and their encoded proteins, encompassing phylogenetic classification, conserved motifs, gene structures, subcellular localization, synteny relationships, cis-acting regulatory elements, and expression profiles. The findings revealed dynamic, tissue-specific expression patterns of BpHSF genes. Significantly, the expression levels of BpHSFA2a, BpHSFA2b, BpHSFB1b, and BpHSFC1a were significantly upregulated under high-temperature stress, with BpHSFA2a exhibiting the highest expression after 24 hours of heat treatment. These results suggest that BpHSFA2a plays a pivotal role in the heat stress response of B. platyphylla. We propose that BpHSFA2a functions as a key regulatory gene, contributing to an effective defense mechanism that enhances heat tolerance. Consequently, BpHSFA2a emerges as a promising candidate for improving high-temperature stress resilience in B. platyphylla.