Characterization of Epithelial-Mesenchymal Transition Identifies a Gene Signature for Predicting Clinical Outcomes and Therapeutic Responses in Bladder Cancer

Abstract

Purpose: The complex etiological variables and high heterogeneity of bladder cancer (BC) make prognostic prediction challenging. We aimed to develop a robust and promising gene signature using advanced machine learning methods for predicting the prognosis and therapy responses of BC patients.

Methods: The single-sample gene set enrichment analysis (ssGSEA) algorithm and univariable Cox regression were used to identify the primary risk hallmark among the various cancer hallmarks. Machine learning methods were then combined with survival and differential gene expression analyses to construct a novel prognostic signature, which would be validated in two additional independent cohorts. Moreover, relationships between this signature and therapy responses were also identified. Functional enrichment analysis and immune cell estimation were also conducted to provide insights into the potential mechanisms of BC.

Results: Epithelial-mesenchymal transition (EMT) was identified as the primary risk factor for the survival of BC patients (HR=1.43, 95% CI: 1.26-1.63). A novel EMT-related gene signature was constructed and validated in three independent cohorts, showing stable and accurate performance in predicting clinical outcomes. Furthermore, high-risk patients had poor prognoses and multivariable Cox regression analysis revealed this to be an independent risk factor for patient survival. CD8+ T cells, Tregs, and M2 macrophages were found abundantly in the tumor microenvironment of high-risk patients. Moreover, it was anticipated that high-risk patients would be more sensitive to chemotherapeutic drugs, while low-risk patients would benefit more from immunotherapy.

Conclusions: We successfully identified and validated a novel EMT-related gene signature for predicting clinical outcomes and therapy responses in BC patients, which may be useful in clinical practice for risk stratification and individualized treatment.

Figure 1

Flowchart of this study. GSEA: gene set enrichment analysis. EMT: epithelial-mesenchymal transition.

Figure 2

Identification of EMT as the leading risk factor for the prognosis of BC patients. The forest plots show that EMT has the highest HR among various cancer hallmarks in the training cohort (a) and multiple cohorts calculated by meta-analysis (b). (c) GSEA plot illustrates that EMT is significantly enriched in BC samples than adjacent normal samples. (d) The heatmap exhibits the distribution of EMT scores and the patient information of grade, M stage, N stage, and T stage in the training cohort. (e) Violin plot displays that patients with higher T stage, N stage, and pathological grade have higher EMT scores. Kaplan–Meier survival curves depict that high-risk patients divided by EMT scores have worse OS (f), DSS (g), and DFI (h) compared with low-risk patients. HR: hazard ratio. BC: bladder cancer. OS: overall survival. DSS: disease-specific survival. DFI: disease-free interval. ∗p < 0.05; ∗∗∗P < 0.001; ns, no significance.

Figure 3

Construction of EMT-related gene signature. Volcano plots show DEGs (a) between BC and adjacent normal samples, and prognostic genes (b) calculated by the univariable Cox regression. Red dots are upregulated genes or risk genes, and blue dots present downregulated genes or protective genes. (c) Venn diagram shows 13 intersected genes between DEGs and prognostic genes. (d) The correlation of 13 EMT-related genes in BC. Upregulated genes and downregulated genes are represented with grey and red colors. Risk genes are described in blue and protective genes are colored in green. The p values of the Cox regression test for 13 genes are represented by the size of circles. Correlation analysis is performed on 13 genes, depicted by the connecting lines between each gene. Red and blue lines present positive and negative correlations. (e) Variable importance plot based on random forest survival analysis for 13 genes. Blue color indicates predictive variables, whereas red color represents nonpredictive variables. (f) Forest plot based on univariable Cox regression analysis shows that four genes in this signature are all significantly associated with overall survival. (g) The coefficient of each gene in the gene signature is depicted by bar plots. DEGs: differentially expressed genes.

Figure 4

The gene signature serves as a robust and promising predictive factor for survival prediction. Kaplan–Meier survival curves illustrate worse survival outcomes in high-risk patients in the TCGA training cohort (a, c, e), validation I cohort (g, i), and validation II cohort (k). ROC curves for 1-year, 3-year, and 5-year survival prediction depict that this gene signature has a promising and stable predictive performance for BC patients in the training cohort (b, d, f), validation I cohort (h, j), and validation II cohort (l).

Figure 5

Functional enrichment and immune cell infiltration analyses based on the gene signature. (a) Bar graph displays significantly enriched pathways in high-risk patients. Violin plots show higher immune score, stromal score, ESTIMATE score, and lower tumor purity in the high-risk patients compared with low-risk patients in the training (b) and merged validation cohorts (c). Box plots depict that CD8+ T cells, Tregs, M1 macrophages, and M2 macrophages are higher infiltrated and B naïve cells are lower infiltrated and in high-risk patients in both training (d) and merged validation (e) cohorts. ∗p < 0.05; ∗∗p < 0.01; ∗∗∗p < 0.001.

Figure 6

EMT-related gene signature predicts chemotherapeutic and immunotherapeutic benefits. Box plots display lower IC₅₀ values of six commonly used chemical drugs in high-risk patients (a-f). Scatter plots illustrate negative correlations between ERS and the estimated IC50 values of cisplatin (g), gemcitabine (h), and doxorubicin (i). Violin plots show that ERS of nonresponders is significantly higher than that of responders in the immunotherapy cohort of IMvigor210 (j) and GSE176307 (l). High-risk patients present significantly lower percentages of responses (CR/PR) and higher percentages of nonresponses (SD/PD) in both IMvigor210 (k) and GSE176307 (m). IC₅₀: half-maximal inhibitory concentration. ERS: EMT-related risk score. CR: complete response. PR: partial response. SD: stable disease. PD: progressive disease. ∗p < 0.05; ∗∗p < 0.01; ∗∗∗p < 0.001.