Development and validation of a novel ubiquitination-related gene prognostic signature based on tumor microenvironment for colon cancer

Abstract

Background: Colon cancer (CC) is one of the most common cancers with high morbidity globally. Ubiquitination is involved in the characterization of multiple biological processes, and some ubiquitinated enzymes are associated with the prognosis of CC. However, the prognostic model associated with ubiquitination-related genes (URGs) for CC is unavailable.

Methods: Gene expression data, somatic mutations, transcriptome profiles, microsatellite instability status (MSI) status, and clinical information for CC were obtained from The Cancer Genome Atlas (TCGA) dataset. Seven URGs were used for establishing a prognostic prediction model, which was constructed and validated in GSE17538. Besides, genomic variance analysis (GSVA) was used to explore further the differences in biological pathway activation status between the high-risk and low-risk groups. Finally, the single-sample gene set enrichment analysis (ssGSEA) and ESTIMATE algorithm analysis were used to characterize the cellular infiltration in the microenvironment.

Results: A seven-URG prognostic signature was established, based on which patients in the training and test groups could be divided into high-risk and low-risk groups. The results demonstrated that the model has a solid ability to predict the prognosis of CC patients.

Conclusions: We established a prognostic prediction model for CC based on ubiquitination. Then we analyzed the genetic characteristics associated with ubiquitination and the tumor microenvironment (TME) cell infiltration in CC. These results are worthy of exploring new clinical treatment strategies for CC.

Keywords: Colon cancer (CC); prognostic; tumor microenvironment (TME); ubiquitination.

Figure 1

Expression and mutation of URGs in the 399 samples. (A) Differences in expression of URGs between normal and tumor samples. (B) Waterfall plot shows each gene’s mutation profile in each CC sample. The left panel shows genes sorted by mutation frequency which is listed in the right panel. (C) Mutation frequencies for the 56 URGs copy numbers. Red represents copy number increases, and green represents copy number deletions. *P<0.05; ***P<0.01; **P<0.01; ***P<0.001, determined by Kruskal-Wallis test. TME, tumor microenvironment; CNV, copy number variations; URGs, ubiquitination-related genes; ns, no significance; CC, colon cancer.

Figure 2

Biological characteristics and prognosis of each ubiquitination-related subtype. (A) Consistent clustering matrix at k=2. (B) Principal component analysis of the transcriptome profiles of the two ubiquitination-related subtypes. (C) Kaplan-Meier analysis of two different groups of patients with ubiquitination subtypes. (D) Abundance of each infiltrating immune cell in two different ubiquitination-associated subtypes. ns, not significant. (E) Heat map showing the degree of enrichment of the biological pathways of the two related isoforms. *P<0.05; ***P<0.01; **P<0.01; ***P<0.001, determined by Kruskal-Wallis test. URGs, ubiquitination-related genes; ns, no significance; PC1, principal component 1; PC2, principal component 2; MDSCs, myeloid-derived suppressor cells; TCGA, The Cancer Genome Atlas; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Figure 3

LASSO Cox regression analysis was used to construct the final prediction model. (A) Lambda in the LASSO model; dashed lines are drawn at the optimal values using the least criterion. (B) LASSO coefficient profiles of the candidate OS-related URGs with nonzero coefficients determined by the optimal lambda. The coefficient profiles are obtained from the logλ sequence. (C) PPI network of the 7-URG. (D) LASSO Cox regression obtained for the seven OS-related URGs. LASSO, least absolute shrinkage and selection operator; Lambda, selection of the optimal parameter; URGs, ubiquitination-related genes; PPI, protein-protein interaction networks; OS, overall survival.

Figure 4

Prognostic value assessment of risk scores and prediction models for CC patients based on the training group and test group of 7-URGs. (A) Risk score distribution of each patient in the training group; (B) survival analysis of the prediction model for the training group, with Kaplan-Meier curves showing the change in survival over time for patients in the high- and low-risk groups; (C) ROC curves of the training group prediction model at 1, 3, and 6 years; (D) risk score distribution of each patient in the test group; (E) the survival analysis of the test group prediction model with Kaplan-Meier curves shows the change in survival over time for patients in the high- and low-risk groups; (F) ROC curves of the Ttest group prediction model at 1, 3, and 6 years. CC, colon cancer; URGs, ubiquitination-related genes; ROC, receiver operating characteristic.

Figure 5

The differences in biological behavior between high and low-risk populations after the merging of training and test groups were explored via GSVA enrichment analysis. The top 20 biological processes with significant differences according to P value from smallest to largest were visualized using heat maps, with red representing activation pathways and blue representing inhibition pathways. GSVA, genomic variance analysis; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Figure 6

Immune cell infiltration and transcriptome characteristics of the TME in the high-risk and low-risk groups. (A) Violin plot showing TME scores for high and low-risk groups. (B) Abundance of each infiltrating immune cell in the high and low-risk groups was determined by Kruskal-Wallis test. (C) Heat map showing the correlation of each infiltrating immune cell with 7-URG, with red representing positive correlation and blue representing negative correlation. (D) Scatter plots showing the correlation between risk scores and Resting Dendritic cells. (E) Scatter plots showing the correlation between risk scores and Neutrophils. (F) Scatter plots showing the correlation between risk scores and CD8⁺ T cells. (G) Scatter plots showing the correlation between risk scores and Gamma Delta T cells. *P<0.05; ***P<0.01; **P<0.01; ***P<0.001. MDSCs, myeloid-derived suppressor cells; TME, the tumor microenvironment; ns, no significance.

Figure 7

Tumor mutation profiles, TMB and MSI status (MSI-H, microsatellite instability-high; MSI-L/MSS, microsatellite instability-low/stable) between the high-risk and low-risk groups. (A) The waterfall plot shows the mutation status of each gene in each CC sample in the low-risk groups. The left panel shows the genes sorted by mutation frequency, and the mutation frequencies are listed in the right panel. (B) The waterfall plot shows the mutation status of each gene in each CC sample in the high-risk groups. The left panel shows the genes sorted by mutation frequency, and the mutation frequencies are listed in the right panel. (C) TMB of TCGA database samples in the high- and low-risk groups. (D) Composition ratios of MSS, MSI-L, and MSI-H in the high-risk and low-risk groups. (E) Risk scores of MSS, MSI-L, and MSI-H in the three groups. TMB, Tumor Mutation Burden; MSI, microsatellite instability status; CC, colon cancer; TCGA, The Cancer Genome Atlas; MSS, microsatellite stable; MSI, microsatellite instability status; MSI-L, MSI-low; MSI-H, MSI-high.