Abstract
Background:
Mitochondrial permeability transition-driven necrosis (MPTDN) has been implicated in a variety of diseases, but its relationship with colorectal cancer (CRC) prognosis is unclear.
Methods:
In this study, TCGA-COAD and TCGA-READ datasets were analyzed to identify differentially expressed genes (DEGs) in CRCs. Differentially expressed MPTDNGRs (DE-MPTDNRGs) were identified by comparing DEGs to MPTDN-related genes (MPTDNRGs). Univariate Cox, Minimum Absolute Contraction and Selection Operator (LASSO) regression and risk scoring analysis divided TCGA-CRC patients into high- and low-risk cohorts.
Results:
The prognostic genes LMNB2, CASP7, PRKCB, GZMB and ENDOG were identified, and the survival rate of high-risk patients was poor. Independent prognostic factors, including risk score, age, and N stage, are effective predictors of survival. Immunoassays revealed that high-risk patients had 9 elevated immune checkpoints, while low-risk patients were more susceptible to pazopanib and temsirolimus. In addition, single-cell analysis showed that PRKCB and GZMB were highly expressed in stem cells, while LMNB2 was more abundantly expressed in mast cells. Real-time PCR (RT-qPCR) confirmed low levels of CASP7, PRKCB, and ENDOG mRNA in CRC tissues, with no significant difference between LMNB2 and GZMB.
Conclusion:
These findings highlight 5 MPTDN-associated prognostic genes in CRC, providing insights for individualized treatment and prognosis.
1 Introduction
Colorectal cancer (CRC) is still a dominant reason for tumor incidence and mortality worldwide (). The global burden of colorectal cancer disease on society is expected to increase significantly by 2040, with new cases rising to 3.2 million and deaths rising to 1.6 million (). Current diagnostic methods for CRC, including colonoscopy and fecal immunochemical tests, are effective but have limitations such as invasiveness, cost, and patient compliance issues. New agents targeting KRAS G12C mutation mark a paradigm shift in terms of managing the influenced patients and may drive progress in developing treatments for the more prevalent KRAS mutations (). Novel innovative combinations with poly (ADP-ribose) polymerase (PARP) suppressors can advance the existing curative landscape (). Despite progress in screening and therapeutic strategies guided by biomarker testing for BRAF mutations, deficient Mismatch Repair (dMMR), microsatellite instability (MSI), and HER2 amplification (), the prognosis for advanced CRC is still compromised by suboptimal responses. This highlights the intricate heterogeneity of CRC and underscores the imperative for the development of innovative diagnostic and therapeutic modalities ().
Mitochondrial Permeability Transition-driven necrosis (MPTDN) represents the regulated form of cell death (RCD) pathway that is induced by oxidative stress and disruptions in the cellular environment. It causes the rapid loss of mitochondrial membrane potential and cell membrane rupture (). MPTDN is implicated in diseases like ischemia-reperfusion injury, neurodegenerative diseases, acute kidney injury, and melanoma (). It involves mitochondrial permeability transition pore (mPTP) opening, causing mitochondrial dysfunction and cell death. For instance, Icaritin, a natural compound, inhibits CRC cell proliferation by activating the JNK signaling pathway and inducing mPTP opening, thereby triggering necrosis (). Other therapeutic agents, such as Curcumin, leverage this pathway (), utilizing mPTP blockers like CsA and SfA to attenuate necrosis. Dynamin-related protein 1 (DRP1), an important effector protein in mitochondrial fission, exerts a dual role in the development of CRC (). On one hand, DRP1-mediated mitochondrial fission is instrumental in metabolic reprogramming within CRC cells, which is essential for the regulation of fatty acids and glucose metabolism. This reprogramming not only enhances cell proliferation, invasion, and migration but also contributes to chemoresistance. On the other hand, the activation of DRP1 also triggers carcinogenic signaling pathways, like the Wnt/β-catenin pathway that fosters cell progress and migration in CRC. These findings suggest that MPTDN-related genes (MPTDNRGs) have the potential for acting as biomarkers and therapeutic targets for CRC, offering possible clinical value. Exploring the roles of other MPTDNRGs in CRC could provide novel perceptions of early diagnosis, therapy, and individualized medicine approaches in such disease.
MPTDN is a well-established mechanism implicated in ischemia-reperfusion injury and neurodegenerative disorders. However, its systematic role in the progression of CRC and its potential as a prognostic marker remain largely unexplored. Currently, there is no comprehensive prognostic classifier based on MPTDN-related genes for CRC, and the influence of these genes on the immune microenvironment and therapeutic responsiveness has not been thoroughly investigated. To address these gaps, we conducted an integrative analysis using data from The Cancer Genome Atlas (TCGA), Gene Expression Omnibus (GEO), and single-cell transcriptomic datasets, complemented by in-house quantitative PCR (qPCR) data. This approach enabled us to develop and externally validate the first MPTDN-based risk signature, explore its immunophenotypic implications, and assess its correlation with drug response. By translating MPTDN biology into a practical stratification tool, we offer a mechanistically informed framework for enhancing prognosis and designing combination therapies in CRC.
2 Materials and methods
2.1 Data extraction
RNA-seq data were obtained from the Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov/) for CRC specimens, specifically, the TCGA-COAD and TCGA-READ datasets were combined to create a unified dataset known as TCGA-CRC. A total of 701 samples were gathered from 627 patients, encompassing 650 samples of tumor tissue and 51 samples of normal tissue. Out of these samples, survival information was accessible for 322 patients. The integration method of TCGA-COAD and TCGA-READ datasets was as follows:GDCquery() was executed separately for TCGA-COAD and TCGA-READ to construct query objects, followed by data downloading via the API using GDCdownload(method=“api”). The downloaded files were organized into SummarizedExperiment objects with GDCprepare().The count matrix (RDrr) was extracted from the SummarizedExperiment using assay() (yielding raw counts in a gene × sample format). Gene annotation information (e.g., gene_type, gene_name) was retrieved via rowRanges() and merged with the count matrix, ultimately generating annotated count matrices for COAD and READ respectively. After aligning at the gene level, samples from the two datasets were combined into a unified expression matrix for downstream analyses.For clinical data, GDCquery(data.category=“Clinical”, data.type=“Clinical Supplement”, data.format=“BCR Biotab”) was applied to TCGA-COAD and TCGA-READ individually to obtain clinical supplement information, which was then downloaded and parsed using GDCdownload() and GDCprepare(…, summarizedExperiment=TRUE) to generate clinical tables directly usable for survival analysis. A total of 268 samples with available clinical information were included.
The GSE39582, GSE282542 and GSE132465 datasets were supplied by the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/gds). The GSE39582 dataset comprised 585 tumor tissue specimens, of which 579 samples had survival information. The GSE282542 dataset contained 8 colorectal cancer samples. The GSE132465 dataset comprises 23 tissue samples from patients with primary colorectal cancer (CRC) and 10 normal mucosal samples. A compilation of genes associated with mitochondrial permeability transition (MPT) originated in the Molecular Signatures Database (MSigDB, https://www.gsea-msigdb.org/gsea/msigdb/index.jsp), with the gene sets M17902, M3873, and M16257, which provided a total of 39 MPTDNRGs ().
2.2 Difference analysis
To detect differentially expressed genes (DEGs) in the TCGA-CRC dataset (tumor vs. normal), the DESeq2 package was employed (p-value< 0.05 & |log2FC (fold change)| > 1) (v 1.34.0) (). The P-value was through the Benjamini-Hochberg (BH). The results were visualized via heatmap and volcano plot.
2.3 Function analysis
The DEGs and MPTDNRGs were intersected to obtain DE-MPTDNRGs. To explore the function of DE-MPTDNRGs, gene ontology (GO) and Kyoto Encyclopedia of Genes and Genes (KEGG) pathway enrichment analyses were conducted through the “clusterprofiler” package (p-value< 0.05) (v 4.0.5) (). The P-value was corrected through the BH. The results were visualized via heatmaps and volcano plots.
2.4 Construction and validation of the risk model
Using 322 tumor samples showing survival data form the TCGA-CRC dataset as the training set, univariate Cox regression (HR≠1, p-value< 0.05) was performed for DE-MPTDNRGs to detect prognostic related genes employing coxph function in survival package (v 3.2-13) () Subsequently, genes with P values greater than 0.05 were screened by PH hypothesis test for subsequent analysis. The proportional hazards (PH) assumption was tested using the Grambsch–Therneau test based on scaled Schoenfeld residuals via survival::cox.zph(), with the reporting of test P-values for each covariate as well as the global test (GLOBAL). Additionally, survminer::ggcoxzph() was utilized to generate smoothed curves of Schoenfeld residuals over time with confidence bands for visual diagnosis. Genes with P values greater than 0.05 were screened by PH hypothesis test for subsequent analysis. Subsequently, the glmnet package (v4.1-2) was employed to perform LASSO regression analysis on the identified prognostic genes (x, y, family=“cox”, maxit=5000, K = 10) (). For the Cox model, cv.glmnet defaults to using the deviance of partial likelihood as the cross-validation evaluation metric (type.measure = “deviance”).The formula was as follows:
Where N is the number of genes, e is the expression level, and w is the LASSO coefficient of the prognostic genes (). The risk score calculation formula provides sparse solutions through L1 regularization, making it suitable for high-dimensional biological data. It ensures good predictive performance and generalizability through cross-validation, while simplifying the model for clinical translation applications. Patients were stratified into high and low-risk cohorts in light of the median risk score. Initially, survival differences between these groups were contrasted through Kaplan-Meier (K-M) curves, which were generated with the survminer package (v 0.4.9) (). Simultaneously, the performance of predicting the survival status of patients was evaluated using receiver operating characteristic (ROC) curves produced with the survivalROC package (v 1.0.3) () for different time points. Finally, the constructed risk model was verified through the validation set GSE39582. Additionally, Based on the optimal cut-off values of biomarker expression, we divided the patient samples into high and low expression groups, and analyzed the survival differences between these two groups using the survminer package (v 0.4.9).
2.5 Screening of independent prognostic factors
Initially, the PH assumption test (p > 0.05) was performed for the risk score and clinical characteristics (age, Clinical Stage, T Stage, N Stage, M Stage) in the TCGA-CRC dataset. Subsequently, univariate Cox regression analysis (HR ≠ 1, P< 0.05) was conducted using the survival package (v 3.2-13) to screen for prognostic factors associated with clinical characteristics and the risk score. Finally, the factors passing the univariate Cox analysis underwent the PH assumption test (p > 0.05) and multivariate Cox regression analysis (HR ≠ 1, P< 0.05), leading to the identification of the final independent prognostic factors (IPF).
2.6 Construction of a nomogram
The rms package (v 6.2-0) () was employed to develop a nomogram incorporating IPFs for predicting the survival probability of patients at different years. Calibration curves and ROC were drawn for the performance evaluation of the nomogram via survivalROC package (v 1.0.3) (AUC > 0.6). Additionally, decision curve analysis (DCA) curves were plotted with the rmda package (v 1.6) to assess the clinical utility of the model. To confirm the unique value of this prognostic model, a direct comparison was conducted with previously published models (–)in the TCGA-CRC dataset, with consistency index (Cindex) > 0.6 and AUC > 0.6 as evaluation criteria. Corresponding ROC curves were also plotted (AUC > 0.6).
2.7 Gene set enrichment analysis
As to the training set, to explore potential differences in function between two risk cohorts, DEGs between two risk cohorts were screened with the DESeq2 package (p-value< 0.05 & |log2FC| > 1) (). The P-value was corrected through the BH. Additionally, to explore their potential functions, correlation coefficients between prognostic genes and all genes in the training set were calculated and sequenced. Secondly, in c5.go.bp.v7.2.symbols. gmt and c2.ccp.kegg.v7.2.symbols.gmt as background gene set, clusterProfiler package was used to carry out GSEA enrichment analysis (p.adjust< 0.05) (v 4.0.5) ().
2.8 Processing and quality control of single cell RNA sequencing data
The single-cell transcriptomic sequencing data (GSE282542) was constructed into Seurat objects through the R package “Seurat,”(v 5.1.0) () with cells having less than 200 genes. Moreover, genes covered by less than three cells excluded. The inclusion criteria were specified as follows: nFeature RNA (the number of genes quantified per cell) was required to fall within the range of 200 to 6,000, nCount RNA (the total gene counts per cell) was restricted to less than 8,000, and percentmt (the proportion of mitochondrial gene expression per cell) was restricted to less than 20%. Doublets were removed using the R package “scDblFinder” (v 1.16.0) (). The data was then normalized using the NormalizeData function from the “Seurat” package (v 5.1.0). The parameters were initialized to “LogNormalize” and scale.factor was initialized to 10000. The FindVariableFeatures function was employed to identify highly variable genes based on the correlation between the mean and variance of gene expression, with the top 2,000 highly variable genes selected for further analysis by default (selection.method = “vst”). Results were observed through the LabelPoints function, with the top 10 most variable genes labeled. To eliminate technical variation among samples, the Harmony algorithm was employed for batch effect correction. Batch effects were evaluated using scatter plots after principal component analysis (PCA) dimensionality reduction, and the DimPlot function was utilized to visualize the distribution of different samples in the PCA space. Subsequently, the Harmony algorithm was run for data integration, with sample origin (orig.ident) as the grouping variable and the theta parameter set to 6 to balance the integration intensity.
To further confirm and analyze the cell populations of different cell groups, the R package “Seurat”(v 5.1.0) was applied to standardize the single-cell transcriptome sequencing data (GSE282542) through the ScaleData function. Subsequently, PCA was executed on the highly variable genes of each sample for dimensionality reduction. The Jackstraw function was then used to generate a Jackstraw plot, and the re-clustering was performed using a permutation test algorithm to select the appropriate principal components (100 simulations were repeated to calculate the p-values of the principal components, with a significance threshold set at p< 0.05). The ElbowPlot function was used to create an elbow plot to identify the usable dimensions. Based on the results showing a plateau in the PCA feature number and strong differentiation capability, combined with the analysis from the Jackstraw plot, suitable principal components were selected for further analysis. The FindNeighbors and FindClusters functions in Seurat (v 5.1.0) were used to identify different cell clusters, and uniform manifold approximation and projection (UMAP) was applied to visualization (resolution = 0.2). For cell type annotation within the clusters, the SingleR method (v 1.831) () and manual labeling were used to annotate each cell cluster. Finally, UMAP plots of the annotated cell clusters were generated, and a DotPlot was formulated to show the distribution of prognostic genes across discrepant groups and cell types. Meanwhile, the expression of prognostic genes across different cell types was validated in the GSE132465 dataset. Since the dataset downloaded from the GEO database had already undergone QC, no additional QC was performed ().
2.9 Immune analysis
Via estimating relative subsets of RNA transcripts (CIBERSORT), cell type recognition was utilized to assess immune cell infiltration by evaluating the relative percentage of 22 immune cell types across all samples in the training set (p< 0.05). The P-value was corrected through the BH. Samples with immune infiltration values exceeding 0.05 and immune cell types with no infiltration were excluded. Additionally, the xCell algorithm within the IOBR package (v2.0.0) () was utilized to calculate the infiltration abundances of 44 immune cell types in CRC samples from the high- and low-risk cohorts (p< 0.05). In addition, similar methods were applied to compare immune checkpoint gene expression between the two risk cohorts. Finally, through the tumor immune dysfunction and exclusion (TIDE) website (http://tide.dfci.harvard.edu/) to figure out the TIDE scores for two risk cohorts of patients and analyze the difference (p< 0.05). To assess the microsatellite instability (MSI) status in the training set, the cBioPortalData package (v2.14.0) () was used to download MSI scores. A score of “MSI_SENSOR_SCORE” greater than 0.3 was used to define MSI, while scores equal to or below 0.3 were defined as microsatellite stable (MSS). The differences in risk scores between the MSI and MSS cohorts were comparedTo predict the response to immunotherapy, the cancer immunome atlas the Cancer Imaging Archive (TCIA, https://tcia.at/) was utilized to obtain the immune prognostic score (IPS) of TCGA-CRC patients, and subsequently, the disparities in IPS between two risk cohorts were compared. Spearman correlation was employed to calculate the correlation (|cor| > 0.3, p< 0.05).
2.10 Drug sensitivity analysis
The pRRophetic package (v 0.5) () was utilized to compute the Half Maximal Inhibitory Concentration (IC50) for each tumor sample in the training set. IC50 values were employed to assess patients’ sensitivity to drug therapy (p< 0.05).
2.11 The protein level analysis of key prognostic genes
In order to further construe the potential relevance between prognostic genes and CRC at the protein level, we adopted the Human Protein Atlas (HPA) online database (http://www.proteinatlas.org/) to analyze the protein expression levels of prognostic genes in CRC tissue specimens.
2.12 Ethical approval and consent to participate
Five CRC patients and five healthy individuals’ tissue samples were gathered by the Fuzhou First General Hospital following the principles of the Helsinki Declaration and approved by the Fuzhou First Hospital Ethics Committee (Ethics Review Number: 202312065).
2.13 RT-qPCR
Total RNA was obtained from the samples using TRIzol reagent as per the protocol. After that, RNA concentration was determined through NanoDrop, and mRNA was converted into cDNA by the SweScript First Strand cDNA Synthesis Kit. Following this, SYBR Green reagent, cDNA, and forward and reverse primers were mixed for the RT-qPCR reaction (Table 1). RT-qPCR amplification was executed more than 40 cycles, with each assay conducted in triplicate. The relative levels of gene mRNA were assessed utilizing the 2–△△Ct method in compliance with GAPDH employed as the internal reference. This has been approved by the relevant ethics committee of Fuzhou First General Hospital. The approval number and date of approval are as follows: [202312010] and [Dec. 27th, 2023].
Table 1
| Primer | Sequence |
|---|---|
| LMNB2 F | 5’-CCTGACTCCTATCACGACGC-3’ |
| LMNB2 R | 5’-GGGAGGGTGTAAACGGTGAG-3’ |
| CASP7 F | 5’-CGCAAAGCAACGTCTAGGAG-3’ |
| CASP7 R | 5’-GAATCCTCAACCCCCTGCTC-3’ |
| PRKCB F | 5’-AAATTCACCGCCCGCTTCTT-3’ |
| PRKCB R | 5’-AAGCAGCAAACTTGGCACTG-3’ |
| GZMB F | 5’-GGCAGATGCAGACTTTTCCT-3’ |
| GZMB R | 5’-CTCGTATCAGGAAGCCACCG-3’ |
| ENDOG F | 5’-GCCATGGACGACACGTTCTA-3’ |
| ENDOG R | 5’-CATTTGCCCACAGCACCAAG-3’ |
| GAPDH F | 5’-CGAAGGTGGAGTCAACGGATTT-3’ |
| GAPDH R | 5’-ATGGGTGGAATCATATTGGAAC-3’ |
RT-qPCR primer sequence list.
2.14 Statistical analysis
Data analytics was implemented through R software. Dissimilarities between cohorts were assessed via the Wilcoxon rank-sum test, with a p-value< 0.05 indicating remarkable differentiation between the two cohorts. For RT-qPCR analysis, the t-test was employed for statistical comparisons (p< 0.05). Notably, **** represented p< 0.0001, *** represented p< 0.001, ** represented p< 0.01, * represented p< 0.05, and ns represented p > 0.05.
3 Results
3.1 Biological functions and pathways of DE-MPTDNRGs
There existed 5,635 DEGs identified in TCGA-CRC datasets (tumor vs normal), with 2,986 up-regulated genes and 2,649 down-regulated genes in the normal cohort (Figures 1A, B). The DEGs were intersected with MPTDNRGs to attain 8 DE-MPTDNRGs, comprising 5 up-regulated genes and 3 down-regulated genes (Figure 1C). To explore the function of DE-MPTDNRGs, enrichment analysis was performed. In results of GO, these genes were mainly concentrated in functions related to mitochondria and neurons, such as neuron death and regulation of mitochondrial membrane permeability, etc. (Figure 1D). In addition, these genes were notably involved in Apoptosis and cancer-related pathways, including apoptosis and p53 signaling pathways, etc. (Figure 1E).
Figure 1
3.2 Functional analysis of prognostic genes
Further, according to 8 DE-MPTDNRGs (LMNB2, CASP7, PRKCB, GZMB, ENDOG, BCL2, ENDOG, and BID), genes related to the prognosis of CRC were further explored. In the tumor cohort, lower expression levels of BCL2, CASP7, and PRKCB compared to the control cohort were observed, while higher expression levels of BID, ENDOG, GZMB, and LMNB2 were noted in the tumor cohort compared to the control cohort (Supplementary Figure 1). Firstly, univariate Cox regression screened 5 prognostic related genes, all of which passed the PH hypothesis test (Figure 2A). Next, these 5 genes were further identified as prognostic genes by LASSO, namely LMNB2, CASP7, PRKCB, GZMB and ENDOG (Figures 2B, C). In addition, we performed survival analysis on the five prognostic genes, and the results showed that patients with low expression levels of GZMB, LMNB2, PRKCB, and CASP7 had better survival outcomes, while no significant difference in prognosis was observed for ENDOG (Supplementary Figure 2). To elucidate the potential function of prognostic genes, enrichment analysis was conducted. The findings revealed that these prognostic genes were similar in functions to a large degree. Several functions might be enriched twice, like ribosome, metabolism of xenobiotics by cytochrome P450, and drug metabolism by cytochrome P450, etc. (Figures 2D–H).
Figure 2
3.3 Moderately predictive risk model performance
The risk score for 5 genes was derived using the formula: riskScore = -0.31130104*LMNB2-0.073557728*CASP7-0.06477476*PRKCB-0.04099177*GZMB + 0.15054822*ENDOG. Afterward, patients in the training dataset were stratified into two risk cohorts using the median risk score (-4.632085) as the cutoff. There was a direct relevance between the increase in risk score and the mortality rate (Figure 3A). As proved by the K-M curve analysis, patients in the low-risk cohort had obviously higher survival rates relative to those in the high-risk cohort during the identical period (p = 0.0069) (Figure 3B). Additionally, the area under the curve (AUC) values at 1-, 3-, and 5- years were all above 0.6, demonstrating that the risk score provided moderately predictive capability for patients survival (Figure 3C). Importantly, similar results were obtained in the GSE39582 validation dataset, highlighting the broad applicability and robustness of the risk model (Figures 3D–F).
Figure 3
3.4 Prognostic value and predictive capability of independent prognostic factors
Independent prognostic factors were further screened. The PH assumption test demonstrated that the effects of the risk score and all clinical characteristics on the survival time of colorectal cancer (CRC) patients did not change over time (p > 0.05) (Supplementary Figure 3A). Subsequently, univariate Cox analysis illustrated that risk score (HR = 3.404), age (HR = 3.309), N2 stage (HR = 3.104), M1 stage (HR = 4.75), and stage IV (HR = 12.322) were significantly associated with the overall survival (OS) prognosis of patients (p< 0.05) (Figure 4A). Due to the extremely small sample size of patients with M stage metastasis (only 1 case), this variable was excluded, and the PH assumption test was re-performed on the aforementioned variables. The results confirmed that their effects on CRC patient survival remained stable (p > 0.05) (Supplementary Figure 3B). Multivariate Cox regression deeply retained age (HR = 4.305), risk score (HR = 3.140), and stage IV (HR = 23.231) as IPFs (Figure 4B). A nomogram integrating IPFs was formulated to predict the 1-, 3-, and 5-years survival probabilities of CRC patients (Figure 4C). The calibration curve revealed that the predicted outcomes for 1-, 3-, and 5- years closely aligned with the diagonal line, indicating a strong predictive capability (Figure 4D). Additionally, the AUC values for different years exceeded 0.65, indicating the model’s effectiveness in prediction (Figure 4E). Finally, the decision curve analysis (DCA) curve indicated that the clinical utility of the nomogram was superior to that of individual prognostic factors and lay above the extreme curves (“None” and “All”) (Figure 4F). In summary, these findings confirm the model’s robust performance.
Figure 4
Furthermore, compared with other published models, the 5-gene prognostic signature in this study exhibited the highest concordance index (Cindex = 0.697) and superior prognostic performance for 1-year, 3-year, and 5-year survival of CRC patients (AUC > 0.7) (Supplementary Figures 4A, B).
3.5 Mechanisms and pathways of prognostic correlation
There existed 565 DEGs between the two risk cohorts, among which 440 were up-regulated and 125 down-regulated (Figure 5A). To further explore the underlying mechanisms of this prognostic relationship, enrichment analysis was implemented (Figures 5B–D). The above-mentioned DEGs were mainly enriched in various biological reactions, like oxidative phosphorylation. Moreover, DEGs primarily contained a variety of neurological diseases, such as Alzheimer’s disease and Parkinson’s disease. Moreover, DEGs were also involved in Cell Adhesion and Cell Junction Organization, etc.
Figure 5
3.6 A series of immune analyses for different risk cohorts
As shown in Figures 6A, the proportion of 22 immune cells in 138 patients with p< 0.05 and cell content not all 0 in the training set. In the high-risk cohort, there existed a comparatively higher percentage of resting memory CD4 T cells. Conversely, the proportions of naive B cells, M1 macrophages, and activated memory CD4 T cells were lower (Figure 6B). Further, GZMB was obviously correlated with 4 different immune cells, and resting memory CD4 T cells was greatly correlated with several other genes except ENDOG (Figure 6C). In the xCell algorithm, compared with the low-risk cohort, the high-risk cohort exhibited higher infiltration abundances of CD4+ effector memory T cells (CD4+ Tem) and megakaryocyte-erythroid progenitors (MEP) (Supplementary Figures 5A, B). In contrast, the high-risk cohort showed lower infiltration abundances of CD4+ memory T cells, CD8+ naive T cells, CD8+ T cells, CD8+ central memory T cells (CD8+ Tcm), mast cells, naive B cells, and T helper 2 cells (Th2 cells) (Supplementary Figures 5A, B). Overall, the correlations between these differentially infiltrated immune cells and the prognostic genes were not strong. Mast cells displayed the strongest negative correlation with ENDOG (cor = -0.42, p< 0.001), while CD4+ memory T cells showed the strongest positive correlation with PRKCB (cor = 0.51, p< 0.001) (Supplementary Figure 5C). A total of 9 immune checkpoints showed significant differences (p< 0.05), namely CD200R1, CD276, CD86, ICOS, LAIR1, PDCD1LG2, TIGIT, TNFRSF9, TNFRSF14, and their expression was high in the high-risk cohort (Figure 6D).
Figure 6
Additionally, no significant difference in TIDE values between two risk cohorts; meanwhile, there was a moderate positive correlation (r = 0.36) observed between the risk score and TIDE values (Figures 6E, F). Furthermore, a noteworthy disparity in risk scores was observed between the MSI and MSS cohorts (p = 0.0014), with no significant correlation found between the MSI sensor score and risk score (Figures 6G, H). Subsequently, substantial disparities were identified between the two risk cohorts regarding the CTLA4-PD1-score and CTLA4+PD1-score, but not in the other two scores (Figure 6I).
3.7 Differential drug sensitivity in risk cohorts
In addition, the sensitivity of different risk cohorts to therapeutic drugs was further compared. The results revealed that among the 8 commonly used chemotherapy and molecular targeted drugs, Pazopanib and Temsirolimus had higher IC50 values in the low-risk cohort (Figures 7A–H).
Figure 7
3.8 Results of single-cell analysis
In the GSE282542 dataset, following stringent QC, 38,037 high-quality cells and 29,655 genes were preserved for further analysis (Supplementary Figures 6A, B). The top 2,000 highly variable genes were selected for subsequent analysis, and 10 of the most variable genes, encompassing IGHA2, IGKC, and IGHG1, were emphasized (Supplementary Figure 7A). Thereafter, the PCA was performed on these genes for dimensionality reduction, and the top 20 principal components, based on statistical significance (P< 0.05), were selected (Supplementary Figures 7B–D). A total of 13 cell populations were ultimately identified (Supplementary Figure 7E). Cell annotation was performed via marker gene analysis, identifying 13 distinct cell types: CD4+ activated T cells, tumor epithelial cells, naive/transitional B cells, CD8+ cytotoxic T cells, activated memory B cells, signaling-competent T cells, CD8+ KLR+ cytotoxic T cells, neutrophils, macrophages, plasma cells, fibroblasts, germinal center B cells, and mast cells (Table 2; Figure 8A). Marker genes exhibited highly specific expression in these cell types, confirming the reliability of the annotations (Figure 8B). Additionally, biomarkers were expressed in multiple cell types: specifically, LMNB2, ENDOG, and CASP7 were expressed in tumor epithelial cells, CD4+ activated T cells, CD8+ cytotoxic T cells, and other cell types; PRKCB exhibited relatively high expression in CD4+ activated T cells, activated memory B cells, naive/transitional B cells, and other cell types; while GZMB showed relatively high expression in tumor epithelial cells and other relevant cell types (Figures 8C, D).
Table 2
| Cell type | Cluster | Cell numbers | Marker gene |
|---|---|---|---|
| Epithelial cell | 0 | 5759 | – |
| T helper 2(Th2) cell | 1 | 4685 | SCGB3A1, PMCH |
| B cell | 2 | 3175 | IGHV3-20 |
| Monocyte | 3 | 1892 | HLA-DRB5 |
| Gamma delta T cell | 4 | 1793 | TRGV9, TRGV4, TRDV2 |
| Oogenesis phase fetal germ cell | 5 | 1621 | TAFA1, ADSSa |
| Mitotic arrest phase fetal germ cell | 6 | 1264 | UTS2B |
| NK cell | 7 | 1235 | GZMA, CCL5, NKG7 |
| Plasma cell | 8 | 1112 | IGHG3 |
| SLC16A7+ cell | 9 | 1045 | MOGAT3 |
| Macrophage | 10 | 638 | S100A8, S100A9, MPEG1, S100A12 |
| Mast cell | 11 | 491 | CPA3, HPGDS, TPSB2, GATA2, TPSAB1 |
| Stem cell | 12 | 76 | MSI1 |
The annotation of marker genes in all cell samples.
Figure 8
In the GSE282542 dataset, cells were clustered into 12 distinct cell clusters (Supplementary Figures 8A–F), and 12 cell types were annotated: Tumor Epithelial cell (Epi_Tumor), Plasma cell-Antibody-Secreting Cell (Plasma_ASC), Enteric glial cell-Schwann cell (Enteric_glia_Schwann), Activated CD4+ T cell-Th17-like T cell (CD4_Th17_like_Act), Fibroblast-SFRP2+ Fibroblast (Fibroblast_SFRP2), Cytotoxic T cell-Natural Killer cell-Gamma delta T cell (Cytotoxic_T_NK_gd), Inflammatory Monocyte-Macrophage-IL1+ Cell (MonoMac_inflam_IL1), Naive B cell (B_Naive), Mast cell (Mast), Mural cell-Pericyte-Vascular Smooth Muscle Cell (Mural_Pericyte_VSMC), Endothelial cell (Endothelial), and plasmacytoid dendritic cell (pDC) (Figure 8E). Bubble plots showed the highly expressed marker genes in each cell type (Figure 8F). The prognostic genes ENDOG, CASP7, and LMNB2 were highly expressed in tumor epithelial cells; GZMB was highly expressed in Cytotoxic_T_NK_gd and pDC cells; and PRKCB was highly expressed in B cells and pDC cells (Figures 8G, H). In summary, the expression patterns of the prognostic genes across cell types were generally consistent between the two single-cell datasets.
3.9 Validation of prognostic genes in clinical samples
To construe the prognostic genes at the protein level, we employed the online database of the Human Protein Atlas (HPA). We construed the protein expression levels of the prognostic genes in CRC tissue specimens. Findings showed that, relative to normal colon tissue samples, LMNB2 was highly expressed in CRC tissue samples; CASP7 and PRKCB possessed lower expressions in CRC tissue specimens; GZMB had low expressions in both normal colon tissue and CRC tissue samples; and ENDOG possessed the elevated expressions in both healthy colon tissue samples and CRC tissue samples (Figures 9A–J). Clinical samples were collected to examine the expression of prognostic genes in CRC tissues and controls. RT-qPCR results implied that, in the CRC cohort, the mRNA expression levels of CASP7, PRKCB, and ENDOG were prominently reduced, relative to the control cohort. Nonetheless, no statistically significant differences in mRNA expression levels of LMNB2 and GZMB between the two cohorts were found (Figures 10A–E).
Figure 9
Figure 10
4 Discussion
CRC represents a major global health concern, standing as the second in terms of mortality and third in incidence among all cancers worldwide (). The debilitating symptoms of CRC and the side effects of its treatment severely compromise patients’ quality of life (). Despite advances in screening and therapy, the prognosis for CRC remains poor, particularly in advanced stages, underscoring the urgent need for improved diagnostic and prognostic tools. In the area of CRC, the focus has primarily been on the importance of mitochondrial dynamics in cancer cell metabolism and survival. Wu et al. highlighted the influence of mitochondrial fusion-fission dynamics on CRC biology (), suggesting a dual role in CRC development. Our study has developed a potentially useful risk model in view of 5 DE-MPTDNRGs, including LMNB2, CASP7, PRKCB, ENDOG, and GZMB, which effectively stratified patients into low- and high-risk cohorts with remarkably different mortality rates and survival outcomes. The risk model exhibited good predictive capability and was further validated in an independent dataset, demonstrating its universal applicability. Moreover, functional and immune analyses revealed potential mechanisms and therapeutic implications of DE-MPTDNRGs in CRC.
CASP7 functions as a regulatory hinge between MPTDN and apoptosis in colorectal cancer. Pan-cancer transcriptomes uniformly show lower CASP7 mRNA in tumors versus matched mucosa (), and reduced expression correlates with a 13% decrease in 5-year survival); this aligns with our findings. At single-cell resolution we mapped CASP7 to malignant epithelia and CD8+ cytotoxic T cells, revealing a lineage-restricted control axis. Post-transcriptionally, TRIM25 docks on the CASP7 3′-UTR and recruits hnRNPH1 to accelerate mRNA decay, blunting chemotherapy-induced apoptosis (). Conversely, SREBP1 blockade restores CASP7 and potentiates gemcitabine lethality (), whereas miR-519a-3p directly targets CASP7 and dampens death signaling (). In ischemia–reperfusion models the MPT couples energy status to cell fate: glycolytic ATP rebound redirects cells from necrosis to caspase-3-mediated apoptosis, a switch abolished by cyclosporin A or caspase-3 inhibitors (). Together, these data establish CASP7 as the cross-point where MPTDN and apoptotic circuits merge to govern chemosensitivity and prognosis in CRC.
Consistent with our findings, PRKCB (Protein Kinase C Beta) is transcriptionally down-regulated in CRC tissues, and its higher expression is significantly associated with better overall survival (). Single-cell RNA-seq data implied a dual role in sculpting the tumor immune micro-environment and maintaining cancer stem-cell properties. A similar inverse correlation between PRKCB promoter hyper-methylation and reduced mRNA abundance has been documented in gastric cancer, where methylation of cg08406370 and cg00735962 CpG islands silences PRKCB and predicts poor prognosis (). Comparable epigenetic repression may therefore underlie the low PRKCB levels observed in our CRC cohort. Extending beyond gastrointestinal malignancies, a large retrospective NSCLC study demonstrated that robust PRKCB expression independently confers favorable prognosis (HR = 0.47, P< 0.001) and negatively correlates with advanced TNM stage (). At the mitochondiral metabolic level, recent work in ischemia-resistant hippocampal neurons has positioned PKCβII as a critical regulator of mitochondrial glutaminase activity and subsequent anaplerotic flux through the TCA cycle (); whether PRKCB-mediated phosphorylation events similarly modulate mitochondrial substrate utilization—glutamine versus pyruvate oxidation—in CRC stem or immune cells remains an attractive, yet untested, mechanism that could link its prognostic impact to metabolic plasticity.
ENDOG is an endonuclease situated in the mitochondrial intermembrane space, is instrumental in orchestrating DNA fragmentation and apoptosis through its translocation to the nucleus, and it also modulates mitochondrial genome cleavage by relocating to the mitochondrial matrix (). In ischemic cardiomyocytes, Bnip3-induced mitochondrial permeability transition drives EndoG release from mitochondria and subsequent caspase-independent DNA fragmentation, establishing a direct mitochondrial-nuclear death axis that operates independently of caspases (). Studies have indicated that ENDOG can be phosphorylated by GSK3B, which strengthens its interaction with YWHAG. This process results in the release of TSC2 and PIK3C3 from YWHAG, subsequent suppression of the MTOR signaling pathway, and initiation of autophagy. Additionally, the endonuclease activity of ENDOG is vital for activating the DNA damage response and thus inducing autophagy (, ). Previous studies have revealed that compounds like amentoflavone () and the LL-37-derived peptide FK-16 () can induce apoptosis, possibly through mechanisms involving ENDOG. Additionally, increased expression of ENDOG in cancers () further supports its potential role in cancer apoptosis. Research has also found that ENDOG deficiency or silencing inhibits proliferation in certain cancer cell lines (e.g., endometrial tumors and glioblastomas), especially when PTEN is missing or p-AKT is highly expressed in these tumor cells, suggesting that low ENDOG expression may limit tumor growth (). However, studies on colorectal cancer indicate that substituted diaryl diselenides show complex cytotoxic and apoptotic effects on colon adenocarcinoma cells (), suggesting that ENDOG’s role in colorectal cancer may not be entirely beneficial and requires further study. Although TCGA data show ENDOG up-regulation in CRC, our RT-qPCR cohort revealed significant down-regulation, and survival analysis detected no stratification between high- and low-ENDOG groups, highlighting the need for expanded samples and functional validation to resolve this discrepancy.
LMNB2 (Lamin B2) is an integral nuclear-membrane protein that orchestrates chromatin remodeling and the rupture–re-assembly of the nuclear envelope during mitosis, thereby governing eukaryotic cell proliferation and tumorigenesis (). It can silence the tumor-suppressor gene p21, creating a cellular environment conducive to colorectal cancer (CRC) advancement (). In hepatocellular carcinoma, high LMNB2 expression is tightly linked to immune-evasive signatures within the tumor microenvironment, implying a broader immunomodulatory role (). Recent work shows that LMNB2 fosters cancer-cell stemness and the Warburg effect via activation of the p38 MAPK cascade—a signaling node that cross-talks with mitochondrial permeability transition and oxidative-stress responses—thus providing an indirect but biologically plausible connection between LMNB2 and mitochondrial functionality (). Additionally, a regulatory circuitry involving long non-coding RNAs and microRNAs (e.g., miR-326) fine-tunes LMNB2 abundance, influencing CRC progression (). Consistent with these findings, our analyses identified significant LMNB2 up-regulation in CRC tissues, and high expression was associated with shorter survival.
Granzyme B (GZMB), a key constituent of cytotoxic granules in natural killer (NK) cells and cytotoxic T lymphocytes (CTLs), is crucial for the immune system’s defense against cancer (). Overexpression of GZMB is often linked to poor prognosis in colorectal cancer (CRC) patients and associated with more aggressive CRC phenotypes (). However, some studies suggest that high GZMB expression may improve survival by enhancing antitumor immune responses ().This dual role indicates that GZMB’s prognostic impact may depend on tumor stage, molecular features, and the immune microenvironment (). Additionally, low GZMB expression may be related to immune evasion in colon cancer, influencing NK cell cytotoxicity and cytokine signaling pathways, therefore promoting tumor progression (). In CRC, mutations in genes like KRAS and BRAF are well-documented. KRAS mutations are closely related to cancer cell proliferation, survival, and metastasis, and can lead to treatment resistance (, , ). In right-sided colon cancer, the BRAF (V600E) mutation is often associated with highly aggressive tumors and poor prognosis (). Similar to KRAS mutations, BRAF mutations also frequently result in treatment resistance, particularly to anti-EGFR monoclonal antibody therapy (, ).
Unlike KRAS and BRAF mutations, GZMB is mainly related to immune responses and has a dual role in the tumor microenvironment. For one thing, the elevated GZMB expression may promote tumor immune evasion and progression. For another, it may enhance antitumor immune responses. Therefore, GZMB’s function is closely related to the tumor immune environment and molecular characteristics. Although GZMB overexpression may serve as a negative prognostic biomarker in CRC, precisely explaining its clinical importance needs a thorough investigation into its interactions with the immune microenvironment.
We conducted enrichment analysis to detect the possible underlying mechanism of MPTDNRGs in CRC. Pathways between different risk groups were enriched in biological processes such as oxidative phosphorylation, which is pivotal for cellular energy metabolism and has been implicated in cancer progression and chemoresistance (). Prior research unveiled that isobavachalcone induced non-apoptotic necrosis among cancer cells by triggering mitochondrial permeability transition, oxidative stress, and the opening of the mitochondrial permeability transition pore, bringing about mitochondrial unbalance and finally cell death (). As claimed by Zhao et al., NPC-26, a new mitochondrion-interfering compound, triggers CRC cell death by activating AMP-activated protein kinase (AMPK) signaling, which is linked to mitochondrial permeability transition pore (mPTP) opening and reactive oxygen species (ROS) production (). Furthermore, the enrichment in neurodegenerative disease pathways is consistent with previous research (). Although these pathways are traditionally associated with neurological disorders, emerging evidence suggests that the molecular mechanisms involved, like protein misfolding and mitochondrial dysfunction, may also play roles in cancer pathogenesis and resistance to therapy (). Such intersection of neurodegenerative and oncogenic pathways could provide novel insights into CRC biology and identify potential biomarkers for disease progression and therapeutic response.
This risk model not only achieved prognostic stratification but also revealed distinctly different tumor immune microenvironments underlying the outcomes. The high-risk group, associated with unfavorable survival, featured a dysfunctional tumor immune microenvironment (TIME): characterized by accumulation of immunosuppressive cells (elevated resting memory CD4 T cells) (), reduction of anti-tumor effector cells (decreased naive B cells, M1 macrophages, and activated memory CD4 T cells), and widespread upregulation of immune checkpoint molecules on both tumor and immune cells (such as elevated PDCD1LG2/PD-L2). These coordinated changes collectively indicated an adaptive immune resistance state (, ). Notably, the alteration in M1 macrophage infiltration was statistically significant only in the CIBERSORT algorithm but not in the xCell algorithm. Nevertheless, M1 macrophages still exhibited a higher proportion in the low-risk group, which may be attributed to sample heterogeneity and functional exhaustion of immune cells. In conjunction with the immunosuppressive tumor microenvironment observed in the high-risk group, the reduced abundance or impaired function of M1 macrophages likely contributes to enhanced tumor immune evasion, consistent with the unfavorable prognosis observed in high-risk patients. Additionally, the higher proportion of MSS patients in the high-risk group aligned with the typically lower immunogenicity and reduced sensitivity to immunotherapy characteristic of MSS tumors, further reinforcing this group’s “cold” or suppressed immune phenotype (). Interestingly, while the TIDE score did not significantly differentiate between groups, the risk score itself positively correlated with TIDE, and the high-risk group showed distinct patterns in IPS predicting differential responses to combination immunotherapy (anti-CTLA-4/PD-1). This suggested that the gene signatures captured by our risk model were intrinsically linked to tumor immune evasion mechanisms. Consequently, the poor prognosis in high-risk patients likely resulted from a microenvironment unfavorable for effective anti-tumor immune responses, collectively constructed by MSS background, suppressive immune cell infiltration, and checkpoint upregulation (). Correspondingly, the low-risk group exhibited favorable immunological characteristics, offering valuable insights for developing personalized immunotherapy strategies in colorectal cancer (CRC). Notably, the high-risk group demonstrated higher predicted sensitivity to Pazopanib. Although currently approved for renal cell carcinoma and advanced soft tissue sarcoma (), this agent has shown potential in CRC by inducing tumor regression when combined with the FOLFOX regimen (), and its safety and tolerability have been established in multiple clinical studies (, ). We hypothesize that the therapeutic efficacy of Pazopanib may be mediated by the prognostic genes identified in this study. However, it is important to emphasize that this recommendation constitutes an in silico prediction rather than a definitive therapeutic guideline. Consequently, its clinical utility requires rigorous validation through functional experiments (e.g., cell lines, organoids, or PDX models) and further clinical investigation.
Our study has several limitations. First, the RT-qPCR validation was constrained by a small sample size. Owing to the restricted sample-use clause in the Ethics Review, we were unable to obtain additional clinical specimens across centers; consequently, the largest feasible cohort was analyzed under the current ethical framework. Although these data provided an initial biological anchor for the in-silico findings, the statistical power is inherently limited, and the observed expression trends require cautious interpretation. Therefore, prior to clinical application, multi-level validation (at both mRNA and protein levels) with an expanded sample size is required. Furthermore, the generalizability of the prognostic gene signature and its potential as a biomarker need to be rigorously validated in large-scale, multi-center prospective cohorts. Second, the study remains predominantly computational; mechanistic insights into how the hub genes modulate MPTDN await functional validation. We propose prioritized experiments—including CRISPR knockout/rescue in CRC organoids, mitochondrial stress tests with high-resolution respirometry, and patient-derived xenograft (PDX) models—to dissect gene-specific effects on oxidative phosphorylation and MPTP opening. Third, drug-sensitivity predictions for pazopanib and temsirolimus were inferred from in-silico IC50 profiling; these findings need prospective verification in 3-D culture systems and PDX models to rule out false-positive biomarker associations. Finally, the current risk model explained a moderate proportion of survival variance; integrating Consensus Molecular Subtypes (CMS) and Immunoscore variables could increment model discrimination and calibration, but this strategy must be formally evaluated through nested cross-validation in external cohorts.
5 Conclusions
This study successfully identified differentially expressed genes associated with mitochondrial permeability transition in CRC and build a strong prognostic risk model. Functional enrichment analysis revealed significant involvement of these genes in various biological processes and signaling pathways. The risk model demonstrated strong predictive performance for patient survival, and its independent prognostic value was confirmed. Furthermore, immune analysis and drug sensitivity assessments highlighted significant differences between high- and low-risk groups, providing insights into potential therapeutic strategies. These findings underscore the potential of mitochondrial permeability transition-related genes as biomarkers and therapeutic targets in CRC. Future research should focus on validating these results through experimental and clinical studies to fully realize their potential in improving patient outcomes.
Statements
Data availability statement
The data presented in the study are deposited in the figshare repository, accession number 10.6084/m9.figshare.31422392.
Ethics statement
The studies involving humans were approved by the Fuzhou First Hospital Ethics Committee (Ethics Review Number: 202312010). The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.
Author contributions
CH: Conceptualization, Data curation, Validation, Visualization, Writing – original draft, Writing – review & editing. DX: Data curation, Validation, Visualization, Writing – review & editing. BL: Visualization, Writing – review & editing. KZ: Conceptualization, Project administration, Supervision, Writing – review & editing.
Funding
The author(s) declared that financial support was received for this work and/or its publication. The research reported in this project was generously supported by Startup Fund for scientific research, Fujian Medical University (Grant number:2023QH1281) and Key Clinical Specialty Discipline Construction Program of Fuzhou, Fujian.P.R.C (Grant number:20220301), and Funding for the Visiting Scholar Program of Fuzhou Municipal Health Commission (Grant number:2025 NO.3). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Acknowledgments
We would like to express our sincere gratitude to all individuals and organizations who supported and assisted us throughout this research. In conclusion, we extend our thanks to everyone who has supported and assisted us along the way. Without your support, this research would not have been possible.
Conflict of interest
The author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declared that generative AI was not used in the creation of this manuscript.
Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2026.1679360/full#supplementary-material
Supplementary Figure 1Expression profiles of 8 de-mptdnrgs.
Supplementary Figure 2Survival curves of the five prognostic genes.
Supplementary Figure 3Proportional Hazards (PH) Assumption Test Results for Independent Prognostic Analysis. (A) PH assumption test results of univariate Cox regression analysis for the risk score and clinical factors. The x-axis (Time) represents the survival time of patients. The y-axis (Beta (t) for riskScore) denotes the time-dependent changes in the regression coefficient corresponding to the risk score, reflecting the temporal dependency of the risk score’s impact intensity on survival risk. The regression coefficient curve in the figure shows no obvious upward or downward trend over the entire follow-up period and fluctuates around the value of 0, which intuitively confirms that the risk effect of the risk score does not change over time, further supporting the statistical conclusion that the proportional hazards assumption is satisfied. The proportional hazards assumption is considered to be met when the P-value is greater than 0.05. (B) PH assumption test of multivariate Cox regression analysis for the factors identified by univariate Cox regression analysis.
Supplementary Figure 4Comparison of the prognostic model with published prognostic models. (A) Comparative analysis results of prognostic performance between the MPTDN-5 model (5 prognostic genes) constructed in this study and the published models. Rows include the prognostic gene model established in this study and four cell death-related subtypes, namely mitophagy, cuproptosis, necroptosis, and pyroptosis. Columns consist of five prognostic performance evaluation metrics, namely Concordance Index (Cindex), 1-year area under the curve (AUC_1y), 3-year area under the curve (AUC_3y), and 5-year area under the curve (AUC_5y). The numbers in the table are raw performance values, and the color of each column corresponds to the z-score normalized values. Published refers to previously reported prognostic models of the same category. (B) Time-dependent receiver operating characteristic (ROC) curves showing the performance of the MPTDN-5 model in predicting 1-year, 3-year, and 5-year survival outcomes of colorectal cancer (CRC) patients across four cell death-related subtypes (cuproptosis, mitophagy, pyroptosis, necroptosis). The x-axis (1−Specificity) represents the false positive rate, indicating the probability that the model incorrectly classifies surviving patients as deceased. The y-axis (Sensitivity) represents the true positive rate, indicating the probability that the model correctly classifies deceased patients as deceased. A larger AUC value indicates a stronger discriminative ability of the model.
Supplementary Figure 5Immune infiltration analysis results via the xCell algorithm. (A) Infiltration proportions of 44 immune cell types in samples from the high-risk and low-risk cohorts. (B) Differences in the proportions of 44 immune cell types between the high-risk and low-risk cohorts. ns: p > 0.05; *: p< 0.05; **: p< 0.01; ***: p< 0.001; ****: p< 0.0001. (C) Correlations between differentially infiltrated immune cells and prognostic genes. The color indicates the correlation coefficient, with yellow representing positive correlation and blue representing negative correlation. Asterisks denote statistical significance: *: p< 0.05; **: p< 0.01; ***: p< 0.001.
Supplementary Figure 6(A, B) nFeature_RNA, nCount_RNA, and percent.mt before and after quality control.
Supplementary Figure 7(A) Acquisition of HVGs. (B-D) Selection of usable PCs. (E) Identification of cell types.
Supplementary Figure 8Dimensionality reduction and clustering of the GSE132465 dataset. (A) Violin plots of the GSE132465 dataset after quality control (QC), including the number of genes detected per cell (Number of Genes), the Unique Molecular Identifier (UMI) count per cell, and the percentage of mitochondrial gene expression per cell (Percentage of MT). The x-axis represents samples. (B) Scatter plot of the top 2,000 highly variable genes (HVGs). The x-axis denotes the average expression level of genes (Average Expression), and the y-axis denotes the standardized variance of genes (Standardized Variance). (C) Principal Component Analysis (PCA) scatter plot of single-cell data before batch effect correction, where each scatter represents a single cell. (D) Jackstraw plot for statistical significance testing of Principal Components (PCs). The x-axis shows the theoretical P-value distribution, and the y-axis shows the actual P-value of each principal component. Each data point in the plot represents one principal component. (E) Elbow plot of PCA. The x-axis represents the serial number of PC sorted by variance explanation ability, and the y-axis represents the standard deviation corresponding to each principal component (Standard Deviation). (F) Uniform Manifold Approximation and Projection (UMAP) plot of cell clustering.
Abbreviations
MPTDN, Mitochondrial permeability transition-driven necrosis; MPTDNRGs, Mitochondrial permeability transition driven necrosis-related genes; CRC, Colorectal cancer; DEGs, Differentially expressed genes; LASSO, Least absolute shrinkage and selection operator; RT-qPCR, Real-time fluorescence quantitative polymerase chain reaction; PARP, Poly (ADP-ribose) polymerase; mPTP, Mitochondrial permeability transition pore; DRP1, Dynamin-related protein 1; DE-MPTDNRGs, Differentially expressed MPTDNGRs; TCGA, The cancer genome atlas; MPT, Mitochondrial permeability transition; GO, Gene ontology; KEGG, Kyoto encyclopedia of genes and genes; K-M, Kaplan-meier; GSEA, Gene set enrichment analysis; MSI, Microsatellite instability; MSS, Microsatellite stable; IPS, Immune prognostic score; lncRNAs, Long non-coding RNAs; GZMB, Granzyme B; CTLs, Cytotoxic T lymphocytes; NK, Natural killer; AMPK, AMP-activated protein kinase; ROS, Reactive oxygen species; mPTP, Mitochondrial permeability transition pore; TIME, Tumor immune microenvironment; DN, Diabetic Nephropathy; YSXZF, Yi-Shen-Xie-Zhuo Formula; TCM, Traditional Chinese Medicine; PDL2, Programmed cell death-ligand 2; PD-1, Programmed cell death-1; PD-L1, Programmed cell death-ligand 1.
References
1
XiYXuP. Global colorectal cancer burden in 2020 and projections to 2040. Transl Oncol. (2021) 14:101174. doi: 10.1016/j.tranon.2021.101174
2
GBD 2019 Colorectal Cancer Collaborators. Global, regional, and national burden of colorectal cancer and its risk factors, 1990-2019: a systematic analysis for the Global Burden of Disease Study 2019. Lancet Gastroenterol Hepatol. (2022) 7:627–47. doi: 10.1016/S2468-1253(22)00044-9
3
HongDSFakihMGStricklerJHDesaiJDurmGAShapiroGIet al. KRAS(G12C) inhibition with sotorasib in advanced solid tumors. N Engl J Med. (2020) 383:1207–17. doi: 10.1056/NEJMoa1917239
4
Salva de TorresCBaraibarISaoudi GonzalezNRosJSalvaFRodriguez-CastellsMet al. Current and emerging treatment paradigms in colorectal cancer: integrating hallmarks of cancer. Int J Mol Sci. (2024) 25:6967. doi: 10.3390/ijms25136967
5
Saoudi GonzalezNSalvaFRosJBaraibarIRodriguez-CastellsMGarciaAet al. Unravelling the complexity of colorectal cancer: heterogeneity, clonal evolution, and clinical implications. Cancers (Basel). (2023) 15:4020. doi: 10.3390/cancers15164020
6
AleseOBWuCChapinWJUlanjaMBZheng-LinBAmankwahMet al. Update on emerging therapies for advanced colorectal cancer. Am Soc Clin Oncol Educ Book. (2023) 43:e389574. doi: 10.1200/EDBK_389574
7
GalluzziLVitaleIAaronsonSAAbramsJMAdamDAgostinisPet al. Molecular mechanisms of cell death: recommendations of the Nomenclature Committee on Cell Death 2018. Cell Death Differ. (2018) 25:486–541. doi: 10.1038/s41418-017-0012-4
8
BonoraMGiorgiCPintonP. Molecular mechanisms and consequences of mitochondrial permeability transition. Nat Rev Mol Cell Biol. (2022) 23:266–85. doi: 10.1038/s41580-021-00433-y
9
ZhouCChenZLuXWuHYangQXuD. Icaritin activates JNK-dependent mPTP necrosis pathway in colorectal cancer cells. Tumour Biol. (2016) 37:3135–44. doi: 10.1007/s13277-015-4134-3
10
QiuYYuTWangWPanKShiDSunH. Curcumin-induced melanoma cell death is associated with mitochondrial permeability transition pore (mPTP) opening. Biochem Biophys Res Commun. (2014) 448:15–21. doi: 10.1016/j.bbrc.2014.04.024
11
WuZXiaoCLiFHuangWYouFLiX. Mitochondrial fusion-fission dynamics and its involvement in colorectal cancer. Mol Oncol. (2024) 18:1058–75. doi: 10.1002/1878-0261.13578
12
LiuJZhangMSunQQinXGaoTXuYet al. Construction of a novel MPT-driven necrosis-related lncRNAs signature for prognosis prediction in laryngeal squamous cell carcinoma. Environ Sci pollut Res Int. (2023) 30:77210–25. doi: 10.1007/s11356-023-26996-1
13
LoveMIHuberWAndersS. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. (2014) 15:550. doi: 10.1186/s13059-014-0550-8
14
YuGWangLGHanYHeQY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. (2012) 16:284–7. doi: 10.1089/omi.2011.0118
15
RamsayISMaSFisherMLoewyRLRaglandJDNiendamTet al. Model selection and prediction of outcomes in recent onset schizophrenia patients who undergo cognitive training. Schizophr Res Cogn. (2018) 11:1–5. doi: 10.1016/j.scog.2017.10.001
16
LiYLuFYinY. Applying logistic LASSO regression for the diagnosis of atypical Crohn’s disease. Sci Rep. (2022) 12:11340. doi: 10.1038/s41598-022-15609-5
17
LiuLChenZShiWLiuHPangW. Breast cancer survival prediction using seven prognostic biomarker genes. Oncol Lett. (2019) 18:2907–16. doi: 10.3892/ol.2019.10635
18
LiuTTLiRHuoCLiJPYaoJJiXLet al. Identification of CDK2-related immune forecast model and ceRNA in lung adenocarcinoma, a pan-cancer analysis. Front Cell Dev Biol. (2021) 9:682002. doi: 10.3389/fcell.2021.682002
19
HeagertyPJLumleyTPepeMS. Time-dependent ROC curves for censored survival data and a diagnostic marker. Biometrics. (2000) 56:337–44. doi: 10.1111/j.0006-341X.2000.00337.x
20
SachsMC. plotROC: A tool for plotting ROC curves. J Stat Softw. (2017) 79:2. doi: 10.18637/jss.v079.c02
21
WangJTaoZWangBXieYWangYLiBet al. Cuproptosis-related risk score predicts prognosis and characterizes the tumor microenvironment in colon adenocarcinoma. Front Oncol. (2023) 13:1152681. doi: 10.3389/fonc.2023.1152681
22
TangYGuoHChenLWangXChenQGouLet al. Development and validation of a prognostic model for mitophagy-related genes in colon adenocarcinoma: A study based on TCGA and GEO databases. PloS One. (2023) 18:e0284089. doi: 10.1371/journal.pone.0284089
23
ZhuLHYangJZhangYFYanLLinWRLiuWQ. Identification and validation of a pyroptosis-related prognostic model for colorectal cancer based on bulk and single-cell RNA sequencing data. World J Clin Oncol. (2024) 15:329–55. doi: 10.5306/wjco.v15.i2.329
24
LuoBLinJCaiWWangM. Identification of the pyroptosis-related gene signature and risk score model for colon adenocarcinoma. Front Genet. (2021) 12:771847. doi: 10.3389/fgene.2021.771847
25
WangYHuangYLiCWangXYangMXuDet al. Establishment of a necroptosis related genes signature to predict prognosis and therapeutic response in colon cancer. Front Cell Dev Biol. (2022) 10:921320. doi: 10.3389/fcell.2022.921320
26
SatijaRFarrellJAGennertDSchierAFRegevA. Spatial reconstruction of single-cell gene expression data. Nat Biotechnol. (2015) 33:495–502. doi: 10.1038/nbt.3192
27
GermainPLLunAGarcia MeixideCMacnairWRobinsonMD. Doublet identification in single-cell sequencing data using scDblFinder. F1000Res. (2021) 10:979. doi: 10.12688/f1000research.73600.1
28
LeeHOHongYEtliogluHEChoYBPomellaVVan den BoschBet al. Lineage-dependent gene expression programs influence the immune landscape of colorectal cancer. Nat Genet. (2020) 52:594–603. doi: 10.1038/s41588-020-0636-z
29
ZengDYeZShenRYuGWuJXiongYet al. IOBR: multi-omics immuno-oncology biological research to decode tumor microenvironment and signatures. Front Immunol. (2021) 12:687975. doi: 10.3389/fimmu.2021.687975
30
CeramiEGaoJDogrusozUGrossBESumerSOAksoyBAet al. The cBio cancer genomics portal: an open platform for exploring multidimensional cancer genomics data. Cancer Discov. (2012) 2:401–4. doi: 10.1158/2159-8290.CD-12-0095
31
JiangHZYangBJiangYLLiuXChenDLLongFXet al. Development and validation of prognostic models for colon adenocarcinoma based on combined immune-and metabolism-related genes. Front Oncol. (2022) 12:1025397. doi: 10.3389/fonc.2022.1025397
32
BrayFLaversanneMSungHFerlayJSiegelRLSoerjomataramIet al. Global cancer statistics 2022: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. (2024) 74:229–63. doi: 10.3322/caac.21834
33
SiegelRLMillerKDJemalA. Cancer statistics, 2020. CA Cancer J Clin. (2020) 70:7–30. doi: 10.3322/caac.21590
34
WuZXiaoCLongJHuangWYouFLiX. Mitochondrial dynamics and colorectal cancer biology: mechanisms and potential targets. Cell Commun Signal. (2024) 22:91. doi: 10.1186/s12964-024-01490-4
35
NasrullahUStankeKRecknagelVBozkurtSWurzelPGauerSet al. The E3 Ligase TRIM25 Impairs Apoptotic Cell Death in Colon Carcinoma Cells via Destabilization of Caspase-7 mRNA: A Possible Role of hnRNPH1. Cells. (2023) 12:201. doi: 10.3390/cells12010201
36
ShenWXuTChenDTanX. Targeting SREBP1 chemosensitizes colorectal cancer cells to gemcitabine by caspase-7 upregulation. Bioengineered. (2019) 10:459–68. doi: 10.1080/21655979.2019.1676485
37
QinYZhuYLuoFChenCChenXWuM. Killing two birds with one stone: dual blockade of integrin and FGF signaling through targeting syndecan-4 in postoperative capsular opacification. Cell Death Dis. (2017) 8:e2920–e. doi: 10.1038/cddis.2017.315
38
KimJ. Mitochondrial permeability transition in the switch from necrotic to apoptotic cell death in ischemic rat hepatocytes. Gastroenterology. (2003) 124:494–503. doi: 10.1053/gast.2003.50059
39
DongZLiuZLiangMPanJLinMLinHet al. Identification of circRNA–miRNA–mRNA networks contributes to explore underlying pathogenesis and therapy strategy of gastric cancer. J Trans Med. (2021) 19:226. doi: 10.1186/s12967-021-02903-5
40
LiLFeiXWangHChenSXuXKeHet al. Genome-wide DNA methylation profiling reveals a novel hypermethylated biomarker PRKCB in gastric cancer. Sci Rep. (2024) 14:26605. doi: 10.1038/s41598-024-78135-6
41
HuGDuJWangBSongPLiuS. Comprehensive analysis of the clinical and prognostic significance of SFRP1 and PRKCB expression in non-small cell lung cancer: a retrospective analysis. Eur J Cancer Prev. (2024) 33:45–52. doi: 10.1097/CEJ.0000000000000832
42
Beręsewicz-HallerMKrupskaOBochomulskiPDudzikDChęcińskaAHilgierWet al. Mitochondrial metabolism behind region-specific resistance to ischemia-reperfusion injury in gerbil hippocampus. Role of PKCβII and phosphate-activated glutaminase. Int J Mol Sci. (2021) 22:8504. doi: 10.3390/ijms22168504
43
WangWLiJZhouQ. The biological function of cytoplasm-translocated ENDOG (endonuclease G). Autophagy. (2024) 20:445–7. doi: 10.1080/15548627.2023.2271750
44
AgoulnikIZhangJYeJAltafajACardonaMBahiNet al. EndoG links bnip3-induced mitochondrial damage and caspase-independent DNA fragmentation in ischemic cardiomyocytes. PloS One. (2011) 6:e17998. doi: 10.1371/journal.pone.0017998
45
LiJWangWZhouQ. Conserved role of ENDOG in promoting autophagy. Autophagy. (2021) 17:1061–2. doi: 10.1080/15548627.2021.1907513
46
LinCHLinKHKuHJLeeKCLinSSHsuFT. Amentoflavone induces caspase-dependent/-independent apoptosis and dysregulates cyclin-dependent kinase-mediated cell cycle in colorectal cancer in vitro and in vivo. Environ Toxicol. (2023) 38:1078–89. doi: 10.1002/tox.23749
47
RenSXShenJChengASLuLChanRLLiZJet al. FK-16 derived from the anticancer peptide LL-37 induces caspase-independent apoptosis and autophagic cell death in colon cancer cells. PloS One. (2013) 8:e63641. doi: 10.1371/journal.pone.0063641
48
YooNJJeongEGKimMSAhnCHKimSSLeeSH. Increased expression of endonuclease G in gastric and colorectal carcinomas. Tumori. (2008) 94:351–5. doi: 10.1177/030089160809400311
49
NedelFCamposVFAlvesDMcBrideAJDellagostinOACollaresTet al. Substituted diaryl diselenides: cytotoxic and apoptotic effect in human colon adenocarcinoma cells. Life Sci. (2012) 91:345–52. doi: 10.1016/j.lfs.2012.07.023
50
DongCHJiangTYinHSongHZhangYGengHet al. LMNB2 promotes the progression of colorectal cancer by silencing p21 expression. Cell Death Dis. (2021) 12:331. doi: 10.1038/s41419-021-03602-1
51
QiuYWangZZhangXHuangPZhangWZhangKet al. A long-acting isomer of Ac-SDKP attenuates pulmonary fibrosis through SRPK1-mediated PI3K/AKT and Smad2 pathway inhibition. IUBMB Life. (2020) 72:2611–26. doi: 10.1002/iub.2389
52
DukelMKalabalikZA. Selective chemical and genetic inhibition of PKCβII strongly sensitizes colon cancer cells to anoikis. Mol Biol Rep. (2025) 52:534. doi: 10.1007/s11033-025-10612-1
53
MuWGuoLLiuYYangHNingSLvG. Long noncoding RNA SNHG1 regulates LMNB2 expression by sponging miR-326 and promotes cancer growth in hepatocellular carcinoma. Front Oncol. (2021) 11:784067. doi: 10.3389/fonc.2021.784067
54
JingZLLiuGLZhouNXuDYFengNLeiYet al. Interferon-γ in the tumor microenvironment promotes the expression of B7H4 in colorectal cancer cells, thereby inhibiting cytotoxic T cells. Sci Rep. (2024) 14:6053. doi: 10.1038/s41598-024-56681-3
55
WangYLiuFDuXShiJYuRLiSet al. Combination of anti-PD-1 and electroacupuncture induces a potent antitumor immune response in microsatellite-stable colorectal cancer. Cancer Immunol Res. (2024) 12:26–35. doi: 10.1158/2326-6066.CIR-23-0309
56
NingYLinKFangJChenXHuXLiuLet al. Pyroptosis-related signature predicts the progression of ulcerative colitis and colitis-associated colorectal cancer as well as the anti-TNF therapeutic response. J Immunol Res. (2023) 2023:7040113. doi: 10.1155/2023/7040113
57
FangYZhanX. Identification of biomarkers associated with the prognoses of colorectal cancer patients. Digestion. (2023) 104:148–62. doi: 10.1159/000528084
58
FengJHuZXiaXLiuXLianZWangHet al. Feedback activation of EGFR/wild-type RAS signaling axis limits KRAS(G12D) inhibitor efficacy in KRAS(G12D)-mutated colorectal cancer. Oncogene. (2023) 42:1620–33. doi: 10.1038/s41388-023-02676-9
59
ZhuQZhangRGuXZhaoZGaoQChenMet al. Honokiol enhances the sensitivity of cetuximab in KRAS(G13D) mutant colorectal cancer through destroying SNX3-retromer complex. Theranostics. (2024) 14:5443–60. doi: 10.7150/thno.97180
60
LiuMLiuQHuKDongYSunXZouZet al. Colorectal cancer with BRAF V600E mutation: Trends in immune checkpoint inhibitor treatment. Crit Rev Oncol Hematol. (2024) 204:104497. doi: 10.1016/j.critrevonc.2024.104497
61
MaZQiZGuCYangZDingYZhouYet al. BRAF(V600E) mutation promoted the growth and chemoresistance of colorectal cancer. Am J Cancer Res. (2023) 13:1486–97.
62
Ruiz-SaenzAAtreyaCEWangCPanBDreyerCABrunenDet al. A reversible SRC-relayed COX2 inflammatory program drives resistance to BRAF and EGFR inhibition in BRAF(V600E) colorectal tumors. Nat Cancer. (2023) 4:240–56. doi: 10.1038/s43018-022-00508-5
63
ZhaoZMeiYWangZHeW. The effect of oxidative phosphorylation on cancer drug resistance. Cancers (Basel). (2022) 15:62. doi: 10.3390/cancers15010062
64
XiaoQZhongBHouYWangMGuoBLinLet al. Fighting cancer by triggering non-canonical mitochondrial permeability transition-driven necrosis through reactive oxygen species induction. Free Radic Biol Med. (2023) 202:35–45. doi: 10.1016/j.freeradbiomed.2023.03.020
65
ZhaoZFengLWangJChengDLiuMLingMet al. NPC-26 kills human colorectal cancer cells via activating AMPK signaling. Oncotarget. (2017) 8:18312–21. doi: 10.18632/oncotarget.15436
66
EfimovaEVTakahashiSShamsiNAWuDLabayEUlanovskayaOAet al. Linking cancer metabolism to DNA repair and accelerated senescence. Mol Cancer Res. (2016) 14:173–84. doi: 10.1158/1541-7786.MCR-15-0263
67
HanYLiuSShiSShuYLuCGuX. Screening of genes associated with immune infiltration of discoid lupus erythematosus based on weighted gene co-expression network analysis. Biochem Genet. (2025) 63:465–82. doi: 10.1007/s10528-023-10603-6
68
RibasA. Adaptive immune resistance: how cancer protects from immune attack. Cancer Discov. (2015) 5:915–9. doi: 10.1158/2159-8290.CD-15-0563
69
KimTKVandsembENHerbstRSChenL. Adaptive immune resistance at the tumour site: mechanisms and therapeutic opportunities. Nat Rev Drug Discov. (2022) 21:529–40. doi: 10.1038/s41573-022-00493-5
70
HuangXSunTWangJHongXChenHYanTet al. Metformin reprograms tryptophan metabolism to stimulate CD8+ T-cell function in colorectal cancer. Cancer Res. (2023) 83:2358–71. doi: 10.1158/0008-5472.CAN-22-3042
71
GaoZLingXShiCWangYLinA. Tumor immune checkpoints and their associated inhibitors. J Zhejiang Univ Sci B. (2022) 23:823–43. doi: 10.1631/jzus.B2200195
72
MiyamotoSKakutaniSSatoYHanashiAKinoshitaYIshikawaA. Drug review: pazopanib. Jpn J Clin Oncol. (2018) 48:503–13. doi: 10.1093/jjco/hyy053
73
ZhuGZhaoMHanQTanYSunYBouvetMet al. Temozolomide and pazopanib combined with FOLFOX regressed a primary colorectal cancer in a patient-derived orthotopic xenograft mouse model. Transl Oncol. (2020) 13:100739. doi: 10.1016/j.tranon.2019.12.011
74
BradyJCorriePChauIDigumartiRAdamsLMBotbylJet al. An open-label study of the safety and tolerability of pazopanib in combination with FOLFOX6 or CapeOx in patients with colorectal cancer. Invest New Drugs. (2013) 31:1228–35. doi: 10.1007/s10637-013-9938-7
75
BennounaJDeslandresMSenellartHde LabareyreCRuiz-SotoRWixonCet al. A phase I open-label study of the safety, tolerability, and pharmacokinetics of pazopanib in combination with irinotecan and cetuximab for relapsed or refractory metastatic colorectal cancer. Invest New Drugs. (2015) 33:138–47. doi: 10.1007/s10637-014-0142-1
Summary
Keywords
colorectal cancer, mitochondrial permeability transition driven necrosis, prognostic genes, risk model, single-cell analysis
Citation
Huang C, Xie D, Lin B and Zhang K (2026) The potential of mitochondrial permeability transition-driven necrosis-related genes in prognostic evaluation of colorectal cancer patients. Front. Oncol. 16:1679360. doi: 10.3389/fonc.2026.1679360
Received
28 October 2025
Revised
07 February 2026
Accepted
16 February 2026
Published
09 March 2026
Volume
16 - 2026
Edited by
Betul Karademir Yilmaz, Marmara University, Türkiye
Reviewed by
Qianxue Wu, First Affiliated Hospital of Chongqing Medical University, China
Ali Sahin, Department of Pediatric Allergy and Immunology, Türkiye
Updates
Copyright
© 2026 Huang, Xie, Lin and Zhang.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Kun Zhang, zhangkun73@aliyun.com
† These authors have contributed equally to this work
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.