Ubiquitination-Related Molecular Subtypes and a Novel Prognostic Index for Bladder Cancer Patients
Hai Cai1†, Hang Chen1†, Qi Huang1†, Jun-Ming Zhu1†, 
Zhi-Bin Ke1†, Yun-Zhi Lin1, 
Qing-Shui Zheng1, 
Yong Wei1, 
Ning Xu1,2* and 
Xue-Yi Xue1*- 1Department of Urology, Urology Research Institute, The First Affiliated Hospital, Fujian Medical University, Fuzhou, China
- 2Fujian Key Laboratory of Precision Medicine for Cancer, The First Affiliated Hospital, Fujian Medical University, Fuzhou, China
Objective: To develop and validate ubiquitination-related molecular subtypes and a novel prognostic index using ubiquitination-related genes (URGs) for patients with bladder cancer (BCa).
Materials and Methods: We downloaded the clinical data and transcriptome data of BCa from TCGA and GEO database. Consensus clustering analysis was conducted to identify ubiquitination-related molecular subtypes for BCa. Besides, we performed univariate and multivariate Cox regression analysis to develop a novel prognostic URGs-related index for BCa. We conducted internal and external verification in TCGA cohort and GEO cohort, respectively. Furthermore, the associations of ubiquitination-related molecular subtypes and prognostic index with tumor immune environment were also investigated.
Results: A total of four ubiquitination-related molecular subtypes of BCa were finally identified. These four molecular subtypes had significantly different clinical characteristics, prognosis, PD-L1 expression level and tumor microenvironment. Besides, we developed a novel prognostic index using six URGs (including HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T and USP5). The difference in OS between high and low-risk group was statistically significant in training cohort, testing cohort, and validating cohort. The area under ROC curve (AUC) for OS prediction was 0.736, 0.723, and 0.683 in training cohort, testing cohort, and validating cohort, respectively. Multivariate survival analysis showed that this index was an independent predictor for OS. This prognostic index was especially suitable for subtype 1 and 3, older, male, high grade, AJCC stage III-IV, stage N0, stage T3-4 BCa patients.
Conclusions: This study identified a total of four ubiquitination-related molecular subtypes with significantly different tumor microenvironment, prognosis, clinical characteristics and PD-L1 expression level. Besides, a novel ubiquitination-related prognostic index for BCa patients was developed and successfully verified, which performed well in predicting prognosis of BCa.
Introduction
The new cases of bladder cancer (BCa) increase by more than 500,000 per year and the deaths caused by BCa increase by approximately 200,000 per year [1,2]. The treatment outcomes are diverse for different BCa patients, especially muscle-invasive BCa (MIBC) [3]. Previously, there were five major subtyping classification systems, including LUND, UROMOL, the University of North Carolina (UNC), The Cancer Genome Atlas (TCGA), and the MD Anderson Cancer Center (MDACC). These five subtyping classification systems not only have evolved independently, but also each taxonomy is different among each other in nomenclature [4–8]. Until now, there is no consistent risk stratification for BCa [9].
Ubiquitination is one of proteins post-translational modification types, and proteins ubiquitination might alter proteins localization, lead to degradation via proteasome, affect proteins activity and interactions [10,11]. Recent studies demonstrated that ubiquitination played vital roles in cancer-related pathways, and that ubiquitinated protein accumulation might be a novel method for treating cancer [12,13]. Akinori Sato et al. [12], reported that ixazomib and ritonavir could suppress the tumorigenesis of BCa by inducing the accumulation of ubiquitinated protein. The development of targetable molecular subtypes and risk stratification tools using ubiquitination-related genes (URGs) is promising. However, as far as we know, there is no previous study exploring the correlations of ubiquitination with prognosis evaluation and molecular subtypes of BCa.
In our study, consensus clustering analysis was used to identify ubiquitination related molecular subtypes. Besides, we developed a novel prognostic index based on URGs by multivariate Cox regression analysis. Next, we conducted internal and external verification, respectively. Finally, we also explored the correlations of the URGs-based prognostic index with molecular subtypes, clinical features, immune function and immune infiltrating cells.
Materials and Methods
Data Acquisition
The transcriptome and clinical data of BCa were obtained from TCGA database (https://portal.gdc.cancer.gov), including 414 BCa cases and 19 normal cases [14]. The clinical features of BCa patients included age, gender, grade, T stage, N stage, M stage, and AJCC stage. A total of 412 BCa cases have clinical data; however, there was only a total of 396 BCa cases with unabridged transcriptome and clinical data simultaneously. Besides, a total of 165 BCa patients with complete mRNA expression profile and clinical data were downloaded from GSE13507 dataset in GEO database (https://www.ncbi.nlm.nih.gov). All these data were processed using Perl language and R language. The clinicopathologic feature of TCGA cohort and GEO cohort was showed in Table 1. We used the Ensemble database (http://asia.ensembl.org/signature.html) to convert Ensemble IDs into gene symbols.
TABLE 1
TABLE 1. clinicopathologic data of TCGA cohort and GEO cohort.
The Molecular Signatures Database (MSigDB) (https://www.gsea-msigdb.org/gsea/msigdb) is a collection of annotated gene sets using GSEA software. We extracted 79 URGs from REACTOME PROTEIN UBIQUITINATION gene set (C2: curated gene sets; systematic name: M27742) from MSigDB database. All these URGs were related to protein ubiquitination. The detailed list of these 79 URGs were presented in Supplementary Table S1.
Identification of Ubiquitination-Related Molecular Subtypes Using Consensus Clustering Analysis
Firstly, the expression matrix of 79 URGs was extracted from BCa transcriptome and merged with overall survival (OS) time using Perl language. There was a total of 396 BCa cases with complete OS data and mRNA expression data in TCGA database. We then performed univariable Cox regression analysis to screen prognostic URGs associated with OS. The cut-off p value was set as 0.05. Next, consensus clustering analysis was performed for identifying ubiquitination-related molecular subtypes of BCa patients using R package “ConsensusClusterPlus”. The association between ubiquitination-related molecular subtypes and OS was explored using R package “survival” and “survminer”. The associations between ubiquitination-related subtypes and clinicopathologic characteristics (including gender, grade, T stage, N stage, stage, and age) were presented by utilizing R package “pheatmap”.
Development and Validation of a Novel Ubiquitination-Based Prognostic Index
First of all, we divided randomly all BCa patients in TCGA database into training cohort and testing cohort (internal verification). All cases in GSE13507 dataset were used as validating cohort (external verification). Then, univariate and multivariate Cox regression analysis were conducted to establish a ubiquitination-based prognostic index for predicting OS of BCa using training cohort.
Next, all cases in training cohort, testing cohort and validating cohort were categorized into low-risk group and high-risk group based on the median risk score. We then performed survival analysis and the time-dependent receiver operating characteristic (ROC) curve to explore the performance of this ubiquitination-based prognostic index in training cohort, testing cohort and validating cohort. The expression heatmap, the distribution of risk score and survival time of training cohort, testing cohort and validating cohort were presented using “pheatmap” R package. Moreover, we performed univariate and multivariate independent prognostic analysis to demonstrate whether this ubiquitination-based prognostic index was an independent predictor of OS in BCa patients.
Exploration of Immune Cells Infiltration, Tumor Microenvironment, PD-L1 Expression Level
The ESTIMATE algorithm was utilized to evaluate tumor microenvironment (TME) scores while the CIBERSORT method was used to calculate the scores of 22 types of immune infiltrating cells [15,16]. Then, we investigated the associations of ubiquitination-related molecular subtypes with BCa tumor microenvironment and immune infiltrating cells. In addition, the association of the ubiquitination-related molecular subtypes with PD-L1 gene expression level [17] was also explored.
Single sample gene set enrichment analysis (ssGSEA) was conducted to calculate the infiltrating score of 16 types of immune infiltrating cell and the activity of 13 types of immune-related function in TCGA cohort using R package “gsva”. We then investigated the associations of the ubiquitination-based prognostic index with immune infiltrating cells and immune function activity. Besides, the associations of ubiquitination-based prognostic index with TME scores were explored.
Validation of Six Risk URGs and Functional Enrichment
UALCAN database (http://ualcan.path.uab.edu/) is a portal for facilitating tumor subgroup gene expression and survival analyses. The mRNA expression levels of risk URGs between normal and tumor tissues were demonstrated using UALCAN database. The prognostic value of risk URGs was validated in GEO cohort using univariate Cox regression analysis. Kyoto Encyclopedia of Genes and Genomes (KEGG) functional enrichment for high-risk group and low-risk group was performed using Gene Set Enrichment Analysis (GSEA) method in the whole TCGA cohort.
Statistical Methods
Statistical analysis was performed utilizing R programming language. Univariate and multivariate Cox regression analysis were performed to establish a ubiquitination-related prognostic index for predicting OS of BCa. Univariate and multivariate independent prognostic analysis were used to demonstrate whether this ubiquitination-related prognostic index was an independent predictor of OS. Survival analysis and the time-dependent receiver operating characteristic (ROC) curve were performed to explore the performance of ubiquitination-related prognostic index. Statistical significance was considered at the level of p value <0.05.
Results
Identification of Four Ubiquitination-Related Molecular Subtypes
The flowchart of this study was presented in Figure 1. Firstly, the expression matrix of 79 URGs was extracted from TCGA database and merged with OS status and time. We then conducted univariable Cox regression analysis in the whole TCGA cohort to screen URGs associated with OS. The results showed that there was a total of six URGs associated with OS, including CDC73, PRKDC, RNF40, TMEM129, UBE2N and USP5. Next, the expression matrix of these six URGs was used to conduct consensus clustering analysis to identify ubiquitination-related molecular subtypes of BCa. Finally, a total of four ubiquitination-related molecular subtypes of BCa were identified, including 116 cases of subtype 1, 90 cases of subtype 2, 132 cases of subtype 3, 58 cases of subtype 4 (Figures 2A,B). The difference of OS among these ubiquitination-related molecular subtypes was statistically significant (p = 0.005, Figure 2C). As indicated by heatmap and histogram, these four molecular subtypes have significantly different grade (p < 0.05, Figures 2D,E).
FIGURE 1
FIGURE 1. The flowchart of this study.
FIGURE 2
FIGURE 2. Identification of four ubiquitination-related molecular subtypes of BCa using consensus clustering analysis (A, B). Comparison of overall survival among these four subtypes (C). The correlation heatmap between this ubiquitination-related molecular subtypes and clinicopathologic features (D). The difference of grade (E) and the expression level of PD-L1 among these four subtypes (F). *p < 0.05; **p < 0.01; ***p < 0.001.
Associations of these Subtypes with PD-L1 Expression and TME
We explored the discrepant expression level of PD-L1 in TCGA database among these four molecular subtypes. The results showed that the expression level of PD-L1 in subtype 3 and subtype 4 was significantly increased compared with that in subtype 1 and subtype 2. However, the PD-L1 expression level was similar between subtype 1 and subtype 2, between subtype 3 and subtype 4, respectively. (Figure 2F).
We used the CIBERSORT method to calculate the scores of 22 types of immune infiltrating cells and the ESTIMATE algorithm to evaluate TME. According to the CIBERSORT method, the activated mast cell was the exclusively discrepant immune infiltrating cell among these four ubiquitination-related subtypes (Figure 3). The ESTIMATE scores, immune scores, stromal scores in subtype 2 and subtype 3 were significantly higher in comparison with that in subtype 1 and subtype 4; the tumor purity scores in subtype 1 and subtype 4 were significantly increased compared with that in subtype 2 and subtype 3. However, the ESTIMATE scores, immune scores, stromal scores and tumor purity scores were similar between subtype 1 and subtype 4, between subtype 2 and subtype 3, respectively (Figure 4).
FIGURE 3
FIGURE 3. Distribution of 22 types of immune cells among these four molecular subtypes (A). Relationship between these ubiquitination-related molecular subtypes and immune cells infiltration (B). *p < 0.05; **p < 0.01; ***p < 0.001.
FIGURE 4
FIGURE 4. Relationship between the ubiquitination-related molecular subtypes and tumor microenvironment. ESTIMATE score (A). Immune score (B). Stromal score (C). Tumor purity (D).
Development and Validation of a Novel Ubiquitination-Based Prognostic Index
There was a total of 200 cases in training cohort, 196 cases in testing cohort and 165 cases in validating cohort. Firstly, univariate and multivariate Cox regression analyses were conducted to develop a novel ubiquitination-based index utilizing six URGs (including HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T and USP5) for prognosis prediction of BCa using training cohort. The calculation formula of risk score is shown as follows: Risk score= (−0.0010853) * HLA-A + (−0.03553995) * TMEM129 + (−0.08324461) * UBE2D1 + (−0.07059857) * UBE2N + (0.01191996) * UBE2T + (0.02752775) * USP5 (Table 2).
TABLE 2
TABLE 2. Univariate and multivariate Cox regression analysis in training cohort to developing ubiquitination-based prognostic index for bladder cancer.
Then, we calculated the risk score of each patient using above calculation formula in training cohort, testing cohort, and validating cohort. All patients were divided into high-risk score group and low-risk score group according to the median risk score. The difference in OS between low-risk and high-risk group was statistically significant in training cohort (p < 0.001), testing cohort (p = 0.003), whole TCGA cohort (p < 0.001) and validating cohort (p = 0.035), respectively. Low-risk score was associated with significantly preferable OS in comparison with high-risk score in training cohort, testing cohort, whole TCGA cohort and validating cohort, respectively. The area under ROC curve (AUC) for OS prediction was 0.736, 0.723, 0.723 and 0.683 in training cohort, testing cohort, whole TCGA cohort and validating cohort, respectively, suggesting the promising value of this novel ubiquitination-based prognostic index for prognosis prediction of BCa (Figure 5). The expression heatmaps, the distributions of risk score and survival time of training cohort, testing cohort, whole TCGA cohort and validating cohort were presented in Figure 6.
FIGURE 5
FIGURE 5. Internal and external verification of a novel ubiquitination-based prognostic index for BCa. The survival analysis and corresponding area under ROC curve in training cohort (A, B), testing cohort (C, D), the whole TCGA cohort (E, F), and validating cohort (G, H).
FIGURE 6
FIGURE 6. The expression heatmap, the distribution of risk score and survival time of training cohort (A), testing cohort (B), the whole TCGA cohort (C), and validating cohort (D).
Univariate and multivariate independent prognostic analysis were used to demonstrate whether this index was an independent predictor of OS in BCa patients. Univariate analysis showed that age (p = 9.97E-05), stage (p = 9.56E-08), T stage (p = 6.08E-05), N stage (p = 1.57E-07) and this prognostic index (p = 1.96E-10) were associated with OS of BCa. Multivariate analysis demonstrated that only age (p = 0.000189) and this prognostic index (p = 3.82E-06) were independent predictors for OS of BCa, indicating the great performance of this prognostic index (Table 3).
TABLE 3
TABLE 3. Univariate and multivariate independent prognostic analysis in whole TCGA cohort.
Associations of this Prognostic Index with Clinicopathologic Features
The difference of clinicopathologic features between high and low risk group were explored. The results demonstrated that ubiquitination-related molecular subtypes (p < 0.001), age (p < 0.001), grade (p = 0.002), T stage (p = 0.003), N stage (p = 0.002), AJCC stage (p < 0.001) and were significantly different between high-risk and low-risk patients. However, gender (p = 0.365) was similar between high-risk and low-risk patients (Figures 7A–G).
FIGURE 7
FIGURE 7. The difference of clinicopathologic features between high and low risk group, including ubiquitination-related molecular subtypes (A), age (B), gender (C), grade (D), T stage (E), N stage (F), and AJCC stage (G). The survival analysis between high and low risk group in patients of high grade (H), age ≤65 (I), age >65 (J), female (K), male (L), AJCC stage I-II (M), AJCC stage III-IV (N), T1-2 stage (O), T3-4 stage (P), N0 stage (Q) and N1-3 stage (R).
Subgroup Survival Analysis
Subgroup survival analysis demonstrated that this prognostic index was especially suitable for older, male, high grade, AJCC stage III-IV, stage N0, stage T3-4 BCa patients. Low-risk score was associated with significantly preferable OS in comparison with high-risk score in patients of age >65 or male or high grade or AJCC stage III-IV or N0 stage or T3-4 stage. However, the difference in OS between low-risk and high-risk group was not statistically significant in patients of age ≤65 or female or AJCC stage I-II or N1-3 stage or T1-2 stage (Figures 7H–R). Besides, considering the extremely low case number of patients with low grade, we did not conduct subgroup survival analysis in this subgroup.
We calculated the risk score for each patient of each ubiquitination-related molecular subtype. The results revealed that this prognostic index was especially applicable in BCa patients of subtype 1 and subtype 3. Low-risk score was associated with significantly preferable OS in comparison with high-risk score in BCa patients of subtype 1 and subtype 3. The area under ROC curve (AUC) for OS prediction was 0.778, 0.713 in BCa patients of subtype 1 and subtype 3, respectively. However, the difference in OS between low-risk and high-risk group was not statistically significant in subtype 2 and subtype 4 patients (Figure 8).
FIGURE 8
FIGURE 8. The survival analysis between high and low risk group, and corresponding area under ROC curve in BCa patients of subtype 1 (A, B), subtype 2 (C, D), subtype 3 (E, F), and subtype 4 (G, H).
Associations of this Prognostic Index with Immune Cell Infiltration and Tumor Microenvironment
The immune functions and immune cells infiltration were quantified using ssGSEA in TCGA cohort. The infiltrating proportions of CD8+ T cells, interdigitating dendritic (iDC) cells, natural killer (NK) cells, T helper 2 (Th2) cells in low-risk group were significantly increased compared with that in high-risk group; patients with low-risk score have higher scores of immune function or pathways than those with high-risk, including cytolytic activity, human leukocyte antigen (HLA), major histocompatibility complex (MHC) class I, para-inflammation, type I interferon (IFN) response (Figures 9A,B). Besides, the associations of ubiquitination-based prognostic index with TME scores were explored. The correlation analysis demonstrated that the risk score was negatively associated with immune score derived from ESTIMATE algorithm. However, the risk score was not associated with ESTIMATE score, stromal score, and tumor purity score (Figures 9C–F).
FIGURE 9
FIGURE 9. Associations of this prognostic index with immune functions (A) and immune infiltrating cells (B) in TCGA cohort. The association between risk score and ESTIMATE score (C), immune score (D), stromal score (E), and tumor purity (F). APC, antigen-presenting cells; CCR, chemokine receptors; HLA, human leukocyte antigen; MHC, major histocompatibility complex; IFN, interferon; DCs, dendritic cells; iDCs, interdigitating dendritic cells; pDCs, plasmacytoid dendritic cells; aDCs, activated dendritic cells; NK, natural killer; Th1, T helper 1; Th2, T helper 2. *p < 0.05; **p < 0.01; ***p < 0.001.
Validation of Six Risk URGs and GSEA Functional Enrichment
As indicated by UALCAN database, the mRNA expression levels of HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T in BCa tissue were significantly increased compared with that in normal tissue; however, the mRNA expression level of USP5 was similar between normal and tumor tissue. We performed univariable Cox regression analysis to validate the prognostic value of risk URGs in validating cohort (GEO cohort). The results showed that TMEM129, UBE2D1, and UBE2T were successfully validated to be associated with OS in GEO cohort. However, the prognostic value of HLA-A, UBE2T, and USP5 was not found to be associated with OS in GEO cohort. KEGG functional enrichment for high-risk group and low-risk group was performed using GSEA method. Toll like receptor signaling pathway, cytosolic DNA sensing pathway, other glycan degradation, and lysosome were major four KEGG enriched pathways in low-risk patients. GAP junction, DNA replication, oocyte meiosis, and RNA polymerase were major four KEGG enriched pathways in high-risk patients (Figure 10).
FIGURE 10
FIGURE 10. Validation of mRNA expression levels of HLA-A (A), TMEM129 (B), UBE2D1 (C), USP5 (F), UBE2T (E), and UBE2N (D) in UALCAN database. Univariable Cox regression analysis in validating cohort (GEO cohort) to validate the prognostic value of risk URGs (G). Kyoto Encyclopedia of Genes and Genomes (KEGG) functional enrichment of high-risk and low-risk group using Gene Set Enrichment Analysis (GSEA) method (H).
Discussions
BCa is the ninth most common cancer worldwide and urothelial carcinoma accounts for 90% of all BCa cases [18–20]. Recently, it has been reported that the inducing accumulation of ubiquitinated proteins was considered as a vital treatment option for BCa [12]. Besides, immune checkpoint inhibitors have been confirmed significantly effective in treating BCa [21]. The identification of ubiquitination-related molecular subtypes and a novel prognostic index for BCa could classify the heterogeneous tumor microenvironment of BCa and contribute to the development of individualized therapy [22,23]. Previously, several gene expression based molecular classification systems with prognostic and predictive relevance have been described for MIBC. However, none of these classifications focused on ubiquitination [8,9].
In this study, we firstly identified a total of four molecular subtypes from the perspective of ubiquitination, which were different from previous studies. These four molecular subtypes were found to exhibit significantly diverse clinical characteristics and tumor microenvironment features. The difference of OS among these ubiquitination-related molecular subtypes was statistically significant. These four molecular subtypes had significantly different grade. The expression level of PD-L1 in subtype 3 and subtype 4 were significantly increased compared with that in subtype 1 and subtype 2. The activated mast cell was the exclusively discrepant immune infiltrating cell among these four ubiquitination-related subtypes. The ESTIMATE scores, immune scores, stromal scores in subtype 2 and subtype 3 were significantly higher in comparison with that in subtype 1 and subtype 4. Therefore, the four ubiquitination molecular subtypes would contribute to stratifying patients for PD-L1 expression level and prognosis. The underlying molecular mechanisms among the four ubiquitination molecular subtypes deserve further investigation.
Furthermore, we used multivariate Cox regression analysis to develop a novel ubiquitination-based index utilizing six URGs (including HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T and USP5) for prognosis prediction of BCa. This novel ubiquitination-based prognostic index was successfully internally tested and externally validated. There were several previous studies developing prognostic index for BCa. Zhang et al. [24] included eight genes (AKAP12, ALDOB, CASP6, DTNA, HS3ST1, JUN, KDELR3, and STC1) to develop a novel hypoxia-related signature for BCa, which could play a precise role in predicting outcome. Liu et al. [25] used four prognostic genes (including ALOX5, FANCD2, HMGCR and FADS2) and multivariate Cox regression analysis to build a ferroptosis-related signature for BCa prognosis prediction. Qu et al. [26] developed a prognostic index for BCa from the perspective of immune-related gene. Luo et al. [27] constructed a poliovirus receptor (CD155)-related risk signature for predicting prognosis and immune infiltration of BCa. Du et al. [28] identified an epithelial to mesenchymal transition related lncRNA signature in order to predict BCa outcome. However, most of them did not conduct both internal and external verification, and this is the first study focusing on ubiquitination-related genes.
In this study, the difference in OS between low-risk and high-risk group was statistically significant in training cohort, testing cohort, and validating cohort, respectively. The AUC for OS prediction was 0.736, 0.723 and 0.683 in training cohort, testing cohort, and validating cohort, respectively, suggesting the promising value of this novel ubiquitination-based prognostic index for prognosis prediction of BCa. Subgroup analysis revealed that this prognostic index was especially suitable for older, male, high grade, AJCC stage III-IV, stage N0, stage T3-4 BCa patients. Moreover, we also found that this prognostic index was especially applicable in subtype 1 and subtype 3 BCa. Low-risk score was associated with significantly preferable OS in comparison with high-risk score in subtype 1 and subtype 3. The prognostic index could only be used to predict OS for subtype 1 and subtype 3 BCa.
The ubiquitination-related prognostic molecular index is composed of six risk URGs (including HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T and USP5). Besides, In the identification of ubiquitination-related molecular subtypes, we revealed that there was a total of six URGs associated with OS, including CDC73, PRKDC, RNF40, TMEM129, UBE2N and USP5. Major histocompatibility complex, class I, A (HLA-A) belongs to the human leukocyte antigen (HLA) class I heavy chain paralogues [29]. Richard Golnik et al. [30] revealed that atypical ubiquitination is of great importance in modulating the processing of MHC class I antigen. Zhou et al. [31] demonstrated in their study that UBE2D1 was crucial in the development and progression of hepatocellular carcinoma. They also reported that high expression of UBE2D1 was attributed to the recurrent genomic copy number gain. However, no research paid attention to the oncogenic role of UBE2D1 in BCa. Jonine D Figueroa et al. [32] examined the data on a prevalent occupational exposure related to increased risk of BCa, and revealed a statistically significant additive interaction for rs798766 (TMEM129-TACC3-FGFR3). Yan et al. [33] demonstrated that miR-934 could bind to UBE2N mRNA and decrease the protein expression of UBE2N, and then lead to the cancer growth of BCa. Gong et al. [34] found that UBE2T might play momentous roles in BCa cells proliferation and cell cycle; hence, UBE2T might be a potential therapeutic targeting point for BCa. Xiao et al. [35] revealed that MAFG-AS1 could directly bind to HuR and then recruit USP5 to prevent HuR from degrading by de-ubiquitination. This process is vital in bladder urothelial carcinoma progression. Michael A Hahn et al. [36] reported that CDC73 is required for maintaining the mono-ubiquitination of H2B-K120 and the reduction of H2B-K120 mono-ubiquitination levels might be a potential mechanism of tumorigenic effect of CDC73. Shang et al. [37] demonstrated that DNA-PKcs (also called PRKDC) could negatively regulate the stability of cyclin B1 via ubiquitination. Fu et al. [38] summarized that E3 ligase RNF40 is vital in cancer development and metastasis. However, indeed, no research has been done yet in the aspect of these URGs in the progressions of BCa.
Toll-like receptors (TLRs) have been reported to drive immune responses in progression as well as treatment of cancer [39]. Recently, various TLRs signaling pathway proteins were proven to be related to the susceptibility of BCa [40]. TLRs agonists were regarded as effective immunostimulants with immunotherapeutic potential against BCa [41]. Besides, Bacillus Calmette-Guerin (BCG) was considered as an agonist of TLR2 and TLR4 and has been widely used in intravesical BCa therapy [41]. Connexins is of great importance for gap junctional intercellular communication and might facilitate the occurrence and metastases of cancer [42]. Ai et al. [43] revealed that GAP junction protein CX43 could be a promising therapeutic target for BCa. H B Grossman et al. [44] found that the alterations of connexin 26 expression may be associated with the malignant phenotype of BCa. In this study, we found that TLRs signaling pathway was a major enriched pathway in low-risk patients while the GAP junction pathway was a main enriched pathway in high-risk patients. However, further investigation into the specific mechanism is required in future.
There were several inescapable limitations of our study. Firstly, this study is retrospective and based on public database. The clinicopathologic information was not comprehensive and the detailed treatment regimens also could not be obtained. Validation using prospective sequencing data is required for generalizability. Secondly, these six risk URGs (including HLA-A, TMEM129, UBE2D1, UBE2N, UBE2T and USP5) might play vital roles in the development and progression of BCa. However, this study was performed via pure bioinformatics methods and the underlying mechanisms were not covered. Further in vitro and in vivo experiments into the underlying mechanism are required. Thirdly, it should be noted that for all parameters, the patient subgroup with the low case number showed statistically unsignificant differences of OS. Hence, one of the reasons for lack of significantly prognostic effect of the prognostic index in some patient subgroups could be the low case number. Validation of this ubiquitination-based prognostic index with larger sample size is required.
Conclusion
In this study, we identified a total of four ubiquitination-related molecular subtypes with significantly different tumor microenvironment, prognosis, clinical characteristics and PD-L1 expression level. Besides, a novel ubiquitination-related prognostic index for BCa patients was developed and successfully verified. The ubiquitination-based molecular subtypes and prognostic index would contribute to further understanding of molecular mechanisms in BCa.
Data Availability Statement
Publicly available datasets were analyzed in this study. This data can be found here: TCGA database (https://portal.gdc.cancer.gov), GEO database (https://www.ncbi.nlm.nih.gov), Molecular Signatures Database (https://www.gsea-msigdb.org/gsea/msigdb) (MSigDB).
Author Contributions
Writing–original draft, HC, HC, QH, and J-MZ. Writing–review and editing, NX and X-YX. Methodology, Y-ZL and Z-BK. Formal analysis, Q-SZ. Data curation, YW. Conceptualization, NX. Visualization, Q-SZ. Project administration, NX and X-YX.
Funding
This study was supported by Natural Science Fundation of Fujian Province (No.2021J01715).
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.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.por-journal.com/articles/10.3389/pore.2021.1609941/full#supplementary-material
References
1. Lenis, A. T., Lec, P. M., Chamie, K., and Mshs, M., Bladder Cancer: A Review. Jama (2020) 324(19):1980–91. doi:10.1001/jama.2020.17598
2. Richters, A., Aben, K. K. H., and Kiemeney, L. A. L. M., The Global burden of Urinary Bladder Cancer: an Update. World J Urol (2020) 38(8):1895–904. doi:10.1007/s00345-019-02984-4
3. Ru, Y., Dancik, G. M., and Theodorescu, D., Biomarkers for Prognosis and Treatment Selection in Advanced Bladder Cancer Patients. Curr Opin Urol (2011) 21(5):420–7. doi:10.1097/MOU.0b013e32834956d6
4. Sjödahl, G., Lauss, M., Lövgren, K., Chebil, G., Gudjonsson, S., Veerla, S., et al. ,A Molecular Taxonomy for Urothelial Carcinoma. Clin Cancer Res (2012) 18(12):3377–86. doi:10.1158/1078-0432.Ccr-12-0077-t
5. Choi, W., Porten, S., Kim, S., Willis, D., Plimack, E. R., Hoffman-Censits, J., et al. ,Identification of Distinct Basal and Luminal Subtypes of Muscle-Invasive Bladder Cancer with Different Sensitivities to Frontline Chemotherapy. Cancer Cell (2014) 25(2):152–65. doi:10.1016/j.ccr.2014.01.009
6. Damrauer, J. S., Hoadley, K. A., Chism, D. D., Fan, C., Tiganelli, C. J., Wobker, S. E., et al. ,Intrinsic Subtypes of High-Grade Bladder Cancer Reflect the Hallmarks of Breast Cancer Biology. Proc Natl Acad Sci (2014) 111(8):3110–5. doi:10.1073/pnas.1318376111
7. Hedegaard, J., Lamy, P., Nordentoft, I., Algaba, F., Høyer, S., Ulhøi, B. P., et al. ,Comprehensive Transcriptional Analysis of Early-Stage Urothelial Carcinoma. Cancer Cell (2016) 30(1):27–42. doi:10.1016/j.ccell.2016.05.004
8. Robertson, A. G, Kim, J., Al-Ahmadie, H., Bellmunt, J., Guo, G., Cherniack, A. D., et al. ,Comprehensive Molecular Characterization of Muscle-Invasive Bladder Cancer. Cell (2017) 171(3):540–556. e525. doi:10.1016/j.cell.2017.09.007
9. Kamoun, A., de Reyniès, A., Allory, Y., Sjödahl, G., Robertson, A. G., Seiler, R., et al. ,A Consensus Molecular Classification of Muscle-Invasive Bladder Cancer. Eur Urol (2020) 77(4):420–33. doi:10.1016/j.eururo.2019.09.006
10. Mansour, M. A., Ubiquitination: Friend and Foe in Cancer. Int J Biochem Cel Biol (2018) 101:80–93. doi:10.1016/j.biocel.2018.06.001
11. Setiawan, B., Rosalina, R., Pangarsa, E. A., Santosa, D., and Suharti, C., Clinical Evaluation for the Role of High-Sensitivity C-Reactive Protein in Combination with D-Dimer and Wells Score Probability Test to Predict the Incidence of Deep Vein Thrombosis Among Cancer Patients. Ijgm (2020) 13:587–94. doi:10.2147/ijgm.S261718
12. Sato, A., Asano, T., Okubo, K., Isono, M., and Asano, T., Ritonavir and Ixazomib Kill Bladder Cancer Cells by Causing Ubiquitinated Protein Accumulation. Cancer Sci (2017) 108(6):1194–202. doi:10.1111/cas.13242
13. Wang, X., Liu, Z., Zhang, L., Yang, Z., Chen, X., Luo, J., et al. ,Targeting Deubiquitinase USP28 for Cancer Therapy. Cell Death Dis (2018) 9(2):186. doi:10.1038/s41419-017-0208-z
14. Xu, N., Ke, Z-B., Lin, X-D., Chen, Y-H., Wu, Y-P., Chen, Y., et al. ,Development and Validation of a Molecular Prognostic index of Bladder Cancer Based on Immunogenomic Landscape Analysis. Cancer Cel Int (2020) 20:302. doi:10.1186/s12935-020-01343-3
15. Ke, Z. B., Wu, Y. P., Huang, P., Hou, J., Chen, Y. H., Dong, R. N., et al. ,Identification of Novel Genes in Testicular Cancer Microenvironment Based on ESTIMATE Algorithm‐derived Immune Scores. J Cel Physiol (2020) 236:706–13. doi:10.1002/jcp.29898
16. Mou, Z., Yang, C., Zhang, Z., Wu, S., Xu, C., Cheng, Z., et al. ,Transcriptomic Analysis of Glycolysis-Related Genes Reveals an Independent Signature of Bladder Carcinoma. Front Genet (2020) 11:566918. doi:10.3389/fgene.2020.566918
17. El Sabbagh, R., Azar, N. S., Eid, A. A., and Azar, S. T., Thyroid Dysfunctions Due to Immune Checkpoint Inhibitors: A Review. Ijgm (2020) 13:1003–9. doi:10.2147/ijgm.S261433
18. Siegel, R., Ma, J., Zou, Z., and Jemal, A., Cancer Statistics, 2014. CA A Cancer J Clinicians (2014) 64(1):9–29. doi:10.3322/caac.21208
20. Cao, J., Li, J., Yang, X., Li, P., Yao, Z., Han, D., et al. ,Integrative Analysis of Immune Molecular Subtypes and Microenvironment Characteristics of Bladder Cancer. Cancer Med (2021) 10(15):5375–91. doi:10.1002/cam4.4071
21. Kim, H. S., and Seo, H. K., Immune Checkpoint Inhibitors for Urothelial Carcinoma. Investig Clin Urol (2018) 59(5):285–96. doi:10.4111/icu.2018.59.5.285
22. Kirkali, Z., Chan, T., Manoharan, M., Algaba, F., Busch, C., Cheng, L., et al. ,Bladder Cancer: Epidemiology, Staging and Grading, and Diagnosis. Urology (2005) 66(6 Suppl. 1):4–34. doi:10.1016/j.urology.2005.07.062
23. Bellmunt, J., Bladder Cancer. Hematology/Oncology Clin North America (2015) 29(2):xiii–xiv. doi:10.1016/j.hoc.2014.12.001
24. Zhang, F., Wang, X., Bai, Y., Hu, H., Yang, Y., Wang, J., et al. ,Development and Validation of a Hypoxia-Related Signature for Predicting Survival Outcomes in Patients with Bladder Cancer. Front Genet (2021) 12:670384. doi:10.3389/fgene.2021.670384
25. Liu, J., Ma, H., Meng, L., Liu, X., Lv, Z., Zhang, Y., et al. ,Construction and External Validation of a Ferroptosis-Related Gene Signature of Predictive Value for the Overall Survival in Bladder Cancer. Front Mol Biosci (2021) 8:675651. doi:10.3389/fmolb.2021.675651
26. Qu, G., Liu, Z., Yang, G., Xu, Y., Xiang, M., and Tang, C., Development of a Prognostic index and Screening of Prognosis Related Genes Based on an Immunogenomic Landscape Analysis of Bladder Cancer. Aging (2021) 13(8):12099–112. doi:10.18632/aging.202917
27. Luo, C., Ye, W., Hu, J., Othmane, B., Li, H., Chen, J., et al. ,A Poliovirus Receptor (CD155)-Related Risk Signature Predicts the Prognosis of Bladder Cancer. Front Oncol (2021) 11:660273. doi:10.3389/fonc.2021.660273
28. Du, Y., Wang, B., Jiang, X., Cao, J., Yu, J., Wang, Y., et al. ,Identification and Validation of a Stromal EMT Related LncRNA Signature as a Potential Marker to Predict Bladder Cancer Outcome. Front Oncol (2021) 11:620674. doi:10.3389/fonc.2021.620674
29. Dard, C., Senoussi, O., Buhler, S., Bardy, B., and Masson, D., Identification of Four Novel HLA‐A Alleles. Hla (2020) 96(2):202–3. doi:10.1111/tan.13899
30. Golnik, R., Lehmann, A., Kloetzel, P-M., and Ebstein, F., Major Histocompatibility Complex (MHC) Class I Processing of the NY-ESO-1 Antigen Is Regulated by Rpn10 and Rpn13 Proteins and Immunoproteasomes Following Non-lysine Ubiquitination. J Biol Chem (2016) 291(16):8805–15. doi:10.1074/jbc.M115.705178
31. Zhou, C., Bi, F., Yuan, J., Yang, F., and Sun, S., Gain of UBE2D1 Facilitates Hepatocellular Carcinoma Progression and Is Associated with DNA Damage Caused by Continuous IL-6. J Exp Clin Cancer Res (2018) 37(1):290. doi:10.1186/s13046-018-0951-8
32. Figueroa, J. D., Koutros, S., Colt, J. S., Kogevinas, M., Garcia-Closas, M., Real, F. X., et al. ,Modification of Occupational Exposures on Bladder Cancer Risk by Common Genetic Polymorphisms. J Natl Cancer Inst (2015) 107(11):djv223. doi:10.1093/jnci/djv223
33. Yan, H., Ren, S., Lin, Q., Yu, Y., Chen, C., Hua, X., et al. ,Inhibition of UBE2N‐dependent CDK6 Protein Degradation by miR‐934 Promotes Human Bladder Cancer Cell Growth. FASEB j. (2019) 33(11):12112–23. doi:10.1096/fj.201900499RR
34. Gong, Y. Q., Peng, D., Ning, X. H., Yang, X. Y., Li, X. S., Zhou, L. Q., et al. ,UBE2T Silencing Suppresses Proliferation and Induces Cell Cycle Arrest and Apoptosis in Bladder Cancer Cells. Oncol Lett (2016) 12(6):4485–92. doi:10.3892/ol.2016.5237
35. Xiao, M., Liu, J., Xiang, L., Zhao, K., He, D., Zeng, Q., et al. ,MAFG‐AS1 Promotes Tumor Progression via Regulation of the HuR/PTBP1 axis in Bladder Urothelial Carcinoma. Clin Translational Med (2020) 10(8):e241. doi:10.1002/ctm2.241
36. Hahn, M. A., Dickson, K-A., Jackson, S., Clarkson, A., Gill, A. J., and Marsh, D. J., The Tumor Suppressor CDC73 Interacts with the Ring finger Proteins RNF20 and RNF40 and Is Required for the Maintenance of Histone 2B Monoubiquitination. Hum Mol Genet (2012) 21(3):559–68. doi:10.1093/hmg/ddr490
37. Shang, Z-F., Tan, W., Liu, X-D., Yu, L., Li, B., Li, M., et al. ,DNA-PKcs Negatively Regulates Cyclin B1 Protein Stability through Facilitating its Ubiquitination Mediated by Cdh1-APC/C Pathway. Int J Biol Sci (2015) 11(9):1026–35. doi:10.7150/ijbs.12443
38. Fu, J., Liao, L., Balaji, K. S., Wei, C., Kim, J., and Peng, J., Epigenetic Modification and a Role for the E3 Ligase RNF40 in Cancer Development and Metastasis. Oncogene (2021) 40(3):465–74. doi:10.1038/s41388-020-01556-w
39. So, E. Y., and Ouchi, T., The Application of Toll like Receptors for Cancer Therapy. Int J Biol Sci (2010) 6(7):675–81. doi:10.7150/ijbs.6.675
40. Cheng, S., Liu, J., Zhang, Y., Lin, Y., Liu, Q., Li, H., et al. ,Association Detection between Genetic Variants in the microRNA Binding Sites of Toll-like Receptors Signaling Pathway Genes and Bladder Cancer Susceptibility. Int J Clin Exp Pathol (2014) 7(11):8118–26.
41. Ohadian Moghadam, S., and Nowroozi, M. R., Toll‐like Receptors: The Role in Bladder Cancer Development, Progression and Immunotherapy. Scand J Immunol (2019) 90(6):e12818. doi:10.1111/sji.12818
42. Tschernig, T., Connexins and Gap Junctions in Cancer of the Urinary Tract. Cancers (2019) 11(5):704. doi:10.3390/cancers11050704
43. Ai, X-l., Chi, Q., Qiu, Y., Li, H-y., Li, D-j., Wang, J-x., et al. ,Gap junction Protein Connexin43 Deregulation Contributes to Bladder Carcinogenesis via Targeting MAPK Pathway. Mol Cel Biochem (2017) 428(1-2):109–18. doi:10.1007/s11010-016-2921-9
Keywords: prognosis, bladder cancer, ubiquitination, molecular subtypes, prognostic index
Citation: Cai H, Chen H, Huang Q, Zhu J-M, Ke Z-B, Lin Y-Z, Zheng Q-S, Wei Y, Xu N and Xue X-Y (2021) Ubiquitination-Related Molecular Subtypes and a Novel Prognostic Index for Bladder Cancer Patients. Pathol. Oncol. Res. 27:1609941. doi: 10.3389/pore.2021.1609941
Received: 08 July 2021; Accepted: 11 October 2021;
Published: 29 October 2021.
Edited by:
Anna Sebestyén, Semmelweis University, HungaryCopyright © 2021 Cai, Chen, Huang, Zhu, Ke, Lin, Zheng, Wei, Xu and Xue. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Xue-Yi Xue, xuexueyi@fjmu.edu.cn; Ning Xu, drxun@fjmu.edu.cn
†These authors have contributed equally to this work
