Skip to main content
Download PDF

Identification of an at-risk subpopulation with high immune infiltration based on the peroxisome pathway and TIM3 in colorectal cancer

BMC Cancer volume 22, Article number: 44 (2022) Cite this article

Abstract

Background

Peroxisomes are pivotal metabolic organelles that exist in almost all eukaryote cells. A reduction in numbers and enzymatic activities of peroxisomes was found in colon adenocarcinomas. However, the role of peroxisomes or the peroxisome pathway in colorectal cancer (CRC) is not defined.

Methods

In the current study, a peroxisome score was calculated to indicate the activity of the peroxisome pathway using gene set variant analysis based on transcriptomic datasets. CIBERSORTx was chosen to infer enriched immune cells for tumors among subgroups. The SubMap algorithm was applied to predict its sensitivity to immunotherapy.

Results

The patients with a relatively low peroxisome score and high level of T-cell immunoglobulin and mucin domain 3 (TIM-3) presented the worse overall survival than others. Moreover, low peroxisome scores were associated with high infiltration of lymphocytes and poor prognosis in those CRC patients. Thus, a PERLowTIM3High CRC risk subpopulation was identified and characterized by high immune infiltration. The results also showed that CD8 T cells and macrophages highly infiltrated tumors of the PERLowTIM3High group, regardless of consortium molecular subtype and microsatellite instability status. This subgroup had the highest tumor mutational burden and overexpression of immune checkpoint genes. Further, the PERLowTIM3High group showed a higher probability of responding to programmed cell death protein-1-based immunotherapy. In addition, genes involved in peroxisomal metabolic processes in CRC were also investigated since peroxisome is a rather pleiotropic and highly metabolic organelle in cell. The results indicated that only those genes involved in fatty acid alpha oxidation could be used to stratify CRC patients as similar as peroxisome pathway genes.

Conclusions

We revealed the favorable prognostic value of the peroxisome pathway in CRC and provided a new CRC stratification based on peroxisomes and TIM3, which might be helpful for CRC diagnostics and personalized treatment.

Peer Review reports

Background

Colorectal cancer (CRC), one of the most prevalent malignancies, is the second leading cause of cancer-related death worldwide [1]. Although the use of chemotherapy and targeted therapy has dramatically prolonged survival time for patients with unresectable CRCs, the efficacy of these therapies can be limited by drug resistance and side effects [2]. In the past few years, a better knowledge of the complex interactions between tumor cells and the immune system has led to the establishment of novel immunotherapies [3, 4]. However, only a small subset of metastatic CRCs with the microsatellite instability (MSI) phenotype benefit from immune checkpoint inhibitors (ICIs) that target co-inhibitory receptors [5, 6]. One major immunotherapeutic obstacle of CRC is the immunosuppressive tumor microenvironment (TME), a highly complex and heterogenous contexture [7, 8]. For CRC tumors with pre-existing strong immune infiltration, restoring the suppressive effects of the TME on immune cells is a potential strategy to improve the efficacy of immunotherapy [9]. Thus, TME-based classifications are valuable for personalized immunotherapy in CRC.

Peroxisomes are membrane-bound organelles first discovered by Johannes Rhodin in 1954. Peroxisomes are involved in multiple metabolic processes, including the metabolism of fatty acid oxidation and D-amino acids, biosynthesis of ether lipids and bile acids, and turnover of reactive oxygen species (ROS) [10]. Apart from being a highly metabolic organelle, the peroxisome has emerged as a pivotal regulator of inflammation and immune responses in these few years [11,12,13]. Furthermore, a recent study suggested that enhanced peroxisome signaling might be involved in reprograming the lipid metabolism of T cells in urological cancers [14]. The attention on the role of peroxisomes or peroxisomal lipid metabolism in cancer, a disease characterized by metabolic reprogramming, is continuously increasing. From earlier publications, an overall reduction in peroxisomal proteins and enzymatic activities has been found in neoplastic tissue, including that of the colon [15,16,17]. Moreover, a more recent study suggested that NUDT7, a peroxisomal gene, is a potent tumor suppressor during KRASG12D-driven CRC development [18]. Although these publications pointed to a potential loss of peroxisome function during CRC tumorigenesis, our current view on the role of peroxisomes or the peroxisome pathway in CRC remains fragmentary and limited.

In the current study, the prognostic role of peroxisome pathway was evaluated in CRC based on the genomic data of public datasets. In addition, an at-risk CRC subpopulation was identified and characterized by massive immune infiltration. Then, we explored the relationship between the newly characterized subpopulation and immunotherapy efficiency. Finally, we further investigated whether the genes specifically involved in peroxisomal metabolic processes can stratify CRC patients individually.

Methods

Data collection and processing

This study used two public patient datasets for discovery and validation. The CRC cohort from the Cancer Genome Atlas (TCGA) was used for discovery. We downloaded the RNA-Seq expression data (FPKM values) of 593 primary colon or rectum tumors with overall survival (OS) data using the GDCquery of the R package TCGAbiolinks [19]. We then transformed the FPKM values into transcripts per kilobase million (TPM) values for subsequent analysis. We used a consensus measurement of tumor purity to remove samples with low tumor purity (< 60%) [20]. We updated the clinical information from TCGA Pan-Cancer Clinical Data Resource [21]. A total of 504 TCGA CRC patients enrolled in this study possessed gene somatic mutation data (MAF files). We analyzed and summarized the mutation data using maftools [22]. A cohort from the Gene Expression Omnibus (GEO), GSE39582 (n = 556), was used for validation. We downloaded the normalized expression data and clinical information for the GSE39582 dataset from the GEO repository. We obtained the MSI status and consortium molecular subtype (CMS) classifications from the Colorectal Cancer Subtyping Consortium synapse data portal (https://www.synapse.org/#!Synapse:syn2623706/files/). We further combined MSI-low and microsatellite stable (MSS) cases and defined them as MSS and defined MSI-high cases as MSI for this analysis.

Gene set expression analyses

Gene set variant analysis (GSVA) is a computational approach that has been utilized to measure the overall activity of biological pathways [23,24,25]. With a pathway-based approach, the expression of pathway member genes can be summarized thus reducing complexity, and yielding more reliable results [26]. We applied it to the transcriptomic data samples to estimate the HALLMARK_PEROXISOME score from the Molecular Signatures Database (MSigDB) Hallmark collection using the GSVA R package (version 3.10) [23]. The median value of the GSVA HALLMARK_PEROXISOME score was used as the cut-off to categorize the CRC patients into either a peroxisome-high (Per-High) or peroxisome-low (Per-Low) group. We used the gene set enrichment analysis (GSEA) method to investigate the relationship between HALLMARK_PEROXISOME and other Hallmark gene sets using the GSEABase package (http://www.bioconductor.org/packages/release/bioc/html/GSEABase.html) [27]. Gene set related to β-oxidation of very-long-chain fatty acids was obtained from Reactome Pathway database (https://reactome.org/) [28]. Gene sets associated with fatty acid alpha oxidation and ROS were from Gene Ontology (http://geneontology.org/) [29, 30]. Gene sets related to biosynthesis of bile acids and ether lipid metabolism were from Kyoto Encyclopedia of Genes and Genomes (KEGG) database (https://www.kegg.jp/kegg/kegg1.html) [31, 32].

Estimation of immune infiltration

The approach utilized for estimating immune infiltration was CIBERSORTx (https://cibersortx.stanford.edu/), which can predict the abundance and proportion of 22 types of tumor-infiltrating lymphocytes based on the expression of 547 genes for each tissue sample [33]. We additionally applied the ESTIMATE (Estimation of STromal and Immune cells in MAlignant Tumours using Expression data) algorithm to assess the enrichment of immune and stromal cells in the TME and to infer purity [34].

Prediction of the immunotherapy response for each subgroup

To predict the sensitivity of each subgroup to immunotherapy, we used Class mapping analysis (SubMap from Gene Pattern, https://www.genepattern.org/) to compare the similarity of the gene expression profiles of our subgroups to those of melanoma patients treated with checkpoint blockades against programmed cell death protein-1 (PD1) or cytotoxic T lymphocyte antigen-4 (CTLA4) [35].

Statistical analysis

We applied the optimal cut-off for T-cell immunoglobulin and mucin domain 3 (TIM3) expression with regard to the associated hazard of death events in a survfit model based on the cutp function of the survMisc R package (version 0.5.5; https://CRAN.R-project.org/package=survMisc). We chose the highest log-rank test score as the optimal cut-off for classifying patients into high- and low-expression groups with different risks. Kaplan–Meier analysis was used to estimate OS. We conducted multivariate Cox proportional hazards regression models to investigate the prognostic value of group III using the survival R package (version 2.42.3; https://CRAN.R-project.org/package=survival). To compare the contingency table variables, we used the chi-square test or Fisher’s exact test. We used the Mann–Whitney test for two-group comparison and the Kruskal–Wallis test for three-group comparisons. All statistical results with a p-value < 0.05 were considered significant.

Results

A low peroxisome pathway score is associated with a worse clinical outcome and high immune cell infiltration in CRC patients

The peroxisome score was calculated by GSVA for CRC patients enrolled in the study and patients were categorized into either a Per-High or Per-Low group by the median of peroxisome scores. Intriguingly, the Per-Low group showed significantly worse OS than the Per-High group in TCGA CRC dataset (p = 0.024; Fig. 1a), findings that were validated in the GSE39582 dataset (p = 0.017; Fig. 1b). We next performed GSEA on 50 hallmark gene sets in TCGA cohort and found lipid-associated gene sets, CHOLESTEROL_HOMEOSTASIS and BILE_ACID_METABOLISM, which were significantly enriched in the Per-High group (Fig. 1c), which was validated in the GSE39582 cohort (Additional file 2: Fig. S1a). Furthermore, we found that immune-related gene sets were highly enriched in the Per-Low group, both in TCGA discovery dataset and the GEO-GSE39582 validation dataset (Fig. 1c; Additional file 2: Fig. S1a; Additional file 1: Table S1–2).

Fig. 1
figure 1

Peroxisome-Low CRC patients had poor clinical outcomes despite high immune infiltrates. a, b Kaplan–Meier analysis of OS in TCGA colorectal cancer cohort (n = 593) and NCBI GEO-GSE39582 cohort (n = 556). P values were calculated by the log-rank test. c The volcano plot shows the enrichment analysis using hallmark gene sets in TCGA cohort. The red gene sets were enriched in the Per-High group, whereas the blue gene sets were enriched in the Per-Low group. d Correlation between the peroxisome score versus immune score, stromal score, and tumor purity in TCGA cohort. e Violin plots of the immune score, stromal score, and tumor purity from ESTIMATE of Per-High and Per-Low groups in TCGA cohort. f Violin plot of total immune infiltrates (sum of absolute scores across 22 immune cell types) among the CRC subgroups in TCGA. g Correlation between peroxisome scores and total immune infiltrates of patients in TCGA cohort. For violin plots, p values in group comparisons with Mann–Whitney U–test are shown. For panels d, Pearson’s rho (r) and statistical difference (p) are indicated. ****P < 0.0001

Full size image

To further investigate the association between the peroxisome pathway and TME in CRC, we applied the ESTIMATE algorithm to TCGA cohort. We observed that the peroxisome pathway activity was negatively correlated with immune and stromal cell infiltration but was positively correlated with tumor purity (Fig. 1d). The Per-Low group had significantly higher immune and stromal scores but a lower level of tumor purity than the Per-High group (Fig. 1e). We also used a deconvolution method, CIBERSORTx, to assess the total immune infiltrates. The Per-Low group demonstrated higher immune cell infiltration than the Per-High group (Fig. 1f). A negative correlation between peroxisome scores and total immune infiltrates was observed (Fig. 1g). These results were validated with the GSE39582 dataset (Additional file 2: Fig. S1b-d). Taken together, our findings revealed that patients with low peroxisome scores had poor clinical outcomes despite having a high degree of immune infiltrates.

Identification and validation of an at-risk CRC subpopulation stratified by peroxisome score and TIM3 expression

To further investigate immunosuppressive mechanisms that could contribute to the poor prognosis observed in the Per-Low group, we analyzed the expression of immune checkpoints (ICs), as these genes are associated with the failure of T cells to efficiently eliminate tumor cells. The Per-Low group showed higher levels of all inhibitory checkpoint genes that we examined (Additional file 2: Fig. S2), including PDCD1 (PD1), PDL1, CTLA4, LAG3, TIGIT, TIM3, OX40, GITR, TNFRSF9 (4-1BB), ICOS, IDO1, CD40, CD70, and CD27. We observed that the peroxisome score was negatively correlated with the respective expression levels of these checkpoint genes in TCGA cohort (Fig. 2a). Among these ICs, TIM3 had the highest negative correlation with the peroxisome score. Additionally, we used an optimal cut-off for TIM3 expression to classify CRC patients into high and low TIM3 groups, followed by OS analysis. We found that for TCGA colorectal tumors, a high level of TIM3 expression was associated with a poor prognosis (Additional file 2: Fig. S3), which is consistent with the findings of previous studies [36, 37]. Intriguingly, Figs. 2b-e show that TIM3 prognostic behaviors were dependent on peroxisome scores. The median peroxisome score was used to classify patients into the Per-High and Per-Low groups. We applied an optimal cut-off based on the log-rank test for each independent group and subsequently defined four groups as follows: group I (PERHighTIM3High), group II (PERHighTIM3Low), group III (PERLowTIM3High), and group IV (PERLowTIM3Low) (Fig. 2b and c). In the Per-Low group, we found that TIM3 served as a prognostic marker of poor outcome; however, this was not seen in the Per-High group (Fig. 2d and e). Group III (PERLowTIM3High), with the highest proportion of MSI-status tumors, showed the most unfavorable outcomes compared with those in group IV and the Per-High group (group I + II), based on an analysis of OS (Fig. 3a and c). In addition, we adopted the same method to validate the existence of risk group III in the GSE39582 validation dataset. Group III consistently had the worst OS, as expected (Fig. 3b and d). Our results showed that group III patients had the worst OS; therefore, it would be of interest to explore differences in clinical and molecular characteristics across three subpopulations.

Fig. 2
figure 2

Identification of an at-risk CRC subpopulation stratified by peroxisome score and TIM3 expression. a Pairwise correlations of the peroxisome score versus expression of individual immune checkpoint genes in TCGA cohort. Pearson’s rho (r) and statistical difference (p) are shown in each graph. b, c Scatter plots of peroxisome score and log2-transformed TIM3 gene expression values are shown for TCGA dataset. MSI (blue diamond) and MSS (black circle) status are labeled for CRC patients. The median value of TIM3 expression is indicated in a gray dashed line. The log-rank test score patterns and candidate cut-offs (blue and red arrows) for high peroxisome and low peroxisome groups are shown in panels b and c, respectively. d, e Kaplan–Meier curves show the OS for the optimal cut-off of TIM3 expression in TCGA CRC subgroups with high peroxisome score and low peroxisome score, respectively. Using the median value of peroxisome score and these two candidate cut-offs, we separated patients into four different groups labeled with their population fraction in percentage (group I-IV). Kaplan–Meier curves are plotted for the risk groups in high peroxisome and low peroxisome groups, respectively. P-value was calculated by the log-rank test and shown for each plot

Full size image
Fig. 3
figure 3

Validation of the risk CRC population using NCBI-GEO GSE39582 dataset. a, b Scatter plots of peroxisome score and log2-transformed TIM3 gene expression values are shown for TCGA and NCBI-GEO GSE39582 cohort. Three CRC subgroups were indicated for each plot (group I + II: high peroxisome group; group III: PERLowTIM3High group; group IV: PERLowTIM3Low group). c, d For OS analysis, Kaplan–Meier curves are plotted for the risk subgroups in TCGA and GSE39582 datasets. The log-rank p values are shown

Full size image

We gathered clinical and molecular data of 1149 CRC patients from TCGA and GSE39582. Characteristics of all patients are shown in Table 1. Advanced age and stages, a higher mutation frequency of BRAF, a CpG-island methylator phenotype (CIMP)-high status, and the MSI phenotype were associated with group III. Furthermore, we investigated the proportions of different CMS subtypes across subgroups and found that group III largely overlapped with CMS1 (MSI immune subtype; strong immune infiltration and activation) and CMS4 (mesenchymal; TGF-β activation; stromal invasion and angiogenesis).

Table 1 Clinical and molecular characteristics of the CRC subgroups stratified by peroxisome score and TIM3 expression
Full size table

To assess the independent poor prognostic value of group III, we performed univariate and multivariate cox regression analysis to investigate the association between traditional clinical and molecular characteristics and OS for the entire cohort. We found three independent indicators for poor OS, stage IV (Stage IV: HR = 6.07, 95% CI: 3.25–11.35, p < 0.001; Stage I as the reference), mutated KRAS (HR = 1.29, 95% CI: 1.00–1.66, p = 0.048) and group III (Group III: HR = 1.41, 95% CI: 1.05–1.90, p = 0.023; Group I + II as the reference). Whereas group III (PERLowTIM3High) had the highest proportions of CMS1 and CMS4 subtypes, of interest, our multivariable OS model showed that group III, but not CMS4, was an independent indicator for predicting poor OS (Table 2).

Table 2 Univariate and multivariate analyses for overall survival in CRC
Full size table

The PERLowTIM3High group has high levels of tumor infiltrated CD8+ T cells and macrophages

To explore the immune heterogeneity across subgroups, we used CIBERSORTx to estimate total immune infiltrates. Group III had the highest enrichment level of immune-infiltrated cells (Fig. 4a; Additional file 2: Fig. S4a). We then used the ESTIMATE-based approach to infer tumor purity; compared with group I + II and IV, group III had the lowest tumor purity (Fig. 4b; Additional file 2: Fig. S4b).

Fig. 4
figure 4

TCGA group III (PERLowTIM3High) tumors were highly infiltrated with CD8 T cells and macrophages. a Violin plot showing the total immune infiltrates of CIBERSORTx for each subgroup in TCGA CRC dataset. b Violin plot showing the ESTIMATE tumor purity for each subgroup in TCGA CRC dataset. c Boxplots showing enrichment levels of CD8 T cells and macrophages for each subgroup in TCGA. d, g Boxplots of enrichment level of CD8 T cells and macrophages for TCGA MSS tumors. e, h Boxplots of enrichment level of CD8 T cells and macrophages for TCGA CMS1 and CMS4 tumors. f, i Boxplots of enrichment level of CD8 T cells and macrophages for TCGA CMS2, CMS3, and indeterminate tumors. ***P < 0.001

Full size image

To further investigate the proportions of immune infiltrates in group III, we aggregated CIBERSORTx results of relative cell-types into six major cell types (B cells, CD4 T cells, CD8 T cells, dendritic cells, macrophages, and neutrophils). We observed that CD8 T cells and macrophages were highly enriched in group III both in TCGA and GSE39582 datasets (Fig. 4c; Additional file 2: Fig. S4c). CRC MSS tumors tend to have poor immunogenicity and lack responses to ICIs [38]; our results showed that group III MSS tumors were also highly infiltrated by CD8 T cells and macrophages (Fig. 4d and g; Additional file 2: Fig. S4d and S4g). As mentioned, group III had high proportions of CMS1 and CMS4, and these types were related to high immune infiltration. Next, we compared tumor infiltrating levels of CD8 T cells and macrophages between CMS1 + CMS4 patients among the three subgroups; CMS1 + CMS4 tumors in group III had the highest enrichment levels (Fig. 4e and h; Additional file 2: Fig. S4e and S4h). Moreover, other CMS subtypes (CMS2, CMS3, and indeterminate) tumors in group III had the highest infiltration of CD8 T cells and macrophages among the three subgroups (Fig. 4f and i; Additional file 2: Fig. S4f and S4i). Together, these results indicate that group III tumors were highly infiltrated by CD8 T cells and macrophages, which was not a simple reflection of high proportions of MSI tumors or CMS1/CMS4 phenotypes.

The PERLowTIM3High group has the highest tumor mutational burden (TMB) and is likely to respond to immunotherapy

To further explore whether these groups have different mutational profiles, we analyzed somatic mutation data from TCGA. Figure 5a shows the 20 most frequently mutated genes for each subgroup. We observed that the mutational load increased significantly in this group (Fig. 5b).

Fig. 5
figure 5

Comparison of tumor mutational burden and possible sensitivity to immunotherapy for CRC subgroups. a Top 20 genes most frequently mutated in TCGA CRC patients in groups I + II, III, and IV. b Comparison of the tumor mutation load across the three CRC subgroups in TCGA. c Plots show immune checkpoint genes of the three CRC subgroups in TCGA. d The same as panel c, but for the GSE39582 cohort. e Heatmaps show the correlation between transcriptomic expression patterns of melanoma patients receiving PD1 or CTLA4 inhibitors and CRC subgroups in TCGA. f The same as panel e, but for the GSE39582 dataset. ***P < 0.001

Full size image

Given that group III was associated with the highest immune infiltration level and TMB, we next sought to assess the possibility of efficient immunotherapy for the three subgroups. We found that group III had the highest expression of ICs in both TCGA and GSE39582 cohorts (Fig. 5c-d), including PDL1, TIM3, LAG3, TIGIT, and VISTA. We then used subclass mapping to measure the similarity of expression profiles between the three CRC subgroups and a public cohort for melanoma patients treated with checkpoint blockades against CTLA4 or PD-1 [35]. When comparing the expression profiles of the TCGA CRC subgroups and this external melanoma cohort, we observed that group III had a significant correlation with the PD-1 response group (Bonferroni corrected p = 0.012; Fig. 5e), a finding that was validated in the GSE39582 dataset (Bonferroni corrected p = 0.012; Fig. 5f). These results suggested that group III patients were likely to benefit from immunotherapy and specifically PD1-based strategies.

Genes involved in fatty acid α-oxidation stratified CRC patients similarly with the peroxisome pathway

Since peroxisomes are highly dynamic and pleiotropic metabolic organelles, we further investigated whether the genes specifically involved in peroxisomal metabolic processes can stratify CRC patients individually. Thus, we attempted to stratify TCGA CRC patients by using gene sets associated with specific peroxisome metabolic functions including REACTOME_BETA_OXIDATION_OF_VERY_LONG_CHAIN_FATTY_ACIDS (BOVLCFA), GOBP_FATTY_ACID_ALPHA_OXIDATION (FAAO), KEGG_PRIMARY_BILE_ACID_BIOSYNTHESIS (PBAB), KEGG_ETHER_LIPID_METABOLISM (ELM) and GOBP_CELLULAR_RESPONSE_TO_REACTIVE_OXYGEN_SPECIES (CRTROS) (Additional file 2: Fig. S5). For TCGA CRC patients with a relatively low BOVLCFA score or FAAO score, high TIM3 expression indicated worse OS (Additional file 2: Fig. S5a and S5b). However, our multivariable OS models showed that FAAOLowTIM3High, but not BOVLCFALowTIM3High, was an independent predictor for poor OS (Additional file 1: Table S3–4). Thus, among these five gene sets, FAAO was found as the only gene set that could stratify CRC patients similarly as HALLMARK_PEROXISOME (Additional file 2: Fig. S5b); TCGA CRC patients with a relatively low FAAO score and high TIM3 expression (FAAOLowTIM3High) had the worst OS, which was validated in GSE39582 (Fig. 6a and b). Similarly, FAAOLowTIM3High group tumors had the lowest tumor purity and were highly infiltrated by immune cells, especially CD8 T cells and macrophages (Fig. 6c-e). We also observed that in TCGA, FAAOLowTIM3High group tumors had significantly increased TMB and expression of ICs (Fig. 6f and g). Then, we applied the SubMap to compare the expression profiles of TCGA CRC subgroups stratified by FAAO and TIM3 with the external melanoma cohort mentioned previously herein. We found that the FAAOLowTIM3High group had a significant correlation with the PD-1 response group, which was validated in GSE39582 (Fig. 6h and i). Taken together, these results showed that genes involved in FAAO could be used to stratify CRC patients as similar as peroxisome pathway genes; moreover, tumors in the PERLowTIM3High group and the FAAOLowTIM3High groups displayed some similar characteristics.

Fig. 6
figure 6

Genes involved in FAAO can stratify CRC patients similarly with genes in the peroxisome pathway. a Kaplan–Meier curves are plotted for the risk subgroups based on FAAO score and TIM3 expression in TCGA. b The same as a, but for the GSE39582 dataset. c, d Violin plots showing the ESTIMATE tumor purity (c) and the enrichment level of total immune infiltrates (d) for each subgroup based on FAAO stratification in TCGA. e Boxplots showing enrichment levels of CD8 T cells and macrophages for the three subgroups based on FAAO stratification in TCGA. f, g Comparisons of the tumor mutation load (f) and expression of immune checkpoint genes (g) across the three subgroups stratified by FAAO score and TIM3 expression in TCGA. h, i Heatmaps show the correlation between transcriptomic expression patterns of melanoma patients receiving PD1 or CTLA4 inhibitors and CRC subgroups stratified by FAAO score and TIM3 in TCGA (h) and GSE39582 (i). *P < 0.05; **P < 0.01; ***P < 0.001

Full size image

Discussion

Although there is a continuously increasing interest in peroxisomes and their functions in health and disease, the prognostic value of the peroxisome pathway has not been investigated. Our work first revealed a prognostic role of the peroxisome pathway in CRC, based on resources of genomic data from public clinical datasets. Then, we identified an at-risk CRC subpopulation based on the peroxisome pathway and TIM3 expression, and patients in this group might exhibit an improved response to immunotherapy. This new stratification was further refined by looking at the specific peroxisome functions, such as beta oxidation of very long-chain fatty acids, FAAO, biosynthesis of bile acids and ether lipids, and ROS turnover. Unexpectedly, we found that only genes involved in FAAO can stratify patients similarly with genes in the peroxisome pathway.

Peroxisomes are critical metabolic organelles associated with lipid metabolism and cellular ROS turnover [10]. Cablé and colleagues found reduced numbers of peroxisomes in colon carcinoma by electron microscopy [16]. However, the role of the peroxisome pathway in colorectal cancer has not been determined. As multiple peroxisomal biogenesis proteins and various oxidases are involved in peroxisome activity, we defined a peroxisome score to summarize the activity of the peroxisome pathway. We found that patients with a low peroxisome score had shorter OS than those with a high peroxisome score. Moreover, the score was associated with peroxisome-related pathways (e.g., the bile acid metabolism pathway and fatty acid metabolism pathway). Within the past few years, the idea that peroxisomes are vital organelles in regulating inflammation and the anti-microbial response has emerged [11,12,13, 39]. Vijayan and colleagues observed that deletion of the peroxisomal biogenesis factors PEX14 or MFP2 leads to a pronounced hyperexpression of COX2 and TNF-α proteins [13]. Deletion of PEX13, another peroxisomal biogenesis factor, causes the activation of Smad-dependent TGF-β signaling and the release of inflammatory cytokines, such as TGF-β and IL-6 [40]. Our results also showed that patients with a low peroxisome score had high immune infiltration and enriched immune-related pathways, despite the poor prognosis.

IC pathways are crucial in preventing autoimmunity and maintaining immune homeostasis; however, their upregulation in the TME of various malignancies can lead to T cell inactivation and immune evasion [41]. By blocking the interactions between ICs and their ligands, ICIs can efficiently initiate anti-tumor immune responses and have shown promise in eliminating tumor cells in some cancers. TIM3 (encoded by Havcr2), as one IC, is an inhibitory receptor expressed by various immune cells and was initially identified as a cell surface marker specific to CD4+ T helper 1 and CD8+ T cells [42, 43]. In cancer, the upregulation of TIM3 expression marks the most terminally exhausted CD8+ T cell subsets [44, 45]. Previous studies have shown that TIM3 expression in the tissues of CRC patients is associated with tumor progression and poor clinical outcomes [36, 37]. Our results showed that for tumors with a low peroxisome score, high TIM3 expression indicated poor OS. Intriguingly, even if it was not statistically significant, high TIM3 expression seemed to be correlated with a better prognosis for tumors with a high peroxisome score, suggesting the possible role of the peroxisome pathway in regulating biological behaviors of TIM3. Further investigations are needed to validate this hypothesis. Our genomic analysis revealed that the prognostic performance of TIM3 was dependent on the peroxisome score, and in agreement with this, we identified an at-risk CRC subpopulation (PERLowTIM3High) with a negative prognostic value.

Beyond the negative prognostic value of group III (PERLowTIM3High), an activated immune phenotype characterized by high levels of pro-tumorous and anti-tumorous immune cells was observed for these tumors. Tumor-infiltrating CD8+ T cells are the most potent tumor-suppressing immune cells. CD8+ T cell abundance was an independent predictor of better disease-free and OS outcomes in early CRC patients, leading to the establishment of the Immunoscore, a reliable independent prognostic marker for CRC [46, 47]. Our data show that tumors in group III had the highest CD8+ T cell abundance. Given that high TIM3 expression is associated with dysfunctional CD8+ T cells, it is not surprising that group III had a poor prognosis in the presence of high CD8+ T cell infiltration. Compared with tumors in other subgroups, group III tumors had the highest level of macrophages, as one of the important types of tumor-supportive immune cells. Tumor-associated macrophages, mostly in the type M2 form, are important for the inhibition of anti-tumor immune responses mediated by T cells within the TME [48, 49]. The risk subgroup largely overlapped with CMS1 and CMS4 subtypes, which are associated with high immune infiltration, but the enrichment of CD8 T cells and macrophages was not a simple reflection of the high proportion of these two CMS subtypes. Moreover, our multivariate cox regression model determined the independent prognostic value of group III for OS. Hence, the peroxisome/TIM3 stratification is distinct from the CMS classification of CRC [50]. We also found an association between group III tumors and MSI status. Our observation might promote refinement of the current clinical dogma that MSI-high tumors are associated with a good prognosis and encourage the further stratification of MSI tumors by peroxisome score and TIM3 expression [51, 52].

Consistent with the concept that MSI-high CRC is associated with high tumor mutations, group III tumors had the highest TMB. Moreover, high TMB was determined to be a predictive biomarker for immunotherapy in several tumor types [53, 54]. Immunotherapy is increasingly being recognized as a major treatment strategy for multiple cancers, including melanoma and a subset of CRC patients [55, 56]. PD1 blockade, as the only immunotherapeutic strategy approved by the FDA, has shown efficacy in metastatic CRC patients with MSI-high status [6, 57, 58]. Further, based on preclinical data, co-blockade of TIM3 and PD1 has shown efficacy for solid tumors, leading to the clinical investigation of anti-TIM3 combined with anti-PD1 in the treatment of various malignancies, including CRC [44, 59,60,61]. The sensitivity to immunotherapy is likely determined by many factors, such as the presence of CD8+ T cells, high TMB, and the expression of ICs [54, 62]. It is well established that metastatic CRC patients with an MSI-high status benefit from anti-PD1/PD-L1 therapy [63]. Our results suggest that patients in group III (enriched with MSI tumors) had poor clinical outcomes but are likely to benefit from treatment targeting PD1. As such, the peroxisome pathway might be an overlooked factor involved in predicting the response to immunotherapy. In addition to anti-PD1 therapy, the co-blockade of PD1 and TIM3 could prove to be a more suitable treatment strategy for group III patients. Although immunotherapy approaches have been unsuccessful for MSS tumors, highly immune-infiltrated MSS tumors in group III could be more likely to respond to ICIs than MSS tumors in other subgroups. Further clinical trials investigating the co-blockade of PD1 and TIM3 in selected group III patients might be warranted.

Since peroxisomes are pleiotropic organelles involved in multiple metabolic processes, the meaning of high peroxisome scores seems obscure. Thus, we also attempted to refine the stratification by genes involved in peroxisomal metabolic processes. We found that CRC patients can be stratified by genes only involved in FAAO, which was similarly with genes in the peroxisome pathway. FAAO is a specific peroxisomal metabolic process in which some unusual fatty acids, such as 2-hydroxy fatty acids and 3-methyl-branched fatty acids (i.e., phytanic acid), are shortened by one carbon atom, which is known to occur in peroxisomes [64, 65]. Peroxisomal FAAO is deficient in patients with peroxisome biogenesis disorders, which is reflected by the accumulation of phytanic acid (PA) in the plasma [66]. However, the role of FAAO or PA in CRC is not clear yet. Our results showed that only for CRC patients with a relatively low peroxisome score or FAAO score was high TIM3 expression indicative of poor OS, suggesting the possible role for normal peroxisomal FAAO in regulating biological functions of TIM3 in CRC, which is a new direction for our future work.

The retrospective nature of this analysis, conducted using public datasets, is the main limitation of our study. Furthermore, although we recapitulated the peroxisome pathway and identified the prognostic role of the peroxisome for CRC, assessing the correlation between the peroxisome score and the abundance of peroxisomes was outside the scope of our study. However, our findings suggested the value of the peroxisome pathway and the peroxisomal FAAO for CRC stratification, as well as the possible role of peroxisomes in TME and immunotherapy.

Conclusions

Collectively, we evaluated the prognostic role of the peroxisome pathway in CRC and identified an at-risk CRC subpopulation (PERLowTIM3High), exhibiting high immune infiltration. Patients with PERLowTIM3High tumors might more likely respond to anti-PD1 based immunotherapy. When looking at individual peroxisomal functions to refine the stratification, only genes involved in FAAO can stratify patients as similar as genes in the peroxisome pathway. Combined evaluation of the expression of TIM3 and genes involved in the peroxisome pathway or FAAO might be used for CRC diagnostics and could be helpful for personalized treatment.

Availability of data and materials

The data from our study were all openly available from the GEO database and TCGA repository (https://portal.gdc.cancer.gov/; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi).

Abbreviations

CRC:

Colorectal cancer

GSVA:

Gene set variant analysis

TIM3:

T-cell immunoglobulin and mucin domain 3

CMS:

Consortium molecular subtype

MSI:

Microsatellite instability

TMB:

Tumor mutational burden

PD-1:

Programmed cell death protein-1

FAAO:

Fatty acid alpha oxidation

ICI:

Immune checkpoint inhibitor

TME:

Tumor microenvironment

ROS:

Reactive oxygen species

TCGA:

The Cancer Genome Atlas

TPM:

Transcripts per kilobase million

GEO:

Gene Expression Omnibus

MSS:

Microsatellite stable

IC:

Immune checkpoint

MSigDB:

The Molecular Signatures Database

GSEA:

Gene set enrichment analysis

KEGG:

Kyoto Encyclopedia of Genes and Genomes

ESTIMATE:

Estimation of STromal and Immune cells in MAlignant Tumours using Expression data

PD1:

programmed cell death protein-1

CTLA4:

cytotoxic T lymphocyte antigen–4

IC:

immune checkpoint

CIMP:

CpG-island methylator phenotype

BOVLCFA:

beta oxidation of very long chain fatty acids

PBAB:

Primary bile acid biosynthesis

ELM:

Ether lipid metabolism

CRTROS:

Cellular response to reactive oxygen species

PA:

Phytanic acid

References

  1. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. 2021;71(3):209–49. https://doi.org/10.3322/caac.21660. Epub 2021 Feb 4.

  2. Biller LH, Schrag D. Diagnosis and treatment of metastatic colorectal Cancer. Jama. 2021;325:669–85.

    Article  PubMed  Google Scholar 

  3. Khalil DN, Smith EL, Brentjens RJ, Wolchok JD. The future of cancer treatment: immunomodulation, CARs and combination immunotherapy. Nat Rev Clin Oncol. 2016;13:273–90.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Hodi FS, O'Day SJ, McDermott DF, Weber RW, Sosman JA, Haanen JB, et al. Improved survival with ipilimumab in patients with metastatic melanoma. N Engl J Med. 2010;363(8):711–23. https://doi.org/10.1056/NEJMoa1003466. Epub 2010 Jun 5. Erratum in: N Engl J Med. 2010;363(13):1290.

  5. Le DT, Durham JN, Smith KN, Wang H, Bartlett BR, Aulakh LK, et al. Mismatch repair deficiency predicts response of solid tumors to PD-1 blockade. Science. 2017;357:409–13.

    Article  PubMed  PubMed Central  Google Scholar 

  6. Le DT, Uram JN, Wang H, Bartlett B, Kemberling H, Eyring A, et al. Programmed death-1 blockade in mismatch repair deficient colorectal cancer. J Clin Oncol. 2016;34:15_suppl:103–103.

    Google Scholar 

  7. Kaymak I, Williams KS, Cantor JR, Jones RG. Immunometabolic interplay in the tumor microenvironment. Cancer Cell. 2020;39:28–37.

    Article  PubMed  PubMed Central  Google Scholar 

  8. Peddareddigari VG, Wang D, DuBois RN. The tumor microenvironment in colorectal carcinogenesis. Cancer Microenviron. 2010;3:149–66.

    Article  PubMed  PubMed Central  Google Scholar 

  9. Kubli SP, Berger T, Araujo DV, Siu LL, Mak TW. Beyond immune checkpoint blockade: emerging immunological strategies. Nat Rev Drug Discov. 2021;20(12):899–919. https://doi.org/10.1038/s41573-021-00155-y. Epub 2021 Mar 8.

  10. Islinger M, Voelkl A, Fahimi HD, Schrader M. The peroxisome: an update on mysteries 2.0. Histochem Cell Biol. 2018;150:443–71.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Cara FD, Sheshachalam A, Braverman NE, Rachubinski RA, Simmonds AJ. Peroxisome-mediated metabolism is required for immune response to microbial infection. Immunity. 2017;47:93–106.e7.

    Article  PubMed  Google Scholar 

  12. Terlecky SR, Terlecky LJ, Giordano CR. Peroxisomes, oxidative stress, and inflammation. World J Biological Chem. 2012;3:93–7.

    Article  Google Scholar 

  13. Vijayan V, Srinu T, Karnati S, Garikapati V, Linke M, Kamalyan L, et al. A new immunomodulatory role for peroxisomes in macrophages activated by the TLR4 ligand lipopolysaccharide. J Immunol. 2017;198:2414–25.

    Article  PubMed  Google Scholar 

  14. Mock A, Zschäbitz S, Kirsten R, Scheffler M, Wolf B, Herold-Mende C, et al. Serum very long-chain fatty acid-containing lipids predict response to immune checkpoint inhibitors in urological cancers. Cancer Immunol Immunother. 2019;68:2005–14.

    Article  PubMed  Google Scholar 

  15. Baur G, Wendel A. The activity of the peroxide-metabolizing system in human colon carcinoma. J Cancer Res Clin. 1980;97:267–73.

    Article  Google Scholar 

  16. Cablé S, Keller JM, Colin S, Haffen K, Kédinger M, Parache RM, et al. Peroxisomes in human colon carcinomas. Virchows Archiv B. 1992;62:221.

    Article  Google Scholar 

  17. Keller J, Cablé S, Bouhtoury F, Heusser S, Scotto C, Armbruster L, et al. Peroxisome through cell differentiation and neoplasia. Biol Cell. 1993;77:77–88.

    Article  PubMed  Google Scholar 

  18. Song J, Park S, Oh J, Kim D, Ryu JH, Park WC, et al. NUDT7 loss promotes KrasG12D CRC development. Cancers. 2020;12:576.

    Article  PubMed Central  Google Scholar 

  19. Mounir M, Lucchetta M, Silva TC, Olsen C, Bontempi G, Chen X, et al. New functionalities in the TCGAbiolinks package for the study and integration of cancer data from GDC and GTEx. PLoS Comput Biol. 2019;15(3):e1006701. https://doi.org/10.1371/journal.pcbi.1006701.

  20. Aran D, Sirota M, Butte AJ. Systematic pan-cancer analysis of tumour purity. Nat Commun. 2015;6:8971.

    Article  PubMed  Google Scholar 

  21. Liu J, Lichtenberg T, Hoadley KA, Poisson LM, Lazar AJ, Cherniack AD, et al. An integrated TCGA Pan-Cancer clinical data resource to drive high-quality survival outcome analytics. Cell. 2018;173:400–416.e11.

    Article  PubMed  PubMed Central  Google Scholar 

  22. Mayakonda A, Lin D-C, Assenov Y, Plass C, Koeffler HP. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 2018;28:1747–56.

    Article  PubMed  PubMed Central  Google Scholar 

  23. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-Seq data. Bmc Bioinformatics. 2013;14:7.

    Article  PubMed  PubMed Central  Google Scholar 

  24. Tokumaru Y, Oshi M, Katsuta E, Yan L, Satyananda V, Matsuhashi N, et al. KRAS signaling enriched triple negative breast cancer is associated with favorable tumor immune microenvironment and better survival. Am J Cancer Res. 2020;10:897–907.

    PubMed  PubMed Central  Google Scholar 

  25. Oshi M, Takahashi H, Tokumaru Y, Yan L, Rashid OM, Nagahashi M, et al. The E2F pathway score as a predictive biomarker of response to neoadjuvant therapy in ER+/HER2− breast Cancer. Cells. 2020;9:1643.

    Article  PubMed Central  Google Scholar 

  26. Khatri P, Sirota M, Butte AJ. Ten years of pathway analysis: current approaches and outstanding challenges. PLoS Comput Biol. 2012;8:e1002375.

    Article  PubMed  PubMed Central  Google Scholar 

  27. Merico D, Isserlin R, Stueker O, Emili A, Bader GD. Enrichment map: a network-based method for gene-set enrichment visualization and interpretation. PLoS One. 2010;5:e13984.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Jassal B, Matthews L, Viteri G, Gong C, Lorente P, Fabregat A, et al. The reactome pathway knowledgebase. Nucleic Acids Res. 2019;48:D498–503.

    PubMed Central  Google Scholar 

  29. Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, et al. Gene ontology: tool for the unification of biology. Nat Genet. 2000;25:25–9.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Consortium TGO, Carbon S, Douglass E, Good BM, Unni DR, Harris NL, et al. The gene ontology resource: enriching a GOld mine. Nucleic Acids Res. 2020;49:D325–34.

    Article  Google Scholar 

  31. Kanehisa M, Sato Y, Kawashima M, Furumichi M, Tanabe M. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Res. 2016;44:D457–62.

    Article  PubMed  Google Scholar 

  32. Kanehisa M, Goto S. KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28:27–30.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol. 2019;37:773–82.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013;4:2612.

    Article  PubMed  Google Scholar 

  35. Roh W, Chen P-L, Reuben A, Spencer CN, Prieto PA, Miller JP, et al. Integrated molecular analysis of tumor biopsies on sequential CTLA-4 and PD-1 blockade reveals markers of response and resistance. Sci Transl Med. 2017;9:eaah3560.

    Article  PubMed  PubMed Central  Google Scholar 

  36. Yu M, Lu B, Liu Y, Me Y, Wang L, Zhang P. Tim-3 is upregulated in human colorectal carcinoma and associated with tumor progression. Mol Med Rep. 2017;15:689–95.

    Article  PubMed  Google Scholar 

  37. Kuai W, Xu X, Yan J, Zhao W, Li Y, Wang B, et al. Prognostic impact of PD-1 and Tim-3 expression in tumor tissue in stage I-III colorectal Cancer. Biomed Res Int. 2020;2020:1–13.

    Article  Google Scholar 

  38. Ghiringhelli F, Fumet J-D. Is there a place for immunotherapy for metastatic microsatellite stable colorectal Cancer? Front Immunol. 2019;10:1816.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Facciotti F, Ramanjaneyulu GS, Lepore M, Sansano S, Cavallari M, Kistowska M, et al. Peroxisome-derived lipids are self antigens that stimulate invariant natural killer T cells in the thymus. Nat Immunol. 2012;13:474–80.

    Article  PubMed  Google Scholar 

  40. Oruqaj G, Karnati S, Vijayan V, Kotarkonda LK, Boateng E, Zhang W, et al. Compromised peroxisomes in idiopathic pulmonary fibrosis, a vicious cycle inducing a higher fibrotic response via TGF-β signaling. Proc National Acad Sci. 2015;112:E2048–57.

    Article  Google Scholar 

  41. Toor SM, Nair VS, Decock J, Elkord E. Immune checkpoints in the tumor microenvironment. Semin Cancer Biol. 2020;65:1–12.

    Article  PubMed  Google Scholar 

  42. Monney L, Sabatos CA, Gaglia JL, Ryu A, Waldner H, Chernova T, et al. Th1-specific cell surface protein Tim-3 regulates macrophage activation and severity of an autoimmune disease. Nature. 2002;415:536–41.

    Article  PubMed  Google Scholar 

  43. Gayden T, Sepulveda FE, Khuong-Quang D-A, Pratt J, Valera ET, Garrigue A, et al. Germline HAVCR2 mutations altering TIM-3 characterize subcutaneous panniculitis-like T cell lymphomas with hemophagocytic lymphohistiocytic syndrome. Nat Genet. 2018;50:1650–7.

    Article  PubMed  Google Scholar 

  44. Sakuishi K, Apetoh L, Sullivan JM, Blazar BR, Kuchroo VK, Anderson AC. Targeting Tim-3 and PD-1 pathways to reverse T cell exhaustion and restore anti-tumor immunity. J Exp Medicine. 2010;207:2187–94.

    Article  Google Scholar 

  45. Fourcade J, Sun Z, Benallaoua M, Guillaume P, Luescher IF, Sander C, et al. Upregulation of Tim-3 and PD-1 expression is associated with tumor antigen–specific CD8+ T cell dysfunction in melanoma patientsNY-ESO-1–specific T cells co-express Tim-3 and PD-1. J Exp Medicine. 2010;207:2175–86.

    Article  Google Scholar 

  46. Galon J, Costes A, Sanchez-Cabo F, Kirilovsky A, Mlecnik B, Lagorce-Pagès C, et al. Type, density, and location of immune cells within human colorectal tumors predict clinical outcome. Science. 2006;313:1960–4.

    Article  PubMed  Google Scholar 

  47. Pagès F, Mlecnik B, Marliot F, Bindea G, Ou F-S, Bifulco C, et al. International validation of the consensus Immunoscore for the classification of colon cancer: a prognostic and accuracy study. Lancet. 2018;391:2128–39.

    Article  PubMed  Google Scholar 

  48. Liu Y, Cao X. The origin and function of tumor-associated macrophages. Cell Mol Immunol. 2015;12:1–4.

    Article  PubMed  Google Scholar 

  49. Qian B-Z, Pollard JW. Macrophage diversity enhances tumor progression and metastasis. Cell. 2010;141:39–51.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Guinney J, Dienstmann R, Wang X, de Reyniès A, Schlicker A, Soneson C, et al. The consensus molecular subtypes of colorectal cancer. Nat Med. 2015;21:1350–6.

    Article  PubMed  PubMed Central  Google Scholar 

  51. Popat S, Hubner R, Houlston RS. Systematic review of microsatellite instability and colorectal Cancer prognosis. J Clin Oncol. 2004;23:609–18.

    Article  Google Scholar 

  52. Ribic CM, Sargent DJ, Moore MJ, Thibodeau SN, French AJ, Goldberg RM, et al. Tumor microsatellite-instability status as a predictor of benefit from fluorouracil-based adjuvant chemotherapy for Colon Cancer. New Engl J Med. 2003;349:247–57.

    Article  PubMed  Google Scholar 

  53. Muzny DM, Bainbridge MN, Chang K, Dinh HH, Drummond JA, Fowler G, et al. Comprehensive molecular characterization of human colon and rectal cancer. Nature. 2012;487:330–7.

    Article  Google Scholar 

  54. Samstein RM, Lee C-H, Shoushtari AN, Hellmann MD, Shen R, Janjigian YY, et al. Tumor mutational load predicts survival after immunotherapy across multiple cancer types. Nat Genet. 2019;51:202–6.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Marin-Acevedo JA, Dholaria B, Soyano AE, Knutson KL, Chumsri S, Lou Y. Next generation of immune checkpoint therapy in cancer: new developments and challenges. J Hematol Oncol. 2018;11:39.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Zappasodi R, Merghoub T, Wolchok JD. Emerging concepts for immune checkpoint blockade-based combination therapies. Cancer Cell. 2018;33:581–98.

    Article  PubMed  PubMed Central  Google Scholar 

  57. Le DT, Uram JN, Wang H, Bartlett BR, Kemberling H, Eyring AD, et al. PD-1 blockade in tumors with mismatch-repair deficiency. New Engl J Med. 2015;372:2509–20.

    Article  PubMed  Google Scholar 

  58. Overman MJ, McDermott R, Leach JL, Lonardi S, Lenz H-J, Morse MA, et al. Nivolumab in patients with metastatic DNA mismatch repair-deficient or microsatellite instability-high colorectal cancer (CheckMate 142): an open-label, multicentre, phase 2 study. Lancet Oncol. 2017;18:1182–91.

    Article  PubMed  PubMed Central  Google Scholar 

  59. Zhou Q, Munger ME, Veenstra RG, Weigel BJ, Hirashima M, Munn DH, et al. Coexpression of Tim-3 and PD-1 identifies a CD8+ T-cell exhaustion phenotype in mice with disseminated acute myelogenous leukemia. Blood. 2011;117:4501–10.

    Article  PubMed  PubMed Central  Google Scholar 

  60. Borate U, Esteve J, Porkka K, Knapper S, Vey N, Scholl S, et al. Phase Ib study of the anti-TIM-3 antibody MBG453 in combination with Decitabine in patients with high-risk myelodysplastic syndrome (MDS) and acute myeloid leukemia (AML). Blood. 2019;134(Supplement_1):570.

    Article  Google Scholar 

  61. Ciardiello D, Vitiello PP, Cardone C, Martini G, Troiani T, Martinelli E, Ciardiello F. Immunotherapy of colorectal cancer: Challenges for therapeutic efficacy. Cancer Treat Rev. 2019;76:22–32. https://doi.org/10.1016/j.ctrv.2019.04.003. Epub 2019 May 4.

  62. Gibney GT, Weiner LM, Atkins MB. Predictive biomarkers for checkpoint inhibitor-based immunotherapy. Lancet Oncol. 2016;17:e542–51.

  63. Overman MJ, Lonardi S, Leone F, McDermott RS, Morse MA, Wong KYM, et al. Nivolumab in patients with DNA mismatch repair deficient/microsatellite instability high metastatic colorectal cancer: Update from CheckMate 142. J Clin Oncol. 2017;35(4_suppl):519.

    Article  Google Scholar 

  64. Casteels M, Sniekers M, Fraccascia P, Mannaerts GP, Van Veldhoven PP. The role of 2-hydroxyacyl-CoA lyase, a thiamin pyrophosphate-dependent enzyme, in the peroxisomal metabolism of 3-methyl-branched fatty acids and 2-hydroxy straight-chain fatty acids. Biochem Soc T. 2007;35:876–80.

    Article  Google Scholar 

  65. Jansen GA, Wanders RJA. Alpha-oxidation. Biochimica Et Biophysica Acta Bba - Mol Cell Res. 2006;1763:1403–12.

    Article  Google Scholar 

  66. Wierzbicki AS. Peroxisomal disorders affecting phytanic acid α-oxidation: a review. Biochem Soc T. 2007;35:881–6.

    Article  Google Scholar 

Download references

Acknowledgements

We would like to thank the members of our laboratory for their suggestions and Yujian Yang for the critical review of the manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (No: 81870390) and the Discipline and Platform Construction Project of Zhongnan hospital of Wuhan University (No: ZLYNXM202017).

Author information

Author notes

  1. Jinwen Yin, Hao Wang and Yuntian Hong contributed equally to this work.

Authors and Affiliations

  1. Department of Gastroenterology, Zhongnan Hospital of Wuhan University, Wuhan, 430000, China

    Jinwen Yin, Hao Wang, Yuntian Hong, Anli Ren, Haizhou Wang, Lan Liu & Qiu Zhao

  2. Hubei Clinical Center and Key Lab of Intestinal and Colorectal Diseases, Wuhan, 430000, China

    Jinwen Yin, Hao Wang, Anli Ren, Haizhou Wang, Lan Liu & Qiu Zhao

  3. Department of Colorectal and Anal Surgery, Zhongnan Hospital of Wuhan University, Wuhan, 430000, China

    Yuntian Hong

Authors
  1. Jinwen Yin
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Hao Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Yuntian Hong
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Anli Ren
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Haizhou Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Lan Liu
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Qiu Zhao
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conception and design: Jinwen Yin, Lan Liu, Qiu Zhao. Acquisition of data: Hao Wang. Analysis and interpretation of data: Yuntian Hong, Hao Wang, Anli Ren. Writing and/or review of the manuscript: Jinwen Yin, Haizhou Wang. Study supervision: Lan Liu, Qiu Zhao. All authors read and proved the final manuscript.

Corresponding authors

Correspondence to Lan Liu or Qiu Zhao.

Ethics declarations

Ethics approval and consent to participate

This research is based on the collection and analysis of publicly available data, and do not involve any human specimens and animal experiments.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Table S1

: Pathway enrichment analysis of peroxisome-high and peroxisome-low groups in TCGA dataset. Table S2: Pathway enrichment analysis of peroxisome-high and peroxisome-low groups in GEO-GSE39582 dataset. Table S3. Univariate and multivariate analyses of the relationship between subgroups stratified by FAAO/TIM3 and clinical characteristics. Table S4. Univariate and multivariate analyses of the relationship between subgroups stratified by BOVLCFA/TIM3 and clinical characteristics.

Additional file 2: Figure S1

. a The volcano plot shows the enrichment analysis using Hallmark gene sets in the GSE39582 cohort. The red gene sets were enriched in the Per-High group, whereas the blue gene sets were enriched in the Per-Low group. b Correlation between the peroxisome score versus immune score, stromal score, and tumor purity in the GSE39582 cohort. c Boxplots of the immune score, stromal score, and tumor purity from ESTIMATE of Per-High and Per-Low groups in the GSE39582 cohort. d Boxplot of total immune infiltrates (sum of absolute scores across 22 immune cell types) and correlation between peroxisome score and total immune infiltrates of patients in TCGA cohort. For Boxplots, p values in group comparison with Mann-Whitney U-test are shown. For panels B, Pearson’s rho (r) and statistical difference (p) are indicated. **P < 0.01; ****P < 0.0001. Figure S2. Boxplots depict the expression of immune checkpoint genes in the TCGA colorectal dataset. Statistical P values between groups were determined by Mann-Whitney U-test. HM_PEROXISOME: Hallmark Peroxisome gene set. **P < 0.01, ***P < 0.001, ****P < 0.0001. Figure S3. a The log-rank test score at the candidate cut-off across the log-transformed TIM3 gene expression values is plotted. b Kaplan-Meier curves are plotted for the TIM3-High group and TIM3-Low group by the optimal cut-off shown in panel a. Figure S4. GSE39582 group III tumors were highly infiltrated with CD8 T cells and macrophages. a Violin plot showing the total immune infiltrates of CIBERSORTx for each subgroup in the GSE39582 CRC dataset. b Violin plot showing the ESTIMATE tumor purity for each subgroup in the GSE39582 CRC dataset. c Boxplots showing enrichment levels of CD8 T cells and macrophages for each subgroup in the GSE39582 dataset. d, g Boxplots of enrichment level of CD8 T cells (d) and macrophages (g) MSS tumors. e, h Boxplots of enrichment level of CD8 T cells (e) and macrophages (h) for GSE39582 CMS1 and CMS4 tumors. f, i Boxplots of enrichment level of CD8 T cells (f) and macrophages (i) for GSE39582 CMS2, CMS3, and indeterminate tumors. *P < 0.05; **P < 0.01; ***P < 0.001. Figure S5. a Scatter plot of BOVLCFA score and log2-transformed TIM3 gene expression values are shown for TCGA cohort. The Pearson’s rho (r) and statistical difference (p) are shown in the scatter plot. MSI (blue diamond) and MSS (black circle) status are labeled for CRC patients. The median value of TIM3 expression is indicated in a gray dashed line. The candidate cut-offs (black lines) for high BOVLCFA and low BOVLCFA groups are shown. Using the median value of BOVLCFA score and these two candidate cut-offs, we separated patients into four different groups labeled. Kaplan-Meier curves show the OS for the optimal cut-off of TIM3 expression in TCGA CRC subgroups with high BOVLCFA score and low BOVLCFA score, respectively. P-value was calculated by the log-rank test and shown for each plot. b The same as panel a, but for FAAO score. c The same as panel a, but for PBAB score. d The same as panel a, but for ELM score. e The same as panel a, but for CRTROS score. BOVLCFA: beta oxidation of very long chain fatty acids; PBAB: primary bile acid biosynthesis; ELM: ether lipid metabolism; CRTROS: cellular response to reactive oxygen species.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Yin, J., Wang, H., Hong, Y. et al. Identification of an at-risk subpopulation with high immune infiltration based on the peroxisome pathway and TIM3 in colorectal cancer. BMC Cancer 22, 44 (2022). https://doi.org/10.1186/s12885-021-09085-9

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12885-021-09085-9

Keywords

BMC Cancer

ISSN: 1471-2407

Contact us