Classification of Muscle-Invasive Bladder Cancer Based on Immunogenomic Profiling

There is a significant heterogeneity in the immunotherapeutic responsiveness of each muscle-invasive bladder cancer (MIBC) patient. In our research, we aimed to identify a novel classification of MIBC based on immunogenomic profiling that may facilitate the reasonable stratification of prognosis and response to immunotherapy. The single-sample gene-set enrichment analysis (ssGSEA) was used to analyze the RNA-seq data of 29 important immune signatures from TCGA. Unsupervised hierarchical clustering was performed to identify an immunogenomic classification of MIBC. Then, we assessed the features of the classification in prognosis, immune infiltration, tumor-infiltration lymphocytes, HLA genes, and PD-L1 expression level. A total of 399 MIBC samples were included and three subtypes named Immunity_High, Immunity_Medium, and Immunity_Low were identified. The Immunity_High had a significant advantage in overall survival over the Immunity_Medium and Immunity_Low (p = 0.046 and p = 0.024). From Immunity_Low to Immunity_High, immune cell infiltration and stromal content showed an upward trend (p < 0.001). Meanwhile, Immunity_High was associated with a significantly higher proportion of TILs including dendritic cells resting, macrophages M1, mast cells resting, T cells CD4 memory activated, and T cells CD8⁺. And the expression levels of all HLA genes and PD-L1 of Immunity_High were the highest, consistent (p < 0.001). Two hundred ninety eight MIBC patients treated with immunotherapy from the IMvigor210 were included to form an independent validation cohort to verify the robustness of immunogenomic classification and the ability to predict the response to immunotherapy. This classification had potential clinical implications for predicting prognosis and immunotherapeutic responsiveness of MIBC patients.

Introduction

Bladder cancer (BC) is one of the most common malignant tumors of the urinary system, and ~25% of BCs invade bladder muscles to become muscle-invasive bladder cancer (MIBC) (1). Patients with MIBC usually have a higher probability of recurrence and a worse survival prognosis (2). In the past 30 years, cisplatin-based combination chemotherapy has become the standard therapeutic option in unresectable, advanced, and metastasis MIBC with median overall survival (OS) of 14–15 months (3). However, due to the toxicity of the drug, more than 50% of patients do not meet the criteria for using cisplatin and mainly use carboplatin-based chemotherapy instead. In this case, the prognosis of these patients often decreases to <12 months (4, 5). Therefore, more effective and safe therapeutic options of MIBC are urgently needed.

Recently, new immunotherapy designed to inhabit several immune checkpoints, such as programmed cell death 1 receptor (PD-1), PD-ligand-1 (PD-L1), and cytotoxic T lymphocyte-associated protein 4 (CTLA-4), have been applied in a variety of tumors including MIBC and got some encouraging results (6). Since 2016, five immune checkpoint inhibitors (ICIs) for MIBC (including atezolizumab, pembrolizumab, nivolumab, durvalumab, and avelumab) have already received approval from the US Food and Drug Administration (US-FDA) for their significant advantages in prolonging patients' survival time (7, 8). However, immunotherapy is not equally effective for all patients. Several recent large clinical studies have demonstrated that only about 20% of patients can benefit significantly from immunotherapy (9, 10). Given the low response rate and high treatment costs of immunotherapy, we need to develop accurate classification methods to identify patients most suitable for immunotherapy.

Nowadays, advances in next-generation sequencing technology have indicated the relevance of the tumor immune microenvironment in terms of tumor behavior and may predict the response to ICIs more comprehensively (11–13). Several recent studies on BC have reported that enrichment in immune cells infiltration of BC tissue was a strong indicator for immunotherapy response (14, 15). Furtherly, the relative proportions of different types of T cells, macrophages, and other tumor-infiltration lymphocytes (TILs) are closely related to the immune characteristics and prognosis of BC and other tumors (16, 17). Besides, the significant correlation between the expressions of the PD-L1 and human leukocyte antigen (HLA) genes and the immunotherapy response has also been demonstrated by previous studies (8, 18, 19). Meanwhile, classification based on a single biomarker (such as PD-L1 expression) is insufficient to describe disease-related immune characteristics and interactions between important cell types comprehensively. A classification method capable of integrating these several immune markers is still an emerging need. Therefore, in order to find a new and promising classification method for disease-related immune characteristics and immunotherapeutic responsiveness, we tried to apply the single-sample gene-set enrichment analysis (ssGSEA) of the transcriptome data of 29 immune gene sets and identify a new classification of MIBC.

Materials and Methods

Data Collection and Processing

This study was based on two cohorts from two independent datasets of BC. The first cohort was the discovery cohort. The gene expression quantification data and clinical data of discovery were downloaded from The Cancer Genome Atlas Urothelial Bladder Carcinoma (TCGA-BLCA) (https://portal.gdc.cancer.gov/) on Jan 19, 2020. From the TCGA-BLCA cohort, we excluded patients who were with non-MIBC (T0 or T1), without complete follow-up records or without complete RNA-seq data. A total of 399 patients were included in the discovery cohort. In order to screen out the matrix data of mRNAs with gene properties, the gene expression profiles were compared with the human genome annotation GTF file, which was downloaded from the GenCode platform (https://www.gencodegenes.org/). Then, the matrix data of the gene expression values were processed by the Perl software.

The second cohort was the validation cohort. The gene expression quantification data and clinical data were obtained from the supplementary data of IMvigor210. IMvigor210 cohort included patients with metastatic urothelial cancer who were treated with immunotherapy (atezolizumab, an anti-PD-L1 agent) (20). We excluded all patients who lacked immunotherapy response information and follow-up information. A total of 298 patients were included in the validation cohort. The validation cohort was used to verify the robustness of the classification of muscle-invasive bladder cancer based on immunogenomic profiling and whether it can distinguish patients with different responses to immunotherapy. The processing of the validation cohort was based on IMvigor210CoreBiologies, a package for the R statistical computing environment.

The detailed clinicopathologic characteristics of the study population are provided in the Supplementary Tables 1, 2.

Clustering of MIBC Samples

In the discovery cohort, ssGSEA was used to quantify the enrichment levels of 29 immune signatures in each MIBC sample in the form of ssGSEA scores. We obtained the signatures for immune cell types from previous publications (21–23). (Supplementary Table 3) With this method, we could calculate the relative abundance of 29 immune gene sets in each tumor sample. Then, we performed unsupervised hierarchical clustering (complete linkage method) of MIBC based on the ssGSEA scores of the 29 immune signatures to establish an immunogenomic classification of MIBC.

After clustering and distinguishing subtypes of MIBC based on immunogenomic profiling, we studied the characteristics and verified the differences for the several MIBC immunogenomic subtypes from multiple aspects of tumor molecular features and phenotypes.

Survival Analyses

We compared the overall survival (OS) of MIBC patients in different immunogenomic subtypes. The log-rank test was performed to test if the differences in survival time were significant. Kaplan-Meier curves were used to visualize the survival of different subtypes and the differences among them. We also performed a univariate Cox regression analysis to determine the hazard ratio (HR) of each subtype.

Evaluating the Correlation of MIBC Immunogenomic Subtypes With Molecular Features

First, by the “Estimation of STromal and Immune cells in MAlignant Tumor tissues using Expression data” (ESTIMATE) method, the level of immune cell infiltration, stromal content, and tumor purity of each MIBC sample were calculated in the form of different scores (24). Then, we compared the scores of immune cell infiltration, stromal content, and tumor purity among different MIBC immunogenomic subtypes to reveal the relationships between MIBC immunogenomic classification and immune infiltration.

Furtherly, tumor-infiltrating immune cells (TILs) in MIBC samples were assessed by applying the “Cell type Identification By Estimating Relative Subsets Of RNA Transcripts CIBERSOFT” (CIBERSORT) deconvolution algorithm (25). CIBERSORT was used to analyze the relative expression levels of 547 genes in samples according to their gene expression profiles, to predict the proportion of 22 types of TILs in each sample (25). The gene expression signature matrix of 22 tumor-infiltrating immune cells was obtained from the CIBERSORT platform (https://cibersortx.stanford.edu/). We set 100 permutations and P < 0.05 as the criteria for the successful deconvolution of a sample. Then, we compared the proportions of the immune cell subsets between MIBC immunogenomic subtypes using the Mann–Whitney U-test.

Then, we calculated the expression levels of the HLA genes and PD-L1 in each subgroup. The differences in expression levels of HLA genes and PD-L1 among MIBC immunogenomic subtypes subgroups were assessed by Kolmogorov-Smirnov test.

Independent Validation of the Existence of Immunogenomic Subtypes

In the validation cohort, in order to verify the robustness of this classification, we calculated the ssGSEA scores of the 29 immune signatures for each patient and performed unsupervised hierarchical clustering to identify the immunogenomic subtypes again. Then, we compared the immune cell infiltration, stromal content, and tumor purity among different immunogenomic subtypes from the validation cohort by the ESTIMATE algorithm.

Evaluating of the Association of Immunotherapeutic Responsiveness With Immunogenomic Subtypes

In the validation cohort, all the patients had been treated with atezolizumab, an anti-PD-L1 agent. The patient's response to immunotherapy was divided into four situations: Complete Response (CR), Partial Response (PR), Stable Disease (SD), Progressive Disease (PD). We compared the response to immunotherapy among different immunogenomic subtypes to verify whether the immunogenomic subtypes have the potential to predict the immunotherapeutic responsiveness. Furtherly, we used Kaplan–meier method to analyze the OS of patients treated with immunotherapy in different immunogenomic subtypes.

Statistical Analysis

The statistical software R (version 3.6.2), Perl (version 5.24.3), and EmpowerStats (version 2.0) were used in the above analyses. A p < 0.05 was considered statistically significant.

Results

Unsupervised Cluster Analysis Identified Three Immunogenomic Subtypes in MIBC Samples

A total of 399 MIBC samples from TCGA were included in the current analyses. Based on the 29 immune-associated gene sets reported before which represented diverse immune cell types, functions, and pathways, we used the ssGSEA scores to quantify the activity or enrichment levels of immune cells, functions, or pathways in the cancer samples (21–23). Unsupervised hierarchical clustering analysis was then conducted to clearly divide the 399 MIBC samples into three clusters (Figure 1A). According to the heatmap of expression levels of the 29 immune-associated gene sets, we define the three clusters as three MIBC immunogenic subtypes named Immunity_High, Immunity_Medium, and Immunity_Low (Figure 1B). The specific score and immune immunogenomic subtype of each sample were shown in Supplementary Table 4.

FIGURE 1

Figure 1. (A) Unsupervised hierarchical clustering of 399 muscle-invasive bladder cancer in the discovery cohort. (B) Enrichment of 29 immune gene sets in the three clusters yielded by the hierarchical clustering. We named the three clusters Immunity_High, Immunity_Medium, and Immunity_Low.

Correlation of MIBC Immunogenomic Subtypes With Prognosis

The results of survival analyses showed that the MIBC immunogenic subtypes had distinct clinical outcomes. In the Kaplan–meier analysis, the Immunity_High subtype had a significant advantage in OS over the Immunity_Medium (p = 0.046) and Immunity_Low subtypes (p = 0.024), while the difference in OS between the Immunity_Medium and the Immunity_Low subtypes was not significant (Figure 2).

FIGURE 2

Figure 2. Kaplan-meier curves for overall survival (OS) of the three MIBC subtypes, including Immunity_High vs. Immunity_Medium vs. Immunity_Low (A), Immunity_High vs. Immunity_Low (B), and Immunity_High vs. Immunity_Medium (C).

Correlation of MIBC Immunogenomic Subtypes With Tumor Microenvironment

Based on the ESTIMATE algorithm, we successfully obtained the scores of immune cell infiltration, stromal content, and tumor purity. After combining the three types of scores with the expression heatmap of the 29 immune-associated gene sets, we could get some interesting findings. From Immunity_Low to Immunity_High, immune cell infiltration and stromal content scores showed an upward trend, indicating that the Immunity_High had a higher level of immune cell infiltration and stromal content. The changing trend of tumor purity was just the opposite. According to the order of Immunity_Low, Immunity_Medium, and Immunity_High, the tumor purity gradually increased (Figure 3A). The Kruskal–Wallis test confirmed that the differences among the three subtypes were all very significant in immune cell infiltration, stromal content, and tumor purity (p < 0.001; Figures 3B–D).

FIGURE 3

Figure 3. (A) The correlation between immunogenomic classification and the level of tumor immune infiltration was showed based immune score, stromal score, tumor purity score. The three types of scores were calculated by the ESTIMATE algorithm. The violin plots specifically showed the differences in the three subtypes in terms of immune cell infiltration level (B), stromal cell content (C), and tumor purity (D). *p < 0.05, **p < 0.01, ***p < 0.001.

Based on CIBERSORT algorithm, we systematically estimated the proportions of various TILs in MIBC samples (Figure 4A). The results showed the TILs subsets with significantly different proportions in different MIBC immunogenomic subtypes. Compared with Immunity_Low subtype, Immunity_High subtype was associated significantly with higher levels of TILs including dendritic cells resting, Macrophages M1, Mast cells resting, T cells CD4 memory activated, T cells CD8 (p < 0.01). We also found that higher proportions of Macrophages M0, NK cells resting, T cells CD4 naive were significantly associated with Immunity_Low (p < 0.001). Besides, we also explored several ratios between specific TIL like T cells CD8 to T cells CD4, T cells CD8 to Treg, and Macrophage M1 to M2 (Figures 4B–D). The results showed that patients in Immunity_H subtype had higher all three ratios significantly higher than those in Immunity_L subtype (p = 0.024, p = 0.018, p = 0.042; Figure 4B).

FIGURE 4

Figure 4. Specific tumor-infiltration lymphocytes with significantly different proportions among the three MIBC subtypes based on the CIBERSORT algorithm (A). The differences of the ratio of T cells CD8 to T cells CD4 (B), T cells CD8 to Treg (C), and Macrophage M1 to M2 (D) among three MIBC subtypes.

Correlation of MIBC Immunogenomic Subtypes With HLA Genes and PD-L1 Gene

After the differential analysis of HLA genes expression among the three MIBC subtypes, we found that, in the expression levels of HLA genes, the rank was Immunity_High > Immunity_Medium > Immunity_Low, consistently. The differences between the three subtypes were very significant (p < 0.001; Figure 5A).

FIGURE 5

Figure 5. (A) Comparison of expression of multiple HLA genes among three MIBC subtypes. (B) Comparison of PD-L1 expression among three MIBC subtypes. *p < 0.05, **p < 0.01, ***p < 0.001.

Similarly, we calculated the expression level of PD-L1 among the three MIBC subtypes. The results also showed that Immunity_High had the highest PD-L1 expression levels, Immunity_Medium had the medium PD-L1 expression levels, and Immunity_Low had the lowest PD-L1 expression levels (p < 0.001; Figure 5B).

Immunogenomic Subtypes in the Independent Validation Cohort

Same as previous analyses in the discovery cohort, 298 patients were successfully divided into three clusters by hierarchical clustering. Similarly, the heatmap of expression levels of the 29 immune-associated gene sets show the obvious differences of the three clusters, we defined the three clusters as three MIBC immunogenomic subtypes named Immunity_High, Immunity_Medium, and Immunity_Low like before (Figure 6). The specific score and immune immunogenomic subtype of each sample were shown in Supplementary Table 5.

FIGURE 6

Figure 6. (A) Unsupervised hierarchical clustering of 298 muscle-invasive bladder cancer in the validation cohort. (B) Enrichment of 29 immune gene sets in the three clusters yielded by the hierarchical clustering. We named the three clusters Immunity_High, Immunity_Medium, and Immunity_Low.

Through the ESTIMATE algorithm, we have calculated the score of immune cell infiltration, stromal content, and tumor purity for three immunogenomic subtypes, respectively. As in the discovery cohort, Immunity_High had the highest score of immune cell infiltration and stromal content, while Immunity_Low had the highest tumor purity (p < 0.001). This showed that the immunogenomic subtypes applied in the independent validation cohort can also reflect the tumor immune microenvironment of patients to a certain extent (Figure 7).

FIGURE 7

Figure 7. The correlation between immunogenomic classification of the validation cohort and the level of tumor immune infiltration was showed based immune score, stromal score, tumor purity score (A). The violin plots specifically showed the differences in the three subtypes in terms of immune cell infiltration level (B), stromal cell content (C), and tumor purity (D). *p < 0.05, **p < 0.01, ***p < 0.001.

Association Between the Immunogenomic Subtypes and Immunotherapeutic Responsiveness

Previous analyses demonstrated that immunogenomic subtypes were closely related to the immune microenvironment and prognosis of MIBC patients. We subsequently tested the association between the immunogenomic subtypes and immunotherapeutic responsiveness in the independent validation cohort with follow-up information on immunotherapy. The results showed that different immunogenomic subtypes had significant differences in immune response to immunotherapy (p = 0.021; Figure 8A) In patients of Immunity_High, 18.8% reached CR, while in patients of Immunity_Medium, and Immunity_Low, only 7.2, and 3% were able to reach CR. And 32% of patients of Immunity_High had a response to immunotherapy, which was also the most among three subtypes. Subsequent survival analysis results also showed that among all patients who were treated with immunotherapy, patients of Immunity_High had significantly better OS (p = 0.032; Figure 8B). Therefore, the immunogenomic subtypes we identified had the potential clinical value to predict the response and effect of immunotherapy.

FIGURE 8

Figure 8. (A) Proportion of the four types of response to immunotherapy among different immunogenomic subtypes. (B) Kaplan-meier survival curves of MIBC patients treated with immunotherapy. CR, complete response; PR, partial response; SD, stable disease; PD, progressive disease.

Discussion

In our study, we used ssGSEA technology to calculate the relative abundance of 29 immune cells in the discovery cohort of 399 MIBC samples. Subsequently, using unsupervised clustering, we can clearly identify three MIBC subtypes based on the ssGSEA scores: Immunity_Low subtype, Immunity_Medium subtype, and Immunity_High subtype. The results of survival analysis showed that Immunity_High subtype was associated with a significantly better prognosis among the three subtypes. We furtherly confirmed with a series of analyses that there was significant immune heterogeneity among the subtypes and Immunity_High subtype showed stronger immune activity in several aspects. In the validation cohort of 298 MIBC patients treated with immunotherapy from the IMvigor210, we validated the robustness of immunogenomic subtypes and the ability to predict the response to immunotherapy.

Researchers have identified some gene biomarkers associated with cancer immunotherapeutic responsiveness such as PD-L1 expression, tumor mutational burden (TMB), and deficient DNA mismatch repair (26–28). But the robustness of these biomarkers to predicting the efficacy of immunotherapy remains controversial. For example, although a clinical trial of BC has reported the median TMB of responders was two times higher than that of non-responders, there are still many responders with low TMB and non-responders with high TMB (9, 29).

There have been a series of studies that identified various subtypes for MIBC based on genomic profiling. These MIBC subtypes have differences in prognosis and biological characteristics, and show clinical significance in many aspects (30–33). In 2017, researchers analyzed the clinical pathological information and RNA-seq data of 412 bladder cancer patients in TCGA, and divided bladder cancer into 5 molecular subtypes (TCGA subtypes 2017), called Luminal papillary, Luminal infiltrated, Basal squamous, and neuronal subtypes. This study provided some new references for treatment options for bladder cancer. For example, the researchers reported that Luminal infiltrated subtype has a relatively higher lymphocyte infiltration and PD-1/PD-L1 expression, which indicates a higher probability of response to immunotherapy. In fact, this study focused on the differences in epithelial-mesenchymal transition status, carcinoma-in-situ scores, histologic features, and survival (30). Another group of researchers explored the characteristics of immune genes related to four distinct molecular subtypes (TCGA subtypes 2014) in MIBC published in 2014 by The Cancer Genome Atlas Network (TCGA) bladder analysis working group (31). The four subtypes are the basis of the TCGA subtypes 2017 mentioned above (32). They studied the immune characteristics of four TCGA subtypes by analyzing the immune, IFN-γ, IFN-α, and cytotoxic genes. They revealed an increased expression of immune-associated genes in Cluster IV and underactive immune environment in Cluster I, which may help increase the applicability of TCGA subtypes in immunotherapy options (31). But the two kinds of TCGA MIBC subtypes did not focus on more comprehensive immune-related genes or gene sets besides the tumor infiltration lymphocytes and immune checkpoints, as well as lacked follow-up data on immunotherapy. Therefore, MIBC subtypes based on more comprehensive immunogenomics characteristics and with immunotherapy follow-up data are still needed.

Another pan-cancer study based on TCGA data identified 6 immune subtypes for 33 non-hematologic tumors (33). The clustering of immune subtypes was based on the main 5 core immunogenomics characteristics including macrophages/monocytes, overall lymphocyte infiltration, TGF-β response, IFN-γ response, and wound healing. The author subsequently reported the differences of the immune characteristics and prognosis among the six immune subtypes, including differences in macrophage or lymphocyte signatures, Th1:Th2 cell ratio, extent of intratumoral heterogeneity, aneuploidy, extent of neoantigen load, overall cell proliferation, expression of immunomodulatory genes, and prognosis. This study comprehensively studied the immunogenomics characteristics of solid tumors and provided important knowledge in this field. And this study confirmed the feasibility of using immunogenomics features for tumor molecular typing. Similar to this study, the immunogenomic subtypes we obtained were also clustering according to a series of specific immunogenomics characteristics. Meanwhile, our immune subtypes also showed significant differences in multiple immune-related features. Our research focused on MIBC, and furtherly used actual clinical immunotherapy follow-up information to verify the predictive value of the immune response of subtypes, which is the strongness of our research compared with the pan-cancer study.

Through the calculation of ESTIMATE algorithm, we demonstrated that Immunity_High subtype had a significantly higher level of immune cell infiltration, stromal content, and lower tumor purity over Immunity_Medium and Imunity_L subtypes. In fact, researchers had demonstrated that immune cell infiltration levels and tumor purity are inversely proportional, and higher immune cell infiltration levels in multiple tumors that were positively correlated with the responsiveness to ICIs (34–36). Additionally, several studies on BC reported that BC could be divided into “hot” (enriched in immune cells infiltration) and “cold” (lack of immune cells infiltration) subtypes as a strong indicator for immunotherapy response (14, 15). Therefore, the obvious heterogeneity of our immunogenomic classification in terms of immune cell infiltration supports its potential clinical value in predicting immunotherapy responses.

After the correlation between the level of immune cell infiltration and patients' prognosis and treatment responsiveness had been confirmed by more and more studies, one emerging question has attracted attention. Was the proportion of immune infiltration cell types different in tumor samples with different prognosis and immunophenotype? Advances in computational methods had allowed us to take advantage of integrating genomic profiles and the state-of-the-art deconvolution accurately analyzed the relative proportions of various immune infiltration cell subsets, such as the CIBORSORT algorithm used in this study. CIBERSORT technology had been applied in a variety of oncology studies such as lung cancer and breast cancer, with an accuracy not less than flow cytology and immunohistochemistry (36–38). When we used CIBERSORT to calculate the proportions of 22 immune cell subsets in MIBC, we found that the three subtypes with different levels of immune cell infiltration also had significant differences in the proportion of immune infiltration cell types. For example, the proportions of activated CD4⁺ T memory cells, resting dendritic cells and CD8⁺ T cells were the highest in Immunity_High. The basic function of ICIs was to exert antitumor effects by enhancing the function of CD8⁺ T cells, while CD4⁺ T cells and dendritic cells could also play a supporting role in this process (16, 17). Previous studies on bladder cancer, prostate cancer, renal cancer, and colorectal cancer had reported that the level of CD8⁺ T cell infiltration was positively correlated with tumor prognosis and responsiveness to immunotherapy (39–41). Therefore, Immunity_High subtype, which was rich in CD8⁺ T cells, CD4⁺ T cells, and dendritic cells, had a higher possibility of immune response with an advantage in survival. We furtherly found that the ratio of CD8⁺/CD4⁺ and the ratio of CD8⁺/Treg were significantly higher in Immunity_High subtype than those in Immunity_Low subtype. Previous studies on breast cancer reported that lower ratio of CD8⁺/CD4⁺ is associated with poor prognosis, which is consistent with our findings (42). Another study believed the balance between the regulatory cells and effector cells of the body's immune system is essential to maintain an adequate immune response to fight diseases (such as cancer) and at the same time prevent damage to healthy tissues. This balance is caused by effector cells (such as CD8⁺ T cells) and regulatory T cells (or Tregs) maintain a mutual regulatory relationship (43). And consistent with our findings, a low CD8⁺/Treg ratio is associated with a poor prognosis.

Meanwhile, we observed that the proportion of macrophage M1 in Imunity H subtype was significantly higher than that in Imunity L subtype, while the proportion of macrophage M0 significantly decreased. The ratio of M1/M2 was also higher in Immunity_High subtype compared with it in Immunity_Low subtype. Tumor-associated macrophages (TAMs) were important components of the tumor immune microenvironment, including two functional states of M1 (anti-tumor) and M2 (tumor-promoting), and macrophage M0 could be transformed to M1 or M2. Previous studies have shown that the level of polarization of M0 to M1 might be associated with improved prognosis of the tumor, which was consistent with our finding that the prognosis of Immunity_High subtype with significantly increased M1/M0 was better (38, 44, 45).

Besides, we also found that Immunity_High subtype very significantly expressed more HLA genes, which indicted stronger immunogenicity and possible stronger immunotherapeutic responsiveness compared to the other subtypes. The HLA genes encode MHC I and MHC II molecules, and both Class I and Class II molecules present pathogen-derived short peptides to T cells to initiate an adaptive immune response. A previous study on metastatic melanoma suggested that HLA expression may affect the response of ICIs. The results showed that high tumor-specific MHC-I expression in response to anti-CTLA-4 therapy was essential, but anti-PD-1 or immune combination therapy was not effective. Meanwhile, tumor-specific MHC-II expression in response to anti-PD-1 therapy was essential, but it is not effective with anti-CTLA-4 or combined immunotherapy (43). Another study involving more than 1,535 cancer patients treated with ICIs found that the presence of more diverse HLA-I molecules was associated with increased survival after receiving immunotherapy. The authors believed that this might be due to the fact that such MHC molecules could bind more diverse peptides and provide a wider range of tumor antigens to T cells (46, 47).

In the Immunity_High subtype, the expression level of PD-L1 was also the highest among the three groups. The positive correlation between PD-L1 expression and immunotherapy response had been widely reported in many tumors including bladder cancer (8, 18, 19). Combining the above analysis, we could find that Immunity_High had very strong immune activity considering the level of immune infiltration, proportion of specific TILs, HLA richness and expression level, and PD-L1 expression level, and had relatively better clinical outcomes. MIBC patients of Immunity_High subtype were more likely to respond to immunotherapy, especially to anti-PD-1 / PD-L1 treatment.

There was an interesting phenomenon in the survival analysis of the discovery cohort, which the difference of OS between Immunity_medium and Immunity_low was not so significant as the difference of their immune characteristics. In a previous study on triple-negative breast cancer, the researchers also achieved similar results, that patients in Immunity_H had the best prognosis, but there was no significant difference in OS between Immunity_M and Immunity_L (23). We speculated that there may be the following reasons. Although our study found that patients in Immunity_M subtypes had significantly higher immune characteristics such as immune cell infiltration levels than those in immunity_L, the degree of this kind of advantages may not be sufficient to bring about significant improvement in survival. It may require sufficient immune cell infiltration and activation of immune-related pathways to reshape the balance between tumor immune surveillance and immune escape. At the same time, the identification of Immunity_High group with significantly better prognosis and immune response was more clinically meaningful and may help to screen out patients suitable for immunotherapy.

Our research identified a new MIBC classification method with potential clinical value in assessing immunotherapeutic responsiveness and prognosis. But we still had some limitations. First, our data was retrospective from TCGA. Although we had used a series of analyses to prove that this classification was distinguishable and clinically significant in many aspects, prospective clinical studies were still needed for further research. Second, although we had studied the potential value of this classification in distinguishing immunotherapeutic responsiveness in four aspects: tumor immunity infiltration, proportion of specific TILs, HLA expression, and PD-L1 expression, real-world data on immunotherapy of the discovery cohort was still missing. But our research demonstrated the feasibility and significance of this classification method, and pointed out the direction for further clinical studies.

Conclusions

We identified a new MIBC classification based on transcriptome differences of 29 immune signatures in tumor samples from patients with MIBC. This classification had potential clinical implications for predicting prognosis and immunotherapeutic responsiveness.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: TCGA data portal (https://portal.gdc.cancer.gov/).

Ethics Statement

Ethical review and approval was not required for the study on human participants in accordance with the local legislation and institutional requirements. Written informed consent for participation was not required for this study in accordance with the national legislation and the institutional requirements.

Author Contributions

QW and LY designed the study. XZho, SQ, and LN partly designed the study, performed the statistical analysis, and drafted the manuscript. DJ and KJ partly performed the statistical analysis. XZhe partly performed the data collection. All authors reviewed and final approval of the manuscript.

Funding

This work was supported by the National key research and development program of China (Grant No. 2017YFC0908003), National Natural Science Foundation of China (Grant Nos. 81902578 and 81974098), China Postdoctoral Science Foundation (2017M612971), Post-doctoral Science Research Foundation of Sichuan University (2020SCU12041), Post-Doctor Research Project, West China Hospital, Sichuan University (2018HXBH085), and National Clinical Research Center for Geriatrics, West China Hospital, Sichuan University (Z2018C01).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank all the staff and researchers who established and maintained The Cancer Genome Atlas. We thank all the researchers and participants of the IMvigor210 study for their contribution to improving the treatment of bladder cancer.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2020.01429/full#supplementary-material

Abbreviations

MIBC, muscle-invasive bladder cancer; ssGSEA, single-sample gene-set enrichment analysis; OS, overall survival; PD-1, programmed cell death 1 receptor; PD-L1, PD-ligand-1; CTLA-4, cytotoxic T lymphocyte-associated protein 4; ICIs, immune checkpoint inhibitors; TILs, tumor-infiltration lymphocytes; HLA, human leukocyte antigen; TAMs, Tumor-associated macrophages; CIBERSORT, Cell type Identification By Estimating Relative Subsets Of RNA Transcripts; ESTIMATE, Estimation of STromal and Immune cells in MAlignant Tumor tissues using Expression data.

References

1. Siegel RL, Miller KD, Jemal A. Cancer statistics, 2019. CA Cancer J Clin. (2019) 69:7–34. doi: 10.3322/caac.21551

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Hussain MHA, MacVicar GR, Petrylak DP, Dunn RL, Vaishampayan U, Primo NLJ, et al. Trastuzumab, paclitaxel, carboplatin, and gemcitabine in advanced human epidermal growth factor receptor-2/neu–positive urothelial carcinoma: results of a multicenter phase II national cancer institute trial. J Clin Oncol. (2007) 25:2218–24. doi: 10.1200/JCO.2006.08.0994

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Clark PE, Spiess PE, Agarwal N, Bangs R, Boorjian SA, Buyyounouski MK, et al. NCCN guidelines insights: bladder cancer, version 2.2016. J Natl Compr Cancer Netw. (2016) 14:1213–24. doi: 10.6004/jnccn.2016.0131

PubMed Abstract | CrossRef Full Text | Google Scholar

4. De Santis M, Bellmunt J, Mead G, Kerst JM, Leahy M, Maroto P, et al. Randomized phase II/III trial assessing gemcitabine/ carboplatin and methotrexate/carboplatin/vinblastine in patients with advanced urothelial cancer “unfit” for cisplatin-based chemotherapy: phase II–results of EORTC study 30986. J Clin Oncol. (2009) 27:5634–9. doi: 10.1200/JCO.2008.21.4924

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Alfred Witjes J, Lebret T, Compérat EM, Cowan NC, De Santis M, Bruins HM, et al. Updated 2016 EAU guidelines on muscle-invasive and metastatic bladder cancer. Eur Urol. (2017) 71:462–75. doi: 10.1016/j.eururo.2016.06.020

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Del Paggio JC. Immunotherapy: cancer immunotherapy and the value of cure. Nat Rev Clin Oncol. (2018) 15:268–70. doi: 10.1038/nrclinonc.2018.27

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Lobo N, Mount C, Omar K, Nair R, Thurairaja R, Khan MS. Landmarks in the treatment of muscle-invasive bladder cancer. Nat Rev Urol. (2017) 14:565–74. doi: 10.1038/nrurol.2017.82

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Kim HS, Seo HK. Immune checkpoint inhibitors for urothelial carcinoma. Invest Clin Urol. (2018) 59:285–96. doi: 10.4111/icu.2018.59.5.285

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Rosenberg JE, Hoffman-Censits J, Powles T, van der Heijden MS, Balar AV, Necchi A, et al. Atezolizumab in patients with locally advanced and metastatic urothelial carcinoma who have progressed following treatment with platinum-based chemotherapy: a single-arm, multicentre, phase 2 trial. Lancet. (2016) 387:1909–20. doi: 10.1016/S0140-6736(16)00561-4

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Sharma P, Retz M, Siefker-Radtke A, Baron A, Necchi A, Bedke J, et al. Nivolumab in metastatic urothelial carcinoma after platinum therapy (CheckMate 275): a multicentre, single-arm, phase 2 trial. Lancet Oncol. (2017) 18:312–22. doi: 10.1016/S1470-2045(17)30065-7

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Liu XS, Mardis ER. Applications of immunogenomics to cancer. Cell. (2017) 168:600–12. doi: 10.1016/j.cell.2017.01.014

PubMed Abstract | CrossRef Full Text | Google Scholar

12. The Cancer Genome Atlas Network. Comprehensive molecular characterization of human colon and rectal cancer. Nature. (2012) 487:330–7. doi: 10.1038/nature11252

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Carosella ED, Ploussard G, LeMaoult J, Desgrandchamps F. A systematic review of immunotherapy in urologic cancer: evolving roles for targeting of CTLA-4, PD-1/PD-L1, and HLA-G. Eur Urol. (2015) 68:267–79. doi: 10.1016/j.eururo.2015.02.032

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Cheng W, Fu D, Xu F, Zhang Z. Unwrapping the genomic characteristics of urothelial bladder cancer and successes with immune checkpoint blockade therapy. Oncogenesis. (2018) 7:2. doi: 10.1038/s41389-017-0013-7

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Sharma P, Allison JP. The future of immune checkpoint therapy. Science. (2015) 348:56–61. doi: 10.1126/science.aaa8172

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Sato Y, Bolzenius JK, Eteleeb AM, Su X, Maher CA, Sehn JK, et al. CD4+ T cells induce rejection of urothelial tumors after immune checkpoint blockade. JCI Insight. (2018) 3:e121062. doi: 10.1172/jci.insight.121062

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Farhood B, Najafi M, Mortezaee K. CD8 cytotoxic T lymphocytes in cancer immunotherapy: a review. J Cell Physiol. (2019) 234:8509–21. doi: 10.1002/jcp.27782

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Patel SP, Kurzrock R. PD-L1 expression as a predictive biomarker in cancer immunotherapy. Mol Cancer Ther. (2015) 14:847–56. doi: 10.1158/1535-7163.MCT-14-0983

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Zhang B, Yu W, Feng X, Zhao Z, Fan Y, Meng Y, et al. Prognostic significance of PD-L1 expression on tumor cells and tumor-infiltrating mononuclear cells in upper tract urothelial carcinoma. Med Oncol. (2017) 34:94. doi: 10.1007/s12032-017-0941-2

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Mariathasan S, Turley SJ, Nickles D, Castiglioni A, Yuen K, Wang Y, et al. TGFβ attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature. (2018) 554:544–8. doi: 10.1038/nature25501

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Bindea G, Mlecnik B, Tosolini M, Kirilovsky A, Waldner M, Obenauf AC, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity. (2013) 39:782–95. doi: 10.1016/j.immuni.2013.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Liu Z, Li M, Jiang Z, Wang X. A comprehensive immunologic portrait of triple-negative breast cancer. Transl Oncol. (2018) 11:311–29. doi: 10.1016/j.tranon.2018.01.011

PubMed Abstract | CrossRef Full Text | Google Scholar

23. He Y, Jiang Z, Chen C, Wang X. Classification of triple-negative breast cancers based on Immunogenomic profiling. J Exp Clin Cancer Res. (2018) 37:327. doi: 10.1186/s13046-018-1002-1

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Yoshihara K, Shahmoradgoli M, Martinez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. (2013) 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text

25. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. (2015) 12:453–7. doi: 10.1038/nmeth.3337

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Fridman WH, Zitvogel L, Sautès-Fridman C, Kroemer G. The immune contexture in cancer prognosis and treatment. Nat Rev Clin Oncol. (2017) 14:717–34. doi: 10.1038/nrclinonc.2017.101

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Brahmer JR, Drake CG, Wollner I, Powderly JD, Picus J, Sharfman WH, et al. Phase I study of single-agent anti-programmed death-1 (MDX-1106) in refractory solid tumors: safety, clinical activity, pharmacodynamics, and immunologic correlates. J Clin Oncol. (2010) 28:3167–75. doi: 10.1200/JCO.2009.26.7609

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Rizvi NA, Hellmann MD, Snyder A, Kvistborg P, Makarov V, Havel JJ, et al. Cancer immunology. Mutational landscape determines sensitivity to PD-1 blockade in non-small cell lung cancer. Science. (2015) 348:124–8. doi: 10.1126/science.aaa1348

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Alexandrov LB, Nik-Zainal S, Wedge DC, S. Aparicio AJR, Behjati S, Biankin AV, et al. Signatures of mutational processes in human cancer. Nature. (2013) 500:415–21. doi: 10.1038/nature12477

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Robertson AG, Kim J, Al-Ahmadie H, Bellmunt J, Guo G, Cherniack AD, et al. Comprehensive molecular characterization of muscle-invasive bladder cancer. Cell. (2017) 171:540–56.e25. doi: 10.1016/j.cell.2017.09.007

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Ren R, Tyryshkin K, Graham CH, Koti M, Siemens DR. Comprehensive immune transcriptomic analysis in bladder cancer reveals subtype specific immune gene expression patterns of prognostic relevance. Oncotarget. (2017) 8:70982–1001. doi: 10.18632/oncotarget.20237

PubMed Abstract | CrossRef Full Text | Google Scholar

32. The Cancer Genome Atlas Research Network. Comprehensive molecular characterization of urothelial bladder carcinoma. Nature. (2014) 507:315–22. doi: 10.1038/nature12965

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang TH, et al. The immune landscape of cancer. Immunity. (2018) 48:812–30.e14. doi: 10.1016/j.immuni.2018.03.023

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Kümpers C, Jokic M, Haase O, Offermann A, Vogel W, Grätz V, et al. Immune cell infiltration of the primary tumor, Not PD-L1 status, is associated with improved response to checkpoint inhibition in metastatic melanoma. Front Med. (2019) 6:27. doi: 10.3389/fmed.2019.00027

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Senbabaoglu Y, Gejman RS, Winer AG, Liu M, van Allen EM, de Velasco G, et al. Tumor immune microenvironment characterization in clear cell renal cell carcinoma identifies prognostic and immunotherapeutically relevant messenger RNA signatures. Genome Biol. (2016) 17:231. doi: 10.1186/s13059-016-1092-z

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Karn T, Jiang T, Hatzis C, Sänger N, El-Balat A, Rody A, et al. Association between genomic metrics and immune infiltration in triple-negative breast cancer. JAMA Oncol. (2017) 3:1707–11. doi: 10.1001/jamaoncol.2017.2140

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Gentles AJ, Newman AM, Liu CL, Bratman SV, Feng W, Kim D, et al. The prognostic landscape of genes and infiltrating immune cells across human cancers. Nat Med. (2015) 21:938–45. doi: 10.1038/nm.3909

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Robert C, Grob JJ, Stroyakovskiy D, Karaszewska B, Hauschild A, Levchenko E, et al. Relevance of tumor-infiltrating immune cell composition and functionality for disease outcome in breast cancer. J Natl Cancer Inst. (2017) 109:djw192. doi: 10.1093/jnci/djw192

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Ali HR, Chlon L, Pharoah PDP, Markowetz F, Caldas C. Patterns of immune infiltration in breast cancer and their clinical implications: a gene-expression-based retrospective study. PLoS Med. (2016) 13:e1002194. doi: 10.1371/journal.pmed.1002194

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Jansen CS, Prokhnevska N, Master VA, Sanda MG, Carlisle JW, Bilen MA, et al. An intra-tumoral niche maintains and differentiates stem-like CD8 T cells. Nature. (2019) 576:465–70. doi: 10.1038/s41586-019-1836-5

PubMed Abstract | CrossRef Full Text | Google Scholar

41. de Andrea CE, Schalper KA, Sanmamed MF, Melero I. Immunodivergence in metastatic colorectal cancer. Cancer Cell. (2018) 34:876-8. doi: 10.1016/j.ccell.2018.11.012

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Wang K, Shen T, Siegal GP, Wei S. The CD4/CD8 ratio of tumor-infiltrating lymphocytes at the tumor-host interface has prognostic value in triple-negative breast cancer. Hum Pathol. (2017) 69:110–7. doi: 10.1016/j.humpath.2017.09.012

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Kareva I. Metabolism and gut microbiota in cancer immunoediting, CD8/treg ratios, immune cell homeostasis, and cancer (immuno)therapy: concise review. Stem Cells. (2019) 37:1273–80. doi: 10.1002/stem.3051

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Zhang C, Zheng J, Lin Z, Lv H, Ye Z, Chen Y, et al. Profiles of immune cell infiltration and immune-related genes in the tumor microenvironment of osteosarcoma. Aging. (2020) 12:3486–501. doi: 10.18632/aging.102824

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Xiong Y, Wang K, Zhou H, Peng L, You W, Fu Z. Profiles of immune infiltration in colorectal cancer and their clinical significant: a gene expression-based study. Cancer Med. (2018) 7:4496–4508. doi: 10.1002/cam4.1745

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Chowell D, Morris LGT, Grigg CM, Weber JK, Samstein RM, Makarov V, et al. Patient HLA class I genotype influences cancer response to checkpoint blockade immunotherapy. Science. (2018) 359:582–7. doi: 10.1126/science.aao4572

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Rodig SJ, Gusenleitner D, Jackson DG, Gjini E, Giobbie-Hurder A, Jin C, et al. MHC proteins confer differential sensitivity to CTLA-4 and PD-1 blockade in untreated metastatic melanoma. Sci Transl Med. (2018) 10:aar3342. doi: 10.1126/scitranslmed.aar3342

PubMed Abstract | CrossRef Full Text | Google Scholar