Multi-Omics Characterization of Tumor Microenvironment Heterogeneity and Immunotherapy Resistance Through Cell States-Based Subtyping in Bladder Cancer

Abstract

Due to the strong heterogeneity of bladder cancer (BC), there is often substantial variation in the prognosis and efficiency of immunotherapy among BC patients. For the precision treatment and assessment of prognosis, the subtyping of BC plays a critical role. Despite various subtyping methods proposed previously, most of them are based on a limited number of molecules, and none of them is developed on the basis of cell states. In this study, we construct a single-cell atlas by integrating single cell RNA-seq, RNA microarray, and bulk RNA-seq data to identify the absolute proportion of 22 different cell states in BC, including immune and nonimmune cell states derived from tumor tissues. To explore the heterogeneity of BC, BC was identified into four different subtypes in multiple cohorts using an improved consensus clustering algorithm based on cell states. Among the four subtypes, C1 had median prognosis and best overall response rate (ORR), which characterized an immunosuppressive tumor microenvironment. C2 was enriched in epithelial-mesenchymal transition/invasion, angiogenesis, immunosuppression, and immune exhaustion. Surely, C2 performed the worst in prognosis and ORR. C3 with worse ORR than C2 was enriched in angiogenesis and almost nonimmune exhaustion. Displaying an immune effective environment, C4 performed the best in prognosis and ORR. We found that patients with just an immunosuppressive environment are suitable for immunotherapy, but patients with an immunosuppressive environment accompanied by immune exhaustion or angiogenesis may resist immunotherapy. Furthermore, we conducted exploration into the heterogeneity of the transcriptome, mutational profiles, and somatic copy-number alterations in four subtypes, which could explain the significant differences related to cell states in prognosis and ORR. We also found that PD-1 in immune and tumor cells could both influence ORR in BC. The level of TGFβ in a cell state can be opposite to the overall level in the tissues, and the level in a specific cell state could predict ORR more accurately. Thus, our work furthers the understanding of heterogeneity and immunotherapy resistance in BC, which is expected to assist clinical practice and serve as a supplement to the current subtyping method from a novel perspective of cell states.

Keywords: bladder cancer; cell states; heterogeneity; immunotherapy resistance; multi-omics; prognosis; subtype; tumor microenvironment.

FIGURE 1

Construction of a single-cell atlas in BC. (A) UAMP plot of a single cell atlas of BC in first identification by Signacx. (B) UAMP plot of a single cell atlas of BC in subclustering. (C) Heatmap plot of the top two markers of cell states in BC.

FIGURE 2

Identification of BC subtypes by clustering cell states. (A) Based on 22 cell states in BC, UROMOL, GSE13507, and GSE48276 cohorts were clustered into four subtypes, and was validated by TCGA-BLCA cohort et al. We only show two heatmaps of UROMOL and TCGA-BLCA cohorts. (B) The line graph demonstrates significant level of cell states among four subtypes in multiple cohorts. *p < .05, **p < .01, and ***p < .001. Abundance means the absolute proportion of cell states.

FIGURE 3

Characterization of clinical features in four subtypes. (A) Kaplan–Meier curves show four subtypes had survival differences in multiple cohorts, and C4 had the best prognosis while C2 has the worst prognosis. The log-rank test p values are shown. (B) The stratified bar chart shows the immune response rates for four subtypes; C1 and C4 both had the highest response rate (30%), and C2 had the lowest immune response rate (20%). C3 had the median immune response rate (23%) (chi-square test, p < .05). (C) In the TCGA-BLCA cohort, the clinical characteristics of C4 were more favorable for survival, whereas C2 was the opposite. C2 had most advanced tumors and most lymph node metastasis (chi-square test, p < .05). (D) In the GSE13507 cohort, the clinical characteristics of C4 were also more favorable for survival, whereas C2 was the opposite (chi-square test, p < .05). (E) In the UROMOL cohort, C4 also had better clinical characteristics than C2 (chi-square test, p < .05). PUNLMP means Papillary Urothelial Neoplasm of Low Malignant Potential.

FIGURE 4

Characterization of TME, inflammation, and immunotherapy heterogeneity by signature in four subtypes. (A) Boxplot shows TME signature of four subtypes; C4 had the highest immune effective and the lowest immunosuppressive and immune exhausted signature, whereas C2 was the opposite. C1 had a relatively high immunosuppressive environment. Kruskal–Wallis (K–W) test was performed among four subtypes. *p < .05, **p < .01, and ***p < .001. (B) Boxplot shows inflammation signature correlated with chemokine, cytokine, interleukin, and TNF family. C2 had the highest tumor-associated inflammation, whereas C4 had the lowest. C1 and C3 had the median level. K–W test was performed among four subtypes. *p < .05, **p < .01, ***p < .001. (C) Boxplot showed immunotherapy signature, immune-resistant related signature and EMT signature were the highest in C4 but the lowest in C2. K–W test was performed among four subtypes. *p < .05, **p < .01, and ***p < .001. (D) The line graph shows the average gene expression of the six immune checkpoints, C2 had the highest level in all six immune checkpoints. (E) Hierarchical bar graph showed the PD1 level in immune cells (ICC) and tumor cells (TCC) in the IMvigor210 cohort (chi-square test, p < .05).

FIGURE 5

Characterization of genomic driver mechanisms by mutation profiles in four subtypes. (A) The vast majority of BC were mutated, and waterfall plots demonstrate significant differences in the mutation rates of the top 20 mutated genes in the overall cohort across the four subtypes as well as significant differences in the mutation rates of top 10 mutated genes in each subtype. (B) Bubble plot shows driver mutations of four subtypes, and C2 had most of the driver mutations, which mostly cause immunosuppression. (C) KRAS, RB1, and EP300 mutations, which are mostly mutated in C2 correlated with worse OS or DFS in TCGA-BLCA cohort. The log-rank test p values are shown. (D) RB1 mutation correlated with bad histologic subtype (chi-squared test, p = .0133), high grade (chi-squared test, p = .0180), high metastasis stage (chi-squared test, p = .0207), and more smoking (chi-squared test, p = .0369).

FIGURE 6

Characterization of genomic driver mechanisms by SCNAs in four subtypes. (A) The four subtypes had significantly different recurrent SCNA regions. Ordinates represent chromosomal regions. (B) The four subtypes had significantly different GISTIC score, mutation frequency, and average amplitude (Kruskal–Wallis test, p < .05). C4 had the highest SCNA level; C2 had the lowest SCNA level though it had the highest mutation frequency and the lowest mutation amplitude. (C) In the C2 subtype, the expression of highly amplified genes tended to correlate with poor prognosis. Cutoff value of high and low expression groups were 75% and 25%.

FIGURE 7

Characterization of associations between the cell-state subtyping and other subtyping methods in BC. Compared with mRNA subtype, basal_squamous was enriched with C1 and C2, and luminal_papillary was enriched with C3 and C4. Compared with miRNA and lncRNA subtypes, C2 had the least miRNA S3 subtype and lncRNA S3 subtype, C1 was the second least, whereas S3 subtype had the best prognosis and low EMT. Conversely C1 and C2 had the most S4 subtype, and the S4 subtype correlated with relatively high EMT and always accompanied the worst prognosis. In comparison with RPPA subtype, C4 was dominated by S1 subtype, which indicated good prognosis and had the least C2, and conversely C4 had the least S4 subtype, which was correlated with the worst prognosis. ND meant that subtype was not determined or unknown.