Palmitoylation-related gene signature for immune subtyping, prognostic stratification, and immunotherapy guidance in hepatocellular carcinoma
Highlight box
Key findings
• A 3-gene palmitoylation signature was constructed to predict prognosis and conduct risk stratification for hepatocellular carcinoma (HCC). PPT1 and ZDHHC18 facilitate HCC cell proliferation and mediate tumor immune suppression. This signature reflects HCC immune evasion features and provides theoretical guidance for HCC immunotherapy.
What is known and what is new?
• Previous studies have confirmed that protein palmitoylation participates in liver tumor progression, while single palmitoyl-related genes cannot accurately evaluate HCC patient prognosis. The crosstalk between palmitoylation modification and HCC tumor immunity remains poorly clarified.
• This work established a novel three-gene palmitoylation risk signature specifically for HCC. We further verified that PPT1 and ZDHHC18 are functional oncogenes driving HCC growth and immune escape, and validated the signature’s value to instruct immunotherapeutic strategies.
What is the implication, and what should change now?
• The palmitoylation signature can assist clinicians in stratifying HCC patients with distinct prognostic risks. This study highlights palmitoylation regulators as promising targets to reverse HCC immune suppression, which may promote personalized immunotherapy design for HCC in future clinical practice.
Introduction
Hepatocellular carcinoma (HCC), which accounts for 75–85% of primary liver cancers, remains a significant global health burden, with over 900,000 new cases and 830,000 deaths reported worldwide in 2020 (1). As the third leading cause of cancer-related mortality (2), HCC incidence is primarily driven by hepatitis B and C infections, as well as exposure to environmental carcinogens such as aflatoxin (3). Moreover, HCC usually develops along the progressive trajectory from metabolic dysfunction-associated steatotic liver disease (MASLD) to liver fibrosis and ultimately HCC (4,5). Despite recent advancements in early diagnostic techniques for liver cancer, the disease’s insidious onset and absence of noticeable early symptoms often result in the majority of cases being identified only during intermediate or advanced stages (6). Current systemic therapies, including multi-kinase inhibitors and immune checkpoint inhibitors, face significant limitations due to the development of resistance and high recurrence rates, with 70% of patients experiencing recurrence within 5 years post-surgery (7). These persistent challenges highlight the critical need to identify novel molecular drivers of HCC progression and therapy resistance, which could pave the way for more effective diagnostic and therapeutic strategies.
Palmitoylation, a reversible form of post-translational modification (PTM), plays a crucial role in modulating protein stability, localization, and functionality, thereby contributing to the development of cancer (8). Research indicates that enzymes responsible for palmitoylation, such as the ZDHHC family of palmitoyltransferases (9), as well as depalmitoylating enzymes like ABHD17C, are key regulators of processes including tumor cell proliferation, invasion, and immune escape mechanisms (10). Emerging evidence highlights the significant involvement of ZDHHC enzymes in tumor progression and the tumor microenvironment. In gliomas, aberrant expression of specific ZDHHC members (e.g., ZDHHC11, ZDHHC12, ZDHHC15, ZDHHC22, and ZDHHC23) regulates tumor cell viability, autophagy, and apoptosis, likely via the PI3K/AKT pathway. ZDHHC inhibition suppresses tumor growth, enhances chemosensitivity, and reduces microglial infiltration, highlighting their therapeutic potential in cancer (11). Furthermore, recent studies have identified ZDHHC enzyme genes, such as ZDHHC3, as prognostic markers and key oncogenic drivers in pancreatic cancer, emphasizing their potential as immunotherapeutic targets. Notably, combination therapy using the ZDHHC inhibitor 2-bromopalmitate (2-BP) and programmed cell death protein 1 (PD-1)/programmed death-ligand 1 (PD-L1) blockade has demonstrated synergistic therapeutic effects in pancreatic cancer mouse models, further supporting the clinical relevance of targeting ZDHHCs in cancer treatment (12). These findings collectively highlight the multifaceted roles of palmitoylation in cancer biology and its promise for advancing therapeutic strategies.
Given the critical role of palmitoylation in cancer progression and its emerging potential as a therapeutic target, alongside other forms of protein acylation such as farnesylation that participates in Ras membrane trafficking, this study aims to explore the diagnostic and prognostic value of palmitoylation-related genes in HCC. Specifically, by integrating transcriptomic and clinical data from HCC patient cohorts, we aim to identify key palmitoylation-related genes that can serve as novel biomarkers for prognosis stratification and risk assessment in HCC. Through this comprehensive approach, we hope to uncover new insights into the molecular mechanisms driven by palmitoylation in HCC and provide a foundation for developing targeted strategies to enhance diagnosis, prognosis, and treatment of this devastating disease. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2026-0745/rc).
Methods
Data acquisition
To develop a predictive model for HCC, transcriptome data along with clinicopathologic information from 377 HCC patients (training cohort) were obtained from The Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov/). These data served as the foundation for constructing the model based on genes associated with palmitoylation modification. For validation purposes, transcriptome data and clinicopathologic details from 231 HCC specimens (validation cohort) were sourced from the International Cancer Genome Consortium (ICGC) database (https://dcc.icgc.org/releases). The ICGC dataset was primarily used to assess the performance and reliability of the model developed with TCGA data. A set of 31 genes implicated in palmitoylation modification was selected based on prior studies (Table S1) (13-15). This study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.
Differential gene expression analysis and palmitoylation-related gene extraction
To identify differentially expressed genes (DEGs) between tumor and adjacent normal tissues, we performed differential expression analysis using the limma package in R software based on HCC transcriptome data. The selection criteria for DEGs were set as |log2FC| >0 and false discovery rate (FDR) <0.05. This ensured that only genes with significant expression differences between tumor and adjacent normal tissues were included for further analysis. Following the identification of these DEGs, we specifically extracted the expression levels of genes associated with palmitoylation modification. The expression levels of these palmitoylation-related genes were then evaluated and incorporated into the subsequent model development, allowing for a deeper exploration of their potential role in HCC.
Univariate Cox and least absolute shrinkage and selection operator (LASSO) regression for model construction
To construct a predictive model based on palmitoylation-related genes, we began with univariate Cox regression analysis, which was applied to the DEGs associated with palmitoylation modification. This analysis aimed to identify genes that exhibited a significant relationship with patient prognosis, allowing us to focus on those genes with the strongest potential as prognostic biomarkers. Afterward, we employed LASSO regression analysis to further refine the set of candidate genes. LASSO regression, which utilizes L1 regularization, effectively reduces the dimensionality of the data by penalizing less relevant predictors, thereby selecting a smaller subset of genes that contribute most significantly to the prediction of outcomes. This dual approach-starting with univariate Cox analysis to identify potential markers, followed by LASSO regression to select the most robust and non-redundant variables-enabled us to construct a reliable and efficient predictive model. The resulting model, based on the expression of key palmitoylation-related genes, provides valuable insights into their prognostic significance in HCC and holds potential for clinical application in predicting patient outcomes.
Model validation and assessment of prognostic independence
To evaluate the predictive performance of the model constructed using TCGA data, we conducted an external validation using the ICGC dataset. This validation step was crucial for assessing the model’s reliability and generalizability across independent cohorts. To further assess the clinical applicability of the model, we incorporated key clinical factors-such as gender, age, tumor stage, and tumor grade-into the analysis. By performing a combined approach of univariate and multivariate Cox regression analysis, we evaluated whether the predictive model, built on palmitoylation-related genes, could serve as an independent prognostic factor for HCC. This approach allowed us to determine whether the model could provide meaningful prognostic information beyond traditional clinical factors, such as tumor stage and grade. Our analysis aimed to establish whether the palmitoylation-related gene-based model could enhance prognostication and potentially offer a more precise tool for predicting patient outcomes in HCC, thus contributing to the development of personalized treatment strategies.
Risk stratification and differential gene expression analysis
Based on the predictive model, we stratified the HCC samples from both the TCGA and ICGC datasets into high-risk and low-risk groups according to the median risk score. To identify genes that were differentially expressed between these two groups, we applied the same filtering criteria used previously: |log2FC| >0 and FDR <0.05. This allowed us to isolate genes that exhibited significant expression differences between the high-risk and low-risk groups in both datasets. Subsequently, we performed functional enrichment analysis on the identified DEGs, which included Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses. Specifically, the GO analysis was divided into three categories: biological process (BP), cellular component (CC), and molecular function (MF), to gain a comprehensive understanding of the biological roles of these genes. The KEGG pathway analysis further helped to reveal the potential pathways that might be involved in HCC progression and prognosis, providing deeper insights into the molecular mechanisms underlying the high- and low-risk stratifications.
In addition to these analyses, we also performed single-sample Gene Set Enrichment Analysis (ssGSEA) to assess the enrichment scores of predefined gene sets in individual HCC samples. This approach allowed us to evaluate the immune and other biological pathways in both the high-risk and low-risk groups, further refining our understanding of the functional differences between the two risk categories and their potential implications in HCC prognosis.
Immune microenvironment analysis based on the predictive model
To further explore the role of the predictive model in assessing the immune microenvironment of HCC, we conducted an analysis using the TCGA database. Based on the stratification of samples into high-risk and low-risk groups according to the predictive model, we evaluated differences between these groups in several key aspects of the immune microenvironment, including ESTIMATE (Estimation of STromal and Immune cells in MAlignant Tumour tissues using Expression data) scores, immune checkpoint expression, and human leukocyte antigen (HLA) gene expression. Specifically, the ESTIMATE algorithm, which is based on ssGSEA, was employed to calculate the stromal score and immune score using the transcriptome profiles of HCC samples, aiming to predict the levels of infiltrating stromal and immune cells in the tumor microenvironment. The comprehensive ESTIMATE score was further generated by combining these two scores (ESTIMATE score = Stromal score + Immune score), which is widely used to analyze tumor purity. We compared the stromal score, immune score, and comprehensive ESTIMATE score between the high- and low-risk groups to systematically assess differences in tumor immune infiltration levels and non-tumor cellular components between the two subgroups. Additionally, we examined the expression levels of immune checkpoint genes, such as PD-1, PD-L1, and CTLA-4, which are crucial in regulating immune responses in cancer. We also evaluated the expression of HLA genes, which are important for antigen presentation and immune system recognition of tumor cells. These analyses provided insights into how the risk groups identified by the predictive model might correlate with immune activity in HCC, potentially revealing the model’s utility in guiding personalized treatment strategies.
Evaluation of prognostic value beyond overall survival (OS) and construction of nomograms
In addition to assessing the predictive value of the model for OS, we extended our evaluation to other important survival metrics, including progression-free survival (PFS), disease-free survival (DFS), and disease-specific survival (DSS). By analyzing these additional endpoints, we aimed to gain a comprehensive understanding of the model’s broader applicability in predicting various aspects of HCC prognosis. To visually represent the model’s prognostic capacity, we constructed nomograms based on the predictive variables identified in the model. These nomograms provide a user-friendly way to quantify individual risk scores and predict patient outcomes. Furthermore, to assess the accuracy of the nomograms, we calculated the Concordance Index (C-index), which quantifies the discriminative ability of the model in predicting patient survival. A higher C-index indicates better predictive accuracy. By evaluating the C-index across OS, PFS, DFS, and DSS, we were able to validate the robustness and reliability of the predictive model, offering insights into its potential utility in clinical decision-making and personalized treatment planning for HCC patients.
Identification of core target genes and exploration of clinical and microenvironmental relevance
To further narrow down the core target genes for HCC, we employed three different machine learning approaches: LASSO regression, random forest, and support vector machine (SVM). All machine learning procedures were implemented within a rigorous cross-validation framework to avoid information leakage. Feature selection was performed independently within each cross-validation fold, and hyperparameter tuning was conducted using grid search with five-fold cross-validation. Genes identified by univariate Cox regression were analyzed using random forest and SVM. The core target genes were defined as the intersection of genes selected by LASSO regression, random forest, and SVM. By analyzing the results from each method, we identified the intersecting genes common to all three analyses, which were considered as the core target genes. Subsequently, we explored the clinical significance of these core target genes using data from the TCGA database. We assessed their potential roles in patient survival and other clinical outcomes to evaluate their value as prognostic biomarkers. To further understand the spatial distribution and functional relevance of these genes within the tumor microenvironment, we leveraged single-cell RNA sequencing (scRNA-seq) data. By analyzing expression patterns across different cell types in the liver cancer microenvironment, we sought to determine how these core target genes are distributed and expressed in various tumor and stromal cells. This analysis provided deeper insights into the molecular mechanisms by which these genes may influence HCC progression and the tumor’s response to therapeutic interventions. We will use the TIDE (Tumor Immune Dysfunction and Exclusion) algorithm to predict potential immunotherapy responses (16); the main purpose of this analysis is to explore the correlation between the expression levels of these core target genes and the response to immunotherapy, so as to clarify whether the expression of core target genes is related to the sensitivity of patients to immunotherapy. Specifically, we will conduct this analysis with the help of the online tool (Home for Researchers).
Cell culture
L02, Hep3B, and Huh7 cells were acquired from the American Type Culture Collection (ATCC) located in Rockville, MD, USA. The human normal liver cell line L02 was cultivated in RPMI-1640 medium that was supplemented with 10% fetal bovine serum (FBS), along with 100 U/mL penicillin and 100 µg/mL streptomycin. Meanwhile, the human HCC cell lines Hep3B and Huh7 were grown in high-glucose Dulbecco’s Modified Eagle Medium (DMEM) containing 10% FBS, 100 U/mL penicillin, and 100 µg/mL streptomycin. All cells were incubated at 37 ℃ in a humidified atmosphere with 5% CO2. The culture medium was refreshed every two days. When the cells reached 80–90% confluence, subculturing was performed using 0.25% trypsin-EDTA (ethylenediaminetetraacetic acid).
RNA extraction and quantitative reverse transcription polymerase chain reaction (RT-qPCR)
RNA was extracted from the following cell lines: L02, Hep3B, and Huh7 cells, using the FastPure Cell/Tissue Total RNA Isolation Kit V2 (No. RC112-01, Vazyme Biotech Co., Ltd., Nanjing, China), following the manufacturer’s instructions. In the process of reverse transcription, 1 µg of total RNA was transformed into complementary DNA (cDNA) with the aid of the HiScript III 1st Strand cDNA Synthesis Kit (+gDNA wiper) (No. R312-01, Vazyme Biotech Co., Ltd., Nanjing, China). The concentration and purity of RNA were determined by measuring the optical density ratio at 260 nm to 280 nm, where an A260/A280 ratio ranging from 1.8 to 2.0 was considered acceptable. For RT-qPCR, the ChamQ SYBR qPCR Master Mix (No. Q311-02, Vazyme Biotech Co., Ltd., Nanjing, China) was employed on a Roche Light-Cycler 480 (Roche, Basel, Switzerland). The relative expression levels of the target genes were computed using the 2−ΔΔCt method, and each experiment was repeated three times to ensure reliability. β-Actin was used as the internal control. The primers for PPT1, ZDHHC18, and β-Actin were designed as follows: PPT1, forward: 5'-TGTTTTTGGACTCCCTCGATG-3', reverse: 5'-CATGCCAGTATTCGGCTTGC-3'; ZDHHC18, forward: 5'-ACCGGCCTCTTCTTCGTCT-3', reverse: 5'-AACTGCCTGTGTTGTCGATCT-3'; β-Actin, forward: 5'-AGCGAGCATCCCCCAAAGTT-3', reverse: 5'-GGGCACGAAGGCTCATCATT-3'.
Statistical analysis
All statistical analyses were performed using R software (version 4.4.0). To compare differences between two groups, the Wilcoxon test was used, while for comparisons among more than two groups, the Kruskal-Wallis test was applied. Immune cell ssGSEA scores and immune functions between low- and high-risk groups were quantified using the Mann-Whitney test. For immune-related analyses and pathway analyses, the Benjamini-Hochberg method was used for multiple-testing correction to control the FDR. Survival differences between low- and high-risk groups were assessed using Kaplan-Meier curves, with significance determined by the log-rank test. The predictive accuracy of the models was evaluated through receiver operating characteristic (ROC) curves and the corresponding area under the ROC curve (AUC). Univariate and multivariate Cox regression analyses were used to identify independent prognostic factors for OS in HCC patients. Additionally, Spearman’s correlation analysis was employed to describe the relationships between non-normally distributed quantitative variables. The threshold for statistical significance was set at a two-tailed P value of ≤0.05.
Results
Identification of differentially expressed palmitoylation-related genes and prognostic analysis
Figure 1 illustrates the overall analysis workflow, outlining the key steps from data collection to the development and validation of the predictive model based on palmitoylation-related genes. Figure 2A shows the 20 differentially expressed palmitoylation-related genes between HCC tissues and adjacent normal tissues, identified through the TCGA dataset. These genes exhibited significant expression differences, highlighting their potential role in HCC progression. To further investigate their prognostic significance, we performed univariate Cox regression analysis, which revealed eight genes that were significantly associated with OS in HCC patients (Figure 2B). These eight genes were selected as key biomarkers and incorporated into the construction of the predictive model, aiming to provide an improved prognostic tool for HCC patients.
Based on the TCGA database, we further visualized the differential expression of these eight genes between HCC tissues and adjacent normal tissues using a heatmap (Figure 2C). This heatmap clearly illustrated the expression patterns of these genes, emphasizing their differential regulation in tumor versus normal tissues. Additionally, we explored the correlations between these eight genes (Figure 2D), revealing potential interactions that may contribute to the biological mechanisms underlying HCC progression.
Construction of predictive model based on LASSO regression
Building upon the eight genes identified through univariate Cox regression analysis, we further utilized LASSO regression to construct a predictive model for HCC prognosis (Figure 2E,2F). This approach enabled us to refine the selection of important genes and assess their combined prognostic value. The risk score for each sample was calculated using the following formula:
Risk score = (0.439×ZDHHC18) + (0.099×ZDHHC23) + (0.268×PPT1)
Where ZDHHC18, ZDHHC23, and PPT1 represent the expression levels of the respective genes in each sample. This risk score serves as a quantitative measure to stratify patients into high-risk and low-risk groups, based on the combined contribution of these genes to OS in HCC patients.
Validation of predictive model performance using TCGA and ICGC datasets
To validate the predictive performance of our model, we utilized both the TCGA and ICGC databases. First, Kaplan-Meier survival curves were generated to assess the OS differences between the high-risk and low-risk groups based on the risk scores derived from the model. In both the TCGA (Figure 3A, P=3.911e-03) and ICGC (Figure 3B, P=2.524e−03) datasets, the high-risk group showed significantly poorer survival rates compared to the low-risk group, confirming that the model effectively stratifies patients by prognosis and validates the risk score as a reliable indicator of patient outcomes.
Additionally, we evaluated the model’s accuracy by constructing ROC curves and calculating the AUC for both datasets. For the TCGA dataset (Figure 3C), the AUC values were 0.700 at 1 year, 0.653 at 2 years, and 0.634 at 3 years, demonstrating moderate predictive performance. In contrast, for the ICGC dataset (Figure 3D), the AUC values were higher, with 0.773 at 1 year, 0.705 at 2 years, and 0.703 at 3 years, indicating a stronger predictive capability.
To further assess the stability and robustness of the model, Principal Component Analysis (PCA) and t-Distributed Stochastic Neighbor Embedding (t-SNE) analyses were performed to visualize the clustering of high- and low-risk patients based on their gene expression profiles and risk scores. Specifically, PCA and t-SNE results for the TCGA dataset are shown in Figure 3E and Figure 3F, respectively, while those for the ICGC dataset are presented in Figure 3G and Figure 3H.
Finally, to explore whether the predictive model could serve as an independent prognostic factor, we incorporated other clinical characteristics, such as age, gender, tumor stage, and grade, into the analysis. We performed both univariate and multivariate Cox regression analyses to evaluate whether the risk score derived from the model remained a significant and independent prognostic factor for HCC patients when adjusted for these clinical traits. Based on the TCGA (Figure 4A) and ICGC (Figure 4B) datasets, we consistently found that the predictive model could serve as an independent prognostic factor, demonstrating its robustness and generalizability across different cohorts.
Functional enrichment analysis of DEGs between high- and low-risk groups in the TCGA and ICGC datasets
To further explore the biological significance of the risk score-based stratification, we compared the DEGs between the high-risk and low-risk groups using both the TCGA and ICGC databases. For each dataset, genes exhibiting significant differential expression between the two risk groups were identified using the threshold of |log2FC| >0 and FDR <0.05. This rigorous filtering ensured that only the most relevant genes, with robust expression differences, were included in the subsequent analyses.
Following DEG identification, we performed GO and KEGG pathway enrichment analyses to further investigate the potential BPs and signaling pathways associated with these DEGs. The GO analysis categorized the DEGs into three main functional categories: BP, CC, and MF. For TCGA (Figure 4C), the top enriched BP terms included immune-related processes such as leukocyte-mediated immunity, immunoglobulin-mediated immune response, and lymphocyte-mediated immunity. Other key BPs included mitotic sister chromatid segregation and nuclear chromosome segregation, which are associated with cell division and cancer cell proliferation. In the CC category, significant enrichment was observed in immunoglobulin complex, outer kinetochore, condensed chromosomes, and the spindle, which are important for immune function and chromosome segregation. For MF, antigen binding, immunoglobulin receptor binding, and single-stranded DNA helicase activity were significantly enriched.
For the ICGC database (Figure 4D), the top enriched BP terms included leukocyte cell-cell adhesion, regulation of T cell activation, and positive regulation of lymphocyte and leukocyte activation. These processes suggest significant involvement of immune responses in the high-risk group. In addition, leukocyte proliferation and lymphocyte proliferation were enriched, reflecting the dynamic immune response in these patients. The CC analysis revealed enrichment in major histocompatibility complex (MHC) class II protein complex, collagen-containing extracellular matrix, and the external side of the plasma membrane, all of which play key roles in immune cell signaling and tissue remodeling. In the MF category, MHC protein complex binding, immune receptor activity, and peptide antigen binding were significantly enriched, indicating a strong immune-related molecular activity.
Additionally, the KEGG pathway analysis for TCGA (Figure 4E) revealed pathways such as the cell cycle, rheumatoid arthritis, and human T-cell leukemia virus 1 infection, which are related to immune regulation and cancer progression. For ICGC (Figure 4F), the top enriched pathways included hematopoietic cell lineage, Staphylococcus aureus infection, phagosome, and rheumatoid arthritis. These pathways are associated with immune cell differentiation, infection response, and immune system modulation.
These functional enrichment analyses not only provide deeper insights into the biological relevance of the DEGs but also validate the robustness of the predictive model by linking it to well-known molecular and cellular processes in HCC. The enrichment of immune-associated pathways and BPs further underscores the critical role of immune modulation in HCC prognosis, thereby validating the potential of our risk model for predicting disease progression and clinical outcomes.
ssGSEA analysis of immune status differences between high- and low-risk groups
To further explore the discrepancies in immune landscape between high- and low-risk subgroups, we implemented ssGSEA. The ssGSEA results revealed significant differences in the enrichment scores of various immune cell subtypes between the high-risk and low-risk groups. Specifically, the high-risk group exhibited altered levels of immune cell infiltration, with increased or decreased activity in several immune cell types (Figure 5A,5B), such as immature dendritic cells (iDCs), follicular helper T cells (Tfh cells), type 2 T helper cells (Th2 cells), regulatory T cells (Tregs), and natural killer cells (NK cells), compared to the low-risk group. These immune cell variations were found to be associated with the poor prognosis observed in the high-risk group.
In addition to characterizing immune cell subtypes, we conducted a comprehensive analysis of immune-related pathways, identifying significant differential activity between high- and low-risk HCC groups in pathways associated with immune evasion and tumor progression (Figure 5C,5D). Specifically, elevated activity in pathways such as antigen-presenting cell (APC) co-inhibition, chemokine receptor (CCR) signaling, immune checkpoint regulation, inflammation-promoting mechanisms, and parainflammation has been linked to poor prognosis in HCC. Conversely, reduced activity in the Type II IFN response pathway is also associated with adverse outcomes. Notably, the high-risk group demonstrated increased activity in these pro-tumorigenic pathways, including enhanced immune checkpoint signaling and reduced antigen presentation capacity, which collectively contribute to immune escape and disease progression. These findings corroborate the model’s ability to reflect critical immune alterations that underpin HCC clinical outcomes, thereby supporting its potential utility in prognostic stratification and therapeutic targeting.
These findings highlight the importance of immune status in the prognosis of HCC and demonstrate that the predictive model not only stratifies patients based on molecular factors but also reflects immune microenvironment differences that may influence tumor behavior and patient survival.
Immune status differences between high- and low-risk groups based on TCGA database: ESTIMATE scores, HLA, and immune checkpoints
To further investigate the immune status differences between the high-risk and low-risk groups, we explored several key immune-related factors using the TCGA database, including ESTIMATE scores, HLA expression levels, and immune checkpoint molecule expression.
First, we compared the ESTIMATE scores (Figure 6A), which reflect the levels of stromal and immune cell infiltration in the tumor microenvironment. We observed that the immune scores and ESTIMATE scores were markedly elevated in the high-risk group relative to the low-risk group, accompanied by a reduction in tumor purity. These findings were statistically significant. However, there was no significant difference in the stromal score between the high-risk group and the low-risk group.
In terms of HLA expression and immune checkpoint molecules, significant differences were observed between the high-risk and low-risk groups (Figure 6B). The high-risk group exhibited altered HLA expression, suggesting an altered immune recognition profile that could contribute to immune evasion. Similarly, higher expression levels of immune checkpoint molecules, such as PD-1, PD-L1, and CTLA-4, were observed in the high-risk group, pointing to the activation of immune escape mechanisms (Figure 6C). These findings indicate that the poor prognosis in the high-risk group is closely associated with immune escape activation and weakened antigen presentation, further supporting the concept that immune evasion plays a critical role in the high-risk group’s adverse outcomes.
These findings underscore significant differences in immune status between the high-risk and low-risk groups, highlighting the potential role of immune-related factors in predicting patient outcomes in HCC. The variations in ESTIMATE scores, immune scores, stromal scores, tumor purity, HLA expression, and immune checkpoint molecule levels further validate the impact of the immune microenvironment on the prognosis of HCC.
Survival differences and nomogram construction based on the predictive model
Based on the TCGA database and the constructed predictive model, we performed an in-depth analysis of the survival differences in PFS, DFS, and DSS among patients in the high-risk and low-risk groups.
The Kaplan-Meier survival curves revealed distinct differences between the two risk groups for all three survival endpoints. For PFS, DFS, and DSS, the survival rates of the high-risk group were consistently lower than those of the low-risk group. Specifically, in the PFS analysis, the P value was P<0.001, indicating a highly significant difference between the two groups (Figure S1A). For DFS, the P value was P=0.01, signifying a statistically significant disparity (Figure S1A). In the DSS analysis, the P value was P=0.004, further highlighting the significant difference in survival outcomes between the high-risk and low-risk groups (Figure S1A). These results suggest that patients in the high-risk group have a poorer prognosis across all these survival measures.
To provide a more intuitive and practical tool for predicting patient survival, we constructed a nomogram that incorporated all four survival endpoints: OS, PFS, DFS, and DSS. This comprehensive nomogram integrated all significant factors identified in the predictive model, enabling individualized risk assessment. Each factor was assigned a specific score, and the total score could be used to calculate the probability of survival at different time points for each of the four survival endpoints. The calibration plots demonstrated good agreement between the predicted probabilities and the observed outcomes, indicating the high accuracy of the nomogram.
To evaluate the discriminatory ability of the nomogram, we conducted internal validation using bootstrap resampling and plotted the C-index curves for 1-, 2-, and 3-year survival. The nomograms showed stable and consistent discrimination across all four endpoints: OS (Figure 7A), PFS (Figure S1B), DFS (Figure S1C), and DSS (Figure S1D). The C-index curves allowed us to visually assess the discriminatory ability of the nomogram across different time points, with the relatively stable trends suggesting that the nomogram could effectively distinguish between high-risk and low-risk patients in terms of their survival prospects.
Furthermore, the results of the internal validation confirmed the robustness of the nomogram, as the C-index curves remained consistent across multiple resamples for all four survival endpoints.
Overall, the construction of this nomogram, encompassing all four survival measures, along with the assessment of its accuracy, provides a valuable tool for clinicians. It enables them to predict the survival outcomes of patients with HCC based on the identified risk factors, facilitating more informed treatment decisions and potentially improving patient management in clinical practice.
Identification and functional characterization of core genes PPT1 and ZDHHC18 in HCC
To screen for core genes, we extracted the intersection genes of the key genes identified by LASSO regression analysis, random forest (Figure 7B), and SVM (Figure 7C), which were designated as PPT1 and ZDHHC18 (Figure 7D). To validate these findings, we performed RT-qPCR experiments using the normal liver cell line L02 and the HCC cell lines Hep3B and Huh7. The results consistently confirmed that the transcriptional levels of PPT1 and ZDHHC18 were significantly higher in the tumor cell lines compared to the normal liver cells (Figure 7E).
To further understand the expression distribution of PPT1 and ZDHHC18 in the immune microenvironment of HCC, we utilized single-cell datasets for analysis. First, based on the single-cell dataset GSE98638, PPT1 was predominantly enriched in tumor proliferating cells, with significantly higher expression in this cell subset relative to other cell types (Figure 7F).
Meanwhile, ZDHHC18 was predominantly enriched in regulatory T cells (Tregs), showing markedly higher expression in Tregs than in other immune cell subsets (Figure 7G).
To validate the reliability of these results, we conducted a similar analysis on another single-cell dataset, GSE140228. The results showed a similar pattern, with PPT1 exhibiting high expression levels in cell proliferation-related cells, further supporting its critical role in the tumor proliferation process (Figure 7H). Meanwhile, ZDHHC18 was primarily expressed in Treg cells (Figure 7I), and its high expression in this cell type was confirmed in this dataset.
After identifying these two core genes, we further explored their expression levels in relation to the staging of HCC patients based on the TCGA database. The results indicated a significant positive correlation between the expression levels of PPT1 and ZDHHC18 and tumor progression. Specifically, as their expression levels increased, the tumor stage tended to be more advanced (Figure 8A). Statistical analysis showed a highly significant positive correlation (P<0.001 for both genes), providing strong evidence for the association between gene expression and tumor staging. This finding suggests that the higher expression of these core genes may promote tumor growth and progression. Additionally, we analyzed the relationship between the expression levels of PPT1 and ZDHHC18 and the response to immunotherapy (Figure 8B). Our results revealed that the expression levels of both genes were significantly correlated with immunotherapy response, indicating their potential roles in modulating the efficacy of immune-based treatments.
We then analyzed the relationship between PPT1/ZDHHC18 expression and biological pathways. PPT1 (Figure 8C) and ZDHHC18 (Figure 8D) were positively correlated with pathways related to cell proliferation and tumor metastasis, including DNA_replication, G2M_checkpoint, Tumor_proliferation_signature, and EMT_markers. These results indicate that they may promote HCC progression by influencing these pathways. Furthermore, we analyzed the correlation between PPT1 and ZDHHC18 expression levels and found a significant positive correlation between the two genes (Figure 8E).
Furthermore, the significant positive correlation between the expression levels of these genes and advanced tumor stages, as well as their association with immunotherapy response, underscores their potential as key drivers of HCC progression and therapeutic targets.
Discussion
The role of palmitoylation in HCC is increasingly recognized as a critical aspect of tumorigenesis and progression. Palmitoylation, a reversible PTM, involves the addition of a palmitate moiety to cysteine residues of proteins, thereby influencing their stability, localization, and function (17). This modification is mediated by a family of enzymes known as ZDHHC20 (18) palmitoyltransferases and is counteracted by depalmitoylating enzymes such as ABHD17C (19). As a classic depalmitoylase, ABHD10 modulates tumor progression by regulating the protein palmitoylation status, maintaining cellular redox homeostasis, and mediating intracellular signal transduction (20). The dysregulation of palmitoylation has been implicated in various cancers, including HCC (21), where it contributes to tumor cell proliferation, invasion, and immune evasion. Additionally, studies have demonstrated that fatty acid synthase (FASN) participates in the palmitoylation of STING via its product palmitic acid, thereby regulating STING’s function in immune responses—a process that plays a crucial role in mitigating liver injury induced by sepsis (22).
In the context of HCC, palmitoylation has been shown to regulate key oncogenic pathways. For instance, the PI3K/AKT pathway, which is frequently hyperactivated in HCC, is influenced by palmitoylation (21). Aberrant palmitoylation of key proteins in this pathway can lead to enhanced tumor cell survival and proliferation. Additionally, palmitoylation can modulate the activity of immune checkpoints, such as TIM-3 (23), thereby affecting the tumor immune microenvironment. The identification of palmitoylation-related genes as potential biomarkers and therapeutic targets underscores the importance of understanding this PTM in HCC.
The identification of palmitoylation-related genes in HCC offers promising avenues for diagnostic and prognostic applications. By integrating transcriptomic profiles and clinical information from HCC patient cohorts, we have identified key palmitoylation-related genes that could serve as novel biomarkers for early diagnosis and prognosis prediction. These genes, such as PPT1 and ZDHHC18, exhibit significant expression differences between tumor and adjacent normal tissues, highlighting their potential as biomarkers for HCC detection.
Moreover, the prognostic significance of PPT1 and ZDHHC18 was tightly correlated with tumor pathological stage. TCGA data analysis showed that the expression of both genes was significantly positively associated with HCC progression, with higher expression corresponding to advanced pathological stage. The utilization of scRNA-seq datasets has provided valuable insights into the expression patterns of palmitoylation-related genes in the HCC immune microenvironment. PPT1 is predominantly enriched in tumor proliferating cells. Of note, a previous study has verified that PPT1 is markedly upregulated in HCC tissues relative to normal tissues, which is consistent with our results and further validates that elevated PPT1 expression correlates with poor clinical prognosis (24). Additionally, PPT1’s expression level has been shown to affect the infiltration level of immune cells within the immune microenvironment of HCC, indicating its broader impact on the tumor immune microenvironment (25). Additionally, studies have shown that targeting PPT1 can enhance the therapeutic effects of anti-PD-1 immunotherapy and sorafenib in HCC: specifically, PPT1 inhibition with ezurpimtrostat promotes lymphocyte infiltration, upregulates MHC-I expression on liver cancer cells, activates cytotoxic CD8+ lymphocytes, enhances anti-PD-1 efficacy, and converts “cold tumors” into “hot tumors” (26); meanwhile, the PPT1-specific inhibitor DC661 inhibits autophagy, induces mitochondrial apoptosis, enhances sorafenib sensitivity, promotes dendritic cell maturation and CD8+ T cell activation, and its combination with sorafenib shows significant anti-tumor effects in animal models (24). These findings further confirm PPT1’s great potential as a novel therapeutic target for HCC. Notably, interferon-alpha (IFN-α) is a common HCC therapeutic agent but has limitations such as low response rates and drug resistance. RORγt agonists can enhance IFN-α efficacy by activating the Akt/Stat3 pathway, promoting Tc17 cell differentiation and CD8+ T cell infiltration via Cxcl10-Cxcr3 interaction (27). Besides, AT7519 (a cyclin-dependent kinase inhibitor) not only suppresses HCC-related genes but also acts on RORγt, regulating CD4+ T cell function and enhancing the effect of anti-PD-1/CTLA-4 therapies (28). These findings further support targeting immune-related molecules to optimize HCC immunotherapy.
In our study, we identified high expression of ZDHHC18 in regulatory T cells (Tregs) within the tumor microenvironment of HCC based on two single-cell datasets (29,30). While the specific role of ZDHHC18 in the HCC microenvironment remains unreported, previous studies have highlighted its involvement in immune regulation. Specifically, ZDHHC18 has been shown to negatively regulate cGAS-mediated innate immunity through palmitoylation (31).
To ensure the robustness and reliability of our findings, we validated our results using an additional single-cell dataset (GSE140228). The similar expression patterns observed in this dataset further support the importance of PPT1 and ZDHHC18 in HCC. The consistency of results across multiple datasets underscores the generalizability of our findings and highlights the potential clinical relevance of these genes as biomarkers and therapeutic targets.
The dysregulation of palmitoylation has been implicated in various aspects of tumorigenesis, including tumor cell proliferation, invasion, and immune evasion (23,32,33). Therefore, targeting palmitoylation enzymes or their downstream effectors could offer a promising strategy for inhibiting tumor growth and enhancing anti-tumor immune responses. As a well-established class of targeted anticancer agents, farnesyltransferase inhibitors (FTIs) have been clinically investigated and applied for tumor therapy by blocking the membrane translocation of N-Ras and suppressing its oncogenic activation (34).
Furthermore, the identification of palmitoylation-related genes as prognostic biomarkers could aid in personalized treatment planning. By stratifying patients based on the expression levels of these genes, clinicians could tailor treatment strategies to individual patient needs. This approach could lead to improved outcomes and reduced side effects compared to traditional one-size-fits-all treatments.
While our study provides valuable insights into the role of palmitoylation-related genes in HCC, it is important to acknowledge its limitations. The sample size of the single-cell datasets used in this study was relatively small, which may limit the generalizability of our findings. Future studies should aim to validate our results using larger and more diverse datasets to ensure robustness.
Additionally, the functional roles of PPT1 and ZDHHC18 in HCC require further investigation. While our analysis suggests their involvement in tumor proliferation and immune regulation, detailed mechanistic studies are needed to elucidate the specific pathways and processes they regulate. This could involve in vitro and in vivo experiments to assess the impact of targeting these genes on tumor growth and immune responses.
Moreover, the clinical utility of palmitoylation-related genes as biomarkers and therapeutic targets needs to be evaluated in larger clinical trials. This would involve assessing their predictive value in a real-world setting and exploring their safety and efficacy in combination with existing treatments. Such studies could provide valuable information for developing personalized treatment strategies for HCC patients.
Conclusions
In conclusion, our study highlights the significant role of palmitoylation-related genes in HCC. The identification of PPT1 and ZDHHC18 as key players in tumor proliferation and immune regulation provides valuable insights into the molecular mechanisms underlying HCC pathogenesis, and the validation of these findings using multiple datasets underscores their robustness and potential clinical relevance. Future research should focus on validating these findings in larger clinical trials and exploring the functional roles of these genes in detail, so as to advance our understanding of HCC and pave the way for more effective and personalized treatment options.
Acknowledgments
None.
Footnote
Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2026-0745/rc
Data Sharing Statement: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2026-0745/dss
Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2026-0745/prf
Funding: None.
Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2026-0745/coif). The authors have no conflicts of interest to declare.
Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. This study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.
Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.
References
- Sung H, Ferlay J, Siegel RL, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin 2021;71:209-49. [Crossref] [PubMed]
- Naidoo S, Daniels A, Habib S, et al. Poly-L-Lysine-Lactobionic Acid-Capped Selenium Nanoparticles for Liver-Targeted Gene Delivery. Int J Mol Sci 2022;23:1492. [Crossref] [PubMed]
- Zi Y, Gao J, Wang C, et al. Pantothenate Kinase 1 Inhibits the Progression of Hepatocellular Carcinoma by Negatively Regulating Wnt/β-catenin Signaling. Int J Biol Sci 2022;18:1539-54. [Crossref] [PubMed]
- Suresh D, Srinivas AN, Kumar DP. Etiology of Hepatocellular Carcinoma: Special Focus on Fatty Liver Disease. Front Oncol 2020;10:601710. [Crossref] [PubMed]
- Singh SP, Madke T, Chand P. Global Epidemiology of Hepatocellular Carcinoma. J Clin Exp Hepatol 2025;15:102446. [Crossref] [PubMed]
- Wang J, Zhou Y, Gu C, et al. LncRNA SAMD12-AS1 Suppresses Proliferation and Migration of Hepatocellular Carcinoma via p53 Signaling Pathway. J Oncol 2022;2022:9096365. [Crossref] [PubMed]
- Sharafi F, Hasani SA, Alesaeidi S, et al. A comprehensive review about the utilization of immune checkpoint inhibitors and combination therapy in hepatocellular carcinoma: an updated review. Cancer Cell Int 2022;22:269. [Crossref] [PubMed]
- Zaballa ME, van der Goot FG. The molecular era of protein S-acylation: spotlight on structure, mechanisms, and dynamics. Crit Rev Biochem Mol Biol 2018;53:420-51. [Crossref] [PubMed]
- Zhao C, Yu H, Fan X, et al. GSK3β palmitoylation mediated by ZDHHC4 promotes tumorigenicity of glioblastoma stem cells in temozolomide-resistant glioblastoma through the EZH2-STAT3 axis. Oncogenesis 2022;11:28. [Crossref] [PubMed]
- Liu Z, Xiao M, Mo Y, et al. Emerging roles of protein palmitoylation and its modifying enzymes in cancer cell signal transduction and cancer therapy. Int J Biol Sci 2022;18:3447-57. [Crossref] [PubMed]
- Tang F, Yang C, Li FP, et al. Palmitoyl transferases act as potential regulators of tumor-infiltrating immune cells and glioma progression. Mol Ther Nucleic Acids 2022;28:716-31. [Crossref] [PubMed]
- Lin Z, Lv Z, Liu X, et al. Palmitoyl transferases act as novel drug targets for pancreatic cancer. J Transl Med 2023;21:249. Retracted Publication.
- Chen Y, Li Y, Wu L. Protein S-palmitoylation modification: implications in tumor and tumor immune microenvironment. Front Immunol 2024;15:1337478. [Crossref] [PubMed]
- Li M, Zhang L, Chen CW. Diverse Roles of Protein Palmitoylation in Cancer Progression, Immunity, Stemness, and Beyond. Cells 2023;12:2209. [Crossref] [PubMed]
- Chamberlain LH, Shipston MJ. The physiology of protein S-acylation. Physiol Rev 2015;95:341-76. [Crossref] [PubMed]
- Jiang P, Gu S, Pan D, et al. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med 2018;24:1550-8. [Crossref] [PubMed]
- Linder ME, Deschenes RJ. Palmitoylation: policing protein stability and traffic. Nat Rev Mol Cell Biol 2007;8:74-84. [Crossref] [PubMed]
- Huang X, Wang M, Zhang D, et al. ZDHHC20 Activates AKT Signaling Pathway to Promote Cell Proliferation in Hepatocellular Carcinoma. J Hepatocell Carcinoma 2024;11:1763-75. [Crossref] [PubMed]
- Wang L, Wang J, Ma X, et al. USP35 promotes HCC development by stabilizing ABHD17C and activating the PI3K/AKT signaling pathway. Cell Death Discov 2023;9:421. [Crossref] [PubMed]
- Cao Y, Qiu T, Kathayat RS, et al. ABHD10 is an S-depalmitoylase affecting redox homeostasis through peroxiredoxin-5. Nat Chem Biol 2019;15:1232-40. [Crossref] [PubMed]
- Bu L, Zhang Z, Chen J, et al. High-fat diet promotes liver tumorigenesis via palmitoylation and activation of AKT. Gut 2024;73:1156-68. [Crossref] [PubMed]
- Kang J, Wu J, Liu Q, et al. FASN regulates STING palmitoylation via malonyl-CoA in macrophages to alleviate sepsis-induced liver injury. Biochim Biophys Acta Mol Basis Dis 2024;1870:167299. [Crossref] [PubMed]
- Zhang Z, Ren C, Xiao R, et al. Palmitoylation of TIM-3 promotes immune exhaustion and restrains antitumor immunity. Sci Immunol 2024;9:eadp7302. [Crossref] [PubMed]
- Xu J, Su Z, Cheng X, et al. High PPT1 expression predicts poor clinical outcome and PPT1 inhibitor DC661 enhances sorafenib sensitivity in hepatocellular carcinoma. Cancer Cell Int 2022;22:115. [Crossref] [PubMed]
- Weng J, Liu S, Zhou Q, et al. Intratumoral PPT1-positive macrophages determine immunosuppressive contexture and immunotherapy response in hepatocellular carcinoma. J Immunother Cancer 2023;11:e006655. [Crossref] [PubMed]
- Bestion E, Rachid M, Tijeras-Raballand A, et al. Ezurpimtrostat, A Palmitoyl-Protein Thioesterase-1 Inhibitor, Combined with PD-1 Inhibition Provides CD8(+) Lymphocyte Repopulation in Hepatocellular Carcinoma. Target Oncol 2024;19:95-106. [Crossref] [PubMed]
- Wang H, Zhang G, Zhao C, et al. RORγt agonist LYC-55716 potentiates IFN-α's efficacy in hepatocellular carcinoma through enhancing cytotoxicity of Tc17 cells and infiltration of CD8(+) T cells. Biochem Pharmacol 2025;238:116963. [Crossref] [PubMed]
- Karaś K, Pastwińska J, Sałkowska A, et al. The cyclin-dependent kinase inhibitor AT7519 is a human RORγt agonist. Immunol Cell Biol 2025;103:317-27. [Crossref] [PubMed]
- Zhang Q, He Y, Luo N, et al. Landscape and Dynamics of Single Immune Cells in Hepatocellular Carcinoma. Cell 2019;179:829-845.e20. [Crossref] [PubMed]
- Zheng C, Zheng L, Yoo JK, et al. Landscape of Infiltrating T Cells in Liver Cancer Revealed by Single-Cell Sequencing. Cell 2017;169:1342-1356.e16. [Crossref] [PubMed]
- Shi C, Yang X, Liu Y, et al. ZDHHC18 negatively regulates cGAS-mediated innate immunity through palmitoylation. EMBO J 2022;41:e109272. [Crossref] [PubMed]
- Zhou B, Wang Y, Zhang L, et al. The palmitoylation of AEG-1 dynamically modulates the progression of hepatocellular carcinoma. Theranostics 2022;12:6898-914. [Crossref] [PubMed]
- Chen QT, Zhang ZY, Huang QL, et al. HK1 from hepatic stellate cell-derived extracellular vesicles promotes progression of hepatocellular carcinoma. Nat Metab 2022;4:1306-21. [Crossref] [PubMed]
- Novotny CJ, Hamilton GL, McCormick F, et al. Farnesyltransferase-Mediated Delivery of a Covalent Inhibitor Overcomes Alternative Prenylation to Mislocalize K-Ras. ACS Chem Biol 2017;12:1956-62. [Crossref] [PubMed]

