Transcriptome data, survival information, somatic mutation data and copy number variation (CNV) data for patients with CRC were retrieved from The Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov/). After excluding patients with a survival time of <30 days, a total of 418 patients were included in the subsequent analyses, and 41 healthy individuals served as negative controls. The GSE17538 dataset, used as an external validation cohort, was obtained from the Gene Expression Omnibus database (https://www.ncbi.nlm.nih.gov/geo/). Additionally, a list of 3,532 palmitoylation-related genes was extracted from the GeneCards database (https://www.genecards.org/).

To identify prognostic genes, univariate Cox regression analysis was first performed, with a significance threshold of P<0.001. To enhance the selection process, Least Absolute Shrinkage and Selection Operator (LASSO) regression was performed using the glmnet package (https://cran.r-project.org/web/packages/glmnet/index.html) in R studio 4.3.1 version software (https://cran.rstudio.com/). LASSO regression analysis was performed employing 10-fold cross-validation to determine the optimal regularization parameter (λ). The λ value selected was the one that maximized model fit while simultaneously minimizing the risk of overfitting. Samples with missing or insufficient clinical data were excluded during univariate Cox analysis, ensuring the accuracy of LASSO regression analysis, which was applied to select the most relevant prognostic genes from the initial set of identified candidates. Following this, a multivariate Cox regression analysis was performed to establish a palmitoylation-related prognostic model. In addition, the differences in the expression profiles of the six model genes between normal and CRC tissues were assessed and were depicted in a heatmap, which was generated using the pheatmap package (https://cran.r-project.org/web/packages/pheatmap/index.html) in R studio 4.3.1 version software (19). Patients in both the training and validation cohorts were stratified into high- and low-risk groups based on the median risk score derived from the model. The receiver operating characteristic (ROC) curve was generated to evaluate the prognostic performance of the model, with the area under the curve (AUC) used to assess sensitivity and specificity. Kaplan-Meier (KM) OS analysis was subsequently performed to assess the association between the palmitoylation-related risk signature and OS. The OS of patients was evaluated by KM analysis using the survival package (https://cran.r-project.org/web/packages/survival/index.html) of R studio 4.3.1 version software (19); log-rank test was used to obtain P-values. Finally, the prognostic impact of the expression levels of individual model genes was evaluated to understand their contribution to patient survival outcomes.

Firstly, the expression levels of the model genes were compared between CRC tissues and adjacent normal tissues to identify potential alterations associated with CRC. Subsequently, the correlation between the expression levels of each model gene and the risk score was assessed to understand their contributions to the risk stratification and their prognostic relevance in patients with CRC with Pearson correlation analysis. A nomogram was constructed using the ‘rms’ package (https://cran.r-project.org/web/packages/rms/index.html) in R studio 4.3.1 version software incorporating clinical features, such as tumor (T) stage, lymph node (N) stage, clinical stage, age, sex and risk score, to predict individualized prognosis. Differences in clinical characteristics, including survival status and tumor stage [T, N and metastasis (M) stages], were then compared between the high- and low-risk groups to evaluate their potential impact on prognosis. Finally, the prognostic accuracy of the model was validated using the independent GSE17538 validation dataset, and survival outcomes for patients in the high- and low-risk groups were compared to assess the robustness of the risk signature.

To explore the molecular mechanisms underlying the palmitoylation-related risk signature, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were performed using the ClusterProfiler package (https://www.bioconductor.org/packages/release/bioc/html/clusterProfiler.html) (20). These analyses aimed to identify key biological processes and signaling pathways associated with the signature. Somatic mutation data for patients with CRC were converted to Mutation Annotation Format using the maftools R package for the subsequent TMB analysis (21). Mutational landscapes were subsequently compared between the high- and low-risk groups to assess potential differences in mutation profiles. Additionally, the TMB was calculated for each sample, and the difference in prognosis between the high-TMB and low-TMB groups was examined to investigate the potential impact of TMB on the prognostic model (22).

CNV was assessed using the Genomic Identification of Significant Targets in Cancer (GISTIC) algorithm (23). This tool identifies genomic regions with recurrent amplifications or deletions across a set of samples, assigning statistical significance to these alterations. The GISTIC analysis aimed to uncover potential oncogenes or tumor suppressor genes associated with CRC. The results of the CNV analysis were visualized through copy number plots, highlighting the genomic regions with significant gains or losses and their potential relevance to cancer development.

Immune cell infiltration and immune function scores were assessed using the single-sample Gene Set Enrichment Analysis algorithm (24), which calculates rank scores for each gene based on gene expression profiles. Differences in immune cell infiltration and immune function scores between the high- and low-risk groups were then compared. A comparative analysis of the differential expression of 50 immune checkpoint genes between the two groups was also performed using the CIBERSORT algorithm (https://github.com/Moonerss/CIBERSORT) (25).

Due to the critical role of immune checkpoint inhibitors in CRC treatment, the correlation between the risk score and the expression levels of programmed cell death protein 1 (PD-1), programmed death-ligand 1 (PD-L1), cytotoxic T-lymphocyte associated protein 4 (CTLA4) and T cell immunoreceptor with Ig and ITIM domains (TIGIT) was examined. The Tumor Immune Dysfunction and Exclusion (TIDE) algorithm (26) was used to predict the response of high- and low-risk groups to immune checkpoint inhibitors. TIDE analyzed gene expression data to assess tumor immune evasion mechanisms, calculate immune evasion scores and predict the likelihood of response to immune therapy for each group (27). Subsequently, SubMap analysis (http://cloud.genepattern.org/gp) was performed to further evaluate the immune therapy responsiveness in the high- and low-risk groups, identifying potential differences in treatment outcomes (28). Finally, drug sensitivity was assessed using the pRRophetic package (29), which calculated the half-maximal inhibitory concentration values for various agents from the Genomics of Drug Sensitivity in Cancer (GDSC) database (http://www.cancerrxgene.org), comparing sensitivity between the high- and low-risk patient groups.

HCT116 (cat. no. CCL-247) and HT29 (cat. no. HTB-38) cell lines were purchased from American Type Culture Collection. Cells were cultured in 90% Dulbecco's Modified Eagle Medium (DMEM; Gibco; Thermo Fisher Scientific, Inc.) supplemented with 10% fetal bovine serum (FBS; Gibco; Thermo Fisher Scientific, Inc.) and 1% penicillin/streptomycin (Shanghai Yeasen Biotechnology Co., Ltd.), and were maintained in a humidified incubator at 37°C with 5% CO₂. Cells were passaged twice, resulting in the use of third-generation cells for subsequent experiments. A lentivirus-based short hairpin RNA (shRNA) approach was employed using three distinct shRNA sequences. These shRNA sequences were cloned into the pLKO.1 lentiviral vector with the aim of downregulating the expression of KRT8P12 in HCT116 cells. In addition, full-length KRT8P12 cDNA was synthesized and cloned into the pLV-CMV-MCS-PGK-Puro lentiviral overexpression vector for the purpose of inducing KRT8P12 overexpression in HT29 cells. Both pLKO.1 and pLV-CMV-MCS-PGK-Puro vectors were procured from Beijing Zhongyuan, Ltd. Lentiviral particles were produced by co-transfecting 293T cells with a packaging plasmid mix (REV:VSVG:PMDL, 2:3:5) in 400 µl serum-free DMEM, 1.5 µg core plasmid and 1.5 µg viral packaging plasmid with 6 µl TurboFect (TurboFect:plasmid, 2:1; cat. no. R0532; Fermentas; Thermo Fisher Scientific, Inc.) at 37°C. The reagent was mixed thoroughly and let sit for 20 min before being added to the cells. Then the cells were harvested 48 h post-transfection. The resultant virus-containing supernatant was collected and filtered. This supernatant was subsequently utilized to infect HCT116 and HT29 cells. Briefly, the HCT116 and HT29 cells were seeded into 6-well plates and cultured until they reached 60–70% confluence, after which, lentiviral particles (multiplicity of infection=10) along with poly-L-lysine (8 µg/ml) were added to enhance viral transduction efficacy. The cells were co-cultured with the virus for 12 h at 37°C, after which they were replaced with fresh culture medium. A total of 48 h post-infection, stable transductants were selected; puromycin was the antibiotic used to confirm successful transduction. The original concentration of puromycin used was 1 mg/ml and 2 µg/ml puromycin was used for selection and maintenance. All culture processes were conducted under conditions of 37°C and 5% CO₂. The sequences for KRT8P12 knockdown and overexpression are shown in Table SI. HCT116, HT29, sh-KRT8P12 and OE-KRT8P12 cells were cultured in DMEM, supplemented with 10% FBS and 1% penicillin-streptomycin. Cells without any treatment served as a blank control group, whereas a pLKO.1 vector containing a nontargeting shRNA and an empty pLV-CMV-MCS-PGK-Puro vector served as negative controls for knockdown and overexpression, respectively. Cells were maintained at 37°C in a 5% CO₂ incubator, with medium changes every 2–3 days. For passaging, cells were detached using 0.25% trypsin-EDTA, centrifuged at 80 × g for 3 min at room temperature and reseeded at a 1:3 to 1:6 dilution. When freezing, cells were suspended in a cryoprotectant medium (10% DMSO, 20% FBS and 70% DMEM) and stored at −80°C or in liquid nitrogen for long-term preservation.

Total RNA was isolated from cell samples utilizing TRIzol^® reagent (Invitrogen; Thermo Fisher Scientific, Inc.), and was subsequently converted into cDNA employing a cDNA synthesis kit (Vazyme Biotech Co., Ltd.) according to the manufacturer's protocol. qPCR was conducted using SYBR Green qRT-PCR Master Mix (Vazyme Biotech Co., Ltd.), with GAPDH serving as the internal control for normalization. qPCR was performed as follows: Pre-denaturation for 2 min at 95°C, followed by 40 cycles of denaturation for 20 sec at 94°C, annealing for 20 sec at 60°C and extension for 20 sec at 72°C. The 2^−ΔΔCq method was used for quantification (30). The primer sequences used are listed in Table SII.

To assess cell proliferation, the CCK-8 assay (Beijing Solarbio Science & Technology Co., Ltd.) was performed according to the manufacturer's instructions. CRC cells were seeded in 96-well plates at a density of 1×10⁴ cells/well and allowed to adhere overnight. 2-Bromopalmitate (2-BP) is a broad-spectrum inhibitor of protein palmitoylation and was purchased from MilliporeSigma. Before CCK-8 assay, the cells underwent different experimental conditions (sh-KRT8P12 or OE-KRT8P12 transduction, and 2-BP treatment). For 2-BP treatment, cells were treated with 50 µM 2-BP for 1 h at 37°C to inhibit palmitoylation in CRC cells. After 24, 48 or 72 h incubation, 10 µl CCK-8 solution was added to each well and the cells were incubated for an additional 2 h at 37°C. The absorbance was measured at 450 nm using a microplate reader. The relative cell proliferation was calculated by comparing the absorbance values of treated cells to the control group. Each experiment was performed in triplicate and repeated at least three times.

To evaluate cell proliferation, the EdU incorporation assay (cat. no. E-CK-A376; Elabscience; Elabscience Bionovation Inc.) was performed. CRC cells were seeded in 24-well plates at a density of 2×10⁴ cells/well and allowed to adhere overnight. After treatment under different experimental conditions (sh-KRT8P12 or OE-KRT8P12 transduction, and 2-BP treatment), the cells were incubated with 10 µM EdU for 2 h at 37°C in a 5% CO₂ incubator. Following incubation, cells were fixed with 4% paraformaldehyde for 15 min at room temperature, and were then permeabilized with 0.3% Triton X-100 and stained with Click-iT EdU reaction cocktail, according to the manufacturer's protocol. The cell nuclei were counterstained with Hoechst to visualize the total cells. Proliferating cells were observed by fluorescence microscopy. Each experiment was performed in triplicate and repeated at least three times.

To assess cell migration, a Transwell assay was performed using 24-well plates with a polycarbonate membrane insert (pore size, 8 µm). CRC cells were seeded into the upper chamber of the Transwell insert at a density of 1×10⁴ cells/well in 200 µl serum-free medium. The lower chamber was filled with 600 µl complete medium containing 10% FBS to act as a chemoattractant. After incubating for 24 h at 37°C in a 5% CO₂ incubator, non-migrated cells on the upper surface of the membrane were gently removed using a cotton swab. Migrated cells on the lower surface of the membrane were fixed with 4% paraformaldehyde for 10 min and stained with 0.1% crystal violet solution for 15 min at room temperature. Images of the stained cells were then captured using an inverted light microscope, and the number of migrated cells was quantified by counting the cells in five random fields per well. Each experiment was performed in triplicate and repeated at least three times.

To assess apoptosis, flow cytometric analysis was performed using an Annexin V-FITC/PI apoptosis detection kit (cat. no. C1062S; Beyotime Institute of Biotechnology). CRC cells were seeded in 6-well plates at a density of 1×10⁵ cells/well and underwent different experimental conditions: sh-KRT8P12 or OE-KRT8P12 transduction, and 2-BP treatment (50 µM 2-BP for 1 h at 37°C in a 5% CO₂ incubator). After treatment, the cells were harvested by trypsinization, washed twice with phosphate-buffered saline and resuspended in 1X binding buffer at a concentration of 1×10⁶ cells/ml. Subsequently, 5 µl Annexin V-FITC and 5 µl PI were added to a 100-µl cell suspension, and the mixture was incubated in the dark at room temperature for 15–20 min. After incubation, 400 µl 1X binding buffer was added to each sample. The stained cells were analyzed using a FACSCalibur flow cytometer (BD Biosciences). Annexin V-positive and PI-negative cells were considered early apoptotic, whereas Annexin V-positive and PI-positive cells were considered late apoptotic or necrotic. The percentage of apoptotic cells was calculated by summing the early and late apoptotic population using FlowJo version 10.8.1 (BD Biosciences). Each experiment was performed in triplicate and repeated at least three times.

All the bioinformatics analyses in this manuscript were conducted using R studio 4.3.1 version software. Data were analyzed using GraphPad Prism 9.0 software (Dotmatics). All experiments were performed in triplicate and results are presented as the mean ± standard deviation. Differences between two groups were analyzed by unpaired Student's t-test. Differences between more than two groups were analyzed using one-way ANOVA with the Bonferroni post hoc test. P<0.05 was considered to indicate a statistically significant difference. All experiments were repeated at least three times to ensure reproducibility.

Univariate Cox regression analysis identified 10 palmitoylation-related genes significantly associated with OS (Fig. 1A). Among these, ZDHHC3, PARGC1A, PMM2, MMAA, LARS2, FDFT1 and CPT2 were protective factors (HR <1), while TIGD1, MPP2 and KRT8P12 were risk factors (HR >1). To refine the model, LASSO regression analysis was performed, narrowing the candidate genes to a manageable subset for multivariate analysis (Fig. 1B). Multivariate Cox regression modeling ultimately identified six critical genes, KRT8P12, ZDHHC3, PCOLCE2, MPP2, LARS2 and MMAA, to construct a palmitoylation-related risk signature (Fig. 1C and D). The risk score formula was defined as follows: RiskScore=(KRT8P12 × 0.105257166581237)-(ZDHHC3 × 0.0240068245440265) + (PCOLCE2 × 0.0952477383419697) + (MPP2 × 0.177664830643707)-(LARS2 × 0.0426016477398653)-(MMAA × 0.1200626992622).

The prognostic efficacy of the model was evaluated using time-dependent ROC curves, which yielded AUC values of 0.733, 0.733 and 0.734 for 1-, 3- and 5-year survival, respectively (Fig. 1E). KM OS analysis demonstrated that patients in the high-risk group exhibited significantly worse survival compared with those in the low-risk group (Fig. 1F). Furthermore, individual KM analyses of the six model genes confirmed their prognostic relevance, with all genes significantly associated with survival outcomes (Fig. 1G). These findings indicated the robust prognostic utility of the palmitoylation-related risk signature in stratifying patients with CRC.

The heatmap of six model genes, highlighted their potential roles in tumor biology (Fig. 2A). Correlation analysis revealed a strong positive association between KRT8P12 expression and the risk score, making it the most influential contributor to the risk signature, as shown in the lollipop plot (Fig. 2B). To facilitate individualized prognostic prediction, a nomogram integrating key clinical parameters, i.e. age, stage, T stage, N stage, sex and the risk score, was constructed, providing a visual tool for risk stratification (Fig. 2C). Analysis of clinical characteristics between high- and low-risk groups demonstrated that patients in the high-risk group had significantly worse prognoses and were more likely to present with advanced disease, including higher stage, T stage, N stage and M stage classifications (Fig. 2D). The prognostic utility of the model was further validated in the independent GSE17538 dataset, where Kaplan-Meier OS analysis confirmed that patients in the high-risk group experienced significantly worse outcomes (Fig. 2E). These findings underscore the robustness of the palmitoylation-related risk signature and its strong association with disease progression and survival outcomes.

GO and KEGG analyses revealed both common and unique pathways in the high- and low-risk groups. Shared pathways between the groups included ‘olfactory transduction’ and ‘olfactory receptor activity’, underscoring certain overlapping molecular mechanisms (Fig. 3A and B). However, distinct pathways were observed in each group: The high-risk group was significantly enriched for ‘intermediate filament’ and ‘sensory perception of smell’, whereas the low-risk group exhibited unique enrichment in ‘odorant binding’ and ‘detection of chemical stimulus involved in sensory perception’. These findings suggested distinct biological underpinnings contributing to the differential prognosis between the two groups. Somatic mutation landscape analyses demonstrated no significant differences in mutation frequencies between the high- and low-risk groups, with the five most frequently mutated genes being APC, TP53, TTN, KRAS and PIK3CA in both groups (Fig. 3C and D). KM OS analysis stratified by TMB revealed that patients with high TMB had significantly worse survival outcomes compared with those with low TMB (Fig. 3E). Furthermore, combining TMB status with risk score stratification provided a more refined prognostic insight: Patients with high-risk/high-TMB demonstrated the worst survival outcomes, highlighting a strong interplay between mutational burden and risk group in determining prognosis (Fig. 3F). These results emphasized the clinical relevance of integrating genetic and molecular signatures for precise prognostic prediction in patients with CRC.

CNV analysis uncovered distinct genomic alterations between high- and low-risk groups. The percentages of copy number gains and losses displayed significant differences, with the high-risk group showing a notably different distribution compared with the low-risk group (Fig. 4A). Additionally, the high-risk group exhibited a higher frequency of amplification and deletion events, reflecting greater genomic instability (Fig. 4B). Further analysis of CNV burden at both focal and arm levels revealed that the high-risk group consistently demonstrated significantly elevated levels of copy number gains and losses compared with the low-risk group (Fig. 4C). These findings suggested that the increased CNV burden in the high-risk group may contribute to the worse prognosis observed in these patients, highlighting the potential role of genomic instability in risk stratification.

Immune cell infiltration analysis demonstrated notable differences in immune cell abundance between the high- and low-risk groups. Differential expression analysis revealed that most immune cell types were significantly more abundant in the low-risk group, suggesting enhanced immune activity (Fig. 5A). Similarly, functional pathway analysis indicated that key immune-related processes were more active in the low-risk group, further highlighting the distinct immune landscape between the two groups (Fig. 5B). Additionally, the expression levels of common immune checkpoint genes exhibited significant differences between the groups, with an overall trend of higher expression in the low-risk group, reflecting its more robust immune environment (Fig. 5C).

Treatment prediction analysis further elucidated the clinical implications of the risk signature. Correlation analysis of the risk score with the expression of immune checkpoint genes, including PD-1, CD274 (PD-L1), CTLA4 and TIGIT, revealed a strong negative association, with higher expression of these genes in the low-risk group; however, although the association was observed, it was not statistically significant (Fig. 6A). TIDE analysis supported these findings, showing that low-risk patients had lower TIDE scores, indicating a greater likelihood of responding to immune checkpoint inhibitors (Fig. 6B). SubMap analysis further demonstrated that low-risk patients might benefit more from immunotherapy (Fig. 6C).

In addition to immunotherapy, drug sensitivity analysis based on the GDSC database identified significant differences in response to several therapeutic agents. The analysis highlighted the top eight drugs with the most pronounced differences in sensitivity between the high- and low-risk groups, with the drug sensitivity being lower in the low risk group, underscoring distinct drug response profiles associated with the risk score (Fig. 6D). Vinblastine is a microtubule inhibitor that prevents cell mitosis by inhibiting microtubule polymerization, which is primarily used for treating Hodgkin's lymphoma, non-Hodgkin's lymphoma and certain solid tumors (31). Afatinib is an irreversible tyrosine kinase inhibitor that targets EGFR mutations, mainly used in non-small cell lung cancer (32). 5-Fluorouracil is a pyrimidine analog that inhibits thymidylate synthase, thus disrupting DNA synthesis, which is widely used in CRC, gastric and breast cancer (33). Gemcitabine is a nucleoside analog that inhibits DNA replication, commonly used in pancreatic, lung and breast cancer (34). YK-4-279 is an EWS-FLI1 inhibitor, primarily targeting EWS-FLI1 fusion gene-associated tumors such as Ewing sarcoma (35). Vincristine is a microtubule inhibitor that blocks mitosis, and is mainly used for leukemia, lymphoma and certain solid tumors (36). Vinorelbine, a vinca alkaloid, inhibits microtubule function and is mainly used for non-small cell lung cancer and breast cancer (37). GSK591 is a protein arginine methyltransferase 5 inhibitor that affects epigenetic regulation and has potential applications in several types of cancer (38). These findings provide a strong foundation for leveraging the risk signature to guide personalized treatment strategies, including both immunotherapy and targeted drug therapies.

KRT8P12 is a gene associated with palmitoylation. To verify the regulatory effect of KRT8P12 on CRC cells, the expression of the KRT8P12 gene was assessed in CRC cell lines HCT116 and HT29 using qPCR. Higher expression was demonstrated in HCT116 cells (Fig. 7A). The expression of KRT8P12 was subsequently knocked down in the HCT116 cell line using shRNA lentiviruses (sh-KRT8P12-3) and overexpression was performed in the HT29 cell line (OE-KRT8P12) (Fig. 7B). Since the target sequence of sh-KRT8P12-2 was absolutely matched with KRT8, it might also knock down KRT8, so it was not chosen for subsequent experiments. Subsequently, the proliferation of the constructed cell lines was assessed and it was demonstrated that sh-KRT8P12 decreased the proliferation of CRC cells, whereas overexpression of KRT8P12 promoted the proliferation of CRC cells. By contrast, 2-BP inhibited the proliferation induced by KRT8P12 overexpression (Figs. 7C and 8A). The effect of KRT8P12 on cell migration was evaluated using Transwell chambers, and it was demonstrated that silencing of KRT8P12 reduced the migration of CRC cells, whereas KRT8P12 overexpression significantly promoted tumor cell migration. Treatment with 2-BP inhibited the migration induced by KRT8P12 overexpression (Fig. 8B). The effect of KRT8P12 on CRC cell apoptosis was assessed using Annexin-V/PI staining, and it was demonstrated that sh-KRT8P12 enhanced apoptosis in CRC cells. Furthermore, 2-BP reversed the effect of KRT8P12 overexpression on inhibiting apoptosis (Fig. 9). In conclusion, the palmitoylation-associated gene KRT8P12 may promote tumor cell proliferation and migration, while inhibiting tumor cell apoptosis.