Integrative Multiomics Nominate GGCT as a Crucial Regulator of Immunosuppression in Colorectal Cancer

TL;DR

This study identifies GGCT as a key gene involved in immune evasion in colorectal cancer, suggesting it could be a new target for treatment.

Contribution

The study introduces GGCT as a novel prognostic biomarker and potential therapeutic target for immunosuppression in CRC.

Findings

GGCT is associated with poor survival in colorectal cancer patients.

GGCT correlates with immunosuppressive tumor microenvironment features.

GGCT shows potential as a predictor of immunotherapy response.

Abstract

Colorectal cancer (CRC) remains a leading cause of cancer‐related mortality worldwide, with tumor microenvironment (TME) heterogeneity playing a critical role in disease progression and therapeutic response. Immune escape (IE) mechanisms facilitate tumor evasion from host immune surveillance, yet their characterization at the single‐cell level in CRC is incomplete. This study integrated single‐cell RNA sequencing (scRNA‐seq) and bulk transcriptomic data from multiple public cohorts to systematically explore IE‐related signatures in CRC. We identified major and minor cell populations within the TME and performed differential gene expression analysis. Using high‐dimensional weighted gene coexpression network analysis (hdWGCNA), we identified gene modules correlated with IE activity. Subsequent survival analysis across six independent cohorts revealed Gamma‐glutamylcyclotransferase (GGCT) as a novel prognostic biomarker associated with poor survival. Functional enrichment analysis indicated GGCT′s involvement in critical oncogenic pathways. Furthermore, GGCT expression correlated with altered immune infiltration profiles and stromal components, suggesting its role in modulating the immunosuppressive TME. Additionally, GGCT demonstrated potential predictive value for response to immunotherapy across multiple datasets. Our findings highlight GGCT as a key player in CRC immune evasion and a promising therapeutic target.

Figures

Identification of DEGs among microenvironment cells at the scRNA‐seq level of CRC.

Functional annotation of GGCT. Bubble plot showing the Metascape‐based pathway enrichment analysis on DEGs between GGCT‐based groups.

Immunotherapy predictive value of GGCT in six cohorts (Hugo cohort, Riaz cohort, Gao cohort, Cho cohort, IMvigor210 cohort, and Nathanson cohort).

Full text

1. Introduction

Colorectal cancer (CRC) is one of the most prevalent malignancies globally, accounting for significant morbidity and mortality [1, 2]. Despite advances in early detection and treatment, metastatic CRC continues to pose therapeutic challenges, with resistance to conventional therapies and immunotherapy being common [3–5]. A deeper understanding of the molecular mechanisms that drive CRC progression is essential for developing novel therapeutic strategies [6–8]. The tumor microenvironment (TME) is a complex ecosystem comprising tumor cells, immune cells, stromal cells, and extracellular matrix components, all of which interact to influence tumor behavior and response to treatment. In particular, the immune contexture of the TME is a critical determinant of patient prognosis and therapeutic efficacy.

Immune escape (IE) is a hallmark of cancer, enabling tumor cells to evade detection and destruction by the immune system [9]. Mechanisms of IE include downregulation of antigen presentation, recruitment of immunosuppressive cells, and expression of immune checkpoint molecules [10]. Although these mechanisms have been studied in various cancers, the specific gene expression programs underlying IE in CRC at the single‐cell resolution remain poorly characterized. A deeper understanding of these programs could identify novel biomarkers and therapeutic targets [11, 12]. Recent advances in single‐cell RNA sequencing (scRNA‐seq) have allowed for unprecedented resolution in dissecting the cellular heterogeneity of the TME [13]. Integrating scRNA‐seq data with bulk transcriptomic data from multiple cohorts can enhance the robustness of identifying consistently prognostic genes.

In this study, we leveraged scRNA‐seq data from public repositories to first identify and characterize the cellular composition of the CRC TME. We then focused on a set of 182 previously identified IE‐related genes to quantify IE activity at the single‐cell level. Using hdWGCNA, we constructed coexpression networks to identify gene modules most associated with this IE signature. Through rigorous survival analysis across six independent CRC cohorts, we pinpointed GGCT as a top prognostic gene. We further explored its biological functions, its correlation with immune infiltration, and its predictive value for immunotherapy response. Our comprehensive analysis is aimed at providing new insights into the mechanisms of immune evasion in CRC and to propose GGCT as a potential target for future therapeutic strategies.

2. Materials and Methods

2.1. Data Collection and Processing

We obtained single‐cell transcriptomic data from GSE132465 and GSE144735 [14]. The R package Seurat was applied to identify major cell types [15]. We obtained the transcriptome data and corresponding clinical information of CRC patients from six publicly available cohorts, including COADREAD_TCGA [16], GSE106584 [17], GSE29621 [18], GSE103479 [19], GSE17538 [20], and GSE72969 [21]. The R package sva was used to remove the batch effect.

2.2. Bioinformatics Analysis

The 182 IE‐related genes were obtained from the previous finding [22]. The IE signature was calculated in each cell using the R package AUCell. The R package hdWGCNA was used to identify the module genes related to the IE signature [23]. The Kaplan–Meier survival curves for overall survival (OS) were plotted and compared using the R package survminer. Drug sensitivity to GGCT was predicted using oncoPredict [24]. The immune infiltration characteristics related to GGCT were predicted using the tumor immune estimation resource (TIMER) [25]. The connection between GGCT and stromal score, immune score, ESTIMATE score, and tumor purity was analyzed using the ESTIMATE algorithm [26]. The Metascape platform was employed for pathway annotation of GGCT [27].

3. Results

3.1. Cellular Heterogeneity of the CRC TME Revealed by scRNA‐seq

We began by analyzing scRNA‐seq data from GEO datasets (GSE132465 and GSE144735) using the Seurat pipeline. This analysis identified major cell types within the CRC TME, including epithelial cells, T cells, B cells, myeloid cells, fibroblasts, and endothelial cells (Figure 1a). Further subclustering revealed minor cell populations, such as distinct CD4+ T cell subsets, plasma cells, and macrophage subtypes (Figure 1b). We also successfully identified malignant epithelial cells based on copy number variation inference and marker gene expression (Figure 1c). This detailed atlas provided the foundation for subsequent cell type–specific analyses.

Figure 1. Identification of microenvironment cells at the scRNA‐seq level of CRC. (a) Major cell types within the microenvironment of CRC. (b) Minor cell types within the microenvironment of CRC. (c) Tumor cells within the microenvironment of CRC.(a)(b)(c)

3.2. Differential Gene Expression Across TME Cell Types

Differential expression analysis between the major cell types identified numerous significantly dysregulated genes. Figure 2 highlights the top differentially expressed genes (DEGs) distinguishing immune cells from stromal and tumor cells. Notably, genes involved in antigen presentation (e.g., HLA genes) were highly expressed in myeloid cells, whereas genes related to matrix remodeling (e.g., COL1A1) were upregulated in fibroblasts.

Identification of DEGs among microenvironment cells at the scRNA‐seq level of CRC.

3.3. Identification of IE‐Related Gene Coexpression Modules via hdWGCNA

To investigate the coordination of IE‐related genes, we performed hdWGCNA on the scRNA‐seq data. We first determined an appropriate soft power threshold based on scale‐free topology fit and connectivity metrics (Figure 3a). The analysis identified several coexpressed gene modules, with the dendrogram illustrating the grouping of genes into distinct modules (Figure 3b). Correlation analysis between module eigengenes and the IE signature score (calculated by AUCell) revealed that the yellow module was most positively correlated with IE activity (Figure 3c), suggesting its strong involvement in immune evasion processes.

Figure 3hdWGCNA on IE‐related gene module. (a) The association between soft power threshold and scale‐free topology model fit and mean/median/max connectivity. (b) hdWGCNA dendrogram showing the distribution of gene modules. (c) Heatmap showing the correlation between IE signature and gene modules.(a)(b)(c)

3.4. GGCT Is a Prognostic Biomarker in Multiple CRC Cohorts

The yellow module from hdWGCNA contained numerous genes. Univariate Cox regression analysis performed on these genes across the TCGA COADREAD cohort identified several significant prognostic genes, with GGCT being among the top hazards (hazard ratio > 1, p < 0.001) (Figure 4a). Examination of GGCT expression at the single‐cell level showed its predominant expression in malignant epithelial cells and a subset of myeloid cells (Figure 4b). We then stratified patients in each of the six bulk transcriptomic cohorts (TCGA, GSE106584, GSE29621, GSE103479, GSE17538, and GSE72969) into high and low GGCT expression groups. Kaplan–Meier survival analysis consistently demonstrated that high GGCT expression was significantly associated with worse OS across all cohorts (Figures 4c, 4d, 4e, 4f, 4g, and 4h, all log–rank p < 0.05), solidifying its role as a robust adverse prognostic marker.

Figure 4. Identification of GGCT as a hazardous marker in CRC. (a) Univariate Cox regression analysis on module genes from the yellow module. (b) Vlnplot showing the GGCT expression in microenvironment cells. (c) Survival curves of GGCT‐based groups in the TCGA cohort. (d) Survival curves of GGCT‐based groups in the GSE103479 cohort. (e) Survival curves of GGCT‐based groups in the GSE106584 cohort. (f) Survival curves of GGCT‐based groups in the GSE17538 cohort. (g) Survival curves of GGCT‐based groups in the GSE29621 cohort. (h) Survival curves of GGCT‐based groups in the GSE72969 cohort.(a)(b)(c)(d)(e)(f)(g)(h)

3.5. Functional Annotation of GGCT

To understand the biological pathways associated with GGCT, we performed differential gene expression analysis between high and low GGCT groups in the TCGA cohort and subjected the DEGs to pathway enrichment analysis using Metascape. The results, presented as a bubble plot in Figure 5, showed that GGCT‐high tumors were significantly enriched in pathways related to cell cycle regulation, DNA replication, immune regulation, and several oncogenic signaling pathways. This suggests that GGCT may promote CRC aggressiveness by driving proliferation and genomic instability.

Functional annotation of GGCT. Bubble plot showing the Metascape‐based pathway enrichment analysis on DEGs between GGCT‐based groups.

3.6. GGCT Is Associated With an Immunosuppressive TME

Using the ESTIMATE algorithm, we found that tumors with high GGCT expression had significantly higher stromal, immune, and ESTIMATE scores, but lower tumor purity (Figure 6a), indicating a greater abundance of nontumor cellular components. Further analysis using TIMER revealed that GGCT‐high tumors had significantly altered immune infiltration, characterized by a decrease in CD8+ T cells and dendritic cells (DCs) and an increase in immunosuppressive macrophages (Figure 6b). This profile is consistent with an immune‐excluded or immunosuppressive TME, which may facilitate immune escape.

Figure 6. Immune features of GGCT. (a) The levels of stromal score, immune score, and ESTIMATE score in GGCT‐based groups. (b) The abundance of B cells, CD4 T cells, CD8 T cells, neutrophils, macrophages, and DCs in GGCT‐based groups.(a)(b)

3.7. GGCT Predicts Immunotherapy Response

Given its association with immune evasion, we evaluated the predictive value of GGCT for immunotherapy response in six publicly available cohorts of patients treated with immune checkpoint inhibitors (Hugo, Riaz, Gao, Cho, IMvigor210, and Nathanson). In most of these cohorts, patients with high GGCT expression had significantly lower rates of clinical response (e.g., complete or partial response) and poorer survival after immunotherapy (Figure 7), suggesting that GGCT could serve as a novel predictive biomarker for resistance to immunotherapy.

Immunotherapy predictive value of GGCT in six cohorts (Hugo cohort, Riaz cohort, Gao cohort, Cho cohort, IMvigor210 cohort, and Nathanson cohort).

4. Discussion

CRC is the world′s fourth most deadly cancer, with almost 900,000 deaths annually [5]. The high mortality rate is predominantly attributable to metastasis [28, 29], which is often facilitated by IE. This study presents a comprehensive analysis of IE‐related genes in CRC, leveraging both single‐cell and bulk transcriptomic data. We constructed a detailed cellular map of the CRC TME, identified a key IE‐related gene module through hdWGCNA, and validated GGCT as a novel and robust prognostic biomarker across multiple independent cohorts. Furthermore, we elucidated its functional roles in promoting cell proliferation and shaping an immunosuppressive TME and demonstrated its potential value in predicting response to immunotherapy.

Our scRNA‐seq analysis confirmed the profound heterogeneity of the CRC TME, consistent with previous studies. The identification of GGCT as predominantly expressed in malignant cells aligns with its known role in other cancers; GGCT is involved in glutathione metabolism [30], which can protect cells from oxidative stress and confer chemoresistance. However, its specific role in immune modulation within CRC is novel. The consistent association between high GGCT expression and poor survival across all six bulk cohorts underscores its clinical relevance as a prognostic factor.

The functional enrichment analysis suggests that GGCT is not merely a passive marker but an active contributor to tumor progression. Its association with cell cycle and DNA replication pathways indicates a role in enhancing tumor cell proliferation. More importantly, the strong correlation between GGCT expression and an immunosuppressive TME landscape provides a mechanistic link to immune evasion. The observed reduction in cytotoxic CD8+ T cells and increase in protumor macrophages are hallmarks of a noninflamed TME, which is typically resistant to immunotherapy [31].

The predictive value of GGCT for immunotherapy resistance is a significant finding. Although biomarkers like PD‐L1 expression and tumor mutational burden are currently used [32, 33], they have limitations. GGCT offers a new, complementary perspective rooted in cellular metabolism and immune contexture. Its ability to predict worse outcomes after immunotherapy in multiple cohorts highlights its potential utility in clinical decision‐making, perhaps for patient stratification or combination therapy design.

Several limitations should be considered. This is a retrospective bioinformatics study based on public data. Experimental validation using in vitro and in vivo models is necessary to confirm the functional role of GGCT in immune evasion. The precise molecular mechanisms by which GGCT influences immune cell infiltration and function require further investigation, potentially involving interactions with specific cytokines or metabolic pathways.

In conclusion, our integrated multiomics approach identifies GGCT as a critical player in CRC immune evasion and prognosis. It impacts the TME by fostering an immunosuppressive milieu and is associated with resistance to immunotherapy. These findings nominate GGCT as a promising therapeutic target. Future work should focus on developing small‐molecule inhibitors against GGCT and testing their efficacy, alone or in combination with immunotherapy, in preclinical models of CRC.

Ethics Statement

Ethics approval was not required for this study, as it analyzed existing genomic datasets with no direct involvement of humans or animals.

Disclosure

S.Z. was responsible for the in vitro experiments.

Conflicts of Interest

The authors declare no conflicts of interest.

Author Contributions

Q.N., S.Y., and C.S. designed the specific details of this study and drafted the manuscript. Y.L. participated in data analysis, whereas K.H. and L.N. contributed to manuscript revisions.

Funding

This study was supported by the Key Research and Development Program Project of Heilongjiang Province (2022ZX06C22), Heart to Heart Foundation (HXXT2021ktyj002), Wu Jieping Medical Foundation (10.13039/100007452; 320.6750.2023‐07‐6, and Beijing Medical Award Foundation (YXJL‐2021‐0763‐0482).

Supporting information

Supporting Information Additional supporting information can be found online in the Supporting Information section.. Figure S1: Functional annotation of GGCT.