Genomic data from The Cancer Genome Atlas (TCGA) were downloaded and analyzed using the NCI Genomic Data Commons (https://portal.gdc.cancer.gov). The corresponding clinical data were retrieved using the University of Santa Cruz (UCSC) Xena platform. We also searched HCC gene expression data from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/), including the GSE6764, GSE19665, GSE112790 and GSE121248 datasets. The box plots are implemented using the R software package ggplot2. For the ICGC-LIRI-JP cohort, RNA-seq and clinical data were downloaded from the International Cancer Genome Consortium (ICGC) portal. All data were downloaded from publicly accessible databases, and no ethics committee approval was required.

UALCAN (http://ualcan.path.uab.edu/) is a publicly accessible and interactive web portal for analyzing transcriptome data from various types of cancer (Liu B. et al., 2021; Fan et al., 2021). HSF2 mRNA levels in HCC and normal liver tissues and different clinicopathological parameters were explored using UALCAN.

The HCCDB database (http://lifeome.net/database/hccdb/home.html), which contains fifteen public HCC transcriptional expression datasets, is a gene expression atlas of HCC that enables visualization of the findings from multiple bioinformatics analyses. We used the HCCDB database to explore the expression patterns and prognostic value of HSF2 in HCC.

GEPIA2 (http://gepia2.cancer-pku.cn), an online cancer microarray database retrieved from the UCSC Xena server, was used to analyze the relationship between HSF2 and other genes.

CBioPortal (www.cbioportal.org) is an online open access database that contains both sequencing and pathological data on 30 different cancers and allows researchers to interactively analyze multidimensional cancer genomics data. According to the online instructions provided by cBioPortal, the genomic profiles, including the genetic alterations, survival curves and correlations of HSF2, were investigated.

The prognostic value of HSF2 mRNA expression was investigated using the Kaplan-Meier plotter database (www.kmplot.com), which contains gene expression data and survival information for clinical patients with cancer. Patients with liver cancer were divided into two groups (high and low expression) according to the median HSF2 expression level. The OS, recurrence-free survival (RFS), progression-free survival (PFS) and disease-specific survival (DSS) of patients with liver cancer were analyzed by constructing a Kaplan-Meier survival plot and calculating hazard ratios (HRs), 95% confidence intervals and log-rank p-values. The relationship between HSF2 expression and clinicopathological parameters, i.e., sex and clinical stage of HCC, was also examined using the Kaplan-Meier plotter database. p < 0.05 was considered to indicate a statistically significant result.

The association between HSF2 expression and the OS of patients with HCC was evaluated by performing univariate and multivariate Cox regression analyses. A forest plot was generated to show the p value, HR and 95% CI of each clinicopathological parameter using the “forestplot” R package. According to the clinical characteristics, such as age, sex, T stage, N stage and M stage, we generated a nomogram to predict the probability of 1-, 3-, and 5-years OS using the R package “rms”(https://www.rdocumentation.org/packages/rms). The concordance index (C-index) was calculated to assess the predictive accuracy.

The correlations between HSF2 expression and the infiltration of each type of immune cell, including B cells, CD4⁺ T cells, CD8⁺ T cells, neutrophils, macrophages and dendritic cells, were investigated using the TIMER web-based platform (https://cistrome.shinyapps.io/timer/).

The CIBERSORT tool (https://cibersort.stanford.edu/), which contains 547 genes and 22 human immune cell subpopulations, was applied to characterize the immune cell composition based on a validated leukocyte gene signature matrix. Our current analysis examined the proportions of tumor-infiltrating immune cells in HCC using the CIBERSORT algorithm and investigated the correlations between HSF2 expression and immune cell subpopulations. A p-value < 0.05 was set as the criterion to select lymphocytes possibly affected by HSF2 expression.

The HSF2 protein levels in normal liver and HCC tissues from the Human Protein Atlas (HPA) database (https://www.proteinatlas.org/) were analyzed based on IHC staining data.

GO analysis, KEGG analysis and GSEA were applied to investigate the biological functions and potential mechanisms of HSF2 in HCC. GO analysis is a powerful bioinformatics method that annotates and categorizes genes according to biological processes (BPs), cellular components (CCs) and molecular functions (MFs) (Chen et al., 2021). GSEA was used to investigate the potential mechanisms of hepcidin. The GO, KEGG and GSEA results were explored using the R package ClusterProfiler (Chen et al., 2021).

GeneMANIA (http://www.genem ania.org), a flexible web interface, generates a list of genes with similar functions to the query gene and constructs a functional association network. In the present study, a gene-gene interaction network for HSF2 was constructed to evaluate the functions of related genes using GeneMANIA.

HCC cell lines, including HepG2 and Hep3B cells, were obtained from the American Type Culture Collection (ATCC, Manassas, VA, United States). The cells were cultured in Dulbecco’s Modified Eagle’s Medium (DMEM, Gibco, United States) supplemented with 10% fetal bovine serum (FBS, ExCell Bio, Shanghai, China) at 37°C with 5% CO₂ (Tan et al., 2015). HepG2 and Hep3B cells were infected with a lentivirus expressing HSF2-shRNA for 72 h and screened with puromycin (Sigma, St. Louis, MO, United States). Puromycin-resistant clones were amplified and screened through western blotting with an anti-HSF2 antibody (1:500, Cat# K003476P, Solarbio, China). The lentivirus was generated by Genechem Co., Ltd. (Shanghai, China).

A CCK-8 kit (MCE, Shanghai, China) was used to examine cell viability according to the manufacturer’s instructions. HepG2 and Hep3B cells were seeded on 96-well plates containing 100 μl of DMEM at a density of 6 × 10³ cells/well. Colony formation assays were performed as previously described (Liu M. et al., 2021). HepG2 and Hep3B cells (8 × 10² cells/well) were seeded into 6-well plates. For the wound healing assay, monolayers of HepG2 and Hep3B cells (5 × 10⁵ cells/well) with or without HSF2 knockdown were scratched using a pipette tip. After 24 h, representative images of HepG2 and Hep3B cells were captured using a microscope. Cell invasion was examined using 24-well Transwell units with polycarbonate filters (pore size, 8.0 μm; Corning, NY, United States) as described in a previous study (Wang et al., 2017). Briefly, the upper Transwell inserts were coated with 100 μl of Matrigel basement membrane (BD Biosciences, San Diego, CA, United States), and then 100 μl of serum-free DMEM medium was added. HepG2 and Hep3B cells (1× 10⁵) were seeded on the upper Transwell inserts. DMEM (600 μl) containing 5% FBS as a chemoattractant was added to the lower chamber. HepG2 and Hep3B cells were allowed to invade from the upper chamber for 24 h; then, 4% paraformaldehyde was used to fix the cells on the filters, and a 0.05% crystal violet solution was used to stain the cells. Invading cells were counted using a microscope.

The survival results obtained from Kaplan-Meier plotter, HCCDB and GEPIA are displayed as the HR and P or Cox p-values calculated using a log-rank test. Pearson’s correlation coefficient, Spearman’s correction coefficient, Kruskal-Wallis test, Wilcox test and statistical significance were used to evaluate the correlation of gene expression, and the strength of the correlation was determined using the absolute values. Potential immune checkpoint blockade (ICB) response was predicted with the TIDE (Tumor Immune Dysfunction and Exclusion) algorithm using R software ggplot2 (v3.3.3) and ggpubr (0.4.0). The results were considered statistically significant at *p < 0.05, **p < 0.01, and ***p < 0.001.

Gene expression was first analyzed in the TIMER database to validate the differential expression of HSF2 between various tumor tissues and adjacent normal tissues (Figure 1A). HSF2 expression was significantly upregulated in cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), esophageal carcinoma (ESCA), head and neck squamous cell carcinoma (HNSC), HCC, lung squamous cell carcinoma (LUSC) and stomach adenocarcinoma (STAD) and downregulated in breast cancer (BRCA), kidney chromophobe renal cell carcinoma (KICH), kidney renal clear cell carcinoma (KIRC), prostate adenocarcinoma (PRAD), thyroid carcinoma (THCA) and uterine corpus endometrial carcinoma (UCEC) (Figure 1A). We further validated the differential expression of HSF2 in the HCCDB database. HSF2 was markedly overexpressed in HCC samples compared with normal liver samples in 12 different datasets (Figure 1B). Consistently, higher expression of HSF2 was observed in HCC tissues than in normal tissues in other GEO datasets, including GSE6764, GSE19665, GSE112790, and GSE121248 (Figure 1C). A further evaluation of HSF2 expression in HCC was performed using the data directly obtained from TCGA database, and HSF2 expression was substantially increased in HCC tissues (Figure 1D). Moreover, marked upregulation of HSF2 expression was detected in 50 paired tumor samples and adjacent normal samples (Figure 1E). We analyzed IHC images from the HPA database to further explore the expression of HSF2 protein in HCC. HSF2 protein expression was significantly upregulated in liver cancer tissues compared with normal liver tissues (Figure 1F). In summary, the aforementioned results indicated that HSF2 was expressed at higher levels in HCC than in healthy controls.

We next analyzed the relationship between HSF2 expression and clinicopathological parameters using the UALCAN database. HSF2 mRNA levels were significantly increased in HCC samples from both males and females compared with the corresponding normal controls (Figure 2A). After stratification according to individual cancer stages, an increase in HSF2 expression was found in patients with stage 1, 2, 3 and 4 HCC compared to normal controls (Figure 2B). Regarding tumor grade, HSF2 expression was significantly elevated in patients with grades 1, 2, 3 and 4 HCC, and the highest levels of HSF2 were observed in patients with grade 4 HCC (Figure 2C). Patients with HCC presenting more advanced tumor grades tended to express higher levels of HSF2 than those with less advanced tumors. After stratification based on histological subtypes, HSF2 was expressed at higher levels in HCC, fibrolamellar carcinoma and hepatocholangiocarcinoma than in normal controls (Figure 2D). In terms of the nodal metastasis status, increased HSF2 expression was observed in patients with HCC classified as N0 and N1 (Figure 2E). In terms of age, HSF2 was dramatically upregulated in patients with HCC compared to normal controls in different age groups (21–40, 41–60, 61–80 and 81–100 years) (Figure 2F). In addition, upregulation of HSF2 expression was observed in patients with HCC of three races, including Caucasian, African-American, and Asian patients (Figure 2G). Moreover, a significant increase in HSF2 expression was detected in patients with TP53-mutant HCC compared to patients with TP53 wild-type HCC and normal controls (Figure 2H).

We investigated the effects of HSF2 expression on the survival rate to explore the correlation between HSF2 expression and the prognosis of patients with HCC. The Kaplan-Meier plotter analysis revealed that patients with HCC presenting upregulated HSF2 levels experienced shorter OS, RFS, PFS and DSS than those without HSF2 upregulation (Figure 3A). Consistent with the aforementioned results, the survival analysis based on TCGA database indicated that patients with HCC presenting high HSF2 expression had a poor prognosis, including OS, PFS, and DSS (Figure 3B). Moreover, the results from the UALCAN database suggested that HSF2 overexpression was significantly correlated with a shorter OS for patients with HCC of different tumor grades and sexes (Supplementary Figure S1A). Then, we further assessed the prognostic value of HSF2 in patients with HCC using the HCCDB database. As expected, patients with HCC presenting high HSF2 expression in the TCGA-LIHC and ICGC-LIRI-JP cohorts had a poor prognosis (Supplementary Figures S1B,C). In contrast, HSF2 had no effect on OS in normal individuals in these analyses (Supplementary Figures S1B,C).

Based on the differential expression of HSF2 and its significant prognostic value in HCC, we investigated the relationship between the expression of HSF2 and different clinicopathological characteristics of HCC using the Kaplan-Meier plotter database. Specifically, upregulated expression of HSF2 was significantly associated with shorter OS of patients with HCC of stages 1, 2, 1+2, 2+3, 3 and 3+4, grades 2 and 3, and AJCC-T stage 1, 2 and 3 (Figure 3C). In addition, high HSF2 expression was significantly related to a shorter PFS for patients with stages 1, 1+2 and 2+3, grades 2 and 3, and AJCC-T stage 1 and 2 HCC (Figure 3C). Moreover, high HSF2 expression was related to poor OS and PFS rates in patients, regardless of race (white and Asian), sex, alcohol abuse, and hepatitis viral infection (Figure 3C). These data suggested that high HSF2 expression may be an independent predictor of a poor prognosis for patients with HCC.

A ROC curve analysis was conducted to analyze the diagnostic accuracy of HSF2 expression for HCC diagnosis. In a comparison between tumor tissues and matched adjacent normal tissues in the TCGA cohort, the areas under the curve (AUCs) for 1-year, 3-years, 5-years and 8-years OS were 0.701, 0.659, 0.621 and 0.669, respectively (Figure 4A). The AUC values for 1-year, 3-years, 5-years and 8-years DFS were 0.672, 0.589, 0.645 and 0.64, respectively (Figure 4A). The AUC values for 1-year, 3-years, 5-years and 8-years PFS were 0.662, 0.6, 0.676 and 0.669, respectively (Figure 4A). The AUC values for 1-year, 3-years, 5-years and 8-years DSS were 0.76, 0.697, 0.608 and 0.686, respectively (Figure 4A).

We performed univariate Cox and multivariate Cox regression analyses to further evaluate whether HSF2 expression was an independent prognostic factor for patients with HCC. In the univariate Cox regression analysis, HSF2 expression, age, T stage and M stage were significantly associated with the OS of patients with HCC (Figure 4B). In the multivariate Cox regression analysis, HSF2 expression and T stage exhibited obvious correlations with the OS of patients with HCC (Figure 4C).

Based on the multivariate regression analysis, we constructed a novel nomogram model to predict the 1-, 3-, and 5-years OS rates of patients with HCC (Figure 4D). The C index of the prognostic nomogram was 0.718 (Figure 4D). The calibration plots of the nomogram showed good agreement between the actual and nomogram-predicted 1-, 3- and 5-years survival rates (Figure 4E).

The frequency of HSF2 alterations in HCC was analyzed using cBioPortal. A total of 372 patients with HCC were analyzed (HCC, TCGA, and PanCancer Atlas). Genetic variations in HSF2 were detected at an incidence rate of 1.71%, and deep depletion was the most common type (Figures 5A,B). However, the results of the Kaplan-Meier plotter analysis indicated that although a statistically significant difference in OS was not observed between patients with HCC presenting with or without alterations in HSF2, patients with HCC presenting alterations in HSF2 experienced shorter DFS and PFS (Figure 5C). We then examined the methylation alteration of HSF2 in HCC. HCC samples showed elevated DNA methylation level of HSF2 compared with normal individuals using the UALCAN database (Supplementary Figure S2A). Next, we constructed the gene-gene interaction network for HSF2 and the altered neighboring genes using GeneMANIA. The 20 most frequently altered genes were closely correlated with HSF2 (Figure 5D). Importantly, other members of the HSF family, including HSF1 and HSF4, were shown to interact with HSF2 (Figure 5D). The correlations between HSF2 and HSF1 and HSF4 were evaluated using the TCGA database and ICGC dataset (Figures 5E–G). HSF2 expression was strongly associated with HSF1 and HSF4 expression in HCC tumors and normal liver tissues (Figures 5E–G).

We performed GO and KEGG pathway analyses with data obtained from the TCGA dataset to understand the role and molecular mechanism of HSF2 in HCC. The top 50 genes that were positively or negatively associated with HSF2 are shown in Figures 6A and B, respectively. GO analysis is a powerful bioinformatics tool to explore the BPs, CCs and MFs of HSF2. The top 5 enriched BP terms were RNA splicing, peptidyl-lysine modification, covalent chromatin modification, regulation of mRNA metabolic process, and RNA splicing via transesterification reactions (Figure 6C). The top 5 enriched MF terms were ubiquitin-like protein transferase activity, ubiquitin-protein transferase activity, tubulin binding, single-stranded DNA binding, and histone binding (Figure 6D). The CC enrichment analysis showed that HSF2 was significantly correlated with chromosomal regions, spindles, nuclear specks, condensed chromosomes and spliceosomal complexes (Figure 6E). In addition, KEGG pathway analysis revealed that HSF2 was involved in signaling pathways related to carcinogenesis, such as ubiquitin-mediated proteolysis, viral carcinogenesis, endocytosis, platinum drug resistance, homologous recombination, cell cycle, p53 signaling pathway, mismatch repair, Hippo signaling pathway and DNA replication (Figure 6F). These findings indicate that HSF2 plays a role in tumor development and progression.

A GSEA was performed to further investigate the molecular mechanisms influenced by HSF2 in HCC. GO analysis revealed that HSF2 was significantly involved in pathways that included mitotic sister chromatid segregation, RNA splicing, DNA damage checkpoints, cell cycle checkpoints, and the ubiquitin ligase complex (Supplementary Figure S2B). Similarly, among the Reactome terms, several pathways related to the immune system, including the Toll-like receptor 4 (TLR4) cascade and the adaptive immune system, were identified (Figure 7A). Moreover, KEGG enrichment analysis showed that HSF2 was related to the cell cycle, ubiquitin-mediated proteolysis, the Fanconi anemia pathway, basal transcription factors, the spliceosome, etc. In addition, herpes simplex virus 1 infection, hepatitis B, shigellosis, Salmonella infection, and human T-cell leukemia virus 1 infection, which are related to the immune response, were also correlated with HSF2 in the KEGG pathways (Figure 7B). These results indicate a potential relationship between HSF2 and immune response regulation.

A growing number of studies have proven an immunoregulatory effect of tumor-infiltrating immune cells on the development and progression of tumors (Mendillo et al., 2012). Therefore, we analyzed the correlation between HSF2 expression and six major types of infiltrating immune cells, including B cells, CD4⁺ T cells, CD8⁺ T cells, neutrophils, macrophages and dendritic cells, in patients with HCC (Figure 8A). HSF2 expression showed a significant positive correlation with the infiltration levels of B cells, CD8⁺ T cells, CD4⁺ T cells, macrophages, neutrophils, and dendritic cells in HCC (Figure 8A). We estimated the association between HSF2 and immune cell infiltration using CIBERSORT to further assess the relationship between HSF2 and the tumor microenvironment (TME). HSF2 expression was significantly and positively correlated with the infiltration levels of native B cells and dendritic cells but negatively correlated with the infiltration levels of activated NK cells and Treg cells in HCC (Figure 8B). After stratification according to HSF2 expression, patients with HCC were separated into high-expression and low-expression groups. The percentage abundances of tumor infiltrating immune cells in each sample and different types of immune cells identified using TIMER are shown in different colors (Figure 8C). Increased percentages of infiltrating B cells, CD4⁺ T cells, neutrophils, macrophages and dendritic cells were observed in the HSF2 high-expression group compared with the HSF2 low-expression group (Figure 8D). These findings indicate that HSF2 may recruit immune cells to the TME of HCC.

We further explored the relationships between HSF2 expression and various tumor-infiltrating immune cells in HCC using the TIMER and GEPIA databases by analyzing the different immune cell gene markers. The correlation was adjusted for tumor purity due to its effect on the immune cell infiltration analysis. HSF2 expression was significantly associated with the expression of most gene markers of T cells, CD8⁺ T cells, B cells, monocytes, tumor-associated macrophages (TAMs), M1 and M2 macrophages, neutrophils, NK cells and dendritic cells in patients with HCC after adjusting for tumor purity (Table 1). Consistent with the results obtained from TIMER, the GEPIA results indicated that HSF2 exhibited positive and significant correlations with gene markers of these immune cells in HCC (Table 2).

Furthermore, HSF2 expression was also significantly correlated with the infiltration of subsets of T cells, including Th1, Th1-like, Th2, Treg, resting Treg cells, effective Treg cells, effective T cells, naïve T cells, effective memory T cells, resistant memory T cells, and exhausted T cells, according to the TIMER database (Table 3). In addition, the expression of many immune checkpoint genes, including CD274, CTLA4, HAVCR2, LAG3, PDCD1, PDCD1LG2 and TIGHT, was much higher in the HSF2 high-expression group than in the HSF2 low-expression group (Figure 8E). Interestingly, HSF2 expression was positively associated with multiple immune checkpoints, including PD-L1, PD-1 and CTLA-4 (Supplementary Figure S3A). The correlations between HSF2 and PD-L1, PD-1 and CTLA-4 expression were further confirmed in both the TIMER and GEPIA databases (Supplementary Figures S3B,C), revealing that increased expression of HSF2 was associated with immunosuppression in HCC. Furthermore, we investigated the effect of HSF2 expression on immune checkpoint blockade therapy in both TCGA and ICGC datasets and found that the HSF2 high-expression group exhibited a higher TIDE (tumor immune dysfunction and exclusion) score, which indicates a worse response to immunotherapy (Figure 8F).

Because a high HSF2 expression level is associated with a poor prognosis and immune cell infiltration, we then investigated the relationship between HSF2 expression and the prognosis of patients with HCC in related immune cell subgroups using the Kaplan-Meier plotter database. Interestingly, patients with HCC presenting high HSF2 expression and either increased or decreased infiltration of B cells, CD4⁺ T cells, CD8⁺ T cells, macrophages, NK cells, Th1, Th2 and Treg cells had a shorter OS (Supplementary Figures S4A-H). Moreover, patients with HCC presenting high HSF2 expression and either increased or decreased infiltration of B cells, macrophages, NK cells, Th1 cells, Th2 cells and Treg cells had a shorter RFS (Supplementary Figures S5A-H).

We genetically knocked down (KD) HSF2 in HepG2 and Hep3B cells to verify the oncogenic role of HSF2 in HCC. The knockdown efficacy was examined using western blotting (Figure 9A). The viability of HSF2-KD HepG2 and Hep3B cells was significantly decreased (Figure 9B). The migration and colony formation abilities of HepG2 and Hep3B cells were also dramatically inhibited when HSF2 was silenced (Figures 9C,D). Matrigel-coated Transwell assays revealed that HSF2 knockdown significantly reduced the invasion ability of both HepG2 and Hep3B cells (Figure 9E). In addition, the correlations between HSF2 expression and EMT-related and invasion-related markers were investigated using the TIEMR database. HSF2 expression was correlated with CDH2 (N-cadherin), VIM (vimentin), SNAI1 (Snail), SNAI2 (Slugc), MMP2 and MMP9 (Figure 9F). Taken together, HSF2 increased the viability, proliferation and metastasis of HCC cells.