Introduction

Hepatocellular carcinoma (HCC) is a prevalent type of liver cancer with a rising incidence of illness and death.¹ Every year, over 500,000 individuals across the globe are newly diagnosed with HCC. According to predictions, HCC is expected to result in 1 million fatalities annually by 2030, making it the third deadliest form of cancer.^2,3 HCC progresses swiftly and surreptitiously, with inconspicuous initial symptoms and a dearth of efficient screening methods, resulting in many patients developing into intermediate or advanced stages upon first diagnosis. At present, alpha-fetoprotein (AFP) serves as the prevailing biomarker for the prognosis and diagnosis of HCC. However, its clinical utility is restricted due to inadequate specificity and sensitivity.^4,5 In addition, HCC is a highly heterogeneous disease that can impact the precision of prognostic forecasting.^6,7 The presence of abnormalities within the tumor immune microenvironment (TIME) has been identified as a significant factor leading to the heterogeneity of HCC.⁸ Given the significant heterogeneity of HCC, the current performance of conventional models in predicting outcomes remains unsatisfactory, thus there is a pressing requirement to develop a new prognostic model to improve the accuracy of prognosis prediction for patients with HCC.

HCC is a hypervascular tumor with significant abnormal angiogenesis, which also promotes tumor growth, invasion and metastasis.⁹ Hypoxia is a critical cause of abnormal angiogenesis.^10,11 Among the many key hub genes of hypoxia-angiogenesis, vascular endothelial growth factor (VEGF)/vascular endothelial growth factor receptor (VEGFR) has been confirmed as the therapeutic target for HCC, targeting to inhibit angiogenesis, of which the typical clinical application is sorafenib in the treatment of advanced HCC.^12,13 Moreover, numerous studies have also identified multiple other genes as biomarkers or regulators of hypoxia-angiogenesis. For instance, hypoxia-inducible factor 1 subunit alpha (HIF-1α), fibroblast growth factor 2 (FGF2) and matrix metalloproteinase 9 (MMP9), which can promote angiogenesis, play crucial roles in HCC invasion and metastasis.^14–16 Conversely, there are also some hypoxia-angiogenesis-related genes, such as angiopoietin-like protein 1 (ANGPTL1) and endoglin (ENG), which can suppress the angiogenesis and metastasis of HCC.^17,18 Notably, hypoxia-angiogenesis is strongly related to the tumor immune microenvironment, and hypoxia-triggered angiogenesis can lead to immunosuppression.¹⁹ Studies have found that hypoxia can facilitate tumor angiogenesis, attract macrophages to chemotactic enter tumor tissues, and subsequently produce a variety of inflammatory cytokines and proangiogenic growth factors, which eventually exacerbate tumor microenvironmental inflammation and accelerate HCC progression.^20,21 Moreover, hypoxia can also induce the production of hypoxia-inducible factor (HIF), which facilitates the maintenance and recruitment of pro-tumor immune cells (such as M2 macrophages and regulatory T cells), suppresses anti-tumor immune cells, and forms an immunosuppressive microenvironment, leading to immune evasion and aggravating HCC.^22–25 Importantly, hypoxia-angiogenesis-related genes also affect the immune microenvironment of HCC, such as HILPDA by promoting the release of inflammatory factors from HCC cells, consequently suppressing the cytotoxicity of natural killer (NK) cells and facilitating HCC metastasis.²⁶ Thus, hypoxia-angiogenesis-related genes might exert an important role in the prognosis of patients with tumors. Nevertheless, the impact of hypoxia-angiogenesis-related genes on the prognosis of HCC patients remains unclear.

In the past few decades, immunotherapy has made tremendous progress in the treatment of HCC, especially immune checkpoint inhibitors (ICIs), with the anticipation of enhancing the therapeutic outcome for patients.^27,28 Nevertheless, as a result of abnormalities within TIME, not all individuals attain the expected result following immunotherapy, thereby heightening the prognostic uncertainty of patients.^29–31 Importantly, numerous studies have shown a strong relationship between hypoxia-angiogenesis and cancer immunotherapy.^19,32 Therefore, we constructed a hypoxia-angiogenesis-related signature HARG-score to forecast the prognosis of HCC patients and reveal its association with immune infiltration and immunotherapy efficacy.

In the current study, we developed a 5-gene prognostic model HARG-score based on differentially expressed HARGs. Experimental data suggested that the key gene HILPDA contributes to HCC by promoting angiogenesis and affecting CTLA4 expression. The present study comprehensively evaluated the connection between HARG-score and prognosis, clinical characteristics, and immune infiltration of patients with HCC, in which HILPDA may be a promising prognostic biomarker and therapeutic target, offering novel insights to improve the prognosis and treatment of HCC.

Materials and Methods

Data Sources

Gene expression information, clinical characteristics, and survival data for patients within 3 HCC cohorts (TCGA-LIHC, n = 365; GSE14520, n = 242; LIRI-JP, n = 208) were obtained from 3 public databases. Subsequent analyses excluded patients with zero survival time in 3 HCC cohorts. Moreover, the R package “IMvigor210CoreBiologies” was employed to download the immunotherapy dataset known as Imvigor210.³³ Meanwhile, 73 HARGs (score >8) used in this study were obtained from the GeneCards database (Table S1). The human data involved in this study were reviewed with the approval of Scientific Research and New Technology of Wannan Medical College Yijishan Hospital IRB (approval number: 2024–102). Figure 1 depicts the workflow diagram of our study.

Figure 1 The workflow diagram of our study.

Identification of Differentially Expressed HARGs

Between 365 HCC samples and 50 non-tumor samples within the TCGA-LIHC cohort, differentially expressed genes (DEGs) were identified utilizing the R package “limma” based on a false discovery rate (FDR) threshold of 0.05, followed by the selection of differentially expressed HARGs. The R package ‘ggplot2’ was used to draw the volcano plot for DEGs. The R package ‘pheatmap’ was utilized to create a heatmap plot of HARGs. Copy number variation (CNV) and somatic mutation data of the TCGA-LIHC cohort were downloaded from the TCGA database. The R packages ‘maftools’ and ‘RCircos’ were used to assess somatic mutations and CNV alterations.

Construction of Hypoxia-Angiogenesis-Related Prognostic Model HARG-Score

In the training cohort (TCGA-LIHC), we identified the differentially expressed HARGs related to overall survival using univariate Cox regression analysis. These HARGs were selected based on a p-value threshold of 0.05. Subsequently, the least absolute shrinkage and selection operator (LASSO) Cox analysis utilized the R package “glmnet” to eliminate overfitting among prognostic genes and conduct penalty parameter tuning through 10-fold cross validation, aiming to reduce the number of prognostic genes. Finally, the hypoxia-angiogenesis-related risk score for every patient with HCC was determined by summing the product of the gene expression levels and the LASSO Cox regression coefficients. Here, we developed a highly effective prognostic prediction model, HARG-score, designed to achieve the best predictive outcome using overall survival within the TCGA-LIHC cohort.

The HARG-score calculation formula was given by the following equation: HARG-score = Σ(Expi * Coefi), where Expi and Coefi represent the expression levels and regression coefficients of prognostic genes, respectively. The 365 patients from the TCGA-LIHC cohort were categorized into 2 groups based on the median HARG-score: high-risk (HARG-score > median value) and low-risk (HARG-score < median value). Subsequently, the predictive ability of HARG-score was evaluated utilizing receiver operating characteristic (ROC) curve as well as Kaplan–Meier survival analysis. Next, we conducted principal component analysis (PCA) utilizing the R package “factoextra”. Likewise, the testing groups (LIRI-JP and GSE14520) also categorized HCC patients into high- and low-risk groups and performed ROC curve and Kaplan–Meier survival analysis.

Clinical Correlation Analysis of the Prognostic Model HARG-Score

To evaluate the relationship between HARG-score and clinical characteristics, the Chi-square test was performed. In 3 cohorts, univariate and multivariate Cox regression analyses were conducted to assess whether HARG-score was a prognostic factor independent of other available clinical characteristics. When the p-value of hypoxia-angiogenesis-related signature HARG-score and other clinical characteristics was less than 0.05, it was identified as an independent prognostic factor.

Development and Validation of the Nomogram

A nomogram model was established utilizing the R packages “survival” and “rms” to forecast the survival time of patients with HCC according to HARG-score and various clinical characteristics. The nomogram assigned a score to each factor, and these scores for every patient were summed to calculate a total score, which predicted the corresponding overall survival time. The C-index was utilized to assess the discriminatory capacity of the nomogram. ROC curve was utilized to evaluate the prognostic predictive ability of the nomogram for 1-, 4-, and 5-year survival. Moreover, we performed calibration plots for the nomogram to illustrate the agreement between the actual observed and predicted survival at 1-, 4-, and 5-year.

Evaluation of Immune Infiltration

The cell-type identification by estimating relative subsets of RNA transcripts (CIBERSORT) algorithm was utilized to assess the proportion of immune cells infiltrating the tumor immune microenvironment, which allowed us to estimate the abundance of 22 infiltrating immune cells within the high- and low-risk groups. Kaplan–Meier survival analysis was utilized to assess the prognostic predictive significance of these infiltrating immune cells. Then, according to the R package “Gene Set Variation Analysis (GSVA)”, we utilized single-sample gene set enrichment analysis (ssGSEA) to calculate the enrichment scores of 13 immune functions as well as 16 immune cells for patients with HCC. Furthermore, the relationships between HARG-score and the enrichment scores of 13 immune functions as well as 16 immune cells were additionally examined.

Immunotherapy Efficacy

The IMvigor210 cohort was utilized to assess the efficacy of immunotherapy between the high- and low-risk groups.

Assessment of the Protein Expression Levels of Prognostic Genes Using the Human Protein Atlas (HPA) Database

To further evaluate the protein expression levels of prognostic genes in normal liver as well as HCC tissues, immunohistochemical (IHC) staining was obtained from the HPA database.

Correlation Analysis Between HILPDA and Immune Checkpoints

The TCGA-LIHC cohort was used to assess the correlation between HILPDA and immune checkpoints, including programmed cell death protein 1 (PD1), CTLA4, lymphocyte-activation gene 3 (LAG3) and T-cell immunoreceptor with Ig and ITIM domains (TIGIT).

Cell Culture and Transfection

Human liver cancer cell line LM3 cells, normal liver cell line LO2 cells, and human umbilical vein endothelial cells (HUVECs) were obtained from Procell Life Science & Technology Co., Ltd. (Wuhan, China), and cultivated in medium containing 9% fetal bovine serum (Cell-Box, Hong Kong, China) and 1% double antibody (UU-Bio, Suzhou, China) in an incubator at 37 °C and 5% CO₂.

Silencing the HILPDA gene using small interfering RNA (siRNA). HILPDA siRNA (siRNA-1: sense: 5’-UGUUGAACCUCUACCUGUUTT-3’, antisense: 5’-AACAGGUAGAGGUUCAACATT-3’; siRNA-2: sense: 5’-CCAGAAGCCAACUAGCCAATT-3’, antisense: 5’-UUGGCUAGUUGGCUUCUGGTT-3’) and control siRNA (siRNA-NC: sense: 5’-GGCUCUAGAAAAGCCUAUGC-3’, antisense: 5’-GCAUAGGCUUUUCUAGAGCC-3’) were obtained from General Biol (Chuzhou, China) and transfected into LM3 cells by Lipofectamine^TM 3000 (Invitrogen, Thermo Fisher Scientific, USA). Specifically, 5 μL of siRNA at a concentration of 20 μM was added to 250 μL of serum-free medium, and 5 μL of Lipofectamine^TM 3000 was added to another 250 μL of serum-free medium, mixed thoroughly and left to stand for 5 minutes, then the two were merged and mixed fully, and left to stand for another 20 minutes. Next, the mixed siRNA and Lipofectamine^TM 3000 were added to 1500 μL of antibiotic-free medium to transfect LM3 cells. After 6 hours, replace with fresh medium and continue to culture for 48 hours.

qPCR Analysis

Trizol reagent (Absin, Shanghai, China) and RNA extract solution (Solarbio, Beijing, China) were utilized to extract total RNA. For cDNA synthesis, total RNA was reverse-transcribed utilizing SweScript All-in-One RT SuperMix for qPCR kit (Servicebio, Wuhan, China). Following the guidelines of SYBR Green Pro Taq HS premixed qPCR kit (Accurate Biology, Hunan, China), primers were added to make the sample system to be tested 20 μL. The 20 μL of reaction mixture was composed of 1 μL of diluted cDNA, 2 μL of primer mix (reverse and forward primers), 7 μL of RNase free water and 10 μL of 2X SYBR^® Green Pro Taq HS Premix (ROX plus). Quantitative real-time PCR was performed using LightCycler^® 480 (Roche Diagnostics GmbH) with an initial step at 95 °C for 30 seconds, followed by 40 cycles of two-step PCR at 95 °C for 5 seconds and 60 °C for 30 seconds, and calculated using 2^−ΔΔCt. All primers used for analysis were commercially synthesized by General Biol (Chuzhou, China) (Table S2). Quantitative analysis was normalized using glyceraldehyde-3-phosphate dehydrogenase (GAPDH).

Western Blot

Protein was obtained from LM3 cells with 100 μL RIPA lysis buffer (Servicebio, Wuhan, China) and transferred to EP tube 10 minutes later, and added with the same ratio (100 μL) of 2×PAGE protein buffer (Servicebio, Wuhan, China), and then boiled at 100 °C for 10 minutes after shaking and mixing. The 5% concentrate glue and 15% separation glue were prepared for protein separation through electrophoresis, followed by the transfer of proteins onto the PVDF membrane (Beyotime, Shanghai, China). The membrane was closed with 5% skim milk for 3 hours, washed 3 times with TBST, and then placed in HILPDA and CTLA4 antibodies (1:1000 diluted) (Santa Cruz, USA), HIF-1α, ERK and P-ERK antibodies (1:1000 diluted) (Zenbio, Chengdu, China) at 4 °C overnight. After 3 washes with TBST, the membrane was closed with a secondary antibody (1:6000 dilution) (Boster Bio, Wuhan, China) for 2 hours at room temperature, and finally visualized utilizing enhanced chemiluminescence (Tanon, Shanghai, China).

Transwell Assay

LM3 cells were suspended in medium without serum, the density of cells was modified to 1–10×10⁵ cells/mL, 200 μL of cell suspension was added to the upper chamber of the 24-well transwell plates with 8 μm-pore (Servicebio, Wuhan, China), and 600 μL of medium containing 5% serum was added to the lower chamber and cultivated in the incubator for 24 hours. Next, the cells were fixed using 600 μL of 4% paraformaldehyde (Biosharp, Guangzhou, China) for 25 minutes, and then stained with 0.1% crystal violet (Servicebio, Wuhan, China) for 20 minutes, and scraped off the cells in the upper layer of the chamber. Finally, photos were taken with microscope.

Wounding Healing Assay

LM3 cells in the phase of logarithmic growth were cultured in 6-well plates at a count of 5×10⁵ cells/well. Upon the cell density increased to 90%–95%, a 200 μL pipette tip was utilized to scratch the cells in every well and recorded with microscope. After incubation at 5% CO₂ and 37 °C for 24 hours, the healing of the scratched cells in 6-well plates was recorded by microscope.

EdU Assay

LM3 cells in the phase of logarithmic growth were cultured in 6-well plates and inoculated on coverslips for 24 hours. 2×Edu working solution (Beyotime, Shanghai, China) was prepared and added to 6-well plates, and continued to incubate for 2 hours. The cells were fixed for 15 minutes with 1mL of 4% fixative solution (Biosharp, Guangzhou, China) per well and then incubated with 1mL of permeabilization solution (PBS containing 0.3% TritonX-100) (Biosharp, Guangzhou, China) for 15 minutes. Next, 250 μL Click reaction solution was introduced into every well and incubated for 30 minutes away from light. Finally, 250 μL of 1×Hoechst 33342 solution prepared with PBS at 1:1000 ratio was introduced into every well and incubated for 10 minutes away from light, and fluorescence detection was performed and analyzed.

Cell Viability Assay (Cell Counting Kit-8, CCK-8)

LM3 cells were collected in the phase of logarithmic growth, inoculated into 96-well plates, 100 μL of cell suspension was introduced into every well, and the density was modified to 3000 cells per well, and cultivated in a 5% CO₂ incubator at 37 °C. 96-well plates were cultivated for 0 hours, 24 hours, 48 hours and 72 hours, and 10 μL of CCK-8 solution (Beyotime, Shanghai, China) was introduced into every well. Finally, after 2 hours of incubation in 96-well plates, the absorbance of the cells in every well was detected at 450 nm by multi-scan spectrum.

Enzyme Linked Immunosorbent Assay (ELISA)

HIF-1α and VEGF in the supernatant of LM3 cells were detected using ELISA kits (Spbio, Wuhan, China) according to the manufacturer’s protocols.

Tube Formation Assay

After transfecting HILPDA siRNA for 24 hours, the LM3 cells supernatant was co-cultured with HUVECs for 24 hours, and then HUVECs (2×10⁵ cells/well) were inoculated into matrix-coated (100 μL/well, Corning, USA) 48-well plates. After 6 hours, we photographed the shape of the tube with a microscope and calculated the grid area to determine the extent of tube formation.

Statistical Analysis

Significant differences between variables were calculated using the Chi-square test as well as Student’s t-test. To assess statistical differences between two and multiple groups of clinical baseline data, the Kruskal–Wallis as well as Wilcoxon tests were conducted. Statistical differences in Kaplan–Meier survival curves were evaluated utilizing the Log rank test. A p-value less than 0.05 denoted a statistically significant difference. R 4.1.3 and Prism 8 software were utilized for all statistical analyses conducted in our study.

Results

Genetic and Expression Alterations of HARGs in HCC

To screen for differentially expressed HARGs in the TCGA-LIHC cohort, we performed differential expression analysis of genes between 365 HCC samples and 50 non-tumor samples, and then determined the differential expression of 73 HARGs (Figure S1A). The DEGs between HCC samples and non-tumor samples were visualized utilizing a volcano plot (Figure S1B). The intersection of DEGs and 73 HARGs was taken, and a total of 54 differentially expressed HARGs were obtained (Table S3). On a genetic perspective, the incidence of somatic mutations in 54 differentially expressed HARGs was examined, and the findings indicated that 154 (42.31%) out of 364 HCC samples exhibited mutations in HARGs (Figure S1C). Catenin beta 1 (CTNNB1) exhibited the top mutation frequency of 24%, trailed by SET-domain-containing 2 (SETD2) and hepatocyte growth factor (HGF). Figure S2 illustrated the relationship network of differentially expressed HARGs. Subsequently, we assessed somatic copy number alterations within 54 differentially expressed HARGs and discovered that the majority of HARGs exhibited copy number alterations. Among differentially expressed HARGs, CNV of angiopoietin 1 (ANGPT1), aryl hydrocarbon receptor nuclear translocator (ARNT), prostaglandin-endoperoxide synthase 2 (PTGS2), Egl-9 family hypoxia inducible factor 1 (EGLN1), vascular endothelial growth factor A (VEGFA) and tumour necrosis factor (TNF) were commonly enhanced, while those of angiopoietin 2 (ANGPT2), caspase-3 (CASP3), vascular endothelial growth factor C (VEGFC), plasminogen (PLG) and Cbp/P300-interacting transactivator with Glu/Asp-rich carboxy-terminal domain 2 (CITED2) was reduced (Figure S1D). Figure S1E exhibited the position of CNV alterations within HARGs on their respective chromosomes.

We performed further analysis to explore the potential connection between CNV alterations and differential expression of HARGs. Notably, HARGs exhibiting CNV gain, for instance, ANGPT1, ARNT, EGLN1 and VEGFA, were markedly increased in HCC samples, whereas HARGs exhibiting CNV loss, such as PLG and CITED2, were remarkably increased in non-tumor samples, indicating that CNVs might modulate HARGs expression. Nevertheless, compared to non-tumor samples, many HARGs exhibiting CNV gain, such as PTGS2 and TNF, were markedly down-regulated in HCC; many HARGs exhibiting CNV loss, such as ANGPT2 and CASP3, were remarkably upregulated in HCC. Hence, while CNV can illuminate the majority of alterations in HARGs expression, it is not a decisive factor in the regulation of gene expression.³⁴ Other important factors that regulate gene expression such as DNA methylation as well as histone modification.^35,36 The above analysis showed that the expression levels as well as genetic alterations of HARGs were significantly different between HCC samples and non-tumor samples, indicating a potential role of HARGs in HCC progression.

Construction and Validation of the Prognostic Model HARG-Score

The prognostic model was established utilizing 365 HCC patients in the TCGA-LIHC cohort. Firstly, a univariate Cox analysis was conducted to determine the correlation between 54 differentially expressed HARGs and overall survival outcomes in HCC patients. Out of the 54 HARGs, 29 prognostic genes were screened (Table S4). Figure 2A exhibited the correlation network of the regulatory relationships of 29 prognostic genes and their prognostic significance in HCC patients. We conducted LASSO Cox analysis utilizing the 29 prognostic genes in the TCGA-LIHC cohort, and to avoid overfitting between genes, we utilized a minimized lambda (Figure 2B). Based on minimized lambda, a hypoxia-angiogenesis-related risk score consisting of HILPDA, EPO, SLC2A1, PSMA7 and CREB1 was constructed, named HARG-score. The HARG-score calculation formula was given by the following equation: Risk score = (0.10833192 * HILPDA) + (0.06310415 * EPO) + (0.10811912 * SLC2A1) + (0.08980190 * PSMA7) + (0.01543502 * CREB1) (Figure S3).

Figure 2 Establishment of the prognostic model HARG-score in the TCGA-LIHC cohort. (A) Correlation network of 29 prognostic genes in the TCGA-LIHC cohort. (B) LASSO coefficients for prognostic genes were obtained using minimized lambda. (C) PCA plot exhibited clustering of the high- and low-risk groups. (D) Kaplan–Meier curves of overall survival in the high- and low-risk groups in the TCGA-LIHC cohort. (E) Distribution of HARG-score as well as survival status for the high- and low-risk groups in the TCGA-LIHC cohort. (F) ROC curves for forecasting 1-, 4-, and 5-year survival according to HARG-score in the TCGA-LIHC cohort.

Using the aforementioned formula, we computed the risk score for each HCC patient in the TCGA-LIHC cohort. Subsequently, according to the median HARG-score, patients with HCC were divided into high- and low-risk groups. The results of PCA suggested that the high- and low-risk groups were clustered dividually (Figure 2C). Patients in the low-risk group demonstrated a remarkably greater overall survival compared with the high-risk group (p < 0.0001, Figure 2D). Figure 2E displays the distribution of HARG-score, survival status as well as survival time for the high- and low-risk groups, with the number of patients in survival status decreasing as HARG-score increases. As illustrated in Figure 2F, in the TCGA-LIHC cohort, the area under curve (AUC) values for HARG-score at 1-, 4-, and 5-year were 0.745, 0.733, and 0.705, respectively. We conducted both univariate and multivariate Cox regression analyses to thoroughly examine whether the prognostic model HARG-score was independent of other clinical characteristics (Table S5). These findings suggested that HARG-score was determined as an independent prognostic factor (HR = 2.80, 95% CI: 1.84–4.30, p < 0.001).

In order to validate the strength and reliability of HARG-score, we additionally assessed HARG-score in both LIRI-JP and GSE14520 cohorts, and achieved comparable results. We utilized an identical formula to categorize the patients in both LIRI-JP and GSE14520 cohorts into high- and low-risk groups. In the LIRI-JP cohort, patients in the low-risk group exhibited markedly prolonged overall survival compared to those in the high-risk group (p = 0.0018, Figure 3A). The distribution of HARG-score, survival status as well as survival time was illustrated in Figure 3B, with the number of patients in survival status decreasing as HARG-score increases. The AUC values for HARG-score at 1-, 4-, and 5-year were 0.811, 0.827, and 0.787, respectively (Figure 3C). Similarly, in the GSE14520 cohort, patients in the low-risk group exhibited a remarkably better overall survival compared to the high-risk group (p = 0.0017, Figure 3D). The distribution of HARG-score, survival status as well as survival time was illustrated in Figure 3E, with the number of patients in death status increasing as HARG-score increases. The AUC values for HARG-score at 1-, 4-, and 5-year were 0.663, 0.667, and 0.670, respectively (Figure 3F). Likewise, both univariate and multivariate Cox regression analyses were conducted on HARG-score as well as clinical characteristics in both LIRI-JP and GSE14520 cohorts (Table S5). The results showed that HARG-score was identified as an independent prognostic factor with superior predictive ability for the prognosis of patients with HCC (HR = 2.48, 95% CI: 1.23–5.02, p = 0.011; HR = 1.59, 95% CI: 1.01–2.50, p = 0.046). Furthermore, we compared the predictive ability of HARG-score with other angiogenesis-related signatures. According to our findings, HARG-score consistently demonstrated better performance in comparison to other published angiogenesis-related signatures (Figure 3G).^37–40 These findings indicated that HARG-score exhibited high reliability as a prognostic prediction model and was superior to other angiogenesis-related signatures.

Figure 3 Validation of the prognostic model HARG-score in the testing cohorts. (A) Kaplan–Meier curves of overall survival in the high- and low-risk groups in the LIRI-JP cohort. (B) Distribution of HARG-score as well as survival status for the high- and low-risk groups in the LIRI-JP cohort. (C) ROC curves for forecasting 1-, 4-, and 5-year survival according to HARG-score in the LIRI-JP cohort. (D) Kaplan–Meier curves of overall survival in the high- and low-risk groups in the GSE14520 cohort. (E) Distribution of HARG-score as well as survival status for the high- and low-risk groups in the GSE14520 cohort. (F) ROC curves for forecasting 1-, 4-, and 5-year survival according to HARG-score in the GSE14520 cohort. (G) The AUC values of HARG-score, other published angiogenesis-related signatures and representative immune checkpoint genes.

Validation of Prognostic Gene Expression Levels for HARG-Score

qPCR was used to measure prognostic gene expression levels in both human normal liver cell line LO2 cells as well as liver cancer cell line LM3 cells. As depicted in Figure 4A, compared to LO2 cells, the expression levels of HILPDA, EPO, SLC2A1, PSMA7 and CREB1 were remarkably elevated in LM3 cells, among which HILPDA was the most significantly increased (p < 0.01). Furthermore, we used the HPA database to evaluate the protein levels of these prognostic genes in normal liver as well as HCC tissues. Figure 4B exhibited the IHC staining markers of prognostic genes in both normal liver as well as HCC tissues in the HPA database. These findings suggested that HILPDA, EPO, SLC2A1, PSMA7 and CREB1 had higher protein levels in HCC tissues, which was consistent with the transcriptome expression trends described above.

Figure 4 Validation of prognostic gene expression levels in HARG-score. (A) Expression of prognostic genes (HILPDA, EPO, SLC2A1, PSMA7 and CREB1) in LO2 and LM3 cells. (B) IHC staining plots illustrating representative markers of protein expression of prognostic genes in both HCC and normal liver tissues sourced from the HPA database. ** p < 0.01; *** p < 0.001.

Associations Between HARG-Score and Clinical Characteristics of HCC Patients

In order to assess the correlation between HARG-score and clinical characteristics, we investigated the clinical presentations of the high- and low-risk groups in the TCGA-LIHC, LIRI-JP and GSE14520cohorts, encompassing age, gender, TNM stage, grade, and status. We plotted heatmaps displaying the clinical characteristics of 3 cohorts (Figure 5A–C). In 3 cohorts, HILPDA, EPO, SLC2A1, PSMA7 and CREB1 were expressed at high levels in the high-risk group. We discovered that the survival status of patients within 3 cohorts was in great agreement with the overall survival benefit depicted in the above studies (Figure 5A-C). The low-risk group in the TCGA-LIHC as well as LIRI-JP cohorts was mainly comprised of patients at early stage (I–II), while the high-risk group consisted mostly of patients at advanced stage (III–IV) (p < 0.01, Figure 5A and B). Similarly, in the TCGA-LIHC as well as GSE14520 cohorts, the low-risk group consisted mainly of patients at early TNM (or T) stage, while the high-risk group consisted predominantly of patients at advanced TNM (or T) stage (p < 0.05, Figure 5A and C). This partly elucidated why the patients in high-risk group were related to a worse overall survival benefit. Significantly, in the GSE14520 cohort, we observed that AFP (<= 300 ng/mL) and tumor size (<= 5 cm) were primarily related to the low-risk group, while tumor AFP (>300 ng/mL) and size (>5 cm) were predominantly related to the high-risk group (p < 0.05, Figure 5C). In the TCGA-LIHC cohort, HARG-score was remarkably correlated with grade (p < 0.01, Figure 5D). In addition, in the TCGA-LIHC cohort, we additionally evaluated whether HARG-score retained good predictive ability in different clinical characteristic subgroups, and the findings demonstrated that HARG-score can be utilized to predict overall survival in different subgroups of clinical characteristics (p < 0.05, Figure S4). The above findings suggested that the high- and low-risk groups of HARG-score fully demonstrated the different clinical characteristics of patients with HCC, and HARG-score retained prominent predictive capability in the subgroups of clinical characteristics.

Figure 5 Correlation between HARG-score and clinical characteristics of HCC patients. (A-C) Correlation between clinical characteristics and the high- and low-risk groups in the TCGA-LIHC, LIRI-JP and GSE14520 cohorts. (D) Correlation analysis of clinical characteristics with high as well as low HARG-score groups in the TCGA-LIHC cohort. * p < 0.05; ** p < 0.01; *** p < 0.001.

Importantly, we analyzed the correlation between HARG-score and HCC vascular invasion within the TCGA-LIHC cohort. Additional evaluation of vascular invasion pathology revealed that patients with HCC in the high-risk group exhibited more vascular invasion (both macrovascular invasion and microvascular invasion) (p < 0.001, Figure 6A). Notably, the HARG-score was remarkably greater in patients with vascular invasion than that of patients with nonvascular invasion, and was highest in patients with macrovascular invasion, followed by patients with microvascular invasion, and lowest in patients with nonvascular invasion (p < 0.001, Figure 6B and C). Survival analysis showed that overall survival was remarkably longer in patients with nonvascular invasion than that of patients with macrovascular invasion (p = 0.028), whereas there were no statistically significant differences in overall survival between patients with nonvascular invasion and microvascular invasion, as well as between patients with microvascular invasion and macrovascular invasion (Figure 6D–F). Not surprisingly, the HARG-score also maintained good performance in the vascular invasion subgroup, which can effectively predict overall survival of patients in both macrovascular invasion and microvascular invasion subgroup (p < 0.05, Figure 6G and H).

Figure 6 The correlation between HARG-score and HCC vascular invasion. (A) Percentage bar plots of vascular invasion distribution in the high- and low-risk groups. (B) Comparison of HARG-score between vascular invasion and nonvascular invasion groups. (C) Comparison of HARG-score between groups with macrovascular invasion, microvascular invasion and nonvascular invasion groups. (D-F) Comparison of overall survival in macrovascular invasion, microvascular invasion and nonvascular invasion groups. (G and H) The Kaplan–Meier curve of overall survival in the high- and low-risk groups according to macrovascular invasion and microvascular invasion.

Development and Assessment of the Predictive Nomogram

Given the significance of HARG-score in forecasting the prognosis of patients with HCC, we subsequently intended to explore its value in clinical practice. In the TCGA-LIHC cohort, we developed a nomogram incorporating HARG-score and clinical characteristics that were easily accessible and known to have an impact on HCC prognosis, to forecast 1-, 4-, and 5-year survival (Figure 7A). The nomogram exhibited an outstanding predictive ability, achieving a C-index of 0.711. Calibration plot for the nomogram showed a robust consistency between the observed and predicted 1-, 4-, and 5-year survival (Figure 7B). The AUC value of the nomogram model at 5-year exhibited better accuracy in forecasting survival (Figure 7C). Herein, the AUC values of the nomogram at 1-, 4-, and 5-year were 0.726, 0.727, and 0.729, respectively. As depicted in Figure 7D, we utilized an alluvial plot to display alterations in the characteristics of patients with HCC.

Figure 7 Construction and evaluation of the nomogram to forecast the overall survival in HCC patients. (A) Nomogram incorporating HARG-score to forecast the overall survival in the TCGA-LIHC cohort. (B) Calibration plot of the nomogram for consistency between the observed and predicted 1-, 4-, and 5-year survival in the TCGA-LIHC cohort. (C) Time-dependent ROC curve for forecasting 1-, 4-, and 5-year survival utilizing the nomogram. (D) Alluvial plot displaying the alterations in HARG-score, age, gender, grade, T stage, M stage as well as N stage in the TCGA-LIHC cohort.

Moreover, in the LIRI-JP as well as GSE14520 cohorts, we additionally developed nomograms to predict 1-, 4-, and 5-year survival (Figure S5A and B; A: LIRI-JP, B: GSE14520). The C-index of the nomograms was 0.768 in the LIRI-JP cohort as well as 0.714 in the GSE14520 cohort, suggesting that both nomograms have excellent predictive ability. Calibration plots indicated that the nomograms were capable of accurately predicting 1-, 4-, and 5-year survival (Figure S5C and D; C: LIRI-JP, D: GSE14520). The nomograms for both LIRI-JP and GSE14520 cohorts showed good accuracy with AUC values of 0.812, 0.738, and 0.739, as well as 0.744, 0.709, and 0.683 at 1-, 4-, and 5-year, respectively (Figure S5E-F; E: LIRI-JP, F: GSE14520). According to these results, we found that HARG-score could exert a crucial role in forecasting the prognosis of patients with HCC.

Evaluation of Immune Infiltration Between Low and High HARG-Score Groups

According to the CIBERSORT algorithm, we analyzed the abundance of 22 infiltrating immune cells in every patient with HCC in the high- and low-risk groups. In the TCGA-LIHC cohort, the infiltration of B cells naive, macrophages M1, natural killer (NK) cells activated as well as mast cells resting was significantly higher in the low-risk group than that in the high-risk group (p < 0.01, Figure 8A), while the infiltration of B cells memory, macrophages M0 as well as T cells regulatory (Tregs) was markedly lower in the low-risk group than that in the high-risk group (p < 0.05, Figure 8A). In the LIRI-JP cohort, the infiltration of macrophages M1 and NK cells resting was markedly elevated in the low-risk group than that in the high-risk group, whereas the infiltration of Tregs and macrophages M0 was significantly decreased (p < 0.01, Figure 8B). In the GSE14520 cohort, the infiltration of macrophages M1 and mast cells resting was remarkably elevated in the low-risk group than that in the high-risk group (p < 0.01, Figure 8C); on the contrary, the infiltration of macrophages M0 and mast cells activated was remarkably elevated in the high-risk group than that in the low-risk group (p < 0.01, Figure 8C). Interestingly, we found that the infiltration of macrophage M1 was significantly elevated in the low-risk group than that in the high-risk group in 3 cohorts (p < 0.01).

Figure 8 Analysis of immune infiltration in the high- and low-risk groups. (A-C) The abundance of infiltration immune cells in the high- and low-risk groups in the TCGA-LIHC, LIRI-JP and GSE14520 cohorts, respectively. (D) Associations between overall survival and T cells CD8, macrophages M0 as well as macrophages M2 in the TCGA-LIHC cohort. (E) Associations between overall survival and macrophages M1, B cells naive as well as macrophages M0 in the GSE14520 cohort. *p < 0.05; **p < 0.01; ***p < 0.001. Abbreviation ns, not significant.

Moreover, we investigated the prognostic significance of 22 infiltrating immune cells. In the TCGA-LIHC cohort, the infiltration of macrophages M0, T cells CD8 as well as macrophages M2 was remarkably associated with overall survival (p < 0.01, Figure 8D). Lower infiltration of macrophages M0 and macrophages M2 was related to better overall survival, while lower infiltration of T cells CD8 was related to poorer overall survival. In the GSE14520 cohort, the infiltration of macrophages M1, B cells naive and macrophages M0 was remarkably correlated with the overall survival of HCC patients (p < 0.05, Figure 8E). The lower infiltration of macrophage M0 exhibited better overall survival, while the lower infiltration of macrophages M1 as well as B cells naive exhibited worse overall survival.

Subsequently, we estimated enrichment scores for 16 immune cells and 13 immune functions in 3 cohorts utilizing ssGSEA to assess the correlation between HARG-score and the scores of immune cells as well as immune functions (Figure S6A and B). Immune cells comprising macrophages, Tregs, as well as dendritic cells (DCs) showed significant positive correlations with HARG-score (p < 0.05). The lower HARG-score was strongly associated with NK cells as well as B cells. In immune functions, major histocompatibility complex (MHC) class I was significantly positively correlated with HARG-score, while interferon (IFN) responses were remarkably negatively related to HARG-score (p < 0.05). These findings indicated that HARG-score was correlated with numerous immune cells and immune functions, suggesting that HARG-score reflects immune infiltration to some extent.

The Role of HARG-Score in Immunotherapy

Considering the association between HARG-score and immune infiltration, we further evaluated the prediction of HARG-score on the efficacy of immunotherapy in the immunotherapy cohort (IMvigor210). Notably, HARG-score was remarkably related to the immune phenotype (desert, excluded and inflamed) of the patients, and HARG-score in patients with immune desert was remarkably higher compared to patients with immune excluded or inflamed (p = 0.0026, Figure 9A). In order to make the IMvigor210 cohort more able to reflect immunotherapy efficacy in HCC patients, we further screened 98 patients with HCC characteristics, that is, these patients presented with liver metastases. Likewise, according to the median HARG-score, 98 patients were grouped into the high- and low-risk groups, and we observed that patients in the low-risk group exhibited superior overall survival in comparison to the high-risk group (p = 0.0044, Figure 9B). The AUC values of HARG-score at 12- and 24-month were 0.675 and 0.963, respectively (Figure 9C). Furthermore, we explored the association between HARG-score and immunotherapy response, tumor cell (IC0, IC1 and IC2+) subgroups, as well as immune cell (TC0, TC1 and TC2+) subgroups. The findings showed that HARG-score was related to immunotherapy response and immune cell (IC0, IC1 and IC2+) subgroups (p < 0.05), but not with tumor cell (TC0, TC1 and TC2+) subgroups (Figures 9D–F). We also found that the lower HARG-score was correlated with immunotherapy efficacy (p < 0.001, Figure 9G). These findings indicated that HARG-score can enhance the accuracy of evaluating immunotherapy efficacy and identify patients who are expected to respond positively to immunotherapy.

Figure 9 The role of HARG-score in immunotherapy. (A) HARG-score in groups with different immune phenotypes. (B) Kaplan–Meier curves of overall survival in the high- and low-risk groups in 98 patients with HCC characteristics. (C) ROC curves for forecasting 12- and 24-month survival according to HARG-score. (D-F) Associations between HARG-score and immunotherapy response, immune cell (IC0, IC1 and IC2+) subgroups as well as tumor cell (TC0, TC1 and TC2+) subgroups. (G) Clinical efficacy of immunotherapy (complete response/partial response [CR/PR] as well as stable disease/progressive disease [SD/PD]) in high and low HARG-score groups.

Inhibition of HILPDA Suppresses the Proliferation and Migration of HCC Cells

The above findings revealed that the prognostic model HAGR-score is of great significance for forecasting the prognosis and immunotherapy of patients with HCC, among which the key gene HILPDA is most significantly upregulated in HCC cells. Therefore, we additionally explored the role of HILPDA in regulating the function of HCC cells. Herein, we utilized siRNA to suppress HILPDA expression in LM3 cells (p < 0.001, Figure 10A). The CCK-8 assay and EdU assay demonstrated that the inhibition of HILPDA remarkably suppressed cell viability as well as cell proliferation (p < 0.01, Figures 10B-C and F). Moreover, the cell scratch and transwell assay suggested that the inhibition of HILPDA remarkably weakened the migration ability of LM3 cells (p < 0.001, Figures 10D-E and G-H). These results indicated that the suppression of HILPDA remarkably inhibited the proliferation as well as migration ability of HCC cells, which is consistent with its HCC-promoting effect in HARG-score.

Figure 10 Inhibition of HILPDA suppresses HCC development. (A) qPCR analysis exhibited the inhibitory efficiency of HILPDA in LM3 cells. (B) The cell viability of LM3 was assessed by CCK8 assay after HILPDA inhibition. (C, F) The cell proliferation capacity of LM3 was detected by EdU assay after HILPDA inhibition. (D-E, G-H) Wounding healing assay and transwell assay were performed to assess the migration ability of LM3 cells after HILPDA inhibition. ** p < 0.01; *** p < 0.001.

HILPDA Inhibition Reduces HCC Cell Hypoxia-Angiogenesis Markers and Affects HUVECs Angiogenesis

To further investigate the effect of HILPDA on hypoxia-angiogenesis markers in HCC cells and angiogenesis in HUVECs. We first used ELISA kits to measure the levels of HIF-1α as well as VEGF in the supernatant of LM3 cells 24 hours after transfection with HILPDA siRNA, and found that suppression of HILPDA reduced the levels of hypoxia-angiogenesis markers HIF-1α as well as VEGF (p < 0.05, Figures 11A–B). Western blot experiment indicated that inhibition of HILPDA alleviated the expression of HIF-1α and key downstream P-ERK of VEGF in LM3 cells (p < 0.001, Figure 11C). In addition, we further validated the effect of HILPDA silencing in LM3 cells on HUVECs angiogenesis. The supernatant of LM3 cells transfected with HILPDA siRNA for 24 hours was co-cultured with HUVECs for 24 hours, and the effect on HUVECs angiogenesis was observed by tube formation assay, and the results suggested that the inhibition of HILPDA could attenuate the angiogenesis of HUVECs (p < 0.001, Figure 11D). These findings indicated that HILPDA inhibition reduces hypoxia-angiogenesis markers in HCC cells and affects HUVECs angiogenesis, thereby participating in the development of HCC.

Figure 11 HILPDA inhibition affects HCC cell hypoxia-angiogenesis markers and HUVECs angiogenesis. (A-B) HIF-1α and VEGF levels in the supernatant of LM3 cells after HILPDA inhibition. (C) HIF-1α and P-ERK protein expression levels in LM3 cells after inhibition of HILPDA. (D) Effect of HILPDA inhibition in LM3 cells on angiogenesis in HUVECs. *p < 0.05; **p < 0.01; ***p < 0.001.

HILPDA Promotes HCC Development by Affecting CTLA4 Expression

Immune checkpoints play an important role in HCC development, and HARG-score has a certain correlation with immune checkpoints, but the relationship between HILPDA and immune checkpoints is still unclear, and whether it will affect the expression of immune checkpoints needs further research to verify. Based on this, we first analyzed the correlation between HILPDA and immune checkpoints, and the findings indicated that HILPDA had the most significant correlation with CTLA4 (p < 0.01, Figure 12A). This was in agreement with the most significant correlation between HARG-score and CTLA4. We then explored whether the inhibition of HILPDA affects the expression of CTLA4. Western blot results showed that the inhibition of HILPDA significantly reduced CTLA4 expression in LM3 cells (p < 0.05, Figure 12B). In conclusion, the findings presented here indicated that HILPDA promotes HCC development by influencing the expression of CTLA4.

Figure 12 Inhibition of HILPDA affects the expression of CTLA4. (A) Correlations between HILPDA and representative immune checkpoint genes. (B) Western blot was conducted to investigate the effect of inhibition or non-inhibition of HILPDA expression on CTLA4 protein expression level. *p < 0.05; ***p < 0.001.

Discussion

Early diagnosis and prognostic prediction of HCC are crucial to patient survival. Currently, the diagnosis and prognosis prediction of HCC primarily rely on pathological evaluation as well as tumor staging, but their sensitivity and accuracy are inadequate.^41,42 Hence, it is imperative to find efficient biomarkers to establish reliable and accurate prognostic model to promote the survival of HCC patients. Multiple studies have shown that hypoxia-mediated angiogenesis not only affects HCC growth, invasion and metastasis but also triggers the infiltration of several distinct immune cells in the HCC microenvironment.^22,43–45 However, most studies have focused only on single hypoxia-angiogenesis-related genes or immune cell infiltration, neglecting to investigate the entire role of multiple HARGs in HCC as well as their impact on immune infiltration.^46–48 Our study provided a comprehensive comparison of genetic as well as expression alterations of HARGs in HCC was performed. We identified differentially expressed HARGs in HCC samples as well as non-tumor samples and developed a prognostic model HARG-score that demonstrated prominent predictive ability. Our previous study found that pyroptosis-related gene prognostic model PRGS could well predict the survival of HCC patients.⁴⁹ Here, compared to PRGS, the HARG-score exhibited better predictive ability. Specifically, in the LIRI-JP cohort, the AUC values of 0.811 and 0.827 for HARG-score at 1- and 4 years were superior to the AUC values of 0.751 and 0.773 for PRGS at 1 and 4 years, respectively. In addition, we investigated the function and mechanism of the key genes that comprise HARG-score, in an attempt to explore its potential as a prognostic biomarker and therapeutic target for HCC. Patients in the high- and low-risk groups determined according to HARG-score demonstrated remarkably distinct prognosis, clinical characteristics, immune infiltration and immunotherapy efficacy. In addition, HARG-score had remarkable survival differences among patients with different subgroups of clinical characteristics. Subsequently, through the integration of HARG-score as well as clinical characteristics, we developed nomograms to additionally augment the predictive ability of HARG-score and promote its clinical practice. Herein, the prognostic ability of the HARG-score-related nomogram was compared with our previous PRGS-related nomogram. The findings showed that in the GSE14520 cohort, the AUC values of 0.744 and 0.709 for HARG-score-related nomogram at 1- and 4 years were better than the AUC values of 0.728 and 0.599 for PRGS-related nomogram at 1- and 4 years, respectively.⁴⁹ Noteworthy, CIBERSORT and ssGSEA analysis revealed that lower and higher HARG-score indicated immune activation and suppression, respectively. In addition, applying HARG-score to the immunotherapy cohort also enabled the anticipation of immunotherapy efficacy. Finally, experimental data showed that the key gene HILPDA contributes to HCC progression by facilitating angiogenesis and affecting the expression of CTLA4.

Currently, hypoxia and angiogenesis have received increasing attention in the study of tumors.^50,51 Notably, hypoxia-induced angiogenesis is a vital characteristic of the tumor microenvironment, which participates in tumor growth and immune evasion.^52,53 Thus, hypoxia-angiogenesis is thought to exert a crucial role in tumor treatment as well as prognostic prediction prognosis prediction. Nevertheless, the impact of hypoxia-angiogenesis on the prognosis of patients with HCC remains unclear. Here, we constructed an independent HCC prognostic model utilizing prognostic HARGs, containing HILPDA, EPO, SLC2A1, PSMA7 and CREB1. EPO was elevated in HCC and negatively associated with the overall prognosis of HCC, which might be attributed to the promotion of HCC-related angiogenesis by EPO.^14,54 SLC2A1 and CREB1 were overexpressed in HCC, regulate HCC cell proliferation, invasion, as well as migration, and can be utilized as prognostic biomarkers for HCC; furthermore, SLC2A1 was also related to immune infiltration of tumors.^55–58 PSMA7 expression was elevated in HCC, and its higher level was related to poor prognosis of HCC.⁵⁹ The initial function of HILPDA was thought to regulate lipid metabolism and lipid droplet abundance, and as research progressed, it was found that it also exerts a critical role in solid tumors.⁶⁰ Studies have demonstrated that HILPDA is upregulated in HCC, and its high expression is related to inferior prognosis in HCC.^26,61 A study showed that HILPDA promotes HCC metastasis by facilitating the secretion of inflammatory factors from HCC cells and suppressing NK cell cytotoxicity.²⁶ However, the role of HILPDA in the progression of HCC as well as its latent mechanism remains incompletely understood. Here, we further conducted experimental studies to verify the role of HILPDA in regulating the function of HCC cells and demonstrate the reliability of our analysis.

In the past few decades, tremendous strides have been made in the immunotherapy of HCC, but there are still heterogeneous in the prognosis of patients, emphasizing the critical role of immune cell infiltration in HCC development as well as prognosis.²⁹ Currently, anti-PD1, anti-programmed death ligand 1 (PD-L1), as well as anti-CTLA4 therapies have become key technologies for current tumor treatment. Nevertheless, there is a subset of patients who do not respond to the above treatments, and numerous studies indicate that PD1 as well as PD-L1 are not effective markers for forecasting treatment with ICIs.^62,63 Hence, it is imperative to develop an effective tool for forecasting the efficacy of immunotherapy in patients. Based on the results of this research, we determined that HARG-score is an effective immune classifier that can be used to predict immunotherapy prognosis. Notably, the key gene HILPDA in HARG-score could facilitate HCC by influencing immune checkpoint genes.

Our study correlated hypoxia-angiogenesis with the prognosis of HCC patients to develop a prognostic model HARG-score, and was superior to other angiogenesis-related signatures. Noteworthy, the key gene HILPDA promoted HCC development by promoting angiogenesis and affecting CTLA4 expression. Nevertheless, the present study has several limitations. Firstly, the data in our study were primarily derived from public databases, and a large amount of clinical data is needed to further verify the predictive ability of the prognostic model, and we have started collecting clinical verification samples. Secondly, the prognostic model was associated with immune infiltration, but its mechanism has not been fully elucidated, and more in vivo as well as in vitro experimental studies are needed to further explore, which is the direction of our future research.

Conclusion

Our comprehensive analysis of HARG-score indicated its association with prognosis, clinical characteristics, as well as immune infiltration in HCC. We also demonstrated the predictive ability of HARG-score in immunotherapy efficacy. These findings emphasized the crucial clinical significance of HARG-score in predicting survival as well as guiding personalized treatment for HCC, in which HILPDA may be a promising prognostic biomarker and therapeutic target, offering new insights to improve the prognosis and treatment of HCC.