Identification and preliminary examination of prognostic genes linked to calcium metabolism in colorectal cancer
Highlight box
Key findings
• A novel prognostic signature consisting of four calcium metabolism-related genes (CAMRGs)—PRKCB, ATP2A3, PLCB4, and SLC25A6—was identified and validated for colorectal cancer (CRC). The risk model effectively stratifies CRC patients into high- and low-risk groups, with the low-risk group demonstrating significantly better overall survival (OS). The clinical relevance of these genes was further confirmed through reverse transcription-quantitative polymerase chain reaction (RT-qPCR) in patient specimens.
What is known and what is new?
• Calcium metabolism plays a critical role in cancer progression and immune signaling, yet its collective genomic impact on CRC prognosis remains poorly understood.
• This study provides the first comprehensive integration of transcriptomic data and clinical validation to establish a CAMRG-based risk model, offering a more robust predictive tool than single-gene markers.
What is the implication, and what should change now?
• The CAMRG signature can serve as a reliable tool for individualized risk assessment and treatment stratification in CRC patients. Clinicians may use this model to identify high-risk patients who require more intensive monitoring or personalized therapeutic interventions targeting calcium signaling pathways.
Introduction
Colorectal cancer (CRC) is a malignancy stemming from the colonic or rectal mucosal epithelium, with over 1.85 million new instances and 850,000 deaths per annum, making it the third leading cause of cancer death globally (1). According to statistics, the incidence of early-onset CRC (CRC diagnosed in patients who are under 50 years old) is on the rise worldwide (2). The high incidence and fatality rates of CRC, coupled with the gradual shift towards younger age groups, have made it a global health issue. The occurrence and development of CRC usually involve a series of processes from adenoma to carcinogenesis, encompassing various genetic mutations, signaling pathway disorders, and microenvironmental remodeling (3,4). Despite the continuous improvement of treatment methods for CRC, delayed diagnosis, treatment resistance, and poor prognosis are still difficult to solve in clinical practice (5) owing to the atypical early symptoms, easy recurrence, and easy metastasis of CRC (6). Therefore, in-depth exploration of CRC molecular mechanisms and identification of novel molecular markers and therapeutic targets are of great significance for improving early diagnostic efficiency and clinical outcomes (7-10).
Abnormal calcium metabolism is intimately linked to the development and progression of various malignant tumors, especially in CRC, where calcium channel proteins, calcium-binding proteins, and their related signaling pathways are significantly associated with tumor proliferation, invasion, drug resistance, and immune regulation (11-13). However, current research mainly focuses on the functions of individual calcium metabolism-related factors, such as calcium channels and calcium-binding proteins, while there is a lack of in-depth analysis of the network of calcium metabolism-related genes (CAMRGs) and their mechanisms in the prognosis of CRC patients (14-16). Therefore, it is necessary to integrate and analyze CAMRGs through comprehensive bioinformatics approaches to uncover potential molecular signatures, thereby providing new biological insights and clinical decision support for precision medicine in CRC.
This study utilized The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) databases to integrate transcriptomic data and clinical information from CRC patients. Through differential expression analysis, functional enrichment analysis, protein-protein interaction (PPI) network construction, and least absolute shrinkage and selection operator (LASSO) and Cox regression methods, prognostic genes related to calcium metabolism were systematically identified and adopted to build a multigene prognostic model. Furthermore, independent prognostic factors for CRC patients were identified through multivariate Cox regression analysis, and an individualized survival prediction column chart model was subsequently constructed. This study provides a novel molecular tool for risk assessment and personalized treatment of CRC patients, while offering theoretical insights into calcium metabolism disorder mechanisms in cancer progression, demonstrating significant research and translational value. We present this article in accordance with the TRIPOD reporting checklist (available at https://jgo.amegroups.com/article/view/10.21037/jgo-2025-1-972/rc).
Methods
Data source
The external validation cohort, GSE17536 (GPL570), was selectively retrieved from the GEO database based on the following inclusion criteria: (I) a substantial sample size (n=177 disease samples); (II) the availability of complete overall survival (OS) and follow-up information; and (III) high data quality on a standardized platform. This systematic selection was intended to ensure the robustness and reproducibility of the model’s external validation. The TCGA-CRC datasets for colon adenocarcinoma (COAD) and rectum adenocarcinoma (READ) were derived from the TCGA database (accessed on January 22, 2025) and used as a training cohort (17). It contained RNA sequencing (RNA-seq) profiles from 624 CRC tumors and 51 normal tissue samples, along with corresponding survival data and clinicopathological parameters including gender, age, TNM stage, and pathological stage. These 624 tumor samples were used for differential expression analysis. For prognostic modeling, samples with OS duration <60 days (n=127) were filtered out to reduce bias from early mortality. As a result, a total of 497 tumor samples, accompanied by complete and detailed survival information, were ultimately retained and utilized for the progression of the risk assessment model (table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-1.xlsx). In addition, a total of 437 CRC samples had complete clinicopathological annotations.
The validation cohort, GSE17536 (GPL570), was derived from the GEO database (17). After filtering samples with survival duration <60 days and incomplete clinical records, 174 CRC tumor samples were retained.
A total of 178 CAMRGs were retrieved from the Molecular Signatures Database (MSigDB, v7.4). Specifically, these genes were identified the C2:KEGG (Kyoto Encyclopedia of Genes and Genomes) curated collection (18). The complete list of these 178 genes is provided in table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-2.xlsx.
Differential expression analysis and identification of candidate genes
Differentially expressed genes (DEGs) were defined between CRC and normal samples via the “DESeq2” package (v 1.40.2) (19). To minimize the risk of false-positive results, P values were adjusted using the Benjamini-Hochberg (BH) method, and an adjusted P value (Padj) <0.05 along with |log2fold change (FC)| >0.5 were set as the significance thresholds. This relatively permissive threshold was intentionally selected to maintain a broad candidate gene pool and avoid the omission of potentially significant calcium metabolism-related factors at the initial screening stage. Any potential noise or redundant features introduced by this threshold were subsequently addressed through rigorous downstream filtering, including LASSO regression with 15-fold cross-validation and independent external validation. Furthermore, this threshold is a widely accepted standard in transcriptomic studies for balancing sensitivity and specificity. Candidate genes were determined by intersecting DEGs with 178 CAMRGs, and a Venn diagram was depicted by the “ggvenn” package (v 0.1.10) (18).
Functional enrichment analysis and PPI network construction
Initially, Gene Ontology (GO) and KEGG analyses were performed on the 125 candidate CAMRGs to characterize the specific functional features of the genes included in the prognostic model. This target-oriented strategy aimed to define the calcium-related biological context of the signature. To avoid pre-selection bias and provide a broader biological landscape of CRC, a complementary enrichment analysis was subsequently conducted using the entire set of 7,511 DEGs. These DEGs were stratified into upregulated and downregulated groups to identify direction-specific biological processes (BPs) and pathways. GO and KEGG analyses were implemented through the “clusterProfiler” package (v 4.8.3) (20) on candidate genes (Padj<0.05). The top 10 most significant GO terms, containing BP, molecular function (MF), cellular component (CC), and KEGG pathways, were visualized by the “ggplot2” package (v 3.5.1) (21). The PPI network was predicted by uploading candidate genes to the STRING database with an interaction confidence threshold of a combined score ≥0.4. Proteins that could not be mapped to the STRING database or were identified as isolated nodes (lacking interactions with any other candidate genes at the specified threshold) were excluded from the final network analysis. Consequently, a refined PPI network consisting of 90 interconnected genes was constructed for further investigation. To identify highly interconnected clusters within the network, the Molecular Complex Detection (MCODE) plugin in Cytoscape (v3.10.3) was employed with the following parameters: degree cutoff =2, node score cutoff =0.2, k-core =2, and max depth =100. Subsequently, GO and KEGG enrichment analyses were performed on the genes within each identified module to clarify their specific biological roles.
Prognostic gene screening
Univariate Cox regression analysis was implemented on candidate genes via the “survival” package (v 3.7.0) (22), with a hazard ratio (HR) ≠1 and P<0.05 as selection criteria. Proportional hazards (PH) assumption testing was executed via the “cox.zph” function from the “survival” package (P>0.05). Subsequently, the LASSO regression analysis was implemented via the “glmnet” package (v 4.1.8) (23). To optimize the penalty parameter (λ) and effectively control the risk of overfitting, a 15-fold cross-validation was performed. Finally, a stepwise multivariate Cox regression analysis was executed to ensure that the genes included in the signature retained independent prognostic value. This hierarchical screening strategy—initial univariate filtering, LASSO-based dimensionality reduction, and multivariate refinement—is a standard approach in transcriptomic prognostic modeling that balances feature selection efficiency with model stability.
Prognostic risk model construction and validation
The risk score formula was defined as:
where β represented regression coefficients derived from multivariate Cox analysis and x denoted gene expression levels. The optimal cutoff value for the risk score was determined in the TCGA training cohort using the “surv_cutpoint” function from the “survminer” package (v0.4.9), which is based on maximally selected rank statistics (minprop =0.15). To rigorously evaluate the stability and generalizability of the model, this fixed cutoff value was subsequently applied to the external validation cohorts (GSE17536 and GSE56699) without recalculation, categorizing patients into high- and low-risk groups for further survival analysis. In addition, Kaplan-Meier survival curves were created to evaluate survival outcomes via the “survival” package. Moreover, time-dependent Receiver Operating Characteristic (ROC) curves, reflecting 1-, 3-, and 5-year survival, were exhibited with the “survivalROC” package (v 1.0.3.1) (24). The area under the curve (AUC) values were computed to estimate the model’s predictive accuracy. Additionally, the generalizability of the model was validated using two independent external datasets: GSE17536 and GSE56699. The consistency of survival outcomes and ROC curves across these diverse cohorts indicates that the risk model effectively captures biological prognostic signals rather than overfitting to the training data.
Clinical characteristic analysis
The Wilcoxon test was employed on CRC samples from the TCGA-CRC cohort to evaluate risk score distributions across clinical subgroups [gender, age, M stage (metastasis), N stage (lymph node involvement), T stage (tumor invasion depth), and pathological stage] (P<0.05).
Independent prognostic analysis and nomogram development
Univariate Cox regression analysis was first executed on risk scores and clinicopathological features (gender, T stage, N stage, M stage, age and pathological stage) via the “survival” package (v 3.7.0) (HR ≠1, P<0.05). Variables that met the PH assumption (P>0.05) were then chosen for multivariate Cox regression analysis to detect independent prognostic factors. Variables meeting the criteria were selected for nomogram construction via the “regplot” package (v 1.1) (25). Calibration curves were created with the “rms” package (v 6.5.0) (26) to assess prediction accuracy. ROC curves for 1-, 3-, and 5-year survival forecasts were plotted via the “survivalROC” package. To evaluate the added prognostic value of the gene-based risk score, a baseline model (comprising age and TNM stage) and a combined model (baseline factors plus risk score) were constructed. The predictive performance was compared using the C-index, time-dependent ROC curves (AUC), Net Reclassification Index (NRI), and Integrated Discrimination Improvement (IDI).
Exploratory gene set enrichment analysis (GSEA)
Differential expression analysis was implemented between risk groups via the “DESeq2” package (v 1.40.2). GSEA was executed via the “clusterProfiler” package (v 4.8.3). An exploratory threshold of false discovery rate (FDR) <0.25 was adopted according to MSigDB recommendations (27). The analysis utilized the “c2.cp.kegg.v7.4.symbols.gmt” gene set, which was derived from the MSigDB. Enriched pathways were visualized via the “GseaVis” package (v 0.1.0) (28) [FDR <0.25, P<0.05, and | normalized enrichment score (NES)| >1].
Descriptive somatic mutation landscape
Somatic mutation data from the TCGA-CRC cohort were downloaded using the “TCGAbiolinks” package (v 2.28.4) (29). To explore the genomic characteristics of different risk subgroups, waterfall plots depicting mutation profiles were generated using the “maftools” package (v 2.16.0) (28). This descriptive analysis focused on identifying the most frequently mutated genes and comparing the mutational architecture between the high- and low-risk groups.
Preliminary tumor microenvironment analysis
Immune cell infiltration was quantified via the CIBERSORT algorithm via the “immunedeconv” package (v 2.1.0) (30). Differential immune cell abundances were evaluated between risk groups by Wilcoxon tests (P<0.05). To minimize the risk of false-positive results from multiple comparisons, P values were adjusted using the BH method, and a Padj<0.05 was considered statistically significant. Moreover, Spearman correlation analysis was implemented between prognostic genes and immune cells on CRC samples from the TCGA-CRC cohort via the “psych” package (v 2.4.6.26) (31).
Preliminary drug sensitivity profiling
Half-maximal inhibitory concentration (IC50) values for 198 chemotherapeutic agents, which were sourced from the GDSC database, were predicted via the “oncoPredict” package (v 1.2) (32). Wilcoxon tests were employed to evaluate the variations in IC50 between risk groups (P<0.05).
Exploratory prediction of the molecular regulatory network
To generate hypotheses regarding potential non-coding RNA-mediated mechanisms, miRNA-mRNA interactions were predicted using the MiRBD (target score >90) and miRNET databases. A predicted lncRNA-miRNA-mRNA regulatory network was subsequently visualized using Cytoscape (v 3.10.3) (33). These analyses were based solely on sequence complementarity and database algorithms, serving as a preliminary exploration of potential regulatory axes.
Validation of prognostic gene expression via reverse transcription-quantitative polymerase chain reaction (RT-qPCR)
Frozen tissue samples (n=5 per group) from CRC and control groups were collected at Guangdong Provincial Hospital of Chinese Medicine (The Second Affiliated Hospital of Guangzhou University of Chinese Medicine). This study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments. This study was approved by the Ethics Committee of Guangdong Provincial Hospital of Chinese Medicine (October 25th 2024/No. ZE2024-361-01). Informed consent was obtained from all individual participants included in the study. Total RNA was extracted from frozen tissue samples using TRIzol reagent (R401-01, Vazyme, Nanjing, China) according to the instructions of the manufacturer. Reverse transcription was conducted using the HP All-in-one qRT Master Mix II (Catalog No. 24Y0124, Yungen Biotech, Kunming, China). cDNA was thinned 5-20-fold with DNase/RNase-free water before qPCR. qPCR was carried out with the Servicebio® 2×Universal Blue SYBR Green Master Mix (G3326-05, Servicebio, Wuhan, China) on a CFX Connect Real-Time PCR System (XLFZ006, BIO-RAD, Hercules, CA, USA). The 10 µL reaction contained 5 µL of 2×Universal Blue SYBR Green qPCR Master Mix, 3 µL of diluted cDNA, and 1 µL of each primer (forward and reverse, 10 µM, Sangon Biotech, Shanghai, China). Primer sequences were provided in Table S1. The thermal cycling protocol was as follows: initial denaturation at 95 °C for 1 minute, followed by 40 cycles of 95 °C for 20 seconds, 55 °C for 20 seconds, and 72 °C for 30 seconds. A melt curve analysis was performed to assess the expression amplification specificity. Reactions were run in technical triplicates, and Ct values were standardized against glyceraldehyde-3-phosphate dehydrogenase (GAPDH) using the 2−ΔΔCt method. Consistency in data was ensured by verifying pipette tip retention and liquid residues during plate loading. Centrifugation was conducted using a H1 16KR High-Speed Centrifuge (XLSY010, Kecheng Instrument, Changsha, China) and a mini centrifuge (SMC5000F, SCILOGEX, Beijing, China).
Statistical analysis
Statistical analysis was performed using Graphpad Prism (v 10.0) (34). Given the small sample size (n=5 per group), differences between the CRC and control groups were assessed using the two-tailed Mann-Whitney U test, as normality assumptions could not be reliably verified. A P value <0.05 was considered the threshold for statistical significance. The R programming language (v 4.3.3) was employed to conduct various bioinformatics analyses. For all analyses, P<0.05 was considered the threshold for statistical significance.
Results
Identification and functional enrichment analysis of candidate genes and construction of PPI
A total of 7,511 DEGs were identified between CRC and normal samples, comprising 3,888 upregulated and 3,623 downregulated genes (Figure 1A,1B). Intersection analysis of these DEGs with 178 CAMRGs revealed 125 candidate genes (Figure 1C). To evaluate the robustness of this gene selection, a sensitivity analysis was performed using a stricter threshold (|log2FC| >1). Under these criteria, 69 CAMRGs remained as overlapping DEGs (Figure S1). Importantly, three of the four core prognostic genes (PRKCB, ATP2A3, and PLCB4) were still captured within this stricter subset. While SLC25A6 did not meet the |log2FC| >1 criterion, it was retained in the primary analysis to ensure a more comprehensive prognostic signature, given its high statistical significance and biological relevance reported in prior literature. In addition, GO analysis identified 954 enriched terms (table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-3.xlsx), with the top 10 BPs including calcium ion transport (e.g., transmembrane transport, cytosolic import, and homeostasis) and calcium-mediated signaling (Figure 1D). CC analysis demonstrated enrichment in calcium channel complexes, ion transporter complexes, and sarcoplasmic reticulum membranes (Figure 1E). MF analysis further revealed calcium binding (e.g., calmodulin binding) and enzymatic functions linked to calcium signaling, such as phosphatidylinositol phospholipase C activity (Figure 1F). KEGG pathway analysis showed 116 markedly enriched pathways (Padj<0.05) (table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-4.xlsx) and confirmed the calcium signaling pathway as the most enriched (Figure 1G), supporting the key role of calcium dysregulation in CRC. Functional enrichment analyses highlighted calcium-related BPs as the most marked pathways. The PPI network analysis was performed on the proteins corresponding to the candidate genes. Among the 125 candidates, 90 genes were successfully mapped to the STRING database and exhibited interactions meeting the predefined threshold, while the remaining 35 genes were excluded as isolated nodes lacking significant connectivity (Figure 1H). To further investigate these interactions, MCODE analysis was performed (Figure S2). Functional enrichment of these modules was revealed. These modular interactions provide candidate interactors for further experimental investigation into the mechanisms of calcium dysregulation in CRC. To further elucidate the systemic alterations in CRC, enrichment analysis was performed on the total 7,511 DEGs stratified by expression direction. Upregulated genes were predominantly enriched in cell cycle-related processes, DNA replication, and the WNT signaling pathway. Conversely, downregulated genes were significantly associated with calcium-mediated signaling, muscle system processes, and metabolic pathways (Figure S3 and table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-3.xlsx). These findings confirm that while the core prognostic signature is deeply rooted in calcium dysregulation, it operates within a broader framework of oncogenic pathway activation and metabolic reprogramming in CRC.
Prognostic genes screening and development and validation of prognostic model
In total, 19 genes related to prognosis, including PRKCB, ATP2A3, and PLCB4, were recognized by univariate Cox regression analysis of 90 candidate genes (HR ≠1, P<0.05, Figure 2A). Subsequently, LASSO regression analysis with 15-fold cross-validation refined 19 to 10 genes (lambda.min =0.01), which were further optimized via multivariate stepwise Cox regression to determine four prognostic genes, including PRKCB, ATP2A3, PLCB4, and SLC25A6 (Figure S4). Notably, PRKCB, ATP2A3, PLCB4, and SLC25A6 exhibited negative coefficients and HR <1, suggesting their elevated expression conferred protective effects against CRC progression. Moreover, the risk score derived from regression coefficients [risk score = (−0.672 × expression level of PRKCB) + (−0.268 × ATP2A3) + (−0.178 × PLCB4) + (−0.456 × SLC25A6)] (Figure 2B) effectively separated TCGA-CRC samples into the high- (n=88) and low-risk groups (n=409) using an optimal cutoff (−7.699613). Kaplan-Meier analysis demonstrated significantly worse OS in high-risk patients (P<0.0001), with the AUC values above 0.6 at 1-, 3-, and 5-year intervals (Figure 2C,2D, Figure S5). External validation using the GSE17536 cohort confirmed model robustness, showing consistent survival stratification (P=0.003) and predictive accuracy (AUC >0.6) (Figure 2E,2F, Figure S5). These findings highlighted PRKCB, ATP2A3, PLCB4, and SLC25A6 as prognostic genes for CRC, providing significant insights for patient risk categorization and individualized therapeutic approaches. Furthermore, to address the potential heterogeneity between tumor sites, subgroup analysis was conducted within the TCGA-CRC training cohort by separating samples into COAD and READ. As illustrated in Figure S6, the risk model demonstrated exceptional prognostic performance in both subgroups (P<0.0001 for both). Notably, in the READ subgroup, the low-risk patients exhibited a remarkably high survival rate with minimal events during the follow-up period, confirming the model’s efficacy across different anatomical locations of CRC.
Relationship between risk score and clinicopathological characteristics
The prognostic risk score demonstrated meaningful associations with advanced tumor stages in CRC patients. Comparative analysis revealed that patients with higher M stage, N stage, T stage, and overall pathological stage exhibited significantly elevated risk scores relative to those with lower-stage disease (P<0.05, Figure 3A-3F). Specifically, a marked variation was noted between M0 and M1 (P<0.05, Figure 3C); for N stage, significant differences were found between N0 and N2 (P<0.001), and between N1 and N2 (P<0.01, Figure 3D); for pathological stage, stage I showed significantly lower risk scores compared to stage III (P<0.01) and stage IV (P<0.01, Figure 3E); and for T stage, a notable variation was detected between T2 and T4 (P<0.05, Figure 3F). This stage-dependent risk stratification pattern suggested that the calcium-related prognostic signature effectively captured tumor progression dynamics, potentially reflecting enhanced metastatic potential and local invasion capabilities in the high-risk group.
Independent prognostic analysis and construction of a nomogram
Univariate Cox regression analysis determined risk score, age, M stage and pathological stage as marked predictors of OS in CRC patients (P<0.05) (Figure 4A). PH assumption testing confirmed the validity of these covariates (P>0.05). Subsequent multivariate Cox analysis demonstrated that risk score, age, and pathological stage retained independent prognostic value (P<0.05), establishing their robustness in survival prediction (Figure 4B). A nomogram integrating these independent factors was constructed to quantify individual survival probabilities (Figure 4C). The calibration curves revealed a close alignment between predicted and observed survival outcomes, demonstrating high model reliability (Figure 4D). Furthermore, time-dependent ROC analysis validated the nomogram’s discriminative capacity, achieving AUC values of 0.912, 0.823, and 0.826 for 1-, 3-, and 5-year survival predictions (Figure 4E). These findings suggested that the nomogram was a dependable and clinically relevant instrument for individualized survival prediction in CRC patients, with potential utility in guiding risk-adapted treatment strategies.
Pathway enrichment analysis between risk groups
GSEA reflected that a total of 62 pathways were markedly enriched (FDR <0.25, P<0.05, |NES| >1), with oxidative phosphorylation, ribosome, hematopoietic cell lineage, chemokine signaling, and cell adhesion molecules ranking as the top five most enriched pathways (Figure 5A and table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-5.xlsx). These findings collectively highlighted metabolic reprogramming and immune-stromal crosstalk as key biological hallmarks distinguishing CRC risk subgroups, providing mechanistic insights into the prognostic signature’s ability to stratify patients based on tumor biology. The combined model demonstrated superior predictive accuracy over the baseline clinical model. The C-index increased from 0.832 (baseline) to 0.870 (combined), representing a 4.57% improvement. Time-dependent ROC analysis showed that the AUCs for 1-, 3-, and 5-year survival improved from 0.898, 0.804, and 0.779 to 0.926, 0.850, and 0.857, respectively. Furthermore, the NRI (13.59%) and IDI (10.58%) were both significantly greater than zero (P<0.05), confirming that the inclusion of the CAMRG signature significantly enhances the reclassification power and discriminative ability of the clinical staging system (Figure S7).
Somatic mutation landscape
Somatic mutation profiling revealed the genomic characteristics of the two risk groups. Among the patients with detectable mutation data, the tumor suppressor gene APC exhibited the highest mutation frequency in both the high- and low-risk groups (77% mutation rate), predominantly characterized by multi-hit alterations (Figure 5B). These findings underscore the central role of APC in CRC pathogenesis and highlight conserved mutational drivers between molecularly distinct prognostic subgroups. This descriptive comparison suggests that risk stratification primarily reflects transcriptional and epigenetic dysregulation rather than fundamental differences in the frequency of common oncogenic mutations. Distinct immune cell infiltration patterns were found between the two risk groups. Among 22 immune cell types analyzed, M0 macrophages exhibited the highest abundance in both groups (Figure 5C). Comparative analysis identified three differentially infiltrated immune cell populations: resting mast cells, activated NK cells, and resting dendritic cells (P<0.05, Figure S8). These associations suggested potential crosstalk between calcium signaling mediators and immune regulatory mechanisms in CRC progression.
Drug sensitivity profiling between risk groups
A comparative analysis of drug response patterns identified 155 compounds with notably varied sensitivities between high- and low-risk groups (P<0.05, table available at https://cdn.amegroups.cn/static/public/jgo-2025-1-972-6.xlsx). The top five differentially sensitive agents included Ribociclib_1632, GSK269962A_1192, BMS.754807_2171, PF.4708671_1129, and PD173074_1049 (P<0.001), with high-risk patients demonstrating enhanced sensitivity to these agents (Figure S8). These findings highlighted the prognostic model’s capacity to identify patient subgroups with distinct molecular dependencies, providing actionable insights for personalized therapy selection. After FDR correction, three immune cell types remained significantly associated with the prognostic genes, including activated NK cells, resting dendritic cells, and resting mast cells (all Padj<0.001, Figure S9). Similarly, the top sensitive compounds, such as Ribociclib and GSK269962A, demonstrated exceptionally strong associations with the high-risk group (Padj<0.001).
Bioinformatic prediction of the potential lncRNA-miRNA-mRNA regulatory network
Integrated analysis of mRNA-miRNA interactions identified three potential prognostic miRNAs (hsa-miR-7-5p, hsa-miR-214-3p, and hsa-miR-6867-5p) through the intersection of predictions from MiRBD and miRNET databases. Subsequent prediction of miRNA-lncRNA interactions using Starbase and miRNET databases revealed 28 prognostic lncRNAs with overlapping binding potential, such as MIR137HG and LINC01003 (Table S2). A mRNA regulatory network was constructed, proposing hypothetical interactions among miRNAs, lncRNAs, and mRNAs. The network highlighted key regulatory relationships, for example, hsa-miR-7-5p-MIR137HG and hsa-miR-214-3p-ATP2A3 (Figure S10). These findings suggested a potential lncRNA-miRNA-mRNA regulatory network underlying CRC prognosis, offering a theoretical basis for future experimental validation.
Expression patterns of PRKCB, ATP2A3, PLCB4, and SLC25A6 in CRC
The expression levels of the four prognostic genes were evaluated in five pairs of CRC and adjacent normal tissues using RT-qPCR. Specifically, PRKCB expression was markedly downregulated in the CRC group relative to the control group (P=0.03, Mann-Whitney U test, Figure 6A). Similarly, ATP2A3 expression was significantly reduced (P=0.03, Figure 6B). In contrast, PLCB4 expression was substantially upregulated in CRC samples (P=0.03, Figure 6C). However, SLC25A6 expression did not exhibit a statistically significant difference between the groups (P=0.20, Figure 6D).
Discussion
Incidence and death rates for CRC are still high, posing a serious threat to human health (35-37). The occurrence and development of CRC represent a multi-step process that gradually progresses from benign adenomas to invasive adenocarcinomas, involving complex molecular mechanisms (e.g., disrupted signaling pathways, gene mutations, and tumor microenvironment remodeling) (38-40). In recent years, calcium metabolism abnormalities have been closely linked to various cancers, but systematic research remains limited in CRC (41,42). This study applied RNA sequencing data and clinical information from CRC patients in TCGA and GEO databases. Through differential expression analysis and regression methods, four prognostic genes (PRKCB, ATP2A3, PLCB4, and SLC25A6) were identified. A risk model was constructed based on these genes and a preliminary examination in an external dataset. The link between risk scores and clinical traits was subsequently assessed, followed by the identification of independent prognostic factors and the construction of a nomogram for individualized survival prediction.
The four identified CAMRGs (PRKCB, ATP2A3, PLCB4, and SLC25A6) consistently exhibited protective roles in CRC progression (HR <1). PRKCB is known to modulate the immune microenvironment by promoting T cell activation and macrophage polarization (43-45). ATP2A3 (SERCA3) maintains intracellular calcium homeostasis, and its downregulation is associated with endoplasmic reticulum stress and genomic instability (46,47). PLCB4 regulates the IP3/Ca2+ signaling pathway, thereby influencing cytoskeletal remodeling and cellular migration. Additionally, SLC25A6 serves as a mitochondrial carrier involved in ATP/ADP exchange, linking calcium signaling to energy metabolism (48,49). Collectively, these genes integrate calcium signaling, metabolic regulation, and immune response to synergistically suppress CRC development, representing promising biomarkers for prognostic prediction.
GSEA analysis revealed disparate pathways between high- and low-risk groups. Several biological pathways closely associated with CRC progression were markedly enriched in the high-risk group, elucidating potential mechanisms by which CAMRGs regulate tumor progression. Among these enriched pathways, the oxidative phosphorylation pathway showed significant upregulation in the high-risk group. This pathway is central to mitochondrial metabolism and is involved in ATP synthesis; excessive activity may lead to the accumulation of reactive oxygen species (ROS), which can induce DNA damage, gene mutations, and metabolic reprogramming of tumor cells (50-52). Calcium ions play a pivotal role in regulating mitochondrial function and are closely tied to oxidative phosphorylation. In this study, SLC25A6, a key factor in mitochondrial energy transport, might influence this pathway’s activity through changes in its expression, thereby participating in the metabolic regulation and progression of CRC. Additionally, although the calcium signaling pathway was not directly ranked among the top results in this GSEA analysis, it is closely related to the candidate genes selected in this study (such as PRKCB, PLCB4). The calcium signaling pathway regulates cell proliferation, apoptosis, adhesion, and migration, acting as a key regulator of CRC cell behavior (53-56). Studies have demonstrated that abnormal activation of this pathway increases CRC cell invasiveness and treatment resistance. These findings suggest that CAMRGs may influence patient risk stratification through regulation of this signaling network. The enrichment of these pathways supports the biological rationale underlying our calcium metabolism-related prognostic model. These findings indicate that calcium homeostasis imbalance and abnormal energy metabolism represent common molecular characteristics of high-risk CRC patients. This knowledge provides valuable insights for future mechanistic investigations and targeted therapeutic interventions.
In terms of the immune microenvironment, the analysis focused on cell types significantly associated with the risk score. Higher infiltration of activated NK cells and resting dendritic cells was observed in the low-risk group, suggesting a more robust anti-tumor immune surveillance (57,58). These findings indicate that the identified calcium-related signature may reflect the intensity of immune clearance, which contributes to the favorable prognosis of low-risk patients.
To assess the expression characteristics of the prognostic genes identified through screening in CRC, this study conducted RT-qPCR experiments on clinical tissue samples. RT-qPCR analysis of clinical tissue samples revealed divergent expression patterns among the four prognostic genes. Consistent with their protective roles in the risk model, PRKCB and ATP2A3 were significantly downregulated in CRC tissues compared to normal controls. In contrast, PLCB4 was significantly upregulated, while SLC25A6 showed no statistically significant difference between the two groups. These findings indicate that the relationship between CAMRG expression and CRC prognosis is not unidirectional and may involve complex post-transcriptional or cell-type-specific regulatory mechanisms. The partial consistency between tissue expression and the model’s coefficients underscores the need for further functional studies and larger-scale investigations stratified by risk groups.
The calcium metabolism-related prognosis model constructed in this study can be used for accurate assessment of survival risk in CRC patients (AUC >0.6), offering effective stratification and guidance for managing the disease. The subsequently developed nomogram, combined with clinical variables, allows for personalized predictions of survival probability, showing high clinical usability and relevance (AUC >0.8). Additionally, the key genes in this model could serve as potential therapeutic targets, and it’s worth exploring their roles in tumor metabolism, immune response, and apoptosis, which could lead to the development of new targeted drugs or combinations of immunotherapy.
While this study demonstrates strengths in integrating multiple datasets and conducting multidimensional analyses, several limitations should be acknowledged. Firstly, the transcriptomic data were primarily retrieved from public databases, and although experimental validation was performed using RT-qPCR on clinical specimens, the sample size was relatively small, with only five pairs of samples included due to the limited availability of clinical tissues. This small sample size may restrict the statistical power and the generalizability of the experimental findings. Consequently, this part of the study should be regarded as a preliminary exploration. Future research involving larger, multicenter cohorts is essential to further validate the expression patterns and clinical utility of these identified prognostic genes. Secondly, although some gene functions have been reported in the literature, this study did not test their specific mechanisms in CRC; further in-depth functional studies using cell and animal models are needed to verify their biological roles and clinical feasibility. Additionally, the inherent heterogeneity between colon and rectal cancers should be acknowledged. While the training cohort integrated both COAD and READ samples, the external validation was performed on a predominantly colon cancer-based dataset (GSE17536). Although internal subgroup analysis demonstrated consistent performance in both rectal and colon cancer patients, the generalizability of the CAMRG signature to independent rectal cancer cohorts warrants further validation in future large-scale studies. Furthermore, the lncRNA-miRNA-mRNA regulatory network constructed in this study is based entirely on computational predictions from established databases. Due to the lack of miRNA expression profiles in the current datasets, these regulatory axes have not been validated through correlation analysis or experimental assays. Therefore, this network should be interpreted as hypothetical and exploratory, providing a framework for future studies to investigate the specific non-coding RNA-mediated mechanisms in CRC. Overall, this study offers valuable insights into how CAMRGs can be used in assessing CRC prognosis and provides new directions for the development of precision medicine.
Conclusions
In summary, this study systematically built and validated a prognostic model for CRC based on CAMRGs, identifying four genes (PRKCB, ATP2A3, PLCB4, SLC25A6) with independent prognostic value. A graphical tool was developed in conjunction with clinical features to achieve individualized survival predictions. The research also uncovered how these genes might be involved in CRC progression through calcium signaling, energy metabolism, and immune regulation, highlighting how crucial calcium metabolism dysregulation is to CRC development. The model demonstrated robust stability and predictive performance in external validation cohorts. This provides a valuable molecular tool for CRC risk assessment and treatment personalization, while establishing a theoretical foundation for investigating relationships between calcium metabolism, tumor immunity, and metabolic reprogramming.
Acknowledgments
We would like to express our sincere gratitude to all individuals and organizations who supported and assisted us throughout this research. Special thanks to the following authors: Tian-ming Jiang, Hai-peng Huang, Ping Tan. In conclusion, we extend our thanks to everyone who has supported and assisted us along the way. Without their supports, this research would not have been possible.
Footnote
Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://jgo.amegroups.com/article/view/10.21037/jgo-2025-1-972/rc
Data Sharing Statement: Available at https://jgo.amegroups.com/article/view/10.21037/jgo-2025-1-972/dss
Peer Review File: Available at https://jgo.amegroups.com/article/view/10.21037/jgo-2025-1-972/prf
Funding: The work was supported by grants from
Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://jgo.amegroups.com/article/view/10.21037/jgo-2025-1-972/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 performed in line with the principles of the Declaration of Helsinki and its subsequent amendments. This study was approved by the Ethics Committee of Guangdong Provincial Hospital of Chinese Medicine (October 25th 2024/No. ZE2024-361-01). Informed consent was obtained from all individual participants included in the study.
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
- Biller LH, Schrag D. Diagnosis and Treatment of Metastatic Colorectal Cancer: A Review. JAMA 2021;325:669-85. [Crossref] [PubMed]
- Patel SG, Karlitz JJ, Yen T, et al. The rising tide of early-onset colorectal cancer: a comprehensive review of epidemiology, clinical features, biology, risk factors, prevention, and early detection. Lancet Gastroenterol Hepatol 2022;7:262-74. [Crossref] [PubMed]
- Zhao H, Ming T, Tang S, et al. Wnt signaling in colorectal cancer: pathogenic role and therapeutic target. Mol Cancer 2022;21:144. [Crossref] [PubMed]
- Shin AE, Giancotti FG, Rustgi AK. Metastatic colorectal cancer: mechanisms and emerging therapeutics. Trends Pharmacol Sci 2023;44:222-36. [Crossref] [PubMed]
- Fan A, Wang B, Wang X, et al. Immunotherapy in colorectal cancer: current achievements and future perspective. Int J Biol Sci 2021;17:3837-49. [Crossref] [PubMed]
- Abedizadeh R, Majidi F, Khorasani HR, et al. Colorectal cancer: a comprehensive review of carcinogenesis, diagnosis, and novel strategies for classified treatments. Cancer Metastasis Rev 2024;43:729-53. [Crossref] [PubMed]
- Walkon LL, Strubbe-Rivera JO, Bazil JN. Calcium Overload and Mitochondrial Metabolism. Biomolecules 2022;12:1891. [Crossref] [PubMed]
- Luan S, Wang C. Calcium Signaling Mechanisms Across Kingdoms. Annu Rev Cell Dev Biol 2021;37:311-40. [Crossref] [PubMed]
- Eisner D, Neher E, Taschenberger H, et al. Physiology of intracellular calcium buffering. Physiol Rev 2023;103:2767-845. [Crossref] [PubMed]
- O’Grady S, Morgan MP. Calcium transport and signalling in breast cancer: Functional and prognostic significance. Semin Cancer Biol 2021;72:19-26. [Crossref] [PubMed]
- Patergnani S, Danese A, Bouhamida E, et al. Various Aspects of Calcium Signaling in the Regulation of Apoptosis, Autophagy, Cell Proliferation, and Cancer. Int J Mol Sci 2020;21:8323. [Crossref] [PubMed]
- Mitaishvili E, Feinsod H, David Z, et al. The Molecular Mechanisms behind Advanced Breast Cancer Metabolism: Warburg Effect, OXPHOS, and Calcium. Front Biosci (Landmark Ed) 2024;29:99. [Crossref] [PubMed]
- Silvestri R, Nicolì V, Gangadharannambiar P, et al. Calcium signalling pathways in prostate cancer initiation and progression. Nat Rev Urol 2023;20:524-43. [Crossref] [PubMed]
- Aslam A, Ahmad J, Baghdadi MA, et al. Chemopreventive effects of vitamin D(3) and its analogue, paricalcitol, in combination with 5-fluorouracil against colorectal cancer: The role of calcium signalling molecules. Biochim Biophys Acta Mol Basis Dis 2021;1867:166040. [Crossref] [PubMed]
- Bertocci LA, Rovatti JR Jr, Wu A, et al. Calcium handling genes are regulated by promoter DNA methylation in colorectal cancer cells. Eur J Pharmacol 2022;915:174698. [Crossref] [PubMed]
- Liu X, Feng C, Yan L, et al. Calcium channels as pharmacological targets for cancer therapy. Clin Exp Med 2025;25:94. [Crossref] [PubMed]
- Yue T, Chen S, Zhu J, et al. The aging-related risk signature in colorectal cancer. Aging (Albany NY) 2021;13:7330-49. [Crossref] [PubMed]
- Mao W, Ding J, Li Y, et al. Inhibition of cell survival and invasion by Tanshinone IIA via FTH1: A key therapeutic target and biomarker in head and neck squamous cell carcinoma. Exp Ther Med 2022;24:521. [Crossref] [PubMed]
- Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014;15:550. [Crossref] [PubMed]
- Wu T, Hu E, Xu S, et al. clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. Innovation (Camb) 2021;2:100141. [Crossref] [PubMed]
- Gustavsson EK, Zhang D, Reynolds RH, et al. ggtranscript: an R package for the visualization and interpretation of transcript isoforms using ggplot2. Bioinformatics 2022;38:3844-6. [Crossref] [PubMed]
- Lei J, Qu T, Cha L, et al. Clinicopathological characteristics of pheochromocytoma/paraganglioma and screening of prognostic markers. J Surg Oncol 2023;128:510-8. [Crossref] [PubMed]
- Sasikumar D, Takano Y, Zhao H, et al. Caging and photo-triggered uncaging of singlet oxygen by excited state engineering of electron donor-acceptor-linked molecular sensors. Sci Rep 2022;12:11371. [Crossref] [PubMed]
- Heagerty PJ, Lumley T, Pepe MS. Time-dependent ROC curves for censored survival data and a diagnostic marker. Biometrics 2000;56:337-44. [Crossref] [PubMed]
- Zhu H, Hu H, Hao B, et al. Insights into a Machine Learning-Based Palmitoylation-Related Gene Model for Predicting the Prognosis and Treatment Response of Breast Cancer Patients. Technol Cancer Res Treat 2024;23:15330338241263434. [Crossref] [PubMed]
- Xu J, Yang T, Wu F, et al. A nomogram for predicting prognosis of patients with cervical cerclage. Heliyon 2023;9:e21147. [Crossref] [PubMed]
- Lu H, Xu Z, Shao L, et al. High infiltration of immune cells with lower immune activity mediated the heterogeneity of gastric adenocarcinoma and promoted metastasis. Heliyon 2024;10:e37092. [Crossref] [PubMed]
- Mayakonda A, Lin DC, Assenov Y, et al. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res 2018;28:1747-56. [Crossref] [PubMed]
- Chen J, Liu Z, Wu Z, et al. Identification of a chemoresistance-related prognostic gene signature by comprehensive analysis and experimental validation in pancreatic cancer. Front Oncol 2023;13:1132424. [Crossref] [PubMed]
- Sturm G, Finotello F, Petitprez F, et al. Comprehensive evaluation of transcriptome-based cell-type quantification methods for immuno-oncology. Bioinformatics 2019;35:i436-45. [Crossref] [PubMed]
- Robles-Jimenez LE, Aranda-Aguirre E, Castelan-Ortega OA, et al. Worldwide Traceability of Antibiotic Residues from Livestock in Wastewater and Soil: A Systematic Review. Animals (Basel) 2021;12:60. [Crossref] [PubMed]
- Maeser D, Gruener RF, Huang RS. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform 2021;22:bbab260. [Crossref] [PubMed]
- Shannon P, Markiel A, Ozier O, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res 2003;13:2498-504. [Crossref] [PubMed]
- Chang J, Wu H, Wu J, et al. Constructing a novel mitochondrial-related gene signature for evaluating the tumor immune microenvironment and predicting survival in stomach adenocarcinoma. J Transl Med 2023;21:191. [Crossref] [PubMed]
- Baidoun F, Elshiwy K, Elkeraie Y, et al. Colorectal Cancer Epidemiology: Recent Trends and Impact on Outcomes. Curr Drug Targets 2021;22:998-1009. [Crossref] [PubMed]
- Eng C, Jácome AA, Agarwal R, et al. A comprehensive framework for early-onset colorectal cancer research. Lancet Oncol 2022;23:e116-28. [Crossref] [PubMed]
- Klimeck L, Heisser T, Hoffmeister M, et al. Colorectal cancer: A health and economic problem. Best Pract Res Clin Gastroenterol 2023;66:101839. [Crossref] [PubMed]
- Li J, Ma X, Chakravarti D, et al. Genetic and biological hallmarks of colorectal cancer. Genes Dev 2021;35:787-820. [Crossref] [PubMed]
- Ionescu VA, Gheorghe G, Bacalbasa N, et al. Colorectal Cancer: From Risk Factors to Oncogenesis. Medicina (Kaunas) 2023;59:1646. [Crossref] [PubMed]
- Cañellas-Socias A, Sancho E, Batlle E. Mechanisms of metastatic colorectal cancer. Nat Rev Gastroenterol Hepatol 2024;21:609-25. [Crossref] [PubMed]
- Zheng S, Wang X, Zhao D, et al. Calcium homeostasis and cancer: insights from endoplasmic reticulum-centered organelle communications. Trends Cell Biol 2023;33:312-23. [Crossref] [PubMed]
- Bartkiewicz P, Kunachowicz D, Filipski M, et al. Hypercalcemia in Cancer: Causes, Effects, and Treatment Strategies. Cells 2024;13:1051. [Crossref] [PubMed]
- Wang J, Shi M, Zhang H, et al. PRKCB is relevant to prognosis of lung adenocarcinoma through methylation and immune infiltration. Thorac Cancer 2022;13:1837-49. [Crossref] [PubMed]
- Liu DZ, Luo XZ, Lu CH, et al. Y4 RNA fragments from cardiosphere-derived cells ameliorate diabetic myocardial ischemia-reperfusion injury by inhibiting protein kinase C β-mediated macrophage polarization. Cardiovasc Diabetol 2024;23:202. [Crossref] [PubMed]
- Zhao L, Li Y, Xu T, et al. Dendritic cell-mediated chronic low-grade inflammation is regulated by the RAGE-TLR4-PKCβ(1) signaling pathway in diabetic atherosclerosis. Mol Med 2022;28:4. [Crossref] [PubMed]
- Ren AJ, Wei C, Liu YJ, et al. ZBTB20 Regulates SERCA2a Activity and Myocardial Contractility Through Phospholamban. Circ Res 2024;134:252-65. [Crossref] [PubMed]
- Li J, Li X, Huang H, et al. Role of SERCA3 in the Prognosis and Immune Function in Pan-Cancer. J Oncol 2022;2022:9359879. [Crossref] [PubMed]
- Ahmed A, Iaconisi GN, Di Molfetta D, et al. The Role of Mitochondrial Solute Carriers SLC25 in Cancer Metabolic Reprogramming: Current Insights and Future Perspectives. Int J Mol Sci 2024;26:92. [Crossref] [PubMed]
- Liu Y, Xu Y, Hao Q, et al. SLC25A21 correlates with the prognosis of adult acute myeloid leukemia through inhibiting the growth of leukemia cells via downregulating CXCL8. Cell Death Dis 2024;15:921. [Crossref] [PubMed]
- Tufail M, Jiang CH, Li N. Altered metabolism in cancer: insights into energy pathways and therapeutic targets. Mol Cancer 2024;23:203. [Crossref] [PubMed]
- Ren L, Meng L, Gao J, et al. PHB2 promotes colorectal cancer cell proliferation and tumorigenesis through NDUFS1-mediated oxidative phosphorylation. Cell Death Dis 2023;14:44. [Crossref] [PubMed]
- Wu Z, Ho WS, Lu R. Targeting Mitochondrial Oxidative Phosphorylation in Glioblastoma Therapy. Neuromolecular Med 2022;24:18-22. [Crossref] [PubMed]
- Li YH, Zheng CR, Liu Y, et al. The role of calcium signaling in organotropic metastasis of cancer. Acta Pharmacol Sin 2025;46:1801-12. [Crossref] [PubMed]
- Cui C, Zhang Y, Liu G, et al. Advances in the study of cancer metastasis and calcium signaling as potential therapeutic targets. Explor Target Antitumor Ther 2021;2:266-91. [Crossref] [PubMed]
- Wi YJ, Na SY. Korean J Gastroenterol 2023;82:47-55. [Calcium, Vitamin D, and Colorectal Cancer]. [Crossref] [PubMed]
- Cruz-Pierard SM, Nestares T, Amaro-Gahete FJ. Vitamin D and Calcium as Key Potential Factors Related to Colorectal Cancer Prevention and Treatment: A Systematic Review. Nutrients 2022;14:4934. [Crossref] [PubMed]
- Li R, Zhang G, Tao Q, et al. Revealing the prognostic potential of natural killer cell-related genes in hepatocellular carcinoma: the key role of NRAS. Discov Oncol 2025;16:807. [Crossref] [PubMed]
- Liu X, Li X, Wei H, et al. Mast cells in colorectal cancer tumour progression, angiogenesis, and lymphangiogenesis. Front Immunol 2023;14:1209056. [Crossref] [PubMed]

