Palmitoylation landscapes across human cancers reveal a role of palmitoylation in tumorigenesis
Journal of Translational Medicine volume 21, Article number: 826 (2023)
Protein palmitoylation, which is catalyzed by palmitoyl-transferase and de-palmitoyl-transferase, plays a crucial role in various biological processes. However, the landscape and dynamics of protein palmitoylation in human cancers are not well understood.
We utilized 23 palmitoyl-acyltransferases and seven de-palmitoyl-acyltransferases as palmitoylation-related genes for protein palmitoylation analysis. Multiple publicly available datasets were employed to conduct pan-cancer analysis, examining the transcriptome, genomic alterations, clinical outcomes, and correlation with c-Myc (Myc) for palmitoylation-related genes. Real-time quantitative PCR and immunoblotting were performed to assess the expression of palmitoylation-related genes and global protein palmitoylation levels in cancer cells treated with Myc depletion or small molecule inhibitors. Protein docking and drug sensitivity analyses were employed to predict small molecules that target palmitoylation-related genes.
We identified associations between palmitoylation and cancer subtype, stage, and patient survival. We discovered that abnormal DNA methylation and oncogenic Myc-driven transcriptional regulation synergistically contribute to the dysregulation of palmitoylation-related genes. This dysregulation of palmitoylation was closely correlated with immune infiltration in the tumor microenvironment and the response to immunotherapy. Importantly, dysregulated palmitoylation was found to modulate canonical cancer-related pathways, thus influencing tumorigenesis. To support our findings, we performed a proof-of-concept experiment showing that depletion of Myc led to reduced expression of most palmitoylation-related genes, resulting in decreased global protein palmitoylation levels. Through mass spectrometry and enrichment analyses, we also identified palmitoyl-acyltransferases ZDHHC7 and ZDHHC23 as significant contributors to mTOR signaling, DNA repair, and immune pathways, highlighting their potential roles in tumorigenesis. Additionally, our study explored the potential of three small molecular (BI-2531, etoposide, and piperlongumine) to modulate palmitoylation by targeting the expression or activity of palmitoylation-related genes or enzymes.
Overall, our findings underscore the critical role of dysregulated palmitoylation in tumorigenesis and the response to immunotherapy, mediated through classical cancer-related pathways and immune cell infiltration. Additionally, we propose that the aforementioned three small molecule hold promise as potential therapeutics for modulating palmitoylation, thereby offering novel avenues for cancer therapy.
Materials and methods
Cell culture and transfection
Human HCT116, SW48, HT29 and HEK293T cells were derived from the American Type Culture Collection (ATCC). HEK293T cells were cultured in DMEM medium (Gibco, NY, USA) supplemented with 10% fetal bovine serum (Gibco, NY, USA) and 1% penicillin–streptomycin (Gibco, CA, USA) at 37 °C in a 5% CO2 incubator. For transfection, after growing to 70% confluence, cells were transfected using Lipofectamine 3000 (Invitrogen, Carlsbad, CA) according to the manufacturer's instructions when they reached 70% confluence.
Reagents and plasmids
Anti-Flag agarose beads (23101) were purchased from Selleck (Houston, USA). Anti-HA agarose beads (KTSM1305) were obtained from Shenzhen KangTi Life Technology (Shenzhen, China). The human ZDHHC7 and ZDHHC23 plasmids were generously provided by Dr. Guojun Zhao from The First Affiliated Hospital of Zhengzhou University, China. Human Myc shRNAs were purchased from Miaoling Biology (Whuhan, China). The shRNA sequences used in the study can be found in Additional file 1: Table S1. Small molecules (BI2536 (T6173), etoposide (T0132), piperlongumine (T6947)) were purchased from Topscience Biotech (Shanghai, China). Hydroxylamine (HAM, 467804), N-Ethylmaleimide (NEM, A600450-0005) and BMCC-biotin ((1-Biotinamido)-4-[4'-(mal-eimidomethyl)cyclohexanecarboxamido]hexane, C100222-0050) were purchased from Sangon Biotech (Shanghai, China).
Data and resources
RNA-Seq, somatic mutation, DNA methylation, copy number variant (CNV), and clinical data were from The Cancer Genome Atlas (TCGA) or Genotype-Tissue Expression (GTEx). Immunogenomic analysis was performed by the Immune Cell Abundance Identifier (ImmuCellAI) [1, 2] algorithm with 24 immune cells. Information on small molecule drugs is from the Cancer Therapeutics Response Portal (CTRP)  and Genomics of Drug Sensitivity in Cancer (GDSC) website. These data were integrated into the Gene Set Cancer Analysis (GSCA) for pan-cancer analysis in this study . All TCGA cancer types included in the study are listed in Additional file 1: Table S2.
Gene Set Variation Analysis (GSVA)
GSVA analysis was performed on nine cancer types (HNSC, LUSC, COAD, STAD, LUAD, GBM, BRCA, KIRC, BLCA) with available subtype information. The analysis was conducted using the R package GSVA, following previously described methods . In this analysis, the GSVA score, which represents the variation of a gene set activity in a specific cancer sample population, was calculated. The GSVA score reflects the integrated level of gene set expression and is positively correlated with the expression of the gene set. By comparing the GSVA scores of tumor tissues with those of adjacent samples, we were able to determine if the overall expression of a particular gene set was higher in the tumor tissues.
Differential expression analysis
Differential expression analysis was conducted for 14 TCGA cancer types (THCA, KIRP, BLCA, LIHC, HNSC, BRCA, LUAD, PRAD, ESCA, KICH, LUSC, KIRC, STAD, COAD) with more than 10 paired tumor and normal samples. This analysis was based on normalized and batch corrected RSEM mRNA expression data. The fold change was calculated by taking the mean expression in tumor samples divided by the mean expression in normal samples. The statistical significance was determined using a t-test, and the resulting P values were adjusted for multiple testing using the false discovery rate (FDR) method.
Differential methylation analysis
Differential methylation analysis was performed using Illumina Human Methylation 450k level data obtained from the TCGA database. This analysis was conducted for 14 different types of cancer, including THCA, KIRP, BLCA, LIHC, HNSC, BRCA, LUAD, PRAD, ESCA, KICH, LUSC, KIRC, STAD, and COAD. Only cancer types with more than 10 paired tumor and normal samples were included in the analysis. Since a region of a gene can have multiple methylation sites, there are multiple tags that store methylation levels at each site. Before conducting the differential methylation analysis, we performed correlation analysis to filter out the sites that were most negatively correlated with gene expression. This was done to ensure that only relevant sites were included in the analysis. To calculate the differential methylation, we used the ratio of the mean methylation levels in tumor samples to the mean methylation levels in normal samples. The statistical significance was determined using a t-test, and the P values were further adjusted using the FDR method.
CNV analysis was performed using data from 11,495 samples obtained from the TCGA database. The data was processed using GISTIC2.0, a program that aims to identify regions of amplification or deletion that are significantly altered in patients. The CNV data is represented by different values: -2 or Deep Deletion indicates a deep loss, possibly a homozygous deletion; -1 or Shallow Deletion indicates a shallow loss, possibly a heterozygous deletion; 0 represents diploid status; 1 or Gain indicates a low-level gain, possibly a heterozygous amplification; 2 or Amplification indicates a high-level amplification, possibly a homozygous amplification.
To analyze the mutations, we obtained mutation data from the TCGA database consisting of 10,234 samples. We focused on seven types of mutations, namely Missense_Mutation, Nonsense_Mutation, Frame_Shift_Ins, Splice_Site, Frame_Shift_Del, In_Frame_Del, and In_Frame_Ins, which are considered deleterious mutations. Additionally, we included non-deleterious mutations such as Silent, Intron, IGR, 3'UTR, 5'UTR, 3'Flank, and 5'Flank in our analysis.
Survival analysis was conducted using clinical data from 33 cancer types from TCGA tumor samples. The survival types analyzed included overall survival (OS), progression-free survival (PFS), disease-specific survival (DSS), and disease-free interval (DFI) . Some uncensored data was excluded, and samples of tumors at risk of competing cancer deaths were filtered out. Regarding data availability, OS data is available for all cancer types, while PFS data is not available for LAML. DSS data is also not available for LAML. and the DFI data are not available for LAML, SKCM, THYM, and UVM. Additionally, samples of tumors at risk of competing cancer deaths were excluded from the analysis.
For expression/methylation and survival analysis, mRNA/methylation data were combined with clinical survival data based on sample barcode. Tumor samples were divided into a high expression/methylation group and a low expression/methylation group based on the median mRNA/methylation value. Cox Proportional-Hazards models and log-rank tests were performed for each gene or gene set in each cancer type.
For mutation and survival analysis, mutation data were combined with clinical survival data based on sample barcode. Tumor samples were classified as the Mutant group when deleterious mutations occurred in specific genes. Cox Proportional-Hazards models and log-rank tests were performed to compare the survival difference between the Wild-Type (WT) and Mutant groups.
Immune infiltration analysis
To assess immune infiltration and mRNA expression analysis, the Spearman correlation was employed to determine the relationship between gene mRNA expression or gene sets and the infiltration of immune cells. The abundance of 24 immune cell types was estimated using gene set signatures, and the infiltrates of these immune cells were evaluated through ImmuCellAI.
Real-time quantitative PCR (RT-qPCR)
The RT-qPCR assays were performed according to the methods described in a previous study . Total RNA was extracted from cells and subjected to reverse transcription. Then, qPCR was conducted using the SYBR Green Supermix (Bio-Rad, Hercules, CA) and standard protocols. The primer sequences used in the study can be found in Additional file 1: Table S1. GAPDH was utilized as the internal control for normalization purposes.
Global palmitoylation detection
To detect global protein palmitoylation, we employed the acyl-biotin exchange (ABE) method . Briefly, cells were lysed using a lysis buffer, and 80 μL of the resulting supernatant was saved as the input sample. The remaining supernatant was utilized for the ABE experiments. The detection of global protein palmitoylation was performed using HRP-conjugated streptavidin (Sangon Biotech, diluted 1:200 in 0.5% BSA). All antibodies used in the study are listed in Additional file 1: Table S3.
Pathway activity analysis
To analyze pathway activity, we utilized RPPA data from the TCPA database. We calculated pathway activity scores for 10 cancer-associated pathways including TSC/mTOR, RTK, RAS/MAPK, PI3K/AKT, Hormone ER, Hormone AR, EMT, DNA Damage Response, Cell Cycle, and Apoptosis pathways. These scores were calculated using 7876 samples from the TCGA database across 32 different cancer types. The RPPA data from RBN were first centered on the median and then normalized by the standard deviation for each component in all samples to obtain relative protein levels. The pathway score for each sample is calculated as the sum of the relative protein levels of all positive regulatory components minus the relative protein levels of negative regulatory components within a given pathway .
To analyze gene set and pathway activity, we initially performed GSVA analysis on the gene set. GSVA scores were calculated using the R package GSVA and represent the variation of gene set activity in a specific population of cancer samples in an unsupervised manner. The GSVA score reflects the overall expression level of the gene set and is positively correlated with gene set expression.
For gene expression and pathway activity analysis, we categorized the samples into high and low groups based on the median gene expression level. The differences in pathway activity score (PAS) between these two groups were determined using a student t-test, and the resulting P values were adjusted using the FDR method. A significance level of FDR ≤0.05 was considered statistically significant. If the PAS (Gene A High expression) > PAS (Gene A Low expression), we inferred that Gene A may have an activating effect on the pathway. Conversely, if the PAS (Gene A High expression) < PAS (Gene A Low expression), Gene A may have an inhibitory effect on the pathway . The global percentage of cancers in which a gene has an effect on the pathway among 32 cancer types is calculated as the percentage (number of cancer types in which the gene is activated or inhibited/32 * 100%).
Immunoprecipitation and mass spectrometry analysis
Immunoprecipitation and mass spectrometry analysis were conducted following previously described methods . HEK293T cells transfected with the ZDHHC23/7 plasmid were lysed in 1 mL of lysis buffer. Anti-Flag/HA agarose beads were washed with 1 mL of lysis buffer three times, and then 0.95 mL of cell lysate was added to the respective groups and incubated overnight at 4 °C. The following day, the agarose beads were centrifuged and washed three times with a wash buffer. The beads were then mixed with a 2 × SDS sample buffer. Lysate samples were boiled for 10 min, followed by SDS-PAGE and Coomassie Blue staining. The bands were excised, subjected to in-gel trypsin digestion, and dried for mass spectrometry analysis, as previously described .
CTRP and GDSC drug sensitivity analysis
To analyze drug sensitivity, we obtained the IC50 values of 481 small molecules in 1001 cell lines from the CTRP database. We also downloaded the corresponding mRNA gene expression data. To assess the correlation between mRNA expression and drug sensitivity, we performed Pearson correlation analysis by combining the mRNA expression and drug sensitivity data. The resulting P values were adjusted for FDR. To validate our findings, we conducted additional drug sensitivity analysis using the GDSC database. The results were presented using the Assistant for Clinical Bioinformatics website (https://www.aclbi.com/static/index.html#/drug_allergy).
Protein structure acquisition
The crystal structures of LYPLA1 (PDB ID: 5SYM), LYPLA2 (PDB ID: 5SYN), PPT1 (PDB ID: 3GRO), and PPT2 (PDB ID: 1PJA) were obtained from the Protein Data Bank (PDB). For proteins without available crystal structures, we utilized the AlphaFold2 tool  to predict their structures. Additionally, the structures of three small molecules (BI2536, etoposide, piperlongumine) were retrieved from PubChem and underwent molecular energy minimization using chem3D 20.0 software. The Vina force field was applied to the proteins, while the small molecules were optimized with polar hydrogen atoms using the Autodock tool.
Determination of protein active pockets and protein–ligand docking
To determine the active pockets in the proteins, we used CavityPlus 2022 . The identified active pockets were then refined based on the enzymatic regions of the proteins. Protein–ligand docking studies were conducted using Autodock Vina 2021 , which allowed us to calculate the binding energy between the protein and small molecules. The resulting conformations from the docking were visualized using Pymol, enabling us to analyze the intermolecular interactions involved in the protein–ligand interactions.
The P value was calculated using the t-test and adjusted for multiple comparisons using the FDR method. Details of the specific statistical methods used can be found in the corresponding methods section or figure legends. For general analysis and visualization, we utilized the GSCA . A significance level of FDR ≤ 0.05 or P < 0.05 was considered statistically significant.
Protein palmitoylation, a reversible post-translational modification, plays a crucial role in altering the localization, stability, and function of proteins [14,15,16,17,18]. The zinc finger DHHC-type containing (ZDHHC) family, consisting of 23 distinct proteins in mammals, catalyzes this reversible modification [14, 15]. On the other hand, de-palmitoylation is typically mediated by three distinct families that include seven genes [16,17,18]. Previous research from our group has demonstrated the promotion of colorectal cancer progression by palmitoylated β-catenin . Moreover, emerging evidence suggests that palmitoyl-acyltransferases and de-palmitoyl-lacyltransferases are associated with various aspects of carcinogenesis, cancer cell growth, survival, and treatment resistance [18, 19]. Consequently, enzymes involved in modulating protein palmitoylation hold great promise as potential targets for cancer therapy. Despite this potential, the landscape of palmitoylation in human cancers has yet to be fully described.
The transcription factor c-Myc (Myc) is a well-known oncogene frequently activated in various human cancers. Numerous studies have demonstrated its role in tumor initiation and progression through the regulation of multiple biological processes [20,21,22]. Recently, our own research has revealed that Myc can contribute to tumorigenesis by modulating GTP metabolic reprogramming . These oncogenic phenotypes of Myc are primarily attributed to its ability to activate or repress gene transcription. Given the substantial influence of Myc on gene expression, it is reasonable to question whether Myc also affects global palmitoylation profiles by regulating the transcription of genes involved in palmitoylation. Investigating whether Myc influences the expression of these palmitoylation-related genes (23 palmitoyl-acyltransferases and seven de-palmitoyl-acyltransferases) in human cancers may provide valuable insights.
Tumor-infiltrating immune cells play a vital role in the tumor microenvironment and have been implicated in tumor progression and response to immunotherapy . Numerous studies have demonstrated that the extent of immune cell infiltration within tumors is associated with clinical outcomes and response to immunotherapy in various cancer types [24,25,26,27,28,29,30,31]. Furthermore, recent investigations have highlighted the significance of palmitoylation in modulating CD8+ T cell infiltration levels in tumor tissues, with palmitoylation of PD-L1 and IFNGR1 emerging as key factors [19, 32]. These findings underscore the critical role of palmitoylation in immunotherapy. To gain a more comprehensive understanding, it is crucial to establish a detailed landscape of palmitoylation and immune cell infiltration in human cancers.
Tumorigenesis is frequently associated with the dysregulation of signaling pathways, including TSC/mTOR, RTK, RAS/MAPK, PI3K/AKT, Hormone ER, Hormone AR, EMT, DNA Damage Response, Cell Cycle, and Apoptosis pathways. These pathways play crucial roles in cancer development and progression. Recent studies have demonstrated that palmitoylation actively regulates several of these cancer-related pathways, impacting tumorigenesis [33,34,35]. To enhance our understanding of tumorigenesis, it is essential to establish a comprehensive landscape that delineates the involvement of palmitoylation in cancer-related pathways across different types of human cancers.
In this study, we present a comprehensive analysis of the palmitoylation landscape in human cancers. We highlight the relationship between Myc and palmitoylation and demonstrate that dysregulated palmitoylation plays a pivotal role in driving tumorigenesis by modulating classical cancer-related pathways and immune infiltration. Moreover, we identify three small molecular compounds, namely BI-2531, etoposide, and piperlongumine, that have the potential to modulate palmitoylation and could serve as promising candidates for cancer therapy.
Dysregulated palmitoylation in human cancers
To investigate the landscape of palmitoylation in human cancers, we analyzed the expression of palmitoylation-related genes using data from the GTEx and TCGA databases, respectively (Additional file 1: Table S4). Initially, we evaluated the expression patterns of these genes across diverse human tissues using the GTEx database and found that most genes, except for a few exceptions such as ZDHHC18, ZDHHC19, and PPT2, exhibited relatively consistent expression levels (Additional file 2: Fig. S1). This suggests that protein palmitoylation may play a comparable role in multiple human tissues.
Further delving into the profiles of protein palmitoylation in human cancers, we performed GSVA in 14 cancer types using more than 10 pairs of tumor-normal samples from the TCGA database. The results revealed dysregulation of protein palmitoylation in these human cancers (Fig. 1A), implying its potential critical involvement in tumorigenesis. Additionally, the expression of palmitoylation-related genes demonstrated associations with cancer subtypes in most cancer types with available subtype information (Fig. 1B). Similarly, the expression patterns of these genes across different stages of cancer were linked to tumor progression in multiple cancer types (Fig. 1C).
Subsequently, we examined the association between the expression of palmitoylation-related genes and patient survival outcomes , including OS, PFS, DSS, and DFI, across the 33 cancer types. Our analysis revealed that dysregulation of palmitoylation profiles correlated with patient survival in several human cancers (Fig. 1D). These findings collectively highlight the dysregulation of protein palmitoylation in human cancers and its potential association with tumor progression.
Dys-expression of palmitoylation-related genes across human cancers
In order to provide a comprehensive description of the expression landscape of palmitoylation in human cancers, we conducted a detailed analysis of the differential expression of each palmitoylation-related gene in normal and tumor tissues across the aforementioned 14 human cancer types. The results demonstrated that a majority of the palmitoylation-related genes exhibited differential expression patterns in various cancer types (Fig. 2A). Moreover, we observed a significant association between all palmitoylation-related genes and cancer subtypes in at least three human cancer types (Additional file 2: Fig. S2A). Notably, two specific genes, ZDHHC9 (a palmitoyl-acyltransferase) and ABHD17C (a de-palmitoyl-acyltransferase), displayed significant dysregulation across more than 10 human cancer types and exhibited correlation with cancer subtypes in seven human cancer types (Fig. 2A, Additional file 2: Fig. S2A). These findings suggest a potential crucial role for ZDHHC9 and ABHD17C in tumorigenesis.
We also assessed the relation between the expression of each individual palmitoylation-related gene and the prognosis of cancer patients, including OS, PFS, DSS, and DFI, across the 33 different human cancer types. The results indicated that all palmitoylation-related genes were associated with at least one survival prognosis in at least one cancer type (Fig. 2B, Additional file 2: Fig. S2B). Collectively, our analysis demonstrates the aberrant expression of palmitoylation-related genes in a variety of human cancers, as well as their association with cancer patient prognosis.
Dys-expression of palmitoylation-related genes is mainly driven by DNA methylation alterations
The analysis of gene mutations is crucial for understanding tumorigenesis [36, 37]. Therefore, we collected mutation data from 10,234 samples across 33 cancer types sourced from the TCGA database. We specifically focused on seven types of deleterious mutations (Missense_Mutation, Nonsense_Mutation, Frame_Shift_Ins, Splice_Site, Frame_Shift_Del, In_Frame_Del, and In_Frame_Ins) to investigate their impact. Intriguingly, pan-cancer analysis revealed that palmitoylation-related genes had generally low mutation rates across most cancer types, with the exception of UCEC (Additional file 2: Fig. S3A). As expected, our survival analysis demonstrated that these mutations in palmitoylation-related genes influenced the prognosis of UCEC patients (Additional file 2: Fig. S3B), emphasizing their role in dysregulated palmitoylation, particularly in UCEC.
Considering the crucial role of CNV in gene expression regulation [38, 39], we also examined the relationship between the expression of palmitoylation-related genes and CNV. Interestingly, in half of the cancer types, none of the palmitoylation-related genes exhibited gains or losses at the CNV level (Additional file 2: Fig. S4A), although the correlations between the expression of palmitoylation-related genes and CNV of these genes were found (Additional file 2: Fig. S4B).
Moreover, given the ubiquitous presence of abnormal DNA methylation in cancer . we analyzed whether dysregulation of palmitoylation-related genes could be attributed to DNA methylation alterations. We found that DNA methylation of all examined palmitoylation-related genes was abnormal in at least half of the analyzed cancer types (Fig. 3A). Moreover, the expression of these genes showed a negative correlation with their methylation levels (Fig. 3B). Importantly, methylation levels of palmitoylation-related genes significantly impacted the survival outcomes of patients in multiple cancer types, particularly LGG (Fig. 3C). In summary, our findings suggest that dysregulation of palmitoylation-related genes, driven by DNA methylation changes, is a notable aspect of tumorigenesis.
Myc drives dysregulated palmitoylation
Transcription factors play a crucial role in regulating gene expression by either activating or suppressing transcription processes [7, 23]. To investigate whether any specific transcription factor is involved in regulating the expression of palmitoylation-related genes, we conducted analysis using the ChEA3 website , which utilizes gene set libraries derived from published Chromatin Immunoprecipitation (ChIP)-seq data from various sources. Interestingly, our analysis revealed a significant positive association between Myc, a transcription factor frequently activated in cancer, and the majority of palmitoylation-related genes (Fig. 4A, B). This suggests that Myc may act as a regulator of dysregulated palmitoylation in human cancers.
To further investigate the regulatory role of Myc in the expression of palmitoylation-related genes, we performed depletion experiments of Myc in colorectal cancer cells, specifically HCT116, SW48 and HT29 cell lines. Consistent with our expectations, heatmap analysis revealed that depletion of Myc led to a significant decrease in the expression of most palmitoylation-related genes (Fig. 4C, Additional file 2: Fig. S5A). As part of our validation process, we discovered that Myc directly binds to the promoter of ZDHHC9 gene (Fig. 4D, Additional file 2: Fig. S5B, C). Additionally, we observed that depletion of Myc resulted in a decrease in global protein palmitoylation levels (Fig. 4D, Additional file 2: Fig. S5D). Collectively, these findings provide evidence that Myc drives the dysregulated profile of palmitoylation in human cancers.
Dysregulated palmitoylation is associated with immune infiltration and immunotherapy response in human cancers
Immune infiltration has emerged as a clinically relevant biomarker, influencing prognosis and therapy response in cancer patients . To explore the relationship between immune infiltration and palmitoylation in human cancers, we performed GSVA analysis focusing on palmitoylation and de-palmitoylation. Through Immune Cell Abundance Identifier , we evaluated the infiltration of 24 immune cell types. Our analysis revealed a significant association between dysregulated palmitoylation and different immune cell infiltrates across 33 human cancers, particularly in PPAD (Additional file 2: Fig. S6A, B). To verify these findings, we performed single-cell analysis, which demonstrated the expression of palmitoylation-related genes in immune cells to some extent in PAAD (Additional file 2: Fig. S6C–F).
In addition, we investigated the relationship between dysregulated palmitoylation and 47 immune regulators in various cancer types. Remarkably, we observed a strong correlation between dysregulated palmitoylation and most immune regulators (Fig. 5A). Furthermore, we explored the predictive value of palmitoylation in assessing the effectiveness of immune checkpoint inhibitors (ICIs) by examining the correlation between dysregulated palmitoylation and tumor mutational burden (TMB) as well as microsatellite instability (MSI). Our results indicated a significant correlation between dysregulated palmitoylation and TMB/MSI across multiple cancer types (Fig. 5B, C). Interestingly, ZDHHC9 and ABHD17C demonstrate similar effects on TMB/MSI (Additional file 2: Fig. S7A, B). These findings suggest that analyzing palmitoylation levels could enable the identification of predictive markers for ICI effectiveness in specific cancer types.
Given the significant advancements in cancer immunotherapy, particularly the use of anti-PD-L1 antibodies, we further investigated the predictive role of the dysregulated palmitoylation in ICI-treated cancer cohorts, specifically focusing on patients with urological tumors. Remarkably, patients with low level of the palmitoylation/depalmitoylation exhibited worse survival rates and shorter survival times compared to those with high level (Fig. 5D). These findings demonstrate a correlation between dysregulated palmitoylation, immune infiltration, and immunotherapy response in human cancers.
Palmitoylation regulates classical cancer-related pathways to affect tumorigenesis
To gain further insights into palmitoylation-related signaling pathways, we utilized RPPA data from the TCPA database. We calculated pathway activity scores for 10 cancer-related pathways using RPPA data from 7,876 samples across 32 cancer types in the TCGA database. Interestingly, we observed associations between dysregulated palmitoylation and at least one cancer-related signaling pathway in each cancer type (Fig. 6A). These findings indicate that palmitoylation may influence classical cancer-related pathways and contribute to tumorigenesis.
To investigate the relationship between individual palmitoylation-related genes and these cancer-related pathways, we generated a heatmap displaying the percentage of cancers in which mRNA expression of a palmitoylation-related gene potentially affects pathway activity (Fig. 6B). For instance, palmitoyl-acyltransferase ZDHHC7 showed a negative association with DNA damage and hormone AR pathways in seven and eight cancer types, respectively (Fig. 6B). Similarly, palmitoyl-acyltransferase ZDHHC23 displayed a positive association with DNA damage and mTOR pathways in eight and three cancer types, respectively. These findings highlight the potential role of ZDHHC7 and ZDHHC23 in modulating cancer-related pathways.
To strengthen our analysis, we performed experimental verification using exogenously expressed ZDHHC7 and ZDHHC23 in HEK293T cells. Immunoprecipitation of the epitope-labeled proteins (Flag or HA-tagged) followed by mass spectrometry analysis allowed us to identify ZDHHC7- and ZDHHC23-associated proteins. Subsequently, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis was conducted on these associations (Fig. 7A). Consistent with our previous analysis, we observed significant correlations between ZDHHC7 or ZDHHC23 and pathways related to mTOR, DNA repair, hormone signaling, and immune response (Fig. 7B, C).Taken together, our data suggest that palmitoylation plays a role in the regulation of classical cancer-related pathways, thereby impacting tumorigenesis.
Potential small molecules regulating palmitoylation are identified
To identify potential drugs that can regulate palmitoylation, we gathered IC50 values for 481 small molecules and their corresponding mRNA gene expressions from the CTRP database. By conducting a correlation analysis between the mRNA expression of each palmitoylation-related gene and the drug IC50 values, we identified a range of small molecules associated with these genes (Additional file 2: Fig. S8A). Specifically, 19 small molecules emerged as potential regulators of palmitoylation (Additional file 2: Fig. S8B).
To validate these findings, we examined data from the GDSC and TCGA databases. Among the 19 small molecules, only three (BI-2536, etoposide, and piperlongumine) matched the available data. We performed pan-cancer susceptibility analysis for these three small molecules and found that they significantly inhibited multiple types of cancer (Fig. 8A, Additional file 2: Fig. S8C, D). Notably, treatment with BI-2536 or piperlongumine resulted in decreased expression of some palmitoylation-related genes in HCT116, SW48 and HT29 cells (Fig. 8B, Additional file 2: Fig. S9), suggesting that these small molecules may regulate palmitoylation by modulating the expression of these genes.
To analyze potential targeted small molecules, we focused on the three aforementioned compounds: BI-2536, etoposide, and piperlongumine. Using AlphaFold2, we conducted protein structure predictions and performed molecular docking calculations with Autodock Vina. Our predictions revealed that all ZDHHC enzymes possess conserved 4 alpha-helical bundles, with the enzymatic regions situated within this region (Additional file 1: Table S5). Through analysis of protein–ligand binding energies and intermolecular interactions, we observed that, except for ZDHHC3, all proteins exhibited favorable binding affinity to BI-2536 (binding energy below -6.0 kcal/mol), with most forming strong hydrogen bond interactions with the compound (Fig. 8C, Additional file 1: Table S5). Etoposide and piperlongumine also displayed favorable binding energies with these proteins (Additional file 2: Fig. S10A).
The binding pockets for the small molecules within the ZDHHC family are located in the hydrophobic interior of the proteins (Fig. 8D, Additional file 2: Fig. S10B–E), while the binding pockets for other proteins (PPT1/2, LYPLA1/2, ABHD17A/B/C) are situated on the protein surface (Fig. 8E, Additional file 2: Fig. S10F, G), as opposed to the enzymatic pockets near the 4 alpha-helical bundles of the ZDHHC enzymes (Fig. 8D, Additional file 2: Fig. S10B–E). Despite belonging to different protein families and exhibiting structural variations, these molecules can effectively regulate the expression of palmitoylation-related genes or bind to palmitoylation-related proteins, making them promising drug candidates for targeting these proteins.
In this study, we have provided a comprehensive overview of the landscape of palmitoylation in human cancers, shedding light on its pivotal role in tumorigenesis. We observed that depletion of Myc, a known transcription factor, leads to alterations in the expression profile of palmitoylation-related genes, further supporting the involvement of Myc in the regulation of palmitoylation. Additionally, we found a correlation between dysregulated palmitoylation and immune infiltration as well as classical cancer-related pathways, providing valuable insights into the mechanisms by which palmitoylation impacts tumorigenesis.
While previous studies have associated palmitoylation with tumorigenesis [18, 42, 43], a comprehensive understanding of the overall landscape of palmitoylation across various human cancers has been lacking. To address this, we utilized palmitoylation-related genes as surrogates for global palmitoylation levels in cancer and conducted pan-cancer analysis. Our results revealed dysregulated palmitoylation patterns in human cancers. Importantly, we also discovered that these dysregulations in palmitoylation are associated with alterations in DNA methylation, highlighting the significant role of epigenetics in modulating palmitoylation.
Oncogenic Myc, a prominent transcription factor, has been extensively studied for its role in promoting tumorigenesis through the regulation of diverse biological processes [44,45,46,47]. However, the specific relationship between Myc and palmitoylation has not been previously reported. In our study, we made a novel discovery that depletion of Myc results in a substantial decrease in global palmitoylation levels. This finding strongly suggests that Myc has the capability to regulate the process of palmitoylation, adding a new dimension to our understanding of Myc’s role in tumorigenesis.
Immune infiltration plays a crucial role in tumor development, prognosis, and response to immunotherapy. Recent research has indicated that the palmitoylation of proteins such as PDL1 and IFNGR1 can impact the levels of CD8+ T cell infiltration and influence the effectiveness of immunotherapy in human cancer [19, 32]. In our study, we conducted a comprehensive pan-cancer analysis to explore the landscape of palmitoylation and its relationship with immune infiltration. Through our analysis, we were able to gain a better understanding of the impact of palmitoylation on immune infiltration in various cancer types. Our findings suggest that modifying the palmitoylation status of certain proteins could be a viable strategy in improving the immune response in immunotolerant patients. This knowledge may offer potential therapeutic targets for enhancing the efficacy of immunotherapy interventions. Overall, our study contributes to the growing understanding of the role of palmitoylation in immune infiltration and its implications for cancer treatment.
The overactivation or inactivation of signaling pathways is responsible for the uncontrolled proliferation and survival of tumor cells [6, 35, 48]. Recent studies have shown that palmitoylation plays a role in regulating cancer-related pathways, influencing tumor initiation and progression [35, 48]. Consistently, our research also found a significant association between palmitoylation and classical cancer-related pathways in human cancers. In addition, we constructed a comprehensive landscape of palmitoylation and its relationship with these pathways, which enhances our understanding of the mechanisms underlying tumorigenesis. Of particular importance, our study identified three drugs—BI-2536, etoposide, and piperlongumine—that have the potential to regulate palmitoylation by targeting specific gene expression or enzymes involved in the process. Notably, these drugs have demonstrated anti-tumor potential in previous studies [49,50,51]. It is worth mentioning that BI-2536, a potent inhibitor of PLK1, has been reported to promote tumor growth by stabilizing Myc. Consequently, it is likely that BI-2536 targets palmitoylation by inhibiting the expression of Myc.
Although our findings in this study have important clinical implications, there are still some limitations that need to be addressed. Firstly, we found that palmitoylation is associated with TMB and MSI, but the relationship remains unclear. Secondly, while we have identified potential small molecule drugs that may regulate palmitoylation, specifically BI-2536, etoposide, and piperamine, it is important to note that further research is necessary to validate these findings. Our future research will be dedicated to conducting additional experiments that specifically address these questions and aim to provide more conclusive evidence.
In summary, our data provide a comprehensive overview of palmitoylation in various human cancers, highlight the role of DNA methylation alterations and Myc in palmitoylation regulation, and elucidate the intricate relationship between palmitoylation and cancer-related pathways. Notably, our identification of three drugs with the potential to regulate palmitoylation offers new possibilities for cancer therapy.
Availability of data and materials
All data that support the findings in this study are available from the corresponding author upon reasonable request.
MYC proto-oncogene, bHLH transcription factor
Zinc finger DHHC-type palmitoyltransferase 7
Zinc finger DHHC-type palmitoyltransferase 23
Mechanistic target of rapamycin kinase
Copy number variant
The cancer genome atlas
Immune cell abundance identifier
Cancer therapeutics response portal
Genomics of drug sensitivity
Gene set cancer analysis
Gene set variation analysis
False discovery rate
Real-time quantitative PCR
Palmitoyl-protein thioesterase 1
Palmitoyl-protein thioesterase 2
Protein data bank
Zinc finger DHHC-type containing
Catenin beta 1
Interferon gamma receptor 1
Zinc finger DHHC-type palmitoyltransferase 18
Zinc finger DHHC-type palmitoyltransferase 19
Zinc finger DHHC-type palmitoyltransferase 9
Abhydrolase domain containing 17C, depalmitoylase
Abhydrolase domain containing 17A, depalmitoylase
Abhydrolase domain containing 17B, depalmitoylase
Immune checkpoint inhibitors
Tumor mutational burden
Miao YR, Zhang Q, Lei Q, Luo M, Xie GY, Wang H, Guo AY. ImmuCellAI: a unique method for comprehensive T-cell subsets abundance prediction and its application in cancer immunotherapy. Adv Sci (Weinh). 2020;7:1902880.
Miao YR, Xia M, Luo M, Luo T, Yang M, Guo AY. ImmuCellAI-mouse: a tool for comprehensive prediction of mouse immune cell abundance and immune microenvironment depiction. Bioinformatics. 2022;38:785–91.
Rees MG, Seashore-Ludlow B, Cheah JH, Adams DJ, Price EV, Gill S, Javaid S, Coletti ME, Jones VL, Bodycombe NE, et al. Correlating chemical sensitivity and basal gene expression reveals mechanism of action. Nat Chem Biol. 2016;12:109–16.
Liu CJ, Hu FF, Xia MX, Han L, Zhang Q, Guo AY. GSCALite: a web server for gene set cancer analysis. Bioinformatics. 2018;34:3771–2.
Liu J, Lichtenberg T, Hoadley KA, Poisson LM, Lazar AJ, Cherniack AD, Kovatich AJ, Benz CC, Levine DA, Lee AV, et al. An integrated TCGA pan-cancer clinical data resource to drive high-quality survival outcome analytics. Cell. 2018;173:400-16.e411.
Zhang Q, Li X, Cui K, Liu C, Wu M, Prochownik EV, Li Y. The MAP3K13-TRIM25-FBXW7alpha axis affects c-Myc protein stability and tumor development. Cell Death Differ. 2020;27:420–33.
Zhang Q, Yang X, Wu J, Ye S, Gong J, Cheng WM, Luo Z, Yu J, Liu Y, Zeng W, et al. Reprogramming of palmitic acid induced by dephosphorylation of ACOX1 promotes beta-catenin palmitoylation to drive colorectal cancer progression. Cell Discov. 2023;9:26.
Akbani R, Ng PK, Werner HM, Shahmoradgoli M, Zhang F, Ju Z, Liu W, Yang JY, Yoshihara K, Li J, et al. A pan-cancer proteomic perspective on The Cancer Genome Atlas. Nat Commun. 2014;5:3887.
Ye Y, Xiang Y, Ozguc FM, Kim Y, Liu CJ, Park PK, Hu Q, Diao L, Lou Y, Lin C, et al. The Genomic landscape and pharmacogenomic interactions of clock genes in cancer chronotherapy. Cell Syst. 2018;6:314-28.e312.
Shevchenko A, Wilm M, Vorm O, Mann M. Mass spectrometric sequencing of proteins silver-stained polyacrylamide gels. Anal Chem. 1996;68:850–8.
Jumper J, Evans R, Pritzel A, Green T, Figurnov M, Ronneberger O, Tunyasuvunakool K, Bates R, Zidek A, Potapenko A, et al. Highly accurate protein structure prediction with AlphaFold. Nature. 2021;596:583–9.
Wang S, Xie J, Pei J, Lai L. CavityPlus 2022 update: an integrated platform for comprehensive protein cavity detection and property analyses with user-friendly tools and cavity databases. J Mol Biol. 2023;435: 168141.
Eberhardt J, Santos-Martins D, Tillack AF, Forli S. AutoDock vina 1.2.0: new docking methods, expanded force field, and python bindings. J Chem Inf Model. 2021;61:3891–8.
Linder ME, Deschenes RJ. Palmitoylation: policing protein stability and traffic. Nat Rev Mol Cell Biol. 2007;8:74–84.
Ko PJ, Dixon SJ. Protein palmitoylation and cancer. EMBO Rep. 2018;19: e46666.
Long JZ, Cravatt BF. The metabolic serine hydrolases and their functions in mammalian physiology and disease. Chem Rev. 2011;111:6022–63.
Lin DTS, Conibear E. ABHD17 proteins are novel protein depalmitoylases that regulate N-Ras palmitate turnover and subcellular localization. Elife. 2015;4: e11306.
Remsberg JR, Suciu RM, Zambetti NA, Hanigan TW, Firestone AJ, Inguva A, Long A, Ngo N, Lum KM, Henry CL, et al. ABHD17 regulation of plasma membrane palmitoylation and N-Ras-dependent cancer growth. Nat Chem Biol. 2021;17:856–64.
Yao H, Lan J, Li CS, Shi HB, Brosseau JP, Wang HB, Lu HJ, Fang CY, Zhang Y, Liang LX, et al. Inhibiting PD-L1 palmitoylation enhances T-cell immune responses against tumours. Nat Biomed Eng. 2019;3:306–17.
Duffy MJ, O’Grady S, Tang MH, Crown J. MYC as a target for cancer treatment. Cancer Treat Rev. 2021;94: 102154.
Chen H, Liu HD, Qing GL. Targeting oncogenic Myc as a strategy for cancer treatment. Signal Transduct Target Ther. 2018;3:5.
Hsieh AL, Walton ZE, Altman BJ, Stine ZE, Dang CV. MYC and metabolism on the path to cancer. Semin Cell Dev Biol. 2015;43:11–21.
Zhang Q, Cui K, Yang X, He Q, Yu J, Yang L, Yao G, Guo W, Luo Z, Liu Y, et al. c-Myc-IMPDH1/2 axis promotes tumourigenesis by regulating GTP metabolic reprogramming. Clin Transl Med. 2023;13: e1164.
Tian J, Cai Y, Li Y, Lu Z, Huang J, Deng Y, Yang N, Wang X, Ying P, Zhang S, et al. CancerImmunityQTL: a database to systematically evaluate the impact of genetic variants on immune infiltration in human cancer. Nucleic Acids Res. 2021;49:D1065-d1073.
Fridman WH, Pagès F, Sautès-Fridman C, Galon J. The immune contexture in human tumours: impact on clinical outcome. Nat Rev Cancer. 2012;12:298–306.
Villegas FR, Coca S, Villarrubia VG, Jiménez R, Chillón MJ, Jareño J, Zuil M, Callol L. Prognostic significance of tumor infiltrating natural killer cells subset CD57 in patients with squamous cell lung cancer. Lung Cancer. 2002;35:23–8.
Coca S, Perez-Piqueras J, Martinez D, Colmenarejo A, Saez MA, Vallejo C, Martos JA, Moreno M. The prognostic significance of intratumoral natural killer cells in patients with colorectal carcinoma. Cancer. 1997;79:2320–8.
Ishigami S, Natsugoe S, Tokuda K, Nakajo A, Che X, Iwashige H, Aridome K, Hokita S, Aikou T. Prognostic value of intratumoral natural killer cells in gastric carcinoma. Cancer. 2000;88:577–83.
Zhu LY, Zhou J, Liu YZ. Pan WD [Prognostic significance of natural killer cell infiltration in hepatocellular carcinoma]. Ai Zheng. 2009;28:1198–202.
Cui K, Yao S, Liu B, Sun S, Gong L, Li Q, Fei B, Huang Z. A novel high-risk subpopulation identified by CTSL and ZBTB7B in gastric cancer. Br J Cancer. 2022;127:1450–60.
Cui K, Yao S, Zhang H, Zhou M, Liu B, Cao Y, Fei B, Huang S, Huang Z. Identification of an immune overdrive high-risk subpopulation with aberrant expression of FOXP3 and CTLA4 in colorectal cancer. Oncogene. 2021;40:2130–45.
Du W, Hua F, Li X, Zhang J, Li S, Wang W, Zhou J, Wang W, Liao P, Yan Y, et al. Loss of optineurin drives cancer immune evasion via palmitoylation-dependent IFNGR1 lysosomal sorting and degradation. Cancer Discov. 2021;11:1826–43.
Hernandez JL, Davda D, Majmudar JD, Won SJ, Prakash A, Choi AI, Martin BR. Correlated S-palmitoylation profiling of Snail-induced epithelial to mesenchymal transition. Mol Biosyst. 2016;12:1799–808.
Hsiao YH, Lin CI, Liao H, Chen YH, Lin SH. Palmitic acid-induced neuron cell cycle G2/M arrest and endoplasmic reticular stress through protein palmitoylation in SH-SY5Y human neuroblastoma cells. Int J Mol Sci. 2014;15:20876–99.
Huang J, Li J, Tang J, Wu Y, Dai F, Yi Z, Wang Y, Li Y, Wu Y, Ren G, Xiang T. ZDHHC22-mediated mTOR palmitoylation restrains breast cancer growth and endocrine therapy resistance. Int J Biol Sci. 2022;18:2833–50.
Helleday T, Eshtad S, Nik-Zainal S. Mechanisms underlying mutational signatures in human cancers. Nat Rev Genet. 2014;15:585–98.
Cui K, Liu C, Li X, Zhang Q, Li Y. Comprehensive characterization of the rRNA metabolism-related genes in human cancer. Oncogene. 2020;39:786–800.
Gamazon ER, Stranger BE. The impact of human copy number variation on gene expression. Brief Funct Genomics. 2015;14:352–7.
Henrichsen CN, Chaignat E, Reymond A. Copy number variants, diseases and gene expression. Hum Mol Genet. 2009;18:R1–8.
Keenan AB, Torre D, Lachmann A, Leong AK, Wojciechowicz ML, Utti V, Jagodnik KM, Kropiwnicki E, Wang Z, Ma’ayan A. ChEA3: transcription factor enrichment analysis by orthogonal omics integration. Nucleic Acids Res. 2019;47:W212–24.
Topalian SL, Taube JM, Pardoll DM. Neoadjuvant checkpoint blockade for cancer immunotherapy. Science. 2020;367:eaax0182.
Zhao C, Yu H, Fan X, Niu W, Fan J, Sun S, Gong M, Zhao B, Fang Z, Chen X. GSK3β palmitoylation mediated by ZDHHC4 promotes tumorigenicity of glioblastoma stem cells in temozolomide-resistant glioblastoma through the EZH2-STAT3 axis. Oncogenesis. 2022;11:28.
Zambetti NA, Firestone AJ, Remsberg JR, Huang BJ, Wong JC, Long AM, Predovic M, Suciu RM, Inguva A, Kogan SC, et al. Genetic disruption of N-RasG12D palmitoylation perturbs hematopoiesis and prevents myeloid transformation in mice. Blood. 2020;135:1772–82.
Stine ZE, Walton ZE, Altman BJ, Hsieh AL, Dang CV. MYC, Metabolism, and Cancer. Cancer Discov. 2015;5:1024–39.
Lourenco C, Resetca D, Redel C, Lin P, MacDonald AS, Ciaccio R, Kenney TMG, Wei Y, Andrews DW, Sunnerhagen M, et al. MYC protein interactors in gene transcription and cancer. Nat Rev Cancer. 2021;21:579–91.
Dhanasekaran R, Deutzmann A, Mahauad-Fernandez WD, Hansen AS, Gouw AM, Felsher DW. The MYC oncogene - the grand orchestrator of cancer growth and immune evasion. Nat Rev Clin Oncol. 2022;19:23–36.
Massó-Vallés D, Beaulieu ME, Soucek L. MYC, MYCL, and MYCN as therapeutic targets in lung cancer. Expert Opin Ther Targets. 2020;24:101–14.
Shi C, Yang X, Liu Y, Li H, Chu H, Li G, Yin H. ZDHHC18 negatively regulates cGAS-mediated innate immunity through palmitoylation. EMBO J. 2022;41: e109272.
Steegmaier M, Hoffmann M, Baum A, Lenart P, Petronczki M, Krssak M, Gurtler U, Garin-Chesa P, Lieb S, Quant J, et al. BI 2536, a potent and selective inhibitor of polo-like kinase 1, inhibits tumor growth in vivo. Curr Biol. 2007;17:316–22.
Saphner T, Weller EA, Tormey DC, Pandya KJ, Falkson CI, Stewart J, Robert NJ. 21-day oral etoposide for metastatic breast cancer—a phase II study and review of the literature. Am J Clin Oncol. 2000;23:258–62.
Mohler H, Pfirrmann RW, Frei K. Redox-directed cancer therapeutics. Taurolidine and Piperlongumine as broadly effective antineoplastic agents. Int J Oncol. 2014;45:1329–36.
We thank Dr. Guojun Zhao (The First Affiliated Hospital of Zhengzhou University, China) for providing the ZDHHC23 and ZDHHC7 plasmids. We acknowledge the TCGA, CTRP, GSCA, GDSC and Assistant for Clinical Bioinformatics projects.
This work was supported by grants from the National Natural Science Foundation of China (82103085; 82074133) and National Science and Technology Major Infectious Diseases Project during the 12th 5-Year Plan Period (2014ZX10003002-003-002).
Ethics approval and consent to participate
Consent for publication
All authors agree to publish.
The authors declare that they have no conflict of interest.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
List of primer and shRNA sequences used in this study. Table S2. List of TCGA cancer types.Table S3. List of antibodies used in this study. Table S4. List of palmitoylation-related genes. Table S5. Docking information of BI-2536 and palmitoylation-associated proteins.
Expression of most palmitoylation-related genes is equal in various normal tissues. Fig. S2. Expression of palmitoylation-related genes is associated with cancer subtype and patient survival. Fig. S3. Mutation landscape of palmitoylation-related genes in human cancers. Fig. S4. CNV landscape of palmitoylation-related genes in human cancers. Fig. S5. Myc regulates palmitoylation. Fig. S6. Dysregulated palmitoylation is associated with immune infiltration in human cancers. Fig. S7. The expression levels of ZDHHC9 and ABHD17C genes are correlated with TMB and MSI. Fig. S8. Etoposide and piperlongumine are potential small molecules for regulating palmitoylation. Fig. S9. BI-2536, etoposide and piperlongumine regulate the expression of palmitoylation-related genes. Fig. S10. BI-2536, etoposide and piperlongumine are potential small molecules for targeting palmitoylation-related proteins.
About this article
Cite this article
Kong, Y., Liu, Y., Li, X. et al. Palmitoylation landscapes across human cancers reveal a role of palmitoylation in tumorigenesis. J Transl Med 21, 826 (2023). https://doi.org/10.1186/s12967-023-04611-8