Genome-Wide Characterization of the sHsp Gene Family in Salix suchowensis Reveals Its Functions under Different Abiotic Stresses

Small heat shock proteins (sHsps) function mainly as molecular chaperones that play vital roles in response to diverse stresses, especially high temperature. However, little is known about the molecular characteristics and evolutionary history of the sHsp family in Salix suchowensis, an important bioenergy woody plant. In this study, 35 non-redundant sHsp genes were identified in S. suchowensis, and they were divided into four subfamilies (C, CP, PX, and MT) based on their phylogenetic relationships and predicted subcellular localization. Though the gene structure and conserved motif were relatively conserved, the sequences of the Hsp20 domain were diversified. Eight paralogous pairs were identified in the Ssu-sHsp family, in which five pairs were generated by tandem duplication events. Ka/Ks analysis indicated that Ssu-sHsps had undergone purifying selection. The expression profiles analysis showed Ssu-Hsps tissue-specific expression patterns, and they were induced by at least one abiotic stress. The expression correlation between two paralogous pairs (Ssu-sHsp22.2-CV/23.0-CV and 23.8-MT/25.6-MT) were less than 0.6, indicating that they were divergent during the evolution. Various cis-acting elements related to stress responses, hormone or development, were detected in the promoter of Ssu-sHsps. Furthermore, the co-expression network revealed the potential mechanism of Ssu-sHsps under stress tolerance and development. These results provide a foundation for further functional research on the Ssu-sHsp gene family in S. suchowensis.


Introduction
Plants are often exposed to complex environmental stimuli, so they inevitably suffer various stresses including extreme temperature, flooding, drought, and salinity, leading to restriction of their growth or even cause death. As a sessile organism, plants have evolved sophisticated mechanisms to accommodate these stresses over long lifespans [1]. Heat shock proteins (Hsps) function as molecular chaperones that are involved in the folding, localization, repairing, and degradation of protein molecules in all living organisms and play vital roles in response to diverse stresses [2].

Genome-Wide Identification and Phylogenetic Analysis of Ssu-sHsp in S. suchowensis
To identify sHsp genes in S. suchowensis, the protein sequences of 19 A. thaliana sHsp (AtsHsp) [9] and 37 P. trichocarpa sHsp (PtsHsp) [13] were used as a query to perform basic local alignment search tool algorithms (BLASTP) and Hidden Markov models (HMMs) search against the S. suchowensis genome database [23]. A total of 35 non-redundant Ssu-sHsp proteins were identified, and the information of open reading frame length, protein length, molecular weight (Mw), and isoelectric point (pI) are shown in Table 1. The length of proteins varied from 116 to 350 amino acids, the predicted pI ranged from 4.60 to 9.65, and the Mw from 13.3 to 38.7 kDa (Table 1 and Table S1).
A phylogenetic tree was constructed using the full-length protein sequences of sHsp from S. suchowensis, A. thaliana, and P. trichocarpa to investigate their evolutionary relationship. In S. suchowensis, the sHsp family was mainly classified into four subfamilies, including C (cytosol), CP (chloroplast), PX (peroxisome), and MT (mitochondria), but lacked the endoplasmic reticulum (ER) subfamily in comparison with AtsHsps and PtsHsps [9,13]. The C subfamily of Ssu-sHsp, which had 26 members, was subdivided into five groups (CI-CV): 16 members in CI, two members in CII and CV, one member in CIII, and five members in CIV, respectively; the CP and MT subfamilies individually contained two members; three members were existed in PX subfamily (Table 1 and Figure 1). On the basis of their phylogenetic relationships, potential subcellular localization/group, and their Mw, we named the Ssu-sHsps the "Ssu-sHsp (Mw)-(group)". For example, the Ssu-sHsp member willow_GLEAN_10019104 with Mw 16.5 kD belonging to the CI group was named Ssu-sHsp16.5-CI (Table 1 and  Multiple alignment was performed using Clustal X2.1, and the phylogenetic tree was constructed using MEGA 6.06 by the neighbor-joining (NJ) method with 1000 bootstrap replicates. sHsp from S. suchowensis, P. trichocarpa, and A. thaliana were marked with green, blue, and red circles, respectively.

Gene Structure and Protein Domain Analysis of Ssu-sHsps
To gain further insight into the structural diversity of different groups of Ssu-sHsp genes, the simplified neighbor-joining (NJ) phylogenetic tree of their Hsp20 domains was constructed, as well as the analysis of gene structures and functional domain ( Figure 2). As shown in Figure 2a, all of the Ssu-sHsp members still remained in the same subfamily as that used by the full-length protein sequences. All the members of CI (except Ssu-sHsp23.5-CI) and CII groups were intronless, four Ssu- Figure 1. Phylogenetic relationships of sHsp family members from Salix suchowensis, Populus trichocarpa, and Arabidopsis thaliana. Multiple alignment was performed using Clustal X2.1, and the phylogenetic tree was constructed using MEGA 6.06 by the neighbor-joining (NJ) method with 1000 bootstrap replicates. sHsp from S. suchowensis, P. trichocarpa, and A. thaliana were marked with green, blue, and red circles, respectively.

Gene Structure and Protein Domain Analysis of Ssu-sHsps
To gain further insight into the structural diversity of different groups of Ssu-sHsp genes, the simplified neighbor-joining (NJ) phylogenetic tree of their Hsp20 domains was constructed, as well as the analysis of gene structures and functional domain ( Figure 2). As shown in Figure 2a, all of the Ssu-sHsp members still remained in the same subfamily as that used by the full-length protein sequences. All the members of CI (except Ssu-sHsp23.5-CI) and CII groups were intronless, four Ssu-sHsps (23.5-CI, 25.6-MT, 38.2-CIV and 38.7-CIV) were comprised of two introns, and the other 14 Ssu-sHsps contained one intron. While the Ssu-sHsps with two-introns showed phase 0 in the first intron, but the second intron of Ssu-sHsp23.5-CI showed phase 1 and others Ssu-sHsp (25.6-CP, 38.2-CIV and 38.7-CIV) showed phase 2 in the second intron. These results indicated that the members of the same subfamily shared similar gene structures and intron phases (Figures 1 and 2). The analysis of the functional domain showed that all the Ssu-sHsp proteins have the Hsp20 domain, but their location was different and there existed two types: the Hsp20 domain of Ssu-sHsp23.1, Ssu-sHsp31.8, and the members in CIV subgroup were located near the N-terminal, while the Hsp20 domain of the others near the C-terminal (Figure 2c). Multiple sequence alignment was conducted to identify the variations in Hsp20 domains. As shown in Figure 3, the Hsp20 domain in most of Ssu-sHsps was composed by CRI, CRII, and β6-loop, which was consistent with the report in switchgrass [24] and pepper [12]. However, Ssu-sHsp17.8-CIII lacked the β2, and ten Ssu-sHsps (21.9-CIV, 23 (Figure 3). After searching with the MEME program, ten conserved motifs were identified in Ssu-sHsps ( Figure 4). Thereinto, the CRI contained the motifs 2, 4, and 5, while CRII contained motifs 1 and 3. However, only 19 of the 35 Ssu-sHsp proteins contained these five motifs, revealing the diversification of Hsp20 domain in S. suchowensis sHsps (Figure 4). The analysis of the functional domain showed that all the Ssu-sHsp proteins have the Hsp20 domain, but their location was different and there existed two types: the Hsp20 domain of Ssu-sHsp23.1, Ssu-sHsp31.8, and the members in CIV subgroup were located near the N-terminal, while the Hsp20 domain of the others near the C-terminal (Figure 2c). Multiple sequence alignment was conducted to identify the variations in Hsp20 domains. As shown in Figure 3, the Hsp20 domain in most of Ssu-sHsps was composed by CRI, CRII, and β6-loop, which was consistent with the report in switchgrass [24] and pepper [12]. However, Ssu-sHsp17.8-CIII lacked the β2, and ten Ssu-sHsps  (Figure 3). After searching with the MEME program, ten conserved motifs were identified in Ssu-sHsps ( Figure 4). Thereinto, the CRI contained the motifs 2, 4, and 5, while CRII contained motifs 1 and 3. However, only 19 of the 35 Ssu-sHsp proteins contained these five motifs, revealing the diversification of Hsp20 domain in S. suchowensis sHsps ( Figure 4).

Expansion of the sHsp Family in S. suchowensis
To analyze the expansion and genomic distribution of Ssu-sHsp genes, the chromosomal location and duplication events of the Ssu-sHsp family were analyzed based on the current genomic information from the S. suchowensis genome database. The 35 Ssu-sHsp genes were unevenly mapped onto 14 chromosomes and six scaffolds in the S. suchowensis genome ( Figure 5). Chr09 contained the largest number of five Ssu-sHsps, followed by Chr04 and Chr19 that each contained four Ssu-sHsps ( Figure 5 and Table 1). On the basis of the chromosomal location, phylogenetic relationships, and

Expansion of the sHsp Family in S. suchowensis
To analyze the expansion and genomic distribution of Ssu-sHsp genes, the chromosomal location and duplication events of the Ssu-sHsp family were analyzed based on the current genomic information from the S. suchowensis genome database. The 35 Ssu-sHsp genes were unevenly mapped onto 14 chromosomes and six scaffolds in the S. suchowensis genome ( Figure 5). Chr09 contained the largest number of five Ssu-sHsps, followed by Chr04 and Chr19 that each contained four Ssu-sHsps ( Figure 5 and Table 1). On the basis of the chromosomal location, phylogenetic relationships, and

Expansion of the sHsp Family in S. suchowensis
To analyze the expansion and genomic distribution of Ssu-sHsp genes, the chromosomal location and duplication events of the Ssu-sHsp family were analyzed based on the current genomic information from the S. suchowensis genome database. The 35 Ssu-sHsp genes were unevenly mapped onto 14 chromosomes and six scaffolds in the S. suchowensis genome ( Figure 5). Chr09 contained the largest number of five Ssu-sHsps, followed by Chr04 and Chr19 that each contained four Ssu-sHsps ( Figure 5 and Table 1). On the basis of the chromosomal location, phylogenetic relationships, and gene structures, eight paralogous pairs were identified in the Ssu-sHsp family ( Table 2). Among of them, three paralogous pairs (Ssu-sHsp18.5-CI/18.6-CI, 22.2-CV/23.0-CV, 25.6-MT/23.8-MT) were generated by whole genome duplication (WGD) events and another five pairs (Ssu-sHsp17.7A-CI/17.7B-CI, 17.8-CI/17.7C-CI, 18.1-CI/17.9-CI, 23.5-CI/13.3-CI, 38.2-CIV/38.7-CIV) were generated by tandem duplication events ( Figure 3 and Table 2). The substitution ratio of non-synonymous (Ka)/synonymous (Ks) was used to measure the selection pressures for duplicated gene pairs [25]. In our study, the Ka/Ks ratios of all the eight Ssu-sHsp paralogous pairs were less than 1, suggesting that these Ssu-sHsp paralogous pairs have evolved from purifying selection (  Figure 3 and Table 2). The substitution ratio of non-synonymous (Ka)/synonymous (Ks) was used to measure the selection pressures for duplicated gene pairs [25]. In our study, the Ka/Ks ratios of all the eight Ssu-sHsp paralogous pairs were less than 1, suggesting that these Ssu-sHsp paralogous pairs have evolved from purifying selection ( Table 2).

Expression Patterns of Ssu-sHsps across Tissues and under Various Stresses
In order to gain the information about the spatial and temporal expression pattern of Ssu-sHsps, qRT-PCR was performed to explore the expression levels of Ssu-sHsps in five different tissues including root (R), stem (S), leaf (L), female catkin (FC), and male catkin (MC). As illustrated in Figure 6a

Expression Patterns of Ssu-sHsps across Tissues and under Various Stresses
In order to gain the information about the spatial and temporal expression pattern of Ssu-sHsps, qRT-PCR was performed to explore the expression levels of Ssu-sHsps in five different tissues including root (R), stem (S), leaf (L), female catkin (FC), and male catkin (MC). As illustrated in Figure  6a Notably, the expression levels of all the Ssu-sHsp genes were relatively low in the leaf across the five analyzed tissues (Figure 6a).  To investigate the responses of Ssu-sHsps to various environmental stresses, the expression profiles of Ssu-sHsps were detected using qRT-PCR under heat, drought, salt, cold, and ABA treatments. Among the detected 30 Ssu-sHsps, 29 Ssu-sHsps were promptly and dramatically up-regulated after heat treatments for 1 h (Figures 6b and 7), and most of them maintained high expression levels until the end (24 h) of the heat stress (Figures 6b and 7). Aside from heat stress, several Ssu-sHsps (e.g., Ssu-sHsp18.0A-CI, 18.1-CI, 18.2A-CI, 22.2-CV, and 15.9A-PX) were responsive to both drought and salt stresses. In contrast to the large number of Ssu-sHsps strongly induced by heat, drought, salt, and ABA stresses, only six, nine, and eight Ssu-sHsps were responsive to cold treatment at 6, 12 and 24 h, respectively. It was noted that Ssu-sHsp18.0A-CI were significantly responsive to all the tested treatments at 1 h and 6 h of treatment (Figures 6 and 7).
Among the eight pairs, the expression of five Ssu-sHsp paralogous pairs ( showed highly similar patterns with R 2 > 0.6 across various tissues and stresses, while the other two detected paralogous pairs showed divergent expression (Figure 6c).
To investigate the responses of Ssu-sHsps to various environmental stresses, the expression profiles of Ssu-sHsps were detected using qRT-PCR under heat, drought, salt, cold, and ABA treatments. Among the detected 30 Ssu-sHsps, 29 Ssu-sHsps were promptly and dramatically upregulated after heat treatments for 1 h (Figures 6b and 7), and most of them maintained high expression levels until the end (24 h) of the heat stress (Figures 6b and 7). Aside from heat stress, several Ssu-sHsps (e.g., Ssu-sHsp18.0A-CI, 18.1-CI, 18.2A-CI, 22.2-CV, and 15.9A-PX) were responsive to both drought and salt stresses. In contrast to the large number of Ssu-sHsps strongly induced by heat, drought, salt, and ABA stresses, only six, nine, and eight Ssu-sHsps were responsive to cold treatment at 6, 12 and 24 h, respectively. It was noted that Ssu-sHsp18.0A-CI were significantly responsive to all the tested treatments at 1 h and 6 h of treatment (Figures 6 and 7).
Among the eight pairs, the expression of five Ssu-sHsp paralogous pairs (

Characterization of Cis-acting Elements in the Promoter Regions of Ssu-sHsps
To further investigate the possible regulation mechanism of Ssu-sHsps in S. suchowensis, the ciselements in the promoter regions of Ssu-sHsps were analyzed using the PlantCARE database. A series of stress-related (e.g., heat, drought, low temperature, and anaerobic induction), hormone-related (e.g., ABA, SA, GA, and MeJA), and development-related (e.g., endosperm expression, circadian control, meristem expression, and light response) cis-acting elements were identified (Figures 8 and  S1). In stress-related cis-acting elements, 77.1% of Ssu-sHsps (27 genes) had the HSE (heat responsive elements) which was the most abundant cis-acting element with its number reaching 56. There were 41 ARE (essential for the anaerobic induction) and 41 TC-rich repeats (stress responsiveness) that were located in 12 and 10 Ssu-sHsp promoters, respectively. In hormone-related cis-acting elements, TCA-element (salicylic acid responsiveness), ABRE (abscisic acid responsiveness), CGTCA-motif (MeJA-responsiveness), and GARE-motif (MeJA-responsiveness) were identified in the promoters of 20, 11, 19, and 19 Ssu-sHsps, respectively. Furthermore, large numbers of cis-acting elements related to endosperm expression (Skn-1) and circadian (GCN4_motif and circadian cis-acting element) were detected, implying that sHsp might be also involved in developmental processes in S. suchowensis.

Characterization of Cis-acting Elements in the Promoter Regions of Ssu-sHsps
To further investigate the possible regulation mechanism of Ssu-sHsps in S. suchowensis, the cis-elements in the promoter regions of Ssu-sHsps were analyzed using the PlantCARE database. A series of stress-related (e.g., heat, drought, low temperature, and anaerobic induction), hormone-related (e.g., ABA, SA, GA, and MeJA), and development-related (e.g., endosperm expression, circadian control, meristem expression, and light response) cis-acting elements were identified ( Figure 8 and Figure S1). In stress-related cis-acting elements, 77.1% of Ssu-sHsps (27 genes) had the HSE (heat responsive elements) which was the most abundant cis-acting element with its number reaching 56. There were 41 ARE (essential for the anaerobic induction) and 41 TC-rich repeats (stress responsiveness) that were located in 12 and 10 Ssu-sHsp promoters, respectively. In hormone-related cis-acting elements, TCA-element (salicylic acid responsiveness), ABRE (abscisic acid responsiveness), CGTCA-motif (MeJA-responsiveness), and GARE-motif (MeJA-responsiveness) were identified in the promoters of 20, 11, 19, and 19 Ssu-sHsps, respectively. Furthermore, large numbers of cis-acting elements related to endosperm expression (Skn-1) and circadian (GCN4_motif and circadian cis-acting element) were detected, implying that sHsp might be also involved in developmental processes in S. suchowensis.

The Co-expression Network of Ssu-sHsps
To further predict the relationship between Ssu-sHsp genes and other genes, a co-expression network of sHsp genes was constructed according to their orthologous in P. trichocarpa (Figure 9). The Ssu-sHsps lacked an ER subfamily compared with the PtsHsp family, and Ssu-sHsp did not correspond one-to-one with P. trichocarpa members (Figure 1). Thus, the ER members of PtsHsps were removed in the co-expression network, and the subfamily was used as a unit to study the co-expression relationship. In this network, 970, 235, 116, and 20 genes were co-expressed with sHsp in the C, CP, MT, and PX subfamily, respectively. Moreover, some of them shared the same co-expressed genes. It was noted that the network contained a transcription factor which related to stress response (e.g., HSFA6B, HSFA2, DREB2C), development (e.g., GRF, Floral homeotic protein HUA1, AGAMOUS-like 62), and hormone (e.g., ARF2, ARR2 and AREB3), and several of them were hub genes ( Figure 9 and Table S2). This indicated that the sHsps in S. suchowensis might be function in stress tolerance or developmental processes through cooperation with specific transcription factors or other functional genes.

The Co-expression Network of Ssu-sHsps
To further predict the relationship between Ssu-sHsp genes and other genes, a co-expression network of sHsp genes was constructed according to their orthologous in P. trichocarpa (Figure 9). The Ssu-sHsps lacked an ER subfamily compared with the PtsHsp family, and Ssu-sHsp did not correspond one-to-one with P. trichocarpa members (Figure 1). Thus, the ER members of PtsHsps were removed in the co-expression network, and the subfamily was used as a unit to study the co-expression relationship. In this network, 970, 235, 116, and 20 genes were co-expressed with sHsp in the C, CP, MT, and PX subfamily, respectively. Moreover, some of them shared the same co-expressed genes. It was noted that the network contained a transcription factor which related to stress response (e.g., HSFA6B, HSFA2, DREB2C), development (e.g., GRF, Floral homeotic protein HUA1, AGAMOUS-like 62), and hormone (e.g., ARF2, ARR2 and AREB3), and several of them were hub genes ( Figure 9 and Table S2). This indicated that the sHsps in S. suchowensis might be function in stress tolerance or developmental processes through cooperation with specific transcription factors or other functional genes.

The Co-expression Network of Ssu-sHsps
To further predict the relationship between Ssu-sHsp genes and other genes, a co-expression network of sHsp genes was constructed according to their orthologous in P. trichocarpa (Figure 9). The Ssu-sHsps lacked an ER subfamily compared with the PtsHsp family, and Ssu-sHsp did not correspond one-to-one with P. trichocarpa members ( Figure 1). Thus, the ER members of PtsHsps were removed in the co-expression network, and the subfamily was used as a unit to study the co-expression relationship. In this network, 970, 235, 116, and 20 genes were co-expressed with sHsp in the C, CP, MT, and PX subfamily, respectively. Moreover, some of them shared the same co-expressed genes. It was noted that the network contained a transcription factor which related to stress response (e.g., HSFA6B, HSFA2, DREB2C), development (e.g., GRF, Floral homeotic protein HUA1, AGAMOUS-like 62), and hormone (e.g., ARF2, ARR2 and AREB3), and several of them were hub genes ( Figure 9 and Table S2). This indicated that the sHsps in S. suchowensis might be function in stress tolerance or developmental processes through cooperation with specific transcription factors or other functional genes.   Table S2.

Discussion
As ubiquitous molecular chaperones, sHsps are present in various organisms and play prominent roles in alleviating injuries from diverse abiotic and biotic stresses, especially in heat stress [26]. The features and functions of the sHsp family have been studied in many model plants (e.g., Arabidopsis, rice, and poplar) [9,10,13,27], however, very little is known regarding sHsp genes in S. suchowensis, which is regarded as a short-rotation woody biomass for bioenergy sources and biofuel.
In this study, 35 Ssu-sHsp genes were identified in the S. suchowensis genome, which is less than that in P. trichocarpa (37 PtsHsps), and this number is consistent with the relative genome size of these two specie (~425 Mb of S. suchowensis and~485 Mb of P. trichocarpa) [13]. However, no ER-family member was identified, possibly because of the incomplete genome assembly of S. suchowensis. Gene duplication including WGD and tandem duplication events play a crucial role in expansions of the gene family during species evolution [28]. In the Ssu-sHsp gene family, eight paralogous pairs were identified, including three and five pairs which were generated by WGD and tandem duplication events, respectively, and contributing to the expansion of the CI (five pairs), CIV (one pair), CV (one pair), and MT (one pair) subgroups. Among of them, CI was the largest subgroup containing most expansion genes, which was consistent with reports in other plants [29]. The duplicated gene pairs shared high amino acid sequence identity and similar gene structures, but the expression patterns of duplicated gene pairs (Ssu-sHsp18.5-CI/18.6-CI, 22.2-CV/23.0-CV) were different across various tissues and under stresses, and they kept a lower correlation with R 2 < 0.6 ( Figure 6) indicating that their function might be divergent during the evolution.
The investigation of gene organization showed that 31 Ssu-sHsp genes (~88.5%) contained one or no intron and confirmed that sHsps have less introns in higher plants. Its reason might be that the length of sHsp itself is too short. It has been reported that the intron number likely affects the efficiency of transcription [30]. When comparing the relationship between the intron numbers and their expression patterns of Ssu-sHsps under stress conditions, we found less introns is helpful for Ssu-sHsps to prompt induction under stresses ( Figure 6). Therefore, the gene structure might also play roles in controlling its transcriptional level under stress condition [31].
In the native state, sHsps always exist as multimeric complexes with variable numbers of subunits [5]. The Hsp20 domain is crucial for the formation of multimeric complexes that are important for their chaperone activity [4]. In our study, all the Ssu-sHsp proteins contained the Hsp20 domain, but multiple sequence alignments and motifs analysis showed that the sequences of Hsp20 were diverse among the Ssu-sHsp family and even parts of them lacked the β2-strand, β6-loop, or the β7-strand (Figures 3 and 4). It is reported that the β2-strand is associated with structure dimerization [11,26]. Additionally, the β6-loop that links the β5-strand and β7-strand is related to dimer formation [32]. Thus, the Ssu-sHsp proteins lacking β2, β6, or β7 might be loss-of-function, or even leading to functional diversity which will need to be fully elucidated.
In the co-expression network, floral development-related genes (AGAMOUS-like 62, MADS-box, and Floral homeotic protein HUA1) were co-expressed with C subfamily members ( Figure 9 and Table S2). In our study, 14 Ssu-sHsps had relatively high expression levels in reproductive organs, and ten of them belonged to the C subfamily (Figure 6a), in which most members were induced by at least one or two stresses (Figures 6b and 7). During the reproductive growth, flowers are always sensitive to high temperature and other environmental stresses [33,34]. Therefore, it can be inferred that Ssu-sHsp with high expression levels in reproductive organs might be cooperating with the flower development-related genes in regulating flower development and maintaining the cellular homeostasis between stress and normal conditions. In addition, 13 Ssu-sHsps were detected as being highly expressed in the stem, a similar result is also reported in poplar [13]. The stem is the product of continuous division and differentiation of the cambium cells in the secondary growth. In addition, the activity of the cambium is susceptible to the external environment [35]. Considering the stem-specific expression patterns of Ssu-sHsps and their conserved function, they might be involved in adjustment stem development by protecting the protein folding which plays roles in maintaining the function of the cambium [16].
Stress response analysis showed that the transcript levels of Ssu-sHsps was generally induced by at least one tested stress, especially heat stress. The similar expression patterns of sHsps from different species in response to stress indicated that they might possess conservation function [9,10,13]. In order to analyze the potential regulatory mechanism, the cis-acting elements and the co-expression network of Ssu-sHsps were analyzed. As the binding sites of transcription factors, cis-acting elements in the promoter of the gene determine its expression patterns. In our study, a series of stress-, hormone-, and development-related cis-acting elements were detected in the promoter of Ssu-sHsps. HSE, the target site of Hsf, was abundantly present in the promoters of Ssu-sHsps (Figure 8). In our study, the Hsf family member (HsfA6B and HsfA2) were co-expressed with the C and CP subfamily sHsps (Figure 9). In Arabidopsis, HsfA6B participate in ABA-mediated salt and drought response and is required for heat stress resistance [36], and HsfA2 is regarded as a heat-inducible transactivator that plays a positive role in enhancing the heat tolerance by activating the expression of Hsp or sHsp [37]. In addition, DREB2C was also detected in our co-expressed network and co-expressed with sHsps from the C and MT subfamily ( Figure 9). In Arabidopsis, DREB2C is responsive to heat stress by activating the expression of HsfA3 [38] and the Thermotolerance-related Phytocystatin 4 (AtCYS4) [39], and DREB2C is also involved in ABA signing by regulating the expression of the ABA responsive gene [40]. In our study, several members in the C subfamily were induced by ABA treatment (Figures 6 and 7), this infers that HsfA6B, HsfA2, and DREB2C may directly or indirectly bind to the promoter region of Ssu-sHsps to activate their expression in the ABA dependent pathway. In addition, some other transcription factors and functional genes that are related to stress response and development were identified in the co-expression network ( Figure 9 and Table S2). This suggests that Ssu-sHsps might be involved in stress tolerance and plant development through cooperating with these genes. Further verification will need to be performed in the future.

Phylogenetic Analysis and Chromosomal Location
The sHsp protein sequences from S. suchowensis, A. thaliana, and P. trichocarpa were used for phylogenetic analysis. Clustal X 2.1 software [46] was used for the multiple alignment of sHsp proteins and unrooted phylogenetic trees were constructed using MEGA 6.06 software [47] with the neighbor-joining method. The information of the chromosomal location of all Ssu-sHsp genes was downloaded from Phytozome 12.1 and the chromosomal locations and duplication of the Ssu-sHsp genes were drawn using Circos software [48]. In the analysis of the duplication, the array of two genes located within a 100 kb region were regarded as duplicated genes. DnaSP 5 software (http: //www.ub.edu/dnasp/) [49] was used to calculate the ratio between nonsynonymous (Ka) and synonymous nucleotide substitutions (Ks).

In Silico Analysis of Regulatory Elements in the Promoter Regions of the Ssu-sHsp Gene Family
The upstream regions (2.0 kb) of the translation initiation sites (ATG) of Ssu-sHsp genes were used as promoter fragments. In addition, the elements were located in the promoter sequences using the program PlantCARE online (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/) [50].

Plant Materials, Growth Conditions, and Treatments
The different tissues of cultivated S. suchowensis, including root (R), stem (S), leaf (L), female catkin (FC), and male catkin (MC), were harvested from triennial S. suchowensis. Four-week-old clones of S. suchowensis were water cultured using Hoagland's nutrient solution in a growth chamber under long-day conditions (16 h light/8 h dark) at 25 • C/22 • C (day/night). In order to investigate the expression profiles of Ssu-sHsps under abiotic stresses and hormone treatment, the clones were treated under 150 mM NaCl (for salt stress), 20% (w/v) polyethylene glycol (PEG 6000, for drought stress), 37 • C (for heat stress), 4 • C (for cold stress), or 100 µM abscisic acid (ABA), respectively. Samples were collected after five time periods (0, 1, 6, 12, and 24 h) of each treatment. Fully matured leaves from five individual clones were immediately harvested at each time period of different treatments, and then rapidly frozen in liquid nitrogen and stored at −80 • C for further analysis. The dosages of the abiotic stresses and hormone treatment were determined based on treatments in S. suchowensis [22]. Three biological replicates were performed in different pots for each treatment.

RNA Isolation and qRT-PCR
The RNA isolation, primer design, qRT-PCR reactions, and program were performed according to Li et al [51]. Because of the high similarity of Ssu-sHsp13.3-CI, Ssu-sHsp17.7A-CI, Ssu-sHsp17.7B-CI, Ssu-sHsp17.7C-CI, and Ssu-sHsp17.8-CI, their expression levels cannot be distinguished by qRT-PCR. On the basis of our previous research, EF1a (Elongation factor-1 alpha) was used as reference gene for different tissue types; HIS (Histone superfamily protein H3) and Actin7 for cold treatment and salt treatment, respectively; UBC (Ubiquitin-conjugating enzyme E2) both for heat and drought treatment [51]. The expression levels of Ssu-sHsps under abiotic stresses and hormone treatment were compared with control (0 h). In addition, the significantly induced genes were selected when the value of Log2 (fold change) ≥1. All the primers used in this study are listed in Table S3.

Integrated Network Analysis
The co-expression data of sHsps was obtained from Phytozome (https://phytozome.jgi.doe.gov/ pz/portal.html). The co-expression network was constructed according to Tian et al [52].

Conclusions
On the basis of the S. suchowensis genome sequence, a total of 35 Ssu-sHsp genes were identified in this study. A comprehensive genome-wide analysis was performed to reveal phylogenetic relationship, gene organization, protein structure, chromosomal localization, co-expression network, and expression patterns across the different tissues and under stress conditions. According to the phylogenetic relationships and the predicted subcellular localization, the 35 Ssu-sHsps were classified into four subfamilies (C, CP, MT, and PX). A total of eight paralogous pairs were identified in the Ssu-sHsp gene family, five pairs were generated by tandem duplication and three pairs were generated by whole genome duplication, and all of them had undergone purifying selection during evolution. Multiple sequence alignment and motif analysis showed that the sequences of the Hsp20 domain were diverse among Ssu-sHsps. Expression profiling and co-expression network analyses revealed that Ssu-sHsps might be involved in plant stress tolerance and development through cooperating with a series of transcription factors and functional genes. The results provide a basis for comprehensively understanding the Ssu-sHsp gene family and will be beneficial to investigate the function of Ssu-sHsps.

Conflicts of Interest:
The authors declare that they have no conflicts of interests.