Multiomics study on the tertiary lymphatic structure signature of colon cancer

Introduction

Colorectal cancer is one of the major types of cancer in the world. According to recent cancer statistics, colorectal cancer, along with prostate cancer and lung cancer, accounts for almost 48% of all incident cancer cases among males in the United States; meanwhile, among females, colorectal cancer, breast cancer, and lung cancer account for 51% of all incident cases (1). In addition, the incidence rate of colorectal cancer varies from region to region and is significantly higher in developed countries than in developing countries (2). However, as the lifestyles and social habits in developing countries continue to change, it is estimated that the incidence rate of colorectal cancer will increase to 2.5 million new cases globally in 2035 (3,4). Moreover, although colon cancer is more common than is colorectal cancer—with an incidence rate about twice that of rectal cancer (5)—and there is evidence suggesting that the widespread use of colonoscopy and colon cancer screening has helped stabilize and reduce the incidence rate, the increasing number of patients under 50 years old with this disease remains a concern (6).

Tertiary lymphoid structures (TLSs), alternatively referred to as ectopic lymphoid aggregates, are immune aggregates that develop due to long-lasting inflammation or infection. These structures mainly comprise B cells, T cells, dendritic cells (DCs), neutrophils, macrophages, lymphatic vessels, and high endothelial venules (HEVs) (7,8). Compared to secondary lymphoid organs (SLOs), TLSs exhibit greater diversity in cellular composition and function. TLSs are capable of antigen processing and presentation, which facilitates the transmission of signals for the proliferation and differentiation of a variety of immune cells (9). In TLSs, loose aggregates of B cells and T cells develop into mature germinal centers, and DCs present antigens near the tumor tissue. After T cells and B cells bind to the antigens and react, they produce adaptive immune responses through outward diffusion of lymphoid tissue, effector T cells, memory B cells, and antibodies, thereby exerting antitumor effects (8,10-12). Studies have shown that TLSs prevent the interaction between programmed cell death protein 1 (PD-1) and its ligands programmed cell death 1 ligand 1 (PD-L1) and programmed cell death 1 ligand 2 (PD-L2) in the tumor microenvironment (TME), helping to restore depleted and cytotoxic-specific T cells. For example, it has been found that for patients with soft tissue sarcoma (13,14), head and neck squamous cell carcinoma (15,16), melanoma (17), intrahepatic cholangiocarcinoma (18), hepatocellular carcinoma (19), or bladder cancer (20), the higher the density of TLSs is after immunotherapy, the better the response and prognosis to the treatment of immune checkpoint block (ICB). The evidence indicates that TLSs can serve as indicators for cancer treatment and prognosis. Hence, evaluating the correlation between TLS and colon cancer may help to clarify the molecular mechanisms underlying the development of colon cancer and provide therapeutic targets for its treatment.

The exact process by which TLSs and colon cancer interact and the manner in which immune cells communicate within TLSs remain unclear. The novelty of this study lies in stratifying colon cancer patients into two TLS-based molecular clusters using TLS-related genes (TRGs), followed by a comprehensive analysis of their correlations with clinical features, tumor immune microenvironment composition, and potential biological pathways. In addition, we constructed two gene clusters based on the differentially expressed genes (DEGs) between the two TLS clusters and established a four-gene prognostic risk model based on these two gene clusters to further improve the accuracy of predicting patient prognosis and treatment sensitivity. We aimed to determine whether the TLS-related signature can stratify prognosis, predict immune microenvironment characteristics, and inform treatment for patients with colon cancer. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1717/rc).

Methods

Collection and classification of multiomics data

A recent study (21) was used to identify 39 relevant TRGs. Data on the sequencing of the transcriptome and clinical details of 458 patients with colon cancer (including stage, age, and gender) were obtained from The Cancer Genome Atlas (TCGA) (https//portal.gdc.cancer.gov). The somatic mutation data for patients with colon cancer were obtained from the TCGA database, while the copy number change matrix was obtained from the University of California Santa Cruz Xena (UCSC Xena) browser (http//xena.ucsc.edu/). For subsequent analysis, we selected 443 colon cancer samples from the GSE39582 dataset in the GEO database. After excluding samples with missing survival data, the gene set file was annotated via Strawberry Perl. To ensure uniform data format, the R packages “limma” and “sva” (The R Foundation for Statistical Computing, Vienna, Austria) were used to merge the TCGA and GEO transcriptome data, thereby eliminating any batch effects. The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Identification and analysis of TRGs

To determine the differential expression of TRGs, we conducted a differential analysis using the Wilcoxon rank-sum test. To further identify genes associated with prognosis-related TRGs, Kaplan-Meier analysis and univariate Cox regression were employed. In addition, we plotted a bar chart and copy number loop chart of the TRGs’ copy number mutation frequency based on the copy number change matrix.

Establishment and analysis of TLS clusters in colon cancer

To identify TLS clusters, we used a nonnegative matrix factorization (NMF) algorithm relying on 39 TRGs. To examine the survival advantage of these TLS clusters, the log-rank test and Kaplan-Meier survival curve were applied to analyze the differential overall survival. We visualized the expression trend of the TRGs for the TLS clusters with heatmaps using the R package “pheatmap”. Pathway differences between various clusters were analyzed via gene set variation analysis (GSVA). Subsequently, we assessed the infiltration of immune cells with various clusters through single-sample gene set enrichment analysis (ssGSEA). In addition, we conducted principal component analysis (PCA) using the R packages “ggplot2” and “limma” to ensure that our clusters were reasonable.

Gene crossover, enrichment analysis, and establishment of gene clusters

We further employed the R software package “limma” to identify the DEGs between TLS clusters A, B, and C with an absolute log₂ fold change (FC) greater than 0.585 and an adjusted P value less than 0.05. Finally, we visualized the results in a Venn diagram. Following this, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed on DEGs with the R packages “clusterProfiler”, “org.Hs.eg.db”, “enrichplot”, and “ggplot2”. We then employed univariate Cox regression analysis to identify the DEGs that had prognostic significance, completed gene clustering accordingly, and conducted prognostic comparisons and other relevant analyses to evaluate the gene clusters.

Generation and verification of the prognostic model

We performed univariate Cox, multivariate Cox, and least absolute shrinkage and selection operator (LASSO) regression analysis on genes that defined gene clusters to select TLS signature genes (TSGs) to construct prognostic models. We then calculated the TSG risk scores for each patient according to the following formula: risk score = sum (expression level of each gene × corresponding coefficient).

Subsequently, patients were divided into high- and low-risk groups based on the median risk score of the training group. We used R packages “ggalluvial”, “ggplot2”, and “dplyr” to draw Sankey diagrams and visualize the connections between different TLS clusters, gene clusters, and patients within the prognostic model. To determine whether different TSG risk scores can effectively distinguish TLS clusters and gene clusters in patients, we conducted the relevant analyses using the R packages “limma” and “ggpubr”. Following this, the R packages “survival” and “survminer” were used to plot the survival curves of the high- and low-risk patients in the training, test, and overall cohorts to evaluate whether there were differences in survival status between the high- and low-risk patients. In addition, we used the R package “timeROC” to plot the receiver operating characteristic (ROC) curves of patients in the training, test, and overall cohorts to evaluate the accuracy of the model in predicting patient prognosis. We further used the R packages “survival”, “regplot”, and “rms” to plot nomograms for predicting the 1-year, 3-year, and 5-year survival rates of patients and tested the accuracy using column charts. Finally, we used the R packages “reshape2”, “tidyverse”, “ggplot2”, “ggpubr” and “ggExtra” to plot the heatmap for the correlation between TSGs and immune cells.

Evaluation of microsatellite status, TME, tumor mutation burden (TMB), and immune checkpoint therapy

The R packages “plyr”, “ggplot2”, and “ggpubr” were used to evaluate the microsatellite status of patients in the high- and low-risk groups and analyzed whether there were differences in microsatellite status between these groups. TME features were quantified using the estimation of stromal and immune cells in malignant tumor tissues using expression data (ESTIMATE) algorithm, with ImmuneScore, StromalScore, and ESTIMATEScore being produced for each sample. Using gene expression signatures, the algorithm estimated the quantities of stromal cells and immune cells present in malignant tumor tissues. The estimation process was performed with the R package “estimate”, after which we conducted somatic mutation analysis to divide all patients into two groups, high- and low-TMB groups, and generated waterfall plots for visualization via the R package “maftools”. To further determine the impact of high and low TMB on patient survival status, as well as the comprehensive effect of TMB and risk scores on patient survival status, we used R packages “survival” and “survminer” to plot Kaplan-Meier curves. Finally, we analyzed the relationship between immune checkpoint and reactivity in patients in the high- and low-risk groups with the R package “ggpubr”.

Analysis of drug sensitivity and correlation analysis of cancer stem cells (CSCs)

The R package “pRRophetic” was used to predict the drug sensitivity and evaluate drugs’ half-maximal inhibitory concentration values in the different risk groups. The correlation analysis of CSCs was conducted with the R packages “limma”, “ggplot2”, “ggpubr”, and “ggExtra”.

Statistical analysis

The software used in this study included R v. 4.4.1 and Strawberry Perl v. 5.30.0.1. A two-sided P value <0.05 (marked with “*” in the figures) was considered statistically significant.

Results

Differential expression and genetic differences in TRGs

First, we identified 39 TRGs reported in previous research (21) (CCL2/3/4/5/8/18/19/20/21, CXCL8/9/11/13, CD200, FBLN7, ICOS, SGPP2, SH2D1A, TIGIT, PDCD1, CD4, CCR5, CXCR3, CSF2, IGSF6, IL2RA, CD38, CD40, CD5, MS4A1, SDC1, GFI1, IL1R1, IL1R2, IL10, IRF4, TRAF6, STAT5A, and TNFRSF17). Subsequently, we examined the distinct expression of TRGs in patients with colon cancer by comparing tumors with normal tissues in the TCGA database. The findings indicated that genes including CCL19, CCL2, and CD4 exhibited elevated expression levels in normal tissues compared to tumor tissues. Conversely, other genes, including CCL20, CCL4, and CXCL8, displayed higher expression levels in tumor tissues than in normal tissues. We further analyzed the frequency of copy number variations (CNVs), and the results showed that all TRGs exhibited a higher frequency of CNVs. Additionally, the exact chromosomal locations at which TRGs undergo CNV alterations were determined (Figure S1).

We used Kaplan-Meier analysis and Cox analysis to evaluate the prognostic impact of TRGs on patients with colon cancer, which indicated that the prognosis of overall survival was linked to 29 genes (Figure S2). It was found that a higher expression of the following genes conferred a survival advantage: CCL19, CCL20, CCR5, CD38, CD40, CXCL8, CXCL9, CXCL11, CXCL13, CXCR3, FBLN7, GFI1, ICOS, IGSF6, IL2RA, IRF4, MS4A1, SGPP2, TIGIT, and TNFRSF17. Moreover, prognosis worsened as the expression of CCL2, CCL3, CCL4, CCL8, CD200, IL1R1, IL10, PDCD1, and STAT5A increased.

Partial subtypes, GSVA, and ssGSEA based on TLS clusters

When k=3, the consistency matrix of subtypes was optimal. Therefore, we classified all patients into three main TLS clusters (cluster A, cluster B, and cluster C) (Figure 1A). The results of PCA showed significant differences in the three clusters (Figure 1B). Next, we conducted survival analysis on the three TLS clusters, and the Kaplan-Meier curve showed significant differences in the prognosis between different clusters. Survival curves indicated that patients in cluster A had the best prognosis, while those in cluster C had the worst prognosis (Figure 1C). A heatmap based on TLS clusters was constructed by combining the stages of patients with colon cancer in the TCGA and GSE39582 cohorts. The results showed that the expression levels of CXCL9/11/13, CD40, and other genes in cluster A were higher than those in cluster B and cluster C (Figure 1D). GSVA identified the 20 biological functions and pathways that differed most significantly between TLS clusters A, B, and C (Figure 1E-1J). The box plot in Figure 1K shows the differences in immune cell infiltration between the three TLS clusters. The infiltration degree of activated B cells, activated CD4/CD8 T cells, immature B cells, myeloid-derived suppressor cells, and regulatory T cells in cluster A was significantly higher than that in cluster B and cluster C.

Figure 1 Establishment and related testing of the TLS cluster. (A) Definition of the consistent matrix heatmap and its related regions for three clusters (k=3). (B) Principal component analysis indicating significant differences between clusters A, B, and C. (C) Survival curve showing a survival difference between TLS clusters A, B, and, C. (D) The difference in TRG expression levels between TLS clusters A, B, and C. (E-J) GSVA heatmap showing the differences in biological functions and pathways between the two clusters of three different subtypes. (K) The difference in the infiltration of immune cells between clusters A, B, and C. *, P<0.05; ***, P<0.001. GSVA, gene set variation analysis; KEGG, Kyoto Encyclopedia of Genes and Genomes; MDSC, myeloid-derived suppressor cell; PC, principal component; TCGA, The Cancer Genome Atlas; TLS, tertiary lymphoid structure; TRG, TLS-related gene.

Establishment and analysis of gene clusters

To further examine the involvement of TLS in colon cancer, we conducted a differential gene analysis on TLS clusters A, B, and C. This analysis yielded a set of 144 genes (Figure 2A) that exhibited differential expression between the TLS clusters. Subsequently, we employed GO and KEGG enrichment analyses to identify potential functions and pathways associated with these genes. The analysis revealed that the differential genes in TLS clusters A, B, and C are primarily enriched in the biological processes (BP) related to the regulation of leukocyte proliferation, regulation of T-cell activation, and immune response-activating signaling pathway; in the cellular components (CC) related to endocytic vesicle, major histocompatibility complex (MHC) class II protein complex, and external side of plasma membrane; and in the molecular functions (MF) related to immune receptor activity, chemokine activity, and MHC protein complex binding (Figure 2B,2C). The KEGG analysis revealed significant enrichment of the intestinal immune network for immunoglobulin A (IgA) production, viral protein interaction with cytokine and cytokine receptor, and T helper 1 and 2 cell differentiation (Figure 2D,2E). Gene clusters (Figure 2F) were subsequently established with these DEGs. According to the genotyping results, the optimal number of gene clusters for our gene cluster typing was determined to be k=3. Subsequently, a survival analysis comparing the three gene clusters was conducted, revealing notable differences in survival rates between the three cohorts (Figure 2G). In addition, a heatmap was created to visualize the expression of genes in these three gene clusters. From the heatmap in Figure 2H, it can be intuitively observed that gene cluster A has a higher level of gene expression compared to gene clusters B and C. According to our findings, most TRGs exhibited a high expression in gene cluster A, including CCL2/3/4/5/8/18/19/21, CCR5, CD4/5/38/40/200, CSF2, CXCL8/9/11/13, CXCR3, GFI1, ICOS, IGSF6, IL10, IL2RA, IRF4, MS4A1, PDCD1, SGPP2, SH2D1A, STAT5A, TIGIT, and TNFRSF17 (Figure 2I).

Figure 2 Correlation analysis of differential genes in TLS clusters and the establishment and related analysis of gene clusters. (A) Venn diagram showing the number of differential genes between TLS clusters A, B, and C. (B,C) The enrichment of BP, CC, and MF items in GO indicating the biological functions of the differentially expressed genes. (D,E) The KEGG pathway enrichment analysis revealed the main enrichment pathways of the differentially expressed genes. (F) Clustering matrix. (G) Kaplan-Meier survival curves of gene clusters A, B, and, C. (H) Heatmap of gene expression in gene clusters A, B, and, C. (I) The expression of TRGs of clusters A, B, and, C. *, P<0.05; **, P<0.01; ***, P<0.001. BP, biological processes; CC, cellular components; GO, Gene Ontology; IgA, immunoglobulin A; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, molecular functions; MHC, major histocompatibility complex; TCGA, The Cancer Genome Atlas; TLS, tertiary lymphoid structure; TRG, TLS-related gene.

Construction of prognostic models

Initially, 898 patients were randomly assigned to two sets at a 1:1 ratio: the training set (n=449) and the test set (n=449). Following this, we analyzed 144 gene cluster–related genes by using univariate Cox regression analysis. Subsequently, we used LASSO regression analysis to minimize gene overfitting in the process of generating a signature, resulting in the identification of six crucial genes (Figure 3A,3B). Next, four TSGs including CXCL10, VSIG4, GZMB, and MMP12 were identified for constructing a prognostic model through multivariate Cox regression analysis. The TSG risk score formula was then established as follows: risk score = exp (CXCL10 × –0.1284) + exp (VSIG4 × 0.2931) + exp (GZMB × –0.0969) + exp (MMP12 × –0.1250). Using the risk score formula, we categorized all patients into low-and high-risk groups based on the median score of the training set.

Figure 3 Construction and detection of the prognostic signature. (A,B) Construction of a risk prediction signature based on the LASSO regression analysis. (C) Sankey diagram showing the connection between the TLS cluster, gene cluster, and prognostic signature. (D) The difference in patient risk scores between TLS clusters A, B, and C. (E) The difference in patient risk scores between gene clusters A, B, and C. (F) Differences in the expression of TRGs between patients in the high- and low-risk groups. (G-I) Kaplan-Meier curves for the overall, test, and training cohorts. *, P<0.05; **, P<0.01; ***, P<0.001. LASSO, least absolute shrinkage and selection operator; TLS, tertiary lymphoid structure; TRG, TLS-related gene.

Predictive significance of the signature

We visualized the relationship between the TLS clusters, gene clusters, and the prognostic model using a Sankey diagram (Figure 3C). We identified variability in the risk scores across different TLS clusters and gene clusters (Figure 3D,3E). In addition, patients in the high- and low-risk groups showed significant differences in TRG expression. The low-risk group had a high expression of CCL3/4/5/8/18/20, CCR5, CD38/40, CSF2, CXCL8/9/11/13, GFI1, ICOS, IGSF6, IL2RA, IRF4, SGPP2, SH2D1A, TIGIT, and TNFRSF17; meanwhile, the high-risk group had a high expression of CCL21, FBLN7, IL1R1, and SDC1 (Figure 3F). We then conducted a prognostic examination of the data for each of the training, test, and overall cohorts. High-risk patients had a worse prognosis compared to low-risk patients (Figure 3G-3I). We further plotted the 1-year, 3-year, and 5-year ROC curves to evaluate the accuracy of TSGs in predicting patient outcomes. The accuracy for predicting the 1-year survival rates for patients in the overall, test, and training cohorts based on TSGs was 65.0%, 63.3%, and 66.9%, respectively; that for the 3-year survival rates was 67.5%, 68.4%, and 66.5%, respectively; and that for the 5-year survival rates was 64.2%, 64.2%, and 63.9%, respectively (Figure 4A-4C). As shown in Figure 4D-4F, the expression of CXCL10, GZMB, and MMP12 in the model decreased with increasing risk score, while the expression of VSIG4 increased with increasing risk score. As shown in Figure 4G-4I, patient survival also decreased with increasing risk score. In addition, we established a nomogram that can predict the 1-, 3-, and 5-year survival rates of patients (Figure 4J), and the results of the calibration curve confirmed the feasibility of our model (Figure 4K). Finally, we demonstrated the correlation between risk score and immune cell infiltration through heatmaps. Risk score was positively correlated with the infiltration of T regulatory cells (Tregs), B cells memory, resting DCs, M2 macrophages, resting mast cells, and neutrophils; however, it was negatively correlated with activated DCs, M1 macrophages, resting natural killer cells, activated memory CD4 T cells, CD8 T cells, and follicular helper T cells (Figure 4L).

Figure 4 The value of the prognostic signature, establishment of the nomogram, and analysis of immune cell correlation. (A-C) The accuracy of the ROC curve in predicting the 1-, 3-, and 5-year survival status of patients. (D-I) The distribution of TSG expression status, risk score, and survival status in the overall, test, and training cohorts. (J) Establishment of nomograms for predicting the prognosis of patients. (K) Calibration curve for validating the accuracy of the nomogram. (L) The correlation between immune cells and four TSGs. *, P<0.05; **, P<0.01; ***, P<0.001. AUC, area under the curve; NK, natural killer; OS, overall survival; ROC, receiver operating characteristic; TLS, tertiary lymphoid structure; TSG, TLS signature gene.

Analysis of microsatellite status, TME, TMB, and immune checkpoints

By analyzing the microsatellite status of patients in the high- and low-risk groups, we found that patients in the high-risk group exhibited higher microsatellite instability (MSI) status (Figure 5A), and the risk score of patients with MSI-high status was different from that of patients with MSI-low or microsatellite-stable status (Figure 5B) (P<0.05). After examining the TME, we observed a significant difference (P<0.05) in stromal and immune cell scores between patients with high- and low-risk scores (Figure 5C). Subsequently, we examined the occurrence rates of somatic mutations in both the high- and low-risk groups. As shown in Figure 5D,5E, APC, TP53, TTN, KRAS, and PIK3CA were the five most mutated genes in both the high-and low-risk groups. In addition, significant differences in patient survival rates were observed between the high-TMB group and the low-TMB group, and the group with a high TMB and high risk score had the worst prognosis (Figure 5F,5G). Patients were then categorized into four following subcategories depending on anti-PD-L1 and anti-CTLA-4 treatment: CTLA-4-negative and PD-1-negative, CTLA-4-negative and PD-1-positive, CTLA positive and PD-1-negative, and CTLA-4-positive and PD-1-positive (Figure 5H-5K). In all four subgroups, it was observed that patients with a low risk score had a greater response to immunotherapy than the high-risk score group.

Figure 5 Analysis of the TME, somatic mutation, immune checkpoint inhibitor treatment, and CSCs. (A) Microsatellite status of patients in the high- and low-risk groups. (B) Significant difference in risk between patients with MSI-H status and those with MSI-L or MSS status. (C) Comparison of ESTIMATE scores, stromal scores, and immune scores between the high-risk and low-risk groups. (D,E) Somatic mutation in the high-risk and low-risk groups. (F,G) Differences in survival between patients with different TMB and risk scores. (H-K) Violin chart showing the differences in immune checkpoint inhibitor treatment between patients with high and low risk scores in the four subgroups. (L) Correlation between CSCs and risk score. *, P<0.05; ***, P<0.001. CSCs, cancer stem cells; ESTIMATE, estimation of stromal and immune cells in malignant tumor tissues using expression data; H, high; L, low; MSI, microsatellite instability; MSS, microsatellite-stable; RNAss, RNA stem score; TMB, tumor mutation burden; TME, tumor microenvironment.

Drug sensitivity analysis and CSC correlation analysis

To evaluate the relationship between TSG scores and CSCs in colon cancer, we conducted a correlation analysis. The correlation between the TSG score and the CSC index is shown in Figure 5L. Our results indicated that the TSG score was negatively correlated with the CSC index (R=–0.41; P<0.05), and thus colon cancer cells with lower TSG scores may have lower levels of cell differentiation. According to the results of the drug sensitivity analysis, 43 drugs showed high sensitivity in the low-risk group (Figure S3), while 33 drugs showed high sensitivity in the high-risk group (Figure S4).

Discussion

TLSs, which are collections of immune cells, play a crucial role in activating T cells and B cells to initiate and sustain immune responses against cancer cells. Recently, Helmink et al. found that TLSs can promote immune checkpoint inhibition response and thus play a crucial role in the antitumor immune response of PD-1/PD-L1 immune checkpoint inhibitors (22). Nevertheless, no research has analyzed the relationship between TLS and the prognosis and targeted treatment of patients with colon cancer. According to our findings, low TSG scores are associated with higher survival rates, while high TSG scores are associated with lower survival rates. In the future, the research and development of TLSs may improve the survival rate of patients with colon cancer and help identify the differences in immune status and survival.

For this study, we initially integrated multiple independent datasets and 39 validated TRGs. Next, we merged and analyzed the TLS data with mutation and copy data. The findings revealed that all TRGs, including CCL21 and CCL8, exhibited higher CNV gain frequencies. CCL8 participates in immune response by interacting with M2 macrophages (23). In cervical cancer, CCL8 induces the recruitment of tumor-associated macrophages by interacting with ZEB1 in hypoxic cancer cells (24). In breast cancer, CCL8 enhances tumor cell activity and promotes tumor metastasis by regulating the TME, and its oncogenic effect is inhibited after macrophage depletion, indicating that CCL8 may promote tumor progression by recruiting macrophages (25). These findings suggest that further exploration of CCL8’s role in colon cancer may be productive. We employed NMF clustering to divide all patients into three TLS clusters according to the correlation between TRGs and the prognosis of patients with colon cancer. We found that patients with TLS cluster A had a better prognosis. On this basis, we analyzed the DEGs between TLS clusters A, B, and C and constructed gene clusters accordingly. Similarly, we also found survival differences between across patients with gene clusters A, B, and C.

Through analyses including univariate, LASSO Cox, and multivariate regression, we identified four TSGs (CXCL10, VSIG4, GZMB, and MMPI2) and their corresponding risk scores. Subsequently, we developed an innovative prognostic risk profile by using four TSGs to enhance the accuracy of prognostic prediction for patients with colon cancer. Furthermore, we evaluated this signature through difference, survival, ROC curve, heatmap, and point map analysis. The findings indicated that TSGs may serve as prognostic risk genes in individuals diagnosed with colon cancer. In the constructed prognostic model, we observed notable variability in TRGs between the high- and low-risk groups. Patients with colon cancer and high risk scores exhibited lower survival rates compared to those with low risk scores and increased expression of TRGs such as SDC1. Chen et al. examined T cells and pancreatic duct cells and found an CCL5-SDC1/4 receptor-ligand interaction occurring between them; they further demonstrated in vitro that CCL5 promotes tumor cell migration through its interaction with SDC1 (26). Liao et al. found that the overexpression of SDC1 in pan-cancer can promote the growth and proliferation of cancer cells and that knocking down the expression of SDC1 can promote the advantageous staging of tumors (27). These findings, which are consistent with our own results, may to some extent indicate a correlation between TRGs and the poor prognosis of patients with high TSG scores.

The correlation between TSGs and risk score and immune cell abundance was assessed via the CIBERSORT algorithm. It was found that the infiltration of naive B cells, activated memory CD4 T cells, and M1 and M2 macrophages was positively correlated with all or a portion of the TSGs and was negatively correlated with risk score. These results indirectly explain why high-risk patients have a poorer prognosis. Chemokines are essential soluble mediators that promote T-cell recruitment into the TME (28). CXCL10 not only promotes the migration of T cells to the tumor site but also drives them to polarize into efficient effector T cells (29,30). Lim et al. found that CXCL10 plays an important role as an ICB-induced antitumor immune mediator in TME; specifically, in a mouse model of non-small cell lung cancer, CXCL10 was found to enhance T-cell infiltration and activation in TME, thereby inhibiting the growth of tumor cells (31). Similarly, Mao et al. found that CXCL10 and Nrf2-β upregulated mesenchymal stem cells rejuvenate T lymphocytes to combat glioblastoma (32). Specifically, overexpression of CLCX10 can promote T-lymphocyte regeneration in the glioblastoma TME, thereby limiting tumor growth. This evidence supports the feasibility of using CXCL10 to establish a prognostic model. In the early stages of malignant tumor development, macrophages play an important role as “scavengers” of malignant tumor cells. However, with the progression and maturation of tumors, the role of macrophages often changes, and their presence can promote a series of processes including tumor cell survival, proliferation, and angiogenesis (33). VSIG4, as a member of the B7 protein family, can be expressed on macrophages and inhibit the proliferation and activation of T cells (34-36). Sazinsky et al. found that the use of anti-VSIG4 antibodies can inhibit M2-like macrophages, ultimately leading to T-cell activation; in their vivo experiments in a mouse model, anti-VSIG4 therapy alone or combined with anti-PD-1 therapy exerted a similar inhibitory effect on tumor growth (37). This constitutes evidence for the relationship between VSIG4 expression, M2-like macrophages, and patient risk score in the prognostic model we constructed. GZMB encodes granzyme B, which is a small-molecule substance that can cleave caspase and cause cell apoptosis (38-40). Yang et al. found that CD4⁺ T cells expressing GZMB have higher cytotoxic activity and can exert stronger killing effects on tumor cells (41). This suggests the auxiliary effect of GZMB in enhancing the killing ability of CD4⁺ T cells. In addition, Wang et al. established an immune-related prognostic model for colorectal cancer and identified MMP12 as a member of the model (42). In summary, our research identified four TRGs as genetic features that can predict prognosis in patients with colon cancer and help guide personalized immunotherapy.

We analyzed the tumor immune microenvironment and found that there was a significant difference in immune score and stromal score. Stromal cells are believed to play an important role in tumor growth, disease progression, and drug resistance (43). Quail et al. reported that stromal cells can release the secretions secreted by fibroblasts to promote the migration of breast cancer cells through Wnt-PCP signal transduction (44). Maller et al. demonstrated that inflammatory stromal cell–mediated collagen cross-linking and matrix sclerosis promote tumor invasion and thus are associated with a poor prognosis (45). Therefore, we compared the TMB between the high- and low-risk groups and found that KRAS and TP53 mutations were more prominent in high- and low-risk patients with colon cancer. KRAS and TP53 oncogenes play a crucial role in the development of cancer. Drosten et al. confirmed a single component of the mitogen-activated protein kinase (MAPK) pathway to be a target for the treatment of KRAS-mutant cancer (46,47). In addition, SDC1 is a protein upregulated by KRAS on the cell surface. The detached extracellular domain HS chain mediates the function of SDC1, which binds to various ligands and affects the growth and reproduction of cancer cells through pathways such as activation of Wnt, flight inhibitory protein long isomers, vascular endothelial growth factor receptor, and MAPK (48-50). We believe that SDC1 guides KRAS in malignant progression, and so the clinical application of SDC1 as an inhibitor may represent an advancement in the treatment of patients with colon cancer.

Immune checkpoint therapy has substantially changed the clinical outcomes of patients with cancer and provided long-lasting clinical benefit. The response rates of different tumor types vary, which requires predictive biomarkers to optimize patient selection to maximize efficacy and minimize toxicity (51). In our study, we discovered that individuals with a low TSG score experienced benefits not just from anti-PD-1/anti-CTLA-4 monotherapy but also from the combination immunotherapy. This finding may provide some assistance for immunotherapy in patients with colon cancer. Different mechanisms lead to heterogeneity within tumors, including gene mutations, microenvironments, and the presence of cancer cell subpopulations, which have enhanced renewal ability and the ability to replicate the heterogeneity found in primary tumors, specifically CSCs (52). Due to the malignant nature of CSCs and cancer cells, as well as the relationship between tumor chemotherapy resistance and metastasis, the eradication of CSCs may be of clinical importance (53). In our study, we combined the TSG score and CSC index values in the context of colon cancer. The findings revealed a negative correlation between the TSG score and CSC index (R=–0.41; P<0.05), suggesting that colon cancer cells with a higher TSG score exhibit more prominent stem cell traits and a reduced level of cell differentiation. This supports heightening the focus on the role of CSCs in the treatment of colon cancer. In addition, we conducted a drug sensitivity analysis and found that 43 drugs had better therapeutic effects on low-risk patients, while 38 drugs had better therapeutic effects on high-risk patients. This finding may provide assistance in guiding clinical medication.

To summarize, we employed TSGs to gain a comprehensive understanding of TLSs role in colon cancer and to assess the collective impact of TLS on the response to immunotherapy in the TME. The correlation between TSG score and immune cell infiltration level is significant, which to some extent indicates the relative maturity of TLSs. The TSG score can also be used to evaluate the prognosis of patients with colon cancer and guide precise medical treatment. However, certain constraints to this study should be noted. To begin, our methods primarily included bioinformatics analysis, with 39 genes being verified as TSGs. Despite this, there remained a deficiency in precise TSG identification, rendering it incomplete and imprecise. The accuracy of the model also needs further improvement. Additionally, the intricate interplay between TLSs and colon cancer, as well as the involvement of diverse immune cells, remains unclear, and thus further investigation into the functional mechanism of TLSs is necessary.

Conclusions

We examined how TLSs participate in colon cancer and developed a novel prognostic model that offers improved accuracy in predicting the prognosis of colon cancer. By using this model, medical professionals will be able to evaluate the overall outlook of patients with colon cancer and gain insights into formulating treatment approaches.

Acknowledgments

This work was supported by Inner Mongolia Natural Science Foundation of China (No. 2021MS08162), and 2025 High-Level Hospital Construction Project of Nanjing Stomatological Hospital, Affiliated Hospital of Medical School, Institute of Stomatology, Nanjing University (No. 0224C055). We gratefully acknowledge the Central Laboratory for providing the necessary equipment for this study. We are also grateful to the participants of the Project Discussion conference for their valuable feedback and discussions.

Footnote

Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1717/rc

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1717/prf

Funding: This work was supported by Inner Mongolia Natural Science Foundation of China (No. 2021MS08162), and 2025 High-Level Hospital Construction Project of Nanjing Stomatological Hospital, Affiliated Hospital of Medical School, Institute of Stomatology, Nanjing University (No. 0224C055).

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1717/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Siegel RL, Kratzer TB, Giaquinto AN, et al. Cancer statistics, 2025. CA Cancer J Clin 2025;75:10-45. [Crossref] [PubMed]
2. Dekker E, Tanis PJ, Vleugels JLA, et al. Colorectal cancer. Lancet 2019;394:1467-80. [Crossref] [PubMed]
3. Bray F, Ferlay J, Soerjomataram I, et al. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin 2018;68:394-424. [Crossref] [PubMed]
4. Arnold M, Sierra MS, Laversanne M, et al. Global patterns and trends in colorectal cancer incidence and mortality. Gut 2017;66:683-91. [Crossref] [PubMed]
5. Labianca R, Beretta GD, Kildani B, et al. Colon cancer. Crit Rev Oncol Hematol 2010;74:106-33. [Crossref] [PubMed]
6. Ait Ouakrim D, Pizot C, Boniol M, et al. Trends in colorectal cancer mortality in Europe: retrospective analysis of the WHO mortality database. BMJ 2015;351:h4970. [Crossref] [PubMed]
7. Lauss M, Donia M, Svane IM, et al. B Cells and Tertiary Lymphoid Structures: Friends or Foes in Cancer Immunotherapy? Clin Cancer Res 2022;28:1751-8. [Crossref] [PubMed]
8. Zou J, Zhang Y, Zeng Y, et al. Tertiary Lymphoid Structures: A Potential Biomarker for Anti-Cancer Therapy. Cancers (Basel) 2022;14:5968. [Crossref] [PubMed]
9. Rossi A, Belmonte B, Carnevale S, et al. Stromal and Immune Cell Dynamics in Tumor Associated Tertiary Lymphoid Structures and Anti-Tumor Immune Responses. Front Cell Dev Biol 2022;10:933113. [Crossref] [PubMed]
10. Vaghjiani RG, Skitzki JJ. Tertiary Lymphoid Structures as Mediators of Immunotherapy Response. Cancers (Basel) 2022;14:3748. [Crossref] [PubMed]
11. Schumacher TN, Thommen DS. Tertiary lymphoid structures in cancer. Science 2022;375:eabf9419. [Crossref] [PubMed]
12. Czepielewski RS, Erlich EC, Onufer EJ, et al. Ileitis-associated tertiary lymphoid organs arise at lymphatic valves and impede mesenteric lymph flow in response to tumor necrosis factor. Immunity 2021;54:2795-2811.e9. [Crossref] [PubMed]
13. Italiano A, Bessede A, Pulido M, et al. Pembrolizumab in soft-tissue sarcomas with tertiary lymphoid structures: a phase 2 PEMBROSARC trial cohort. Nat Med 2022;28:1199-206. [Crossref] [PubMed]
14. Petitprez F, de Reyniès A, Keung EZ, et al. B cells are associated with survival and immunotherapy response in sarcoma. Nature 2020;577:556-60. [Crossref] [PubMed]
15. Liu Z, Meng X, Tang X, et al. Intratumoral tertiary lymphoid structures promote patient survival and immunotherapy response in head neck squamous cell carcinoma. Cancer Immunol Immunother 2023;72:1505-21. [Crossref] [PubMed]
16. Ruffin AT, Cillo AR, Tabib T, et al. B cell signatures and tertiary lymphoid structures contribute to outcome in head and neck squamous cell carcinoma. Nat Commun 2021;12:3349. [Crossref] [PubMed]
17. Chelvanambi M, Fecek RJ, Taylor JL, et al. STING agonist-based treatment promotes vascular normalization and tertiary lymphoid structure formation in the therapeutic melanoma microenvironment. J Immunother Cancer 2021;9:e001906. [Crossref] [PubMed]
18. Ding GY, Ma JQ, Yun JP, et al. Distribution and density of tertiary lymphoid structures predict clinical outcome in intrahepatic cholangiocarcinoma. J Hepatol 2022;76:608-18. [Crossref] [PubMed]
19. Lu Y, Yang A, Quan C, et al. A single-cell atlas of the multicellular ecosystem of primary and metastatic hepatocellular carcinoma. Nat Commun 2022;13:4594. [Crossref] [PubMed]
20. An Y, Sun JX, Xu MY, et al. Tertiary lymphoid structure patterns aid in identification of tumor microenvironment infiltration and selection of therapeutic agents in bladder cancer. Front Immunol 2022;13:1049884. [Crossref] [PubMed]
21. Sautès-Fridman C, Petitprez F, Calderaro J, et al. Tertiary lymphoid structures in the era of cancer immunotherapy. Nat Rev Cancer 2019;19:307-25. [Crossref] [PubMed]
22. Helmink BA, Reddy SM, Gao J, et al. B cells and tertiary lymphoid structures promote immunotherapy response. Nature 2020;577:549-55. [Crossref] [PubMed]
23. Lou X, Zhao K, Xu J, et al. CCL8 as a promising prognostic factor in diffuse large B-cell lymphoma via M2 macrophage interactions: A bioinformatic analysis of the tumor microenvironment. Front Immunol 2022;13:950213. [Crossref] [PubMed]
24. Chen XJ, Deng YR, Wang ZC, et al. Hypoxia-induced ZEB1 promotes cervical cancer progression via CCL8-dependent tumour-associated macrophage recruitment. Cell Death Dis 2019;10:508. [Crossref] [PubMed]
25. Farmaki E, Kaza V, Chatzistamou I, et al. CCL8 Promotes Postpartum Breast Cancer by Recruiting M2 Macrophages. iScience 2020;23:101217. [Crossref] [PubMed]
26. Chen K, Wang Y, Hou Y, et al. Single cell RNA-seq reveals the CCL5/SDC1 receptor-ligand interaction between T cells and tumor cells in pancreatic cancer. Cancer Lett 2022;545:215834. [Crossref] [PubMed]
27. Liao S, Liu C, Zhu G, et al. Relationship between SDC1 and cadherin signalling activation in cancer. Pathol Res Pract 2020;216:152756. [Crossref] [PubMed]
28. Kohli K, Pillarisetty VG, Kim TS. Key chemokines direct migration of immune cells in solid tumors. Cancer Gene Ther 2022;29:10-21. [Crossref] [PubMed]
29. Chow MT, Ozga AJ, Servis RL, et al. Intratumoral Activity of the CXCR3 Chemokine System Is Required for the Efficacy of Anti-PD-1 Therapy. Immunity 2019;50:1498-1512.e5. [Crossref] [PubMed]
30. Mikucki ME, Fisher DT, Matsuzaki J, et al. Non-redundant requirement for CXCR3 signalling during tumoricidal T-cell trafficking across tumour vascular checkpoints. Nat Commun 2015;6:7458. [Crossref] [PubMed]
31. Lim RJ, Salehi-Rad R, Tran LM, et al. CXCL9/10-engineered dendritic cells promote T cell activation and enhance immune checkpoint blockade for lung cancer. Cell Rep Med 2024;5:101479. [Crossref] [PubMed]
32. Mao J, Li J, Chen J, et al. CXCL10 and Nrf2-upregulated mesenchymal stem cells reinvigorate T lymphocytes for combating glioblastoma. J Immunother Cancer 2023;11:e007481. [Crossref] [PubMed]
33. Mantovani A, Allavena P, Marchesi F, et al. Macrophages as tools and targets in cancer therapy. Nat Rev Drug Discov 2022;21:799-820. [Crossref] [PubMed]
34. Greisen SR, Aspari M, Deleuran B. Co-Inhibitory Molecules - Their Role in Health and Autoimmunity; Highlighted by Immune Related Adverse Events. Front Immunol 2022;13:883733. [Crossref] [PubMed]
35. Jung K, Kang M, Park C, et al. Protective role of V-set and immunoglobulin domain-containing 4 expressed on kupffer cells during immune-mediated liver injury by inducing tolerance of liver T- and natural killer T-cells. Hepatology 2012;56:1838-48. [Crossref] [PubMed]
36. Zhou X, Khan S, Huang D, et al. V-Set and immunoglobulin domain containing (VSIG) proteins as emerging immune checkpoint targets for cancer immunotherapy. Front Immunol 2022;13:938470. [Crossref] [PubMed]
37. Sazinsky S, Zafari M, Klebanov B, et al. Antibodies Targeting Human or Mouse VSIG4 Repolarize Tumor-Associated Macrophages Providing the Potential of Potent and Specific Clinical Anti-Tumor Response Induced across Multiple Cancer Types. Int J Mol Sci 2024;25:6160. [Crossref] [PubMed]
38. Hameed A, Lowrey DM, Lichtenheld M, et al. Characterization of three serine esterases isolated from human IL-2 activated killer cells. J Immunol 1988;141:3142-7.
39. Krähenbühl O, Rey C, Jenne D, et al. Characterization of granzymes A and B isolated from granules of cloned human cytotoxic T lymphocytes. J Immunol 1988;141:3471-7.
40. Poe M, Blake JT, Boulton DA, et al. Human cytotoxic lymphocyte granzyme B. Its purification from granules and the characterization of substrate and inhibitor specificity. J Biol Chem 1991;266:98-103.
41. Yang X, Wu J, Fan L, et al. Single-Cell Analysis Identifies Distinct Populations of Cytotoxic CD4(+) T Cells Linked to the Therapeutic Efficacy of Immune Checkpoint Inhibitors in Metastatic Renal Cell Carcinoma. J Inflamm Res 2024;17:4505-23. [Crossref] [PubMed]
42. Wang L, Dong X, Yu M, et al. Association between immune-related hub genes CD36, CXCL13, FGFR4, GABBR1, LAMP3, MMP12, and PPM1H and colorectal cancer prognosis. Am J Transl Res 2024;16:216-23. [Crossref] [PubMed]
43. Yoshihara K, Shahmoradgoli M, Martínez E, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun 2013;4:2612. [Crossref] [PubMed]
44. Quail DF, Joyce JA. Microenvironmental regulation of tumor progression and metastasis. Nat Med 2013;19:1423-37. [Crossref] [PubMed]
45. Maller O, Drain AP, Barrett AS, et al. Tumour-associated macrophages drive stromal cell-dependent collagen crosslinking and stiffening to promote breast cancer aggression. Nat Mater 2021;20:548-59. [Crossref] [PubMed]
46. Luo J. KRAS mutation in pancreatic cancer. Semin Oncol 2021;48:10-8. [Crossref] [PubMed]
47. Drosten M, Barbacid M. Targeting the MAPK Pathway in KRAS-Driven Tumors. Cancer Cell 2020;37:543-50. [Crossref] [PubMed]
48. Yu L, Xu H, Zhang S, et al. SDC1 promotes cisplatin resistance in hepatic carcinoma cells via PI3K-AKT pathway. Hum Cell 2020;33:721-9. [Crossref] [PubMed]
49. Zeng L, Zheng W, Liu X, et al. SDC1-TGM2-FLOT1-BHMT complex determines radiosensitivity of glioblastoma by influencing the fusion of autophagosomes with lysosomes. Theranostics 2023;13:3725-43. [Crossref] [PubMed]
50. Yao W, Rose JL, Wang W, et al. Syndecan 1 is a critical mediator of macropinocytosis in pancreatic cancer. Nature 2019;568:410-4. [Crossref] [PubMed]
51. Sharma P, Goswami S, Raychaudhuri D, et al. Immune checkpoint therapy-current perspectives and future directions. Cell 2023;186:1652-69. [Crossref] [PubMed]
52. Nassar D, Blanpain C. Cancer Stem Cells: Basic Concepts and Therapeutic Implications. Annu Rev Pathol 2016;11:47-76. [Crossref] [PubMed]
53. Ishiguro T, Ohata H, Sato A, et al. Tumor-derived spheroids: Relevance to cancer stem cells and clinical applications. Cancer Sci 2017;108:283-9. [Crossref] [PubMed]

(English Language Editor: J. Gray)