The Expression of Ferroptosis-Related Genes in Hepatocellular Carcinoma and Their Relationships With Prognosis

Introduction

Global Cancer Statistics 2024 shows that primary liver cancer is the 12^th most common cancer in the world in 2024, with the sixth highest mortality rate, and it is an extremely malignant tumor.¹ Among them, hepatocellular carcinoma (HCC) accounts for 75%–85% of primary liver cancers. The screening and diagnosis of HCC are mainly carried out through several methods nowadays, including screening of high-risk groups, imaging tests, hematological molecular markers, puncture biopsy, and pathological diagnosis, depending on the different conditions of patients. The field of HCC treatment is characterized by multidisciplinary participation and coexistence of multiple treatment methods, and general treatment methods include hepatectomy, liver transplantation, ablation therapy, precision radiotherapy, and systemic antitumor therapy. Choosing reasonable treatments for HCC patients with different stages can maximize the therapeutic efficacy.² Therefore, multidisciplinary and multicenter joint efforts are needed to achieve breakthrough progress in liver cancer treatment, and more high-quality studies should be carried out to improve the level of liver cancer diagnosis and treatment.³ Although HCC treatment has advanced rapidly, further high-level research evidence is needed to explore and supplement current knowledge.

Cell death is the irreversible cessation of life. To date, multiple programmed cell death patterns have been identified, such as apoptosis, necrosis, autophagy, and pyroptosis. Cell death is essential for normal development, maintaining homeostasis in the body, and preventing cancers and other hyperproliferative diseases.⁴ Ferroptosis, a type of cell death discovered in recent years, is a regulated cell death caused by the accumulation of lipid peroxidation products and reactive oxygen species, and is distinct from other programmed cell deaths such as apoptosis. Ferroptosis is not characterized by typical apoptosis and necrosis, but mainly by cell membrane vesiculation and rupture, mitochondrial atrophy, increased membrane density, reduction or even disappearance of mitochondrial cristae, and chromatin condensation. Dixon et al named this unique iron-dependent nonapoptotic cell death ferroptosis in 2012.⁵

Research has been conducted on the regulation of ferroptosis, and some regulatory signaling pathways related to ferroptosis have been explored. Among them, Glutathione Peroxidase 4 (GPX4) is an antioxidant enzyme, that is a key regulator of ferroptosis, and it uses glutathione (GSH) as a cofactor to catalyze the reduction of lipid peroxides. Meanwhile, GSH is also affected by Solute Carrier Family 7 Member 11 (SLC7A11), which leads to the reduction of GSH synthesis when its activity is inhibited and triggers oxidative damage that ultimately leads to the occurrence of ferroptosis. The SLC7A11-GSH-GPX4 pathway is the classical pathway of ferroptosis.^6,7 Ferroptosis suppressor protein 1 (FSP1), previously known as apoptosis-inducing factor mitochondrial 2, is a key protein that resists ferroptosis.⁸ Beyond the GPX4-centered ferroptosis pathway, the NAD(P)H-FSP1-CoQ10 axis represents a parallel antioxidant mechanism.

FerrDb,⁹ a database of ferroptosis-related genes, has been created to help researchers acquire insights into ferroptosis. FerrDb is the first manually organized ferroptosis database to manage and characterize ferroptosis-related markers and regulators, as well as ferroptosis-associated diseases. FerrDb has downloaded 784 articles on ferroptosis from the PubMed database, and extracted and organized 259 regulatory genes, 111 markers, and 95 ferroptosis-related diseases.

Correlation between ferroptosis and HCC has already been demonstrated in several studies,¹⁰ and inducing or promoting cellular ferroptosis may become a promising tumor therapy, but ferroptosis-related genes with prognostic and clinical diagnostic significance in HCC have not yet been more comprehensively researched. Therefore, in this study, based on the sample data in The Cancer Genome Atlas (TCGA) database, a bioinformatics approach was utilized to screen for ferroptosis-related genes in the FerrDb database, to explore the differential expression of ferroptosis-related genes in HCC, whether the genes were at risk factors for their survival, and the correlation between the gene expression and the clinical traits, and validated with the Gene Expression Profiling Interactive Analysis 2 (GEPIA2) database. This study utilized Gene Set Enrichment Analysis (GSEA) to explore functionally related pathways, aiming to identify potential biomarkers and therapeutic targets for early diagnosis and prognosis of HCC, thereby contributing to improved clinical strategies.

Materials and Methods

Data Extraction and Arrangement

The liver cancer-related data that needed to be added to the cart for downloading was selected through the TCGA Genomic Data Commons (GDC) website (https://portal.gdc.cancer.gov/repository/). The data, including transcriptome analysis and gene expression quantitative data of HCC and normal samples, have been downloaded. A total of 465 files were collected that contained 407 tumor samples and 58 normal samples. Clinically relevant information on 418 HCC patients was obtained through download. The clinical data is comprised of 377 HCC samples and 41 metastatic HCC samples (Cholangiocellular carcinoma, CHOL). Perl (version 5.32.1) has been used to organize the downloaded data for further analysis. The Ferroptosis-related gene clusters “Driver”, “Suppressor”, and “Marker” were downloaded from the FerrDb database (http://www.zhounan.org/ferrdb/legacy/). After removing the intersections for further selection, 259 genes were obtained after sorting out 108 Driver genes, 69 Suppressor genes, and 111 Marker genes. The study utilized the anonymized secondary data that contained no personally identifiable information. Therefore, ethical approval is exempted. The requirement for consent was waived because this study was based on publicly available genomic data from established databases and did not involve direct interaction with patients or any invasive procedures. According to relevant ethical regulations and guidelines for using publicly available, this study meets the conditions for exemption from ethical approval and informed consent.

Preliminary Identification of Candidate Ferroptosis-Related Genes by Significance of Scatter Differential Expression Analysis and Kaplan–Meier Survival Analysis

Candidate ferroptosis-related genes were preliminarily selected by the significance of the scatter difference analysis and Kaplan–Meier survival analysis, and genes with significant differences with P-value < 0.05 were screened, and then the genes derived from the screening were subjected to the following analyses.

Differential expression analysis and survival analyses were conducted with R software (version 3.6.3). Scatter plot differential expression analysis was performed using Wilcoxon signed-rank test, and the gene expression differences between normal and HCC tissues were plotted and analyzed by using the R packages “limma” and “beeswarm”. Ferroptosis-related genes with a P-value < 0.05 were selected. Paired difference analysis was performed by pairing normal samples with HCC samples, comparing the gene expression of the two samples, and analyzing the expression differences between HCC tissues and normal tissues of the same samples by the Wilcoxon signed-rank test for paired difference analysis. The gene expression information related to the paired normal and HCC samples was organized into the format required for plotting heatmaps, as well as the gene-related expression information of all samples in the same way, and the organized gene expression data was uploaded to the website for plotting heatmap. The website (https://www.bioinformatics.com.cn) is used to perform bioinformatics analysis for free.

To carry out the survival analysis, the organized clinical information files were first used to delete the samples that lacked information on survival status and survival time. Survival curves were plotted using the R package “survival”, and the Kaplan–Meier survival curves were evaluated using the Log rank test. The ferroptosis-related genes with P-value < 0.05 were screened based on the output of the survival curves.

Clinical Correlation Analysis With Univariate and Multivariate Cox Regression Analysis

Clinical correlation analysis was performed using R software. The arranged clinical information was deleted from the samples with unknown information and other data, and only the data on the clinicopathological factors to be analyzed were kept. The five clinicopathological factors of the clinical stage (stage), histological grading (G), tumor stage (T), distant metastasis (M), and lymph node metastasis (N) have been analyzed. The analysis was conducted with the Wilcoxon signed rank test between two groups and the Kruskal–Wallis test between multiple groups, and the images were plotted.

Univariate Cox regression analysis and multivariate Cox regression analysis were performed using R software. Unknown information was removed from the organized clinical information data to ensure data integrity. Univariate Cox regression analysis was conducted using the R package “survival”, and logistic regression was used to analyze the relationship between different levels of clinicopathological factors and gene expression. The objective of the univariate Cox regression analysis was to confirm whether a single clinical factor or gene expression was linked to survival. If the analysis produced a P-value of 0.05, the factor was deemed significant, and further multivariate Cox regression analysis was conducted.

The multivariate Cox regression analysis was performed using the R packages “survival” and “survminer”. The clinical data were organized, and a multivariate Cox regression analysis was conducted to analyze simultaneously whether several factors, including clinical traits and gene expression, were correlated with survival. P-value < 0.05 was considered significant, and independent prognostic markers were confirmed. A forest plot was also plotted to show the hazard ratio (HR) for each clinical trait and gene expression. For a more aesthetically pleasing graph, gene expression was logarithmically formulated as log₂(TPM+1)(TPM, Transcripts Per Million). In general, HR > 1 indicates that the gene is a hazard factor and HR < 1 indicates that the gene is a protective factor.

Validation of the Screened Genes in GEPIA2 and HPA Databases

The GEPIA2 database¹¹ and the Human Protein Atlas (HPA) database¹² were used to validate the results. The results of the screened ferroptosis-related genes in the previous steps were verified by plotting gene expression difference box plots, clinical staging plots, and survival analysis plots in the GEPIA2 database (http://gepia2.cancer-pku.cn). To use the GEPIA2 database for plotting, enter the gene name or ID in “Gene” and select HCC (Liver cancer, LIHC). To create a box plot, set the P-value cutoff to 0.05, select ‘Multiple Datasets’ and “Match TCGA normal and Genotype-Tissue Expression Project (GTEx) data”, and use log₂(TPM+1). More data from normal liver samples were added by selecting data from the GTEx database. In the process of plotting clinical stages, log₂(TPM+1) was employed and major stages were chosen for plotting. The plotting of survival analysis involved selecting “Overall Survival (OS)” and “Median”, and using default values for other options. The GEPIA2 database was used to validate the results of bioinformatics analyses.¹³

The HPA database portal (https://www.proteinatlas.org/) verified the protein expression of the screened ferroptosis-related genes in HCC tissues and normal liver tissues.

Selection of Gene Function-Related Pathways by Gene Set Enrichment Analysis (GSEA)

Gene expression datasets were obtained from The Cancer Genome Atlas (TCGA) database, processed using Perl scripts, and analyzed for functional enrichment via Gene Set Enrichment Analysis (GSEA) software (version 4.3.2) based on pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG). The genomic parameters and parameters for running the enrichment test were set. “c2.cp.kegg.v2023.1.Hs.symbols.gmt” was selected as the gene set database. When filtering pathways and plotting multi-GSEA enrichment maps, signaling pathways with false discovery rate (FDR) q-values < 0.05 were usually selected to have significant enrichment. When all the q-values > 0.05, then the P-value < 0.05 was used at this point to select enriched signaling pathways.

Cell Lines and Culture

HepG2 cells and Huh-7 cells were purchased from the Cell Culture Center of the Institute of Basic Medical Sciences, Chinese Academy of Sciences (Beijing, China). HCC-LM3 cell was obtained from the China Center for Type Culture Collection of Wuhan University (Wuhan, China). MHCC97-L cells were purchased from Meisen Chinese Tissue Culture Collections (Zhejiang, China). The normal hepatocyte cell line L02 was presented by the Wuhan Churuike Pharmaceutical Technology Co., Ltd (CRK Pharma, Wuhan, China) and derived from the National Collection of Authenticated Cell Cultures. RPMI-1640 (Gibco, USA) medium having 10% fetal bovine serum (FBS, Hyclone) was used to culture the L02 cells. MHCC97-L, HepG2, HCC-LM3, and Huh-7 cells were cultured in Dulbecco’s modified Eagle’s medium (DMEM, Hyclone, USA) supplemented with 10% FBS, penicillin, and streptomycin. All the cells were cultured at 37 °C with a moistened environment of 5% CO₂.

Isolation of RNA and Quantitative Reverse Transcription PCR (RT-qPCR)

Total RNA was extracted from L02, MHCC97-L, HepG2, HCC-LM3, and Huh-7 cell lines by the RNAiso Plus (Takara, Japan) and the RNA was reverse transcribed with a PrimeScript RT reagent Kit (Takara). The paired primers (Tsingke, Beijing, China) used for amplification were illustrated (Table S1). The quantitative reverse transcription PCR (RT-qPCR) of various genes was performed using the TransStart Tip Green qPCR SuperMix (TransGen Biotech Co., Beijing, China). All the procedures were conducted following the manufacturer’s instructions. Eventually, the relative expression level mRNA was calculated and quantified with the 2^−ΔΔCt method after normalization regarding the expression of Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH).

The Establishment of a Protein–Protein Interaction (PPI) Network

Identify significantly differentially expressed genes and obtain the protein interaction information through the STRING database (Search Tool for the Retrieval of Interacting Genes/Proteins), which can be accessed at https://string-db.org/. Next, import the protein interaction data into the Cytoscape software (version 3.10.3) and enhance the visualization to generate a protein interaction network diagram.

Results

Screening of Ferroptosis-Related Genes

After differential expression analysis and survival analysis of 259 ferroptosis-related genes inside the FerrDb database (Table S2), a total of 82 genes with P-values < 0.05 were obtained (Table 1). By searching for publications, genes that had already been published in the field of liver cancer using the same bioinformatics methods as this study were eliminated. The remaining genes were filtered by the standard of having a significant univariate Cox analysis result, retaining genes with significant results. Finally, eight ferroptosis-related genes associated with the prognosis of HCC were obtained. The eight screened genes include Acyl-CoA synthetase long-chain family member 3 (ACSL3), Asparagine synthetase (ASNS), Charged multivesicular body protein 5 (CHMP5), MYB proto-oncogene, transcription factor (MYB), Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene (PCK2), Phosphogluconate Dehydrogenase (PGD), Solute Carrier Family 38 Member 1 (SLC38A1), YY1 Associated Protein 1 (YY1AP1).

Table 1 The Preliminary Screening Results of Ferroptosis-Related Genes

Differential Expression of Genes in Hepatocellular Carcinoma (HCC) Tumor Tissue and Normal Tissue

The TCGA HCC dataset has a total of 465 patients with HCC, including 407 tumor tissue samples and 58 normal tissue samples with mRNA sequencing data. The differential expression plot demonstrated that the mRNA expression of the PCK2 gene in HCC tumor tissues was significantly lower than in normal tissues, with P-value < 0.001 showing a very significant difference. Except for PCK2, the mRNA expression of all other genes in HCC tumor tissues was significantly higher than that in the normal tissues, and the P-values obtained were < 0.001. Matching normal tissues with HCC tumor tissues, 58 pairs of normal tissues and cancer tissues were obtained from the same patients for analyzing the differences in the expression of ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 genes. As shown by the paired analysis graph, the mRNA expression levels of ACSL3, ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1 showed an overall up-regulation trend in HCC tissues, with P-values of < 0.001. Meanwhile, mRNA expression of PCK2 showed a down-regulation trend in most of the tissues, with significant P-values of < 0.001 (Figure 1). The expression of the screened ferroptosis-related genes in 58 paired samples (Figure 2A) and the expression of the screened ferroptosis-related genes in all 465 samples were organized and plotted as heatmaps (Figure 2B). The heatmap used colors to represent the data, which can reflect a large amount of data more intuitively.

Figure 1 The mRNA Expression of screened ferroptosis-related gene in the adjacent normal and hepatocellular carcinoma(HCC) tissues, including ACSL3 (A), ASNS (B), CHMP5 (C), MYB (D), PCK2 (E), PGD (F), SLC38A1 (G), YY1AP1 (H). Abbreviations: ACSL3, Acyl-CoA synthetase long-chain family member 3; ASNS, Asparagine synthetase; CHMP5, Charged multivesicular body protein 5; MYB, MYB proto-oncogene, transcription factor; PCK2, Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene; PGD, Phosphogluconate Dehydrogenase; SLC38A1, Solute Carrier Family 38 Member 1; YY1AP1, YY1 Associated Protein 1.

Figure 2 The Heatmap of screened ferroptosis-related gene expression in the adjacent normal and hepatocellular carcinoma (HCC) tissues. (A) The expression of 58 paired samples; (B) The expression of all 465 samples.

Kaplan–Meier Survival Analysis of Patients With Hepatocellular Carcinoma (HCC) From TCGA Database

Kaplan–Meier survival curves were plotted based on survival status and survival time of clinical data of HCC in the TCGA database. The difference of ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 was statistically significant (P < 0.05, Figure 3). There were 201 clinical samples in the high-expression group and 202 clinical samples in the low-expression group. The five-year survival rate of the PCK2 gene was lower for the low-expression group (P < 0.001, Figure 3E), which showed that the low expression of PCK2 was associated with the low survival rate of HCC. The five-year survival rate of the rest of genes was lower for the high-expression group, and the high expression of ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 was highly correlated with the low survival rate of HCC.

Figure 3 The Kaplan-Meier survival curves of screened ferroptosis-related gene expression in hepatocellular carcinoma (HCC), including ACSL3 (A), ASNS (B), CHMP5 (C), MYB (D), PCK2 (E), PGD (F), SLC38A1 (G), YY1AP1 (H). Abbreviations: ACSL3, Acyl-CoA synthetase long-chain family member 3; ASNS, Asparagine synthetase; CHMP5, Charged multivesicular body protein 5; MYB, MYB proto-oncogene, transcription factor; PCK2, Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene; PGD, Phosphogluconate Dehydrogenase; SLC38A1, Solute Carrier Family 38 Member 1; YY1AP1, YY1 Associated Protein 1.

Clinical Correlation Analysis of Ferroptosis-Related Genes in Hepatocellular Carcinoma (HCC)

ASNS (P = 0.023), MYB (P = 0.045), PCK2 (P < 0.001), PGD (P < 0.001), SLC38A1 (P = 0.022), and YY1AP1 (P = 0.048) expression was significantly correlated with the clinical stage (STAGE), ACSL3 (P = 0.009), ASNS (P = 0.01), MYB (P = 0.01), PCK2 (P < 0.001), PGD (P = 0.017), SLC38A1 (P < 0.001), and YY1AP1 (P < 0.001) were significantly correlated with histologic grading (G), and ACSL3 (P = 0.024), MYB (P = 0.011), PCK2 (P < 0.001), PGD (P = 0.019), SLC38A1 (P = 0.008), YY1AP1 (P = 0.047) expression was significantly correlated with tumor staging (T), ASNS (P = 0.007), CHMP5 (P < 0.001), PCK2 (P = 0.001), and SLC38A1 (P = 0.007) were significantly associated with lymph node metastasis (N). However, no genes significantly correlated with distant metastasis (M) were present in this analysis (Figure 4).

Figure 4 The result of the correlation between screened ferroptosis-related gene expression levels and various clinicopathological features in hepatocellular carcinoma (HCC), including ACSL3 (A), ASNS (B), CHMP5 (C), MYB (D), PCK2 (E), PGD (F), SLC38A1 (G), YY1AP1 (H). Abbreviations: ACSL3, Acyl-CoA synthetase long-chain family member 3; ASNS, Asparagine synthetase; CHMP5, Charged multivesicular body protein 5; MYB, MYB proto-oncogene, transcription factor; PCK2, Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene; PGD, Phosphogluconate Dehydrogenase; SLC38A1, Solute Carrier Family 38 Member 1; YY1AP1, YY1 Associated Protein 1 stage, clinical stage; grade, histologic grade; T, tumor stage; M, distant metastasis; N, lymph node metastasis.

Univariate and Multivariate Cox Regression Analysis

Univariate Cox regression analysis showed that the increased clinical stage, increased tumor stage (T), and high expression of ACSL3, MYB, PGD, and SLC38A1 were significantly associated with an increased risk of survival in HCC patients (P < 0.001). The development of distant metastases (M), low expression of PCK2, and high expression of ASNS, CHMP5, and YY1AP1 were significantly associated with poorer survival in HCC patients (P < 0.05). However, no correlation was found between age (Age), gender (Gender), histologic grading (G), and lymph node metastasis (N), and survival of HCC patients (P > 0.05) (Table 2).

Table 2 The Univariate Cox Regression Analysis for Prognosis Risk Factors of Hepatocellular Carcinoma Patients

Multivariate Cox regression analysis indicated that the expression of ACSL3, ASNS, MYB, PGD, SLC38A1, and YY1AP1 could serve as an independent prognostic factor for OS in patients with HCC (P < 0.05) (Table 3 and Figure 5).

Table 3 The Multivariate Cox Regression Analysis for Prognosis Risk Factors of Hepatocellular Carcinoma Patients

Figure 5 The forest plots of multivariate Cox analysis about screened ferroptosis-related gene expression levels among hepatocellular carcinoma(HCC), including ACSL3 (A), ASNS (B), CHMP5 (C), MYB (D), PCK2 (E), PGD (F), SLC38A1 (G), YY1AP1 (H). Abbreviations: ACSL3, Acyl-CoA synthetase long-chain family member 3; ASNS, Asparagine synthetase; CHMP5, Charged multivesicular body protein 5; MYB, MYB proto-oncogene, transcription factor; PCK2, Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene; PGD, Phosphogluconate Dehydrogenase; SLC38A1, Solute Carrier Family 38 Member 1; YY1AP1, YY1 Associated Protein 1. Notes: “*”, P < 0.05; “**”, P < 0.01; “***”, P < 0.001. “# Events”, the number of clinical sample. The italicized text in the lower left corner provides statistical information about the model.

Validation of GEPIA2 and HPA Databases

GEPIA2 analysis showsed the expression trends of ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 in liver tumor tissues were consistent with the analysis in this research. In the survival analysis (Figure S1), the OS time of ACSL3, ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1 high-expression group was significantly higher than that of the low-expression group (P < 0.05), and the survival time of the PCK2 low-expression group was significantly higher than that of the high-expression group (P < 0.05) under the data obtained from this research (Figure S2). Gene expression by pathological stage was plotted in Stage plot, and in stage plot ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1 were expressed significantly differently in different clinical stages (P < 0.05), and there was no significant difference in ACSL3 (P > 0.05) (Figure S3). Except for the results of CHMP5, they were consistent with the conclusions obtained in this research (Figure S3).

Immunohistochemistry (IHC) analysis in the HPA database showed that the screened ferroptosis-related genes ACSL3, ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1 had high expression in HCC tissues (Figure S4), and PCK2 had low expression in HCC tissues, which could be used as supportive evidence for the conclusions of this research. According to other analyses of the HPA database, a comprehensive analysis showed that the screened ferroptosis-related genes ACSL3, ASNS, PCK2, PGD, and SLC38A1 were associated with the prognosis of HCC. CHMP5, MYB, and YY1AP1 were not associated with the prognosis.

Gene Set Enrichment Analysis (GSEA) Screening for Ferroptosis-Related Genes Functionally Related Signaling Pathways

To identify signaling pathways specifically regulated in HCC, the screened ferroptosis-related genes were subjected to GSEA gene enrichment analysis, and classified into high-expression and low-expression gene datasets based on median values. The signaling pathways that were the most significantly enriched with a q-value of 0.05 or P-value of 0.05 were examined.

Different degrees of high expression of ACSL3, ASNS, CHMP5, MYB, SLC38A1, and YY1AP1 were associated with apoptosis, endocytosis, pathways in cancer, ubiquitin-mediated proteolysis, and Wnt signaling pathways. Base excision repair, DNA replication, Cell cycle Glutathione metabolism, and Glutathione metabolism are enriched in the high PGD expression. Low expression of ASNS, MYB, and SLC38A1 was enriched in the PPAR signaling pathway, primary bile acid biosynthesis, fatty acid metabolism pathway, peroxisome, and drug metabolism cytochrome P450. Low expression of YY1AP1 was associated with all four pathways except peroxisome, ACSL3, and PGD were in the low expression associated with primary bile acid biosynthesis, PGD in the low expression was also associated with Linoleic acid metabolism and Complement and coagulation cascades (Figure 6).

Figure 6 The enrichment plots of screened ferroptosis-related gene expression levels in hepatocellular carcinoma (HCC) patients from Gene Set Enrichment Analysis (GSEA), including ACSL3 (A), ASNS (B), CHMP5 (C), MYB (D), PCK2 (E), PGD (F), SLC38A1 (G), YY1AP1 (H). Abbreviations: ACSL3, Acyl-CoA synthetase long-chain family member 3; ASNS, Asparagine synthetase; CHMP5, Charged multivesicular body protein 5; MYB, MYB proto-oncogene, transcription factor; PCK2, Phosphoenolpyruvate carboxykinase 2, mitochondrial Gene; PGD, Phosphogluconate Dehydrogenase; SLC38A1, Solute Carrier Family 38 Member 1; YY1AP1, YY1 Associated Protein 1.

PCK2 distinguishes itself from the other seven genes in that it is associated with apoptosis, endocytosis, pathways in cancer, ubiquitin-mediated proteolysis, and Wnt signaling pathways when expressed at low levels. While the PPAR signaling pathway, primary bile acid biosynthesis, fatty acid metabolism pathway, peroxisome, and drug metabolism cytochrome P450 were enriched in the PCK2 high expression (Figure 6E).

The mRNA Expression Validation of MYB, ASNS, SLC38A1, ACSL3, CHMP5, and PGD by Quantitative Reverse Transcription PCR (RT-qPCR)

The mRNA expression levels of MYB, ASNS, SLC38A1, ACSL3, CHMP5, and PGD were further validated by RT-qPCR in HCC cell lines and normal hepatocyte cell lines. Similar to the results of the TCGA database bioinformatics method for HCC samples, MYB, ASNS, SLC38A1, ACSL3, CHMP5, and PGD were checked in the MHCC97-L, HepG2, HCC-LM3, and Huh-7 cell lines compared to the L02 cell lines (Figure 7).

Figure 7 The validation of mRNA expressions of screened ferroptosis-related genes in L02, MHCC97-L, HepG2, and HCC-LM3 cell lines by RT-qPCR analysis, including MYB (A), ASNS(B), ASCL3(C), SLC38A1(D), CHMP5(E), PGD(F). Abbreviations: MYB, MYB proto-oncogene, transcription factor; ASNS, Asparagine synthetase; ACSL3, Acyl-CoA synthetase long-chain family member 3; SLC38A1, Solute Carrier Family 38 Member 1; CHMP5, Charged multivesicular body protein 5; PGD, Phosphogluconate Dehydrogenase. Notes: All data were presented as means ± standard deviation (SD), “*”: p<0.05, “**”: p<0.01, “***”: p<0.001, “ns”, p>0.05.

Differentially Expressed Gene Protein–Protein Interaction (PPI) Network

A protein–protein interaction (PPI) network comprising 165 protein nodes and 695 interaction edges has been established. The size and shading of the nodes denote the significance of the proteins in biological processes, with TP53 having the highest number of links, totaling 66. CHMP5 is connected to CHMP4 but has no association with other genes (Figure 8). This PPI network provides a comprehensive view of the complex interactions among ferroptosis-related genes differentially expressed in HCC. Within this network, core hub genes such as GPX4, SLC7A11, Acyl-CoA Synthetase Long-Chain Family Member 4 (ACSL4), Ferritin Heavy Chain 1 (FTH1), and Nuclear Receptor Coactivator 4 (NCOA4) emerge as central nodes exhibiting high connectivity. These genes are intimately involved in the regulation of iron metabolism, lipid peroxidation, and oxidative stress—key processes that influence ferroptotic cell death and, ultimately, tumor behavior. Their central positioning suggests they may serve as pivotal regulators of ferroptosis pathways, influencing cellular proliferation, apoptosis, and immune surveillance in the tumor microenvironment.

Figure 8 The Protein–Protein Interaction (PPI) network plot of ferroptosis-related hepatocellular carcinoma differential expression genes.

Notably, TP53, a well-characterized tumor suppressor gene, also occupies a critical position in the network. Its high degree of connectivity and interaction with other ferroptosis-associated genes underscores TP53’s multifaceted role in tumor biology. Through its regulatory effects on genes involved in iron homeostasis, redox balance, and lipid peroxidation, TP53 may modulate ferroptotic sensitivity in HCC cells, thereby impacting tumor progression, response to treatment, and overall patient prognosis.

Peripheral genes, connected indirectly to these central hubs, contribute to a multi-layered and dynamic signaling framework. Though less connected, these genes still participate in modulating the ferroptosis pathway at different levels, potentially influencing how tumor cells respond to external cues or therapeutic interventions. The resulting hierarchical and interactive structure underscores the complexity of ferroptosis regulation, as multiple genetic components interact to shape the phenotypic landscape of HCC.

Discussion

HCC is one of the most common malignant tumors worldwide. In the diagnosis and treatment of HCC, only a minority of HCC patients are detected in the early stage, while most of them are diagnosed in the middle or late stage, when they have lost the opportunity for surgical treatment, which is an important reason for the high mortality rate of HCC patients.^14,15 Ferroptosis has been shown to play a significant role in the pathophysiological process of many liver diseases and the development and progression of HCC.^16–19 Therefore, the research of ferroptosis-related genes is expected to provide new targets and prognostic markers for the diagnosis and treatment of HCC, which is of great clinical significance in guiding the prognosis of patients and prolonging survival.

In this research, based on the sample data of 407 HCCs and 58 normal cases in the TCGA database, 259 genes were screened in the FerrDb database for genes related to HCC. Based on the results of the literature search and the analysis of the data, eight ferroptosis-related genes were finally screened, which were ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1. Differential expression analysis of these eight genes revealed that ACSL3, ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1 were more expressed in HCC tumors than in normal tissues. The expression of PCK2 was lower in tumor tissues than in normal tissues. The expression of the eight genes, except PCK2, was significantly correlated with a lower survival rate of HCC, while the expression of PCK2 showed a correlation with a higher survival rate of HCC. Clinical correlation analyses yielded that the eight screened genes were all associated with multiple clinical factors. Multivariate Cox analysis revealed that mRNA expression of ACSL3, ASNS, MYB, PGD, SLC38A1, and YY1AP1 could serve as independent prognostic factors for HCC. To ensure the accuracy of the results, the expression difference analysis, clinical stage correlation analysis, and survival analysis of these eight genes were performed again with the use of GEPIA2, and the results were consistent with our results. The results also confirmed the protein expression differences of these eight genes between HCC tissues and the normal tissue using the HPA database. Finally, GSEA gene enrichment analysis was employed to uncover the potential signaling pathways associated with the screened eight genes. The analysis revealed differential enrichment of these eight genes in critical pathways, including apoptosis, endocytosis, pathways in cancer, Wnt signaling pathway, ubiquitin-mediated proteolysis, drug metabolism cytochrome P450, PPAR signaling pathway, peroxisome, primary bile acid biosynthesis, and fatty acid metabolism.

Multiple lines of evidence link the ferroptosis genes ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 to HCC and other cancers. ACSL3 is crucial for ferroptosis,²⁰ it is also a suppressor gene that inhibits the generation of ferroptosis in cancer cells leading to tumor growth, and it has been shown that oleic acid protects melanoma cells from ferroptosis in an ACSL3-dependent mode and increases the ability to develop metastatic tumors.^21,22 Another research showed that steatocytes can induce ferroptosis resistance in breast cancer cells by secreting specific fatty acids and that the process is dependent on the fatty acid synthase ACSL3.²³ It is promising to reverse the protection of adipocytes against ferroptosis by inhibiting fatty acid metabolism. In determining the resemblance between HCC cell lines and human liver tumors, the result indicated that HUH7 cells showed a higher expression of ASNS upregulated in tumors.²⁴ There is evidence that ASNS expression levels may also be inversely correlated with asparaginase efficacy in certain solid tumors as well.²⁵ Arginine is a second messenger-like molecule that promotes liver tumorigenesis by binding to RNA Binding Motif Protein 39 (RBM39) and controlling metabolic gene expression. Arginine-dependent ASNS expression further enhances arginine uptake, and promotes liver tumorigenesis.²⁶ In 2018, a finding showed that CHMP5 is a direct target of MIR429 in human colon cancer cell lines and suggested that CHMP5 up-regulation as a result of reduced MIR429 expression in DSS-induced mice colitis tissues and human ulcerative colitis (UC) tissues may restrict apoptosis and promote cell proliferation.²⁷ Li’s team bioinformatically determined that MYB expressions were associated with intensive infiltration of B cell, CD4⁺ T cell, CD8⁺ T cell, and macrophage, which was then validated by the correlation in clinical samples by IHC.²⁸ The CD36-BATF2/MYB signature serves as a robust predictor of anti-PD-1 immunotherapy response in Gastric Cancer (GC).²⁹ Phosphoenolpyruvate carboxykinase (PEPCK or PCK) catalyzes the first rate-limiting step in the hepatic gluconeogenesis pathway to maintain blood glucose levels. A study has reported that the expressions of both PCK1 and PCK2 genes are downregulated in primary HCC and low PCK2 expression was associated with poor prognosis in patients with HCC. Forced expression of either PCK1 or PCK2 in liver cancer cell lines results in severe apoptosis under the condition of glucose deprivation and suppressed liver tumorigenesis in mice.³⁰ Treatment of HepG2 cells with ribulose-5-phosphate, a catalytic product of PGD, gave rise to a concentration-dependent upregulation of Nrf2, and then Nrf2 promotes hepatoma cell growth and progression.³¹ Wnt/β-catenin works in association with SLC38A1 to induce the development of hepatoblastoma.³² Studies also revealed that YY1AP1 silencing eliminates oncogene addiction by altering the chromatin landscape and triggering massive apoptosis in vitro and tumor suppression in vivo. YY1AP1 expression promotes HCC proliferation.³³ YY1AP1 expression plays an important role in the tumorigenesis and progression of colon adenocarcinoma (COAD) and could serve as a clinical prognostic indicator for COAD.³⁴

The Wnt/β-catenin signaling pathway, which is involved in biological processes that regulate development, differentiation, and adult tissue homeostasis, and its aberrant activation is highly correlated with tumorigenesis and metastasis.^35,36 It was found that Keratin 15 (KRT15) may have promoted breast cancer progression through mechanisms such as regulating the Wnt/β-catenin signaling pathway and inhibiting ferroptosis.³⁷ The up-regulation of ferroptosis-related genes ACSL3, ASNS, CHMP5, MYB, PCK2, PGD, SLC38A1, and YY1AP1 screened in this research were enriched in the Wnt signaling pathway, and it was hypothesized that ferroptosis might be correlated with the Wnt signaling pathway.

In the constructed differential gene PPI network, the TP53 gene was observed to interact extensively with multiple genes. TP53 is an important tumor suppressor gene within cells, and mutations in the TP53 gene have been found in various cancers.³⁸ The interaction map indicates that the TP53 gene plays a significant role in the mechanism by which ferroptosis affects the pathology of HCC.³⁹

The treatment of cancer by utilizing ferroptosis is still in the exploratory stage, as the currently used apoptosis-inducing tumor therapeutic pathways stimulate cellular anti-apoptotic mechanisms, thus making cancer cells resistant to anticancer drugs. A study showed that ferroptosis heterogeneity in triple-negative breast cancer reveals an innovative immunotherapy combination strategy.⁴⁰ Intriguingly, ferroptosis has been linked to cancer therapy resistance, and induction of ferroptosis can potentially reverse this resistance. In recent years, certain drugs and compounds have been found to have the ability to induce ferroptosis and demonstrate anti-tumor activity.⁴¹ Therefore, inducing new forms of regulated cell death like ferroptosis, and conducting research on ferroptosis-associated genes to search for new targets and prognostic markers, are of great value for the development of new therapeutic approaches for tumors.

Conclusions

In this research, it was found that the high expression of ferroptosis-related genes ACSL3, ASNS, CHMP5, MYB, PGD, SLC38A1, and YY1AP1, as well as the low expression of PCK2, were intimately associated with the prognosis of HCC. These genes may become markers and novel targets for early diagnosis, precise treatment, and prognosis assessment of HCC.