Integrating Stemness Features and Immune Microenvironment in a Prognostic Model for Bladder Cancer Treatment

Integrating Stemness Features and Immune Microenvironment in a Prognostic Model for Bladder Cancer Treatment | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Integrating Stemness Features and Immune Microenvironment in a Prognostic Model for Bladder Cancer Treatment Zhiwei Li, Qiqi Zhu, Yidong Cheng, Junjie Li, Zhiyang Xiao, Dong Li, and 1 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-9059350/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract Background The high heterogeneity and treatment resistance of bladder cancer (BLCA) are closely associated with cancer stemness, yet an integrated prognostic-therapeutic model incorporating stemness features and the immune microenvironment remains lacking. Methods We integrated data from the TCGA-BLCA (n = 453), GSE13507 (n = 167), and GSE32894 (n = 308) cohorts. Stemness indices were calculated via ssGSEA/PCA, and consensus clustering was employed for subtyping. A prognostic model was constructed using LASSO-Cox regression, while immune microenvironment analysis was performed via CIBERSORT/TIDE. Drug sensitivity was predicted using oncoPredict. Results Based on 35 stemness-related genes, BLCA was classified into Cluster 1 (immune-responsive, n = 173) and Cluster 2 (immune-desert, n = 230). The latter exhibited significantly elevated stemness indices, TIDE scores, and immunosuppressive cell infiltration (MDSC/CAF/TAM-M2, P < 0.01). A 7-gene prognostic model (POLE2/UTP6/XPOT/RRAS2/PLAA/CENPH/DIAPH3) was established, demonstrating that high-risk patients had a threefold shorter median survival (720 vs. 2880 days, P 1.5, P 0.85, P < 10⁻¹¹), whereas low PLAA expression enhanced sensitivity to SB-743921 (ρ=-0.33). Risk genes were found to drive poor prognosis by modulating the "M2/Treg suppression axis" (XPOT/POLE2 positively correlated, ρ > 0.37) and the "CD8⁺/NK effector axis" (PLAA/CENPH negatively correlated, ρ<-0.32). Conclusion This study establishes the first stemness-immune classification system for BLCA and a cross-subtype prognostic model, offering novel targeted therapeutic strategies (e.g.Aurora inhibitor combination therapy) for immune checkpoint inhibitor (ICI)-resistant patients. Bladder cancer (BLCA) Cancer stemness Molecular subtyping Prognostic model Tumor immune microenvironment (TME) Targeted therapy Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Introduction Bladder cancer (BLCA) is one of the most common malignancies of the urinary system, with approximately 570,000 new cases and 210,000 deaths annually worldwide, and its occurrence rate continues to rise [ 1 ]. In China, the age-standardized incidence rate of BLCA has reached 6.8 per 100,000[ 2 ], with muscle-invasive bladder cancer (MIBC) accounting for ~ 25% of new cases but contributing to over 50% of cancer-specific deaths [ 3 ]. Although cisplatin-based neoadjuvant chemotherapy combined with radical cystectomy remains the standard treatment for MIBC, ~ 50% of patients develop distant metastases within two years post-surgery, with a 5-year overall survival rate below 50% [ 4 , 5 ]. In recent years, immune checkpoint blockers (ICBs) have significantly transformed the therapeutic landscape for advanced BLCA, with programmed death-1/programmed death-ligand 1 (PD-1/PD-L1) monoclonal antibodies approved for second- and first-line treatment after chemotherapy failure [ 6 , 7 , 8 ]. However, real-world data indicate that only 20%–30% of patients exhibit durable responses to ICBs, and no reliable biomarkers currently predict long-term benefits [9.10]. Thus, elucidating the molecular basis of BLCA heterogeneity, uncovering mechanisms driving recurrence and metastasis, and developing integrated prognostic-therapeutic models remain critical challenges [ 11 , 12 ]. The concept of "cancer stemness" provides a novel perspective for understanding tumor heterogeneity and treatment failure [ 13 , 14 ]. Stemness does not refer to a single cell population but rather encompasses a phenotypic spectrum—including self-renewal, unlimited proliferation, multipotent differentiation, and stress resistance—regulated by specific genetic-epigenetic networks [ 15 , 16 ]. Single-cell transcriptomic studies confirm the existence of rare stemness-associated subpopulations in BLCA, characterized by high expression of core transcription factors (e.g., SOX2, NANOG, OCT4) and enrichment in PI3K/AKT/mTOR, Hippo-YAP/TAZ, and epithelial-mesenchymal transition (EMT) pathways [ 17 , 18 ]. Functionally, these cells exhibit marked resistance to conventional chemo- and radiotherapy and evade immune surveillance via PD-L1/CD47 upregulation and TGF-β secretion [ 19 ]. Clinically, large-scale TCGA-BLCA analyses demonstrate that stemness indices (e.g., mRNAsi, ssGSEA stemness score) correlate significantly with pathological grade, T-stage, lymph node metastasis, and poor prognosis, suggesting stemness as an independent risk factor for BLCA progression [ 20 ]. However, prior studies have focused on isolated stemness-related genes or pathways, lacking systematic integration of stemness features with tumor microenvironment (TME) interactions or validation in prospective cohorts or drug screening platforms [ 21 ]. Building on recent advances in BLCA and stemness research, we hypothesize that BLCA stemness features may serve as independent prognostic-therapeutic biomarkers beyond conventional clinical parameters. To test this, we integrated TCGA-BLCA, GEO, applying machine learning algorithms (consensus clustering, LASSO-Cox, Random Forest) to construct and validate a stemness-related prognostic index (SRPI) comprising nine core genes [ 22 ]. Mechanistically, the SRPI high-risk group exhibited an "immune-desert" phenotype with elevated tumor mutational burden (TMB) but reduced CD8 + T-cell infiltration, indicative of ICB resistance [ 23 ]. Collectively, this study translates BLCA stemness into a quantifiable prognostic-therapeutic decision tool, offering novel insights to overcome current treatment limitations . Materials and Methods 2.1. Data Collection and Integration In this study, transcriptomic data from multiple public databases were comprehensively analyzed. The primary dataset was obtained from The Cancer Genome Atlas Program (TCGA) bladder cancer (BLCA) project, comprising 453 bladder cancer tumor tissue samples and 51 adjacent normal tissue samples. These data provided transcriptomic expression profiles and corresponding clinical information for bladder cancer tissues. To validate our findings, two independent datasets were used as external validation sets. The first validation set, GSE13507, was downloaded from the Gene Expression Omnibus (GEO) database and included gene expression profiles of bladder cancer tumor tissues and normal tissues. The second validation set, GSE32894, also sourced from GEO, contained gene expression data from bladder cancer and normal tissue samples. 2.2. Stemness Feature Calculation and Gene Screening In the TCGA-BLCA transcriptomic FPKM dataset, pathway activity scores were computed using the GSVA package via single-sample gene set enrichment analysis (ssGSEA) based on the stemness gene set defined by Miranda et al. [ 22 ]. Principal component analysis (PCA) was employed to integrate multi-pathway information, with the first principal component extracted as the comprehensive stemness index. Subsequently, gene screening was performed by calculating Spearman correlations between all gene expressions and the stemness index. Genes with significant correlations (|cor| > 0.5 and P < 0.01) were defined as stemness-associated genes. Samples were stratified into tumor (01) and normal (11) groups based on the 14th–15th digits of TCGA barcodes. A two-sided Wilcoxon rank-sum test was applied to compare stemness index differences between groups. Statistical analyses were conducted using the ggpubr package, with significance thresholds set at P < 0.001, P < 0.01, and P < 0.05 (uncorrected). All analyses were performed in R 4.3.2. 2.3. Subtype Identification and Cluster Analysis To elucidate heterogeneity among urothelial carcinoma (UC) patients, consensus clustering analysis was performed to stratify UC patients into distinct subtypes. In the development set (GSE92415 UC samples), the ConsensusClusterPlus R package and K-means algorithm were applied for consensus clustering based on key genes [ 23 ]. The optimal number of clusters (ranging from 2 to 10) was determined, with two clusters generally exhibiting the best separation. Heatmaps generated using the ggplot2 package visualized expression differences of key genes across subtypes (P < 0.05). 2.4. Tumor Immune Microenvironment Assessment The Tumor Immune Dysfunction and Exclusion (TIDE) platform ( http://tide.dfci.harvard.edu ) was utilized to evaluate the immune status of TCGA-BLCA samples [ 24 ]. The TIDE algorithm computed a composite TIDE score for each sample, along with infiltration abundances of immune cells (e.g., CD8 + T cells, regulatory T cells), expression levels of immune checkpoint molecules (e.g., CD274/PD-L1, PDCD1/PD-1, CTLA4), and feature scores of immunosuppressive factors (e.g., myeloid-derived suppressor cells [MDSCs], cancer-associated fibroblasts [CAFs], M2 tumor-associated macrophages [TAMs]). Spearman or Pearson correlation analyses assessed associations between immune features (including TIDE scores) and stemness markers, visualized via heatmaps. Subtype-specific comparisons employed the Mann-Whitney U test to evaluate differences in TIDE scores, immune effector cell infiltration (e.g., CD8 + T cells), immune checkpoint expression, immunosuppressive factor scores, and immune activity markers (e.g., interferon-gamma [IFNG]-related signaling). Statistical significance was set at P < 0.05. This analysis systematically untangled immune microenvironment composition, activation states, and suppression levels across BLCA subtypes, offering insights into immune evasion mechanisms. 2.5. Differential Gene Expression Analysis The DESeq2 R package (v3.5.2) was used to analyze differential expression in the TCGA colorectal cancer dataset (453 tumors vs. 51 normal tissues). Differentially expressed genes (DEGs) were identified with thresholds of |log2FoldChange| ≥ 1 and P < 0.05. Volcano plots (ggplot2) illustrated DEG distributions, while standardized (Z-score) expression profiles of significant DEGs were visualized via heatmaps (ComplexHeatmap). 2.6. Functional Pathway and Gene Set Enrichment Analysis Gene set enrichment analysis (GSEA) was performed using the clusterProfiler R package with the MSigDB human Canonical Pathways dataset (3,917 gene sets) [ 25 ]. Results were clustered by similarity (aPEAR) and presented as network plots. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses annotated candidate genes across biological processes (BP), cellular components (CC), and molecular functions (MF). KEGG elucidated metabolic pathways and signal transduction mechanisms. 2.7. Risk Score Model Construction and Validation In the integrated TCGA-BLCA expression and clinical matrix (453 tumors, 51 normals; OS > 0, missing values excluded), LASSO-Cox regression (glmnet v4.1-8) with 10-fold cross-validation (nfolds = 10, α = 1, nlambda = 1000) selected seven non-zero-coefficient genes (lambda.min). Risk scores (Σ(βi × gene expressioni)) stratified samples into high/low groups for prognostic nomograms, Kaplan-Meier (KM) curves, and subtype-specific survival validation. 2.8. Survival Prognosis and Nomogram Development Prognostic performance was validated in GSE13507 and GSE32894 using KM survival curves. A nomogram (rms package) integrating key genes and clinical data predicted BLCA risk probabilities. Variables were mapped to "Points," summed to "Total Points," and translated to risk probabilities. 2.9. Immune Infiltration and Correlation Analysis The CIBERSORT algorithm estimated 22 immune cell infiltration levels in BLCA subtypes (samples with P > 0.05 excluded). Wilcoxon tests identified subtype-specific differences (P 0.3, P < 0.05), plotted as heatmaps (pheatmap). 2.10. Drug Sensitivity Prediction and Gene Association Analysis The oncoPredict package (v0.2) leveraged GDSC2022 data to predict IC50 values for 265 compounds in BLCA subtypes [ 26 ]. limma (v3.54) identified differentially sensitive drugs (adj. P 1). Spearman correlations (|cor| > 0.4, P < 0.001) linked key genes to drug IC50 values, highlighting candidate therapeutics. Results 3.1 Identification of Stemness Features and Molecular Subtypes in Bladder Cancer Based on the TCGA-BLCA transcriptomic dataset, this study employed the stemness gene set defined by Miranda et al. to calculate pathway activity scores via ssGSEA and integrated principal component analysis (PCA) to derive a composite stemness index. The analysis revealed that the stemness index was significantly higher in the tumor group (TCGA barcode positions 14–15: "01") compared to the normal group ("11") (Wilcoxon rank-sum test, P 0.5 and P < 0.01), which were enriched in multiple biological functions (Fig. 1 B) and subsequently used for subtype classification. Functional enrichment analysis of the stemness signature genes using ClusterProfiler revealed significant enrichment in biological processes such as cell movement, cell proliferation, and signal transduction, all of which are closely associated with the maintenance of tumor stemness and progression. In the molecular function category, enriched terms included protein binding, enzyme activity, and transcription factor activity, which are critical for gene expression regulation and cell fate determination. Cellular component analysis showed predominant enrichment in the nucleus, cytoplasm, and cell membrane, reflecting the broad functional roles of stemness genes in cellular architecture and function. Consensus clustering analysis of TCGA-BLCA tumor samples based on these stemness signature genes (Fig. 1 C) revealed two major subtypes (Cluster 1: n = 173; Cluster 2: n = 230; Fig. 1 D). Heatmap analysis demonstrated significant heterogeneity in stemness gene expression between the two subtypes, with Cluster 1 exhibiting markedly higher expression than Cluster 2, suggesting that Cluster 1 may possess stronger stemness characteristics and greater tumorigenic potential. This differential expression pattern may correlate with tumor aggressiveness, metastatic capacity, and treatment resistance, indicating that Cluster 1 represents a more aggressive and therapeutically challenging subtype. Conversely, Cluster 2 may correspond to a less stem-like subtype with potentially better treatment responsiveness. These findings provide novel insights into the role of tumor stemness in bladder cancer progression and may facilitate the development of personalized therapeutic strategies for distinct molecular subtypes. 3.2 Immune Microenvironment Profiling of Bladder Cancer Subtypes via TIDE In the TCGA-BLCA cohort (n = 453), TIDE platform analysis revealed distinct immune signatures between the two subtypes. Spearman correlation coefficients were calculated for eight core immune metrics (Fig. 2 A), showing a strong positive correlation between CD8⁺ T-cell infiltration and IFN-γ signaling (ρ = 0.813), both of which exhibited moderate positive correlations with PD-L1 (CD274) (ρ = 0.611 and 0.755, respectively), suggesting IFN-γ-driven adaptive immune upregulation alongside PD-L1-mediated feedback inhibition. Cancer-associated fibroblasts (CAFs) positively correlated with TIDE score (ρ = 0.405), whereas TAM-M2 showed a strong negative correlation with CD274 (ρ = − 0.763), indicating that pro-fibrotic CAFs and M2-polarized macrophages may cooperatively suppress antitumor immunity through distinct mechanisms. MSI.Expr.Sig exhibited weak correlations with all immune metrics (|ρ| ≤ 0.40), only showing a negative correlation with CAFs (ρ = − 0.395), implying that tumor stemness exerts limited direct regulatory effects on the immune microenvironment and requires further investigation in conjunction with other stemness dimensions. The composite TIDE score was significantly higher in Cluster 2 (n = 230) than in Cluster 1 (n = 173) (Mann-Whitney U test, P = 4.3 × 10⁻⁷; Fig. 2 B), indicating stronger immune dysfunction and immune exclusion in Cluster 2, which may confer immunotherapy resistance. All immunosuppressive components exhibited significant differential distribution between the two subtypes (Fig. 2 C). Cluster 2 displayed higher myeloid-derived suppressor cell (MDSC) infiltration (median: 0.62 vs. 0.33; P < 0.001), elevated CAF signature score (0.35 vs. 0.28, P = 0.009), and increased TAM-M2 score (0.63 vs. 0.45, P < 0.001). Combined with the positive correlation network in Fig. 2 A, these findings suggest that Cluster 2 establishes a more robust immune barrier through enrichment of the MDSC-CAF-TAM-M2 inhibitory axis, likely driving its high TIDE score. Further comparison of key immune effector molecules (Fig. 2 D) revealed significantly upregulated PD-L1 (CD274) in Cluster 2 (median: 0.61 vs. 0.46, P < 0.001), consistent with its positive correlation with TIDE score, indicating an amplified immune-suppressive feedback loop in high-TIDE subtypes. CD8⁺ T-cell infiltration was significantly higher in Cluster 1 (median: 0.40 vs. 0.25, P < 0.01) but exhibited a concomitant decline in IFN-γ signaling (Cluster 1 vs. Cluster 2: 0.81 vs. 0.75, P = 0.028), suggesting that although Cluster 1 harbors more effector T cells, their functionality may be partially suppressed by the microenvironment. MSI.Expr.Sig (stemness index) showed no significant inter-subtype difference (P = 0.17), further supporting that stemness per se is not the primary driver of immune divergence but may act indirectly via immunosuppressive pathways. In summary, our multi-dimensional immune profiling delineated two polarized BLCA subtypes: - Cluster 1 (Immune-Responsive Subtype): Low TIDE score, high CD8⁺ T-cell infiltration, active IFN-γ signaling, and potential antitumor immunity. - Cluster 2 (Immune-Desert Subtype): High TIDE score, enrichment of the MDSC-CAF-TAM-M2 inhibitory axis, PD-L1-mediated T-cell exhaustion, and predicted immunotherapy resistance. This stark contrast provides a molecular basis for subtype-specific therapeutic strategies. 3.3 Differential Gene Expression and Functional Characteristics of Tumor Subgroups To elucidate molecular differences between the two stemness subtypes, we performed DESeq2 differential expression analysis on TCGA-BLCA samples (Cluster 2: n = 230 vs. Cluster 1: n = 173), applying thresholds of |log₂FC| ≥ 1 and P < 0.05. A total of 3,872 differentially expressed genes (DEGs) were identified, including 2,097 upregulated and 1,775 downregulated genes (Fig. 3 A). Integrated GSEA, functional, and pathway enrichment analyses revealed that upregulated genes were significantly enriched in keratinization, extracellular matrix remodeling, and lipid metabolic pathways, whereas downregulated genes were associated with neuroactive ligand-receptor interactions and inhibition of cellular differentiation. Collectively, Cluster 2 exhibited a "high-keratinization, high-metabolism, low-differentiation" stemness-enhanced phenotype, whereas Cluster 1 retained more neuroendocrine and stromal homeostasis features, providing a gene expression foundation for subtype-specific therapeutic targeting. 3.4 Construction and Cross-Subtype Validation of a 7-Gene Prognostic Model for Bladder Cancer Based on LASSO-Cox Regression Using the integrated TCGA-BLCA expression-clinical dataset (n = 395 tumor samples, excluding cases with missing survival data), candidate genes were initially screened via univariate Cox regression, followed by 10-fold LASSO cross-validation (α = 1, nλ = 1,000) to construct the prognostic model. At λ = λ min, seven genes were identified as key prognostic factors: POLE2, UTP6, and XPOT (risk coefficients > 0) as well as RRAS2, PLAA, CENPH, and DIAPH3 (risk coefficients < 0) (Fig. 4 A–B). The risk score formula was defined as: risk_score = exp(0.468×POLE2 + 0.380×UTP6–0.183×RRAS2–0.439×PLAA + 0.237×DIAPH3–0.417×CENPH + 0.429×XPOT). Multivariate Cox regression confirmed the model’s significant independent prognostic value (global test P = 1.4×10 − 8, C-index = 0.690), with POLE2 (HR = 1.63, P < 0.001), XPOT (HR = 1.30, P < 0.001), and PLAA (HR = 0.64, P = 0.002) exhibiting the strongest associations (Fig. 4 C). A nomogram integrating risk scores and clinical variables further demonstrated robust predictive performance (C-index = 0.69; Fig. 4 D). Using the median risk score (0.07) as the cutoff, Kaplan-Meier analysis revealed significantly shorter median overall survival in the high-risk group (720 vs. 2,880 days, log-rank P < 0.001; Fig. 4 E). Subtype-stratified validation confirmed the model’s prognostic discrimination in both Cluster 1 (P < 0.001) and Cluster 2 (P < 0.001) (Fig. 4 F), demonstrating its robustness across stemness subtypes. In summary, the 7-gene LASSO-Cox risk model provides a reliable prognostic tool for bladder cancer patients and lays the foundation for targeted intervention strategies. 3.5 Cross-Cohort Validation and Development of an Individualized Survival Prediction Tool The 7-gene risk model (risk_score = Σβ i ×gene expression i ) was validated in two independent external cohorts: GSE13507 (n = 167) and GSE32894 (n = 308). Survival analysis showed that in GSE13507, the high-risk group (n = 83) had a median survival of 24 months, significantly lower than the low-risk group (n = 84; 48 months, log-rank P = 0.032; Fig. 5 A). Similarly, in GSE32894, the high-risk group (n = 154) exhibited further reduced median survival (48 vs. 96 months, P = 0.0017; Fig. 5 B), with both cohorts showing hazard ratios (HR) > 1.5 and sustained survival curve separation, confirming the model’s cross-platform robustness. An individualized nomogram integrating risk score, age, and pathological stage (Fig. 5 C) demonstrated high concordance between predicted and observed 1–5-year survival rates (mean absolute error = 0.014; Fig. 5 D). Decision curve analysis revealed that the model significantly improved clinical net benefit over traditional TNM staging within the 0.2–0.8 threshold probability range (AUC = 0.751; Fig. 5 E). The "Points→Total Points→Survival Probability" mapping (Fig. 5 C) enables rapid quantification of individual mortality risk (e.g., a 60-year-old Stage III patient with risk_score = 0.5 has an estimated 5-year mortality probability of ≈ 60%), facilitating tailored follow-up or combination therapy. Thus, cross-cohort validation and nomogram implementation enable precise risk stratification and clinical translation of bladder cancer prognosis. 3.6 Immune Microenvironment Characterization and Risk Gene–Immune Cell Interaction Network in Bladder Cancer CIBERSORT-based immune infiltration analysis of the TCGA-BLCA development set (n = 395) revealed significant immune polarization between Cluster 2 (n = 230) and Cluster 1 (n = 173). Cluster 2 exhibited effector cell exhaustion (CD8⁺ T cells: 1.5% vs. 2.8%, P = 0.002; activated NK cells: 0.4% vs. 1.1%, P = 0.007) and immunosuppressive cell enrichment (M2 macrophages: 15.6% vs. 8.3%, P < 0.001; Tregs: 3.4% vs. 2.0%, P = 0.004), whereas Cluster 1 maintained antigen-presenting capacity (resting DCs: 1.1% vs. 0.5%, P = 0.011), aligning with TIDE-classified "immune-desert" (Cluster 2) and "immune-responsive" (Cluster 1) phenotypes (Fig. 6 A). Spearman correlation analysis (|ρ| > 0.3, P < 0.05) further untangled the regulatory network between the seven risk genes and immune cells: - Prognosis-promoting genes (XPOT, UTP6, POLE2) correlated positively with M2 macrophages (ρ = 0.37–0.42) and suppressed CD8⁺ T cells (ρ = −0.31 to − 0.29). - Prognosis-suppressing genes (RRAS2, PLAA, CENPH) enhanced activated NK cell activity (ρ = 0.33–0.36) and attenuated Treg suppression (ρ = −0.32 to − 0.28). - DIAPH3 promoted antigen presentation via resting DC activation (ρ = 0.34) (Fig. 6 B). The risk score positively correlated with M2 macrophage abundance (ρ = 0.46, P < 0.001) and negatively with CD8⁺ T cells (ρ = −0.39, P < 0.001). After adjusting for immune variables, multivariate Cox regression retained its independent prognostic value (HR = 2.27, P < 0.001), indicating that risk genes drive poor outcomes by reshaping the "M2/Treg immunosuppressive axis" while suppressing the "CD8⁺/NK effector axis", concurrently capturing stemness-related regulatory signals beyond the immune microenvironment (Fig. 6 C). 3.7 Stemness Subtype-Specific Drug Resistance Mechanisms and Risk Gene–Drug Interaction Profiles Using the oncoPredict algorithm and GDSC2022 database, drug sensitivity to 265 compounds was predicted in the TCGA-BLCA cohort (n = 453). Differential analysis (adj. P < 0.01, |log 2 FC| > 1) revealed that the high-stemness Cluster 2 exhibited marked resistance to three chemotherapeutics: Aurora B inhibitor SB-743921 (log 2 FC = 1.12, P = 2.7×10 − 11), microtubule stabilizer docetaxel (log 2 FC = 0.98, P = 4.8×10 − 14), and paclitaxel (log 2 FC = 0.85, P = 3.2×10 − 13), with systematic rightward shifts in IC 50 values (Fig. 7 A–B), suggesting cell cycle–microtubule dynamics and Aurora kinase activation jointly mediate chemotherapy resistance. Spearman correlation (|ρ| > 0.4, P < 0.001) between risk genes (PLAA, CENPH, XPOT, RRAS2) and drug sensitivity revealed: - PLAA overexpression reduced SB-743921 IC 50 (ρ = −0.33), indicating enhanced Aurora inhibitor sensitivity. - CENPH overexpression inversely correlated with docetaxel resistance (ρ = −0.35), implicating cell cycle regulation in overcoming microtubule drug evasion. - XPOT overexpression reversed paclitaxel resistance (ρ = −0.37), suggesting nuclear-cytoplasmic transport dysfunction exacerbates insensitivity (Fig. 7 C). These findings demonstrate that Cluster 2 establishes cross-resistance via synergistic stemness pathway activation and risk gene dysregulation, whereas Cluster 1 may benefit from targeted interventions, providing a computational framework for dual stemness-risk gene biomarker-guided precision therapy. Discussion This study systematically elucidates the core drivers of bladder cancer heterogeneity from a "stemness–immune–prognosis–drug" integrative perspective, establishing a clinically actionable risk stratification and targeting system to address ICIs resistance and chemotherapy failure. Consensus clustering of 35 stemness signature genes divided TCGA-BLCA into Cluster 1 (immune-responsive) and Cluster 2 (immune-desert). Cluster 2 displayed a "high-keratinization–high-metabolism–low-differentiation" transcriptional program with upregulated EMT, lipid metabolic reprogramming, and keratinization pathways, whereas Cluster 1 retained neuroendocrine and stromal homeostasis features. This classification was validated not only transcriptionally but also immunologically: Cluster 2 was enriched with M2 macrophages, CAFs, and Tregs but depleted of CD8⁺ T and activated NK cells, forming an immunosuppressive loop. Notably, limited correlation between stemness indices and immune metrics suggests stemness indirectly shapes the immune landscape via extracellular matrix remodeling and inhibitory cell recruitment, expanding the traditional "stemness–immune" unidirectional model to emphasize microenvironment-mediated mechanisms. The 7-gene LASSO-Cox model (POLE2, UTP6, XPOT [prognosis-promoting]; RRAS2, PLAA, CENPH, DIAPH3 [prognosis-suppressing]) robustly discriminated survival in the TCGA training set (C-index = 0.69) and external cohorts (HR > 1.5, P < 0.05), confirming cross-platform generalizability. Its independence from stemness subtypes, age, or stage implies capture of biological features beyond conventional clinical variables. Mechanistically, XPOT/POLE2 activate the "M2/Treg axis" to impair antitumor immunity, whereas PLAA/CENPH sustain the "CD8⁺/NK axis", forming a dual-axis regulatory network. Correlation with immune cells further reveals the model’s integration of stemness-cycling and microenvironmental signals, ensuring prognostic precision. Cluster 2’s resistance to microtubule inhibitors (paclitaxel, docetaxel) and SB-743921 (elevated IC 50) implicates cell cycle–microtubule dynamics and Aurora kinase activation in therapeutic escape. Key risk genes negatively correlated with drug resistance: PLAA enhanced SB-743921 sensitivity (ρ = −0.33), CENPH reversed docetaxel resistance (ρ = −0.35), and XPOT reduced paclitaxel IC 50 (ρ = −0.37). These findings validate risk gene functionality and propose a dual "stemness-risk gene" biomarker framework: - Cluster 2: Aurora inhibitor + microtubule stabilizer combinations to overcome resistance. - Low-risk patients: Mono-/immunotherapy strategies. Additionally, risk gene overexpression in immune-desert subtypes suggests their targeting may synergize with ICIs to disrupt immunosuppression, enabling "stemness–immune" dual intervention. Limitations and Future Directions The retrospective nature of public data warrants prospective trials to validate the model’s predictive power and targeting safety. Single-cell multi-omics could delineate spatiotemporal dynamics of risk gene–immune cell crosstalk and optimize combination therapy timing/dosing. Conclusion This study establishes the first "stemness subtyping–prognostic modeling–targeted therapy" decision framework for bladder cancer. The 7-gene risk score overcomes ICIs response prediction barriers, while Cluster 2’s drug sensitivity profile unveils novel targets to reverse chemotherapy resistance, advancing precision therapy from "ICI monotherapy" to "stemness–immune–targeted" tripartite strategies. Declarations Ethics approval and consent to participate Not applicable. Competing interests The authors declare that they have no competing interests. Funding Not applicable. Author Contribution Z L: Writing–original draft, Writing–review & editing, Software, Methodology, Conceptualization, Data curation, Formal analysis. Q Z: Validation, Writing–original draft, Writing–review & editing, Methodology, Visualization. D C: Resources, Writing–review & editing, Software, Visualization, Formal analysis. J L: Resources, Writing–review & editing, Software, Visualization, Formal analysis. Z X:Writing–review & editing, Resources. D L: Writing–review & editing, Resources. Y X :Supervision, Conceptualization, Writing–review & editing. Acknowledgements Not applicable. Data Availability All data were sourced from public databases: TCGA-BLCA (via GDC Portal), GSE13507, and GSE32894 (via GEO). Drug sensitivity data were derived from GDSC2022 (oncoPredict R package). TIDE scores were computed via http://tide.dfci.harvard.edu. Data access complies with original database policies, ensuring transparency and reproducibility. The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request. References Siegel RL, Miller KD, Fuchs HE, Jemal A, Cancer statistics. 2023. CA: A Cancer Journal for Clinicians. 2023;73(1):17–48. Sun L, Zhao K, Liu X, Meng X. Global, regional, and national burden of esophageal cancer using the 2019 global burden of disease study. Sci Rep. 2025;15(1):3284. Boegemann M, Krabbe LM. Prognostic Implications of Immunohistochemical Biomarkers in Non-muscle-invasive Blad Cancer and Muscle-invasive Bladder Cancer. Mini Rev Med Chem. 2020;20(12):1133–52. Funt SA, Rosenberg JE. Systemic, perioperative management of muscle-invasive bladder cancer and future horizons. Nat Rev Clin Oncol. 2017;14(4):221–34. Kim KH, Lee HW, Ha HK, Seo HK. Perioperative systemic therapy in muscle invasive bladder cancer: Current standard method, biomarkers and emerging strategies. Investig Clin Urol. 2023;64(3):202–18. Funt SA, Lattanzi M, Whiting K, et al. Neoadjuvant atezolizumab with gemcitabine and cisplatin in patients with muscle-invasive bladder cancer: a multicenter, single-arm, phase II trial. J Clin Oncol. 2022;40(12):1312–22. Chen M, Jiang M, Wang X, Shen L, Li J. Envafolimab - first PD-1/PD-L1 antibody to be administered by subcutaneous injection for microsatellite instability-high or deficient mismatch repair advanced solid tumors. Expert Opin Biol Ther. 2022;22(10):1227–32. Wu Y, Mou H. Complete remission in a patient with advanced renal pelvis carcinoma with lung metastasis treated with durvalumab immunotherapy: a case report. Ann Palliat Med. 2020;9(4):2408–13. Powles T, Rosenberg JE, Sonpavde CP, et al. Enfortumab vedotin in previously treated advanced urothelial carcinoma. N Engl J Med. 2021;384(12):1125–35. Fan B, Zheng X, Wang Y, Hu X. Predicting prognosis and clinical efficacy of immune checkpoint blockade therapy via interferon-alpha response in muscle-invasive bladder cancer. Pathol Oncol Res. 2023;29:1611117. Mariathasan S, Turley SJ, Nickles D, et al. TGF-β attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature. 2018;554(7693):544–8. Chen Y, Yan Z, Li L, Liang Y, Wei X, Zhao Y, Cao Y, Zhang H, Tang L. Identification and validation of a 9-RBPs-related gene signature associated with prognosis and immune infiltration in bladder cancer based on bioinformatics analysis and machine learning. Transl Androl Urol. 2025;14(4):1066–81. Loh JJ, Ma S. Hallmarks of cancer stemness. Cell Stem Cell. 2024;31(5):617–39. Sato R, Semba T, Saya H, Arima Y. Concise Review: Stem Cells and Epithelial-Mesenchymal Transition in Cancer: Biological Implications and Therapeutic Targets. Stem Cells. 2016;34(8):1997–2007. Lin G, Lin J, Wang X, et al. Characterization of the stem cell landscape and identification of a stemness-associated prognostic signature in bladder cancer. Cancer Cell Int. 2024;24:45. Hass R, von der Ohe J, Ungefroren H. Impact of the Tumor Microenvironment on Tumor Heterogeneity and Consequences for Cancer Cell Plasticity and Stemness. Cancers (Basel). 2020;12(12):3716. Zhang Y, Zhang X, Li H, et al. Tumor stemness score to estimate epithelial-to-mesenchymal transition and cancer stem cells characterization and to predict prognosis and immunotherapy response in bladder urothelial carcinoma. Stem Cell Res Ther. 2023;14(1):15. Sathe A, Guerth F, Cronauer MV, Heck MM, Thalgott M, Gschwend JE, Retz M, Nawroth R. Mutant PIK3CA controls DUSP1-dependent ERK 1/2 activity to confer response to AKT target therapy. Br J Cancer. 2014;111(11):2103–13. Fu S, Tan Z, Shi H, et al. Development of a stemness-related prognostic index to provide therapeutic strategies for bladder cancer. NPJ Precision Oncol. 2024;8(1):14. Miranda A, Hamilton PT, Zhang AW, et al. Cancer stemness, intratumoral heterogeneity, and immune response across cancers. Proc Natl Acad Sci U S A. 2019;116(18):9020–9. Lai Y, Lin P, Lin F, et al. Identification of immune microenvironment subtypes and signature genes for Alzheimer's disease diagnosis and risk prediction based on explainable machine learning. Front Immunol. 2022;13:1046410. Fu S, Tan Z, Shi H, Chen J, Zhang Y, Guo C, Feng W, Xu H, Wang J, Wang H. Development of a stemness-related prognostic index to provide therapeutic strategies for bladder cancer. NPJ Precis Oncol. 2024;8(1):14. 10.1038/s41698-024-00510-3 . Mellman I, Chen DS, Powles T, Turley SJ. The cancer-immunity cycle: Indication, genotype, and immunotype. Immunity. 2023;56(10):2188–205. Fu J, Li K, Zhang W, et al. Large-scale public data reuse to model immunotherapy response and resistance. Genome Med. 2020;12(1):21. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst. 2015;1(6):417–25. Maeser D, Gruener RF, Huang RS. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform. 2021;22(6):bbab260. Additional Declarations No competing interests reported. Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-9059350","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":629599552,"identity":"2848cec3-def8-4c1b-b126-6a7cc8f7f214","order_by":0,"name":"Zhiwei Li","email":"","orcid":"","institution":"Nanjing Integrated Traditional Chinese and Western Medicine Hospital","correspondingAuthor":false,"prefix":"","firstName":"Zhiwei","middleName":"","lastName":"Li","suffix":""},{"id":629599553,"identity":"0332469e-a64e-43ea-a323-a2ce317d3ea7","order_by":1,"name":"Qiqi Zhu","email":"","orcid":"","institution":"The First Affiliated Hospital of Nanchang University","correspondingAuthor":false,"prefix":"","firstName":"Qiqi","middleName":"","lastName":"Zhu","suffix":""},{"id":629599554,"identity":"bafae5d9-d822-4b37-8983-b662c48924bb","order_by":2,"name":"Yidong Cheng","email":"","orcid":"","institution":"Nanjing Hospital of Chinese Medicine, Nanjing University of Chinese Medicine","correspondingAuthor":false,"prefix":"","firstName":"Yidong","middleName":"","lastName":"Cheng","suffix":""},{"id":629599555,"identity":"7311eade-6772-4987-adac-86c5feffa25e","order_by":3,"name":"Junjie Li","email":"","orcid":"","institution":"Jinggangshan University","correspondingAuthor":false,"prefix":"","firstName":"Junjie","middleName":"","lastName":"Li","suffix":""},{"id":629599556,"identity":"957b64a3-012b-4496-b699-96a1f08ddf95","order_by":4,"name":"Zhiyang Xiao","email":"","orcid":"","institution":"Ganzhou Hospital of Traditional Chinese Medicine","correspondingAuthor":false,"prefix":"","firstName":"Zhiyang","middleName":"","lastName":"Xiao","suffix":""},{"id":629599557,"identity":"f8418c64-99e7-4af1-98db-27ecdacd4a7d","order_by":5,"name":"Dong Li","email":"","orcid":"","institution":"Nanjing Integrated Traditional Chinese and Western Medicine Hospital","correspondingAuthor":false,"prefix":"","firstName":"Dong","middleName":"","lastName":"Li","suffix":""},{"id":629599558,"identity":"94543026-30c5-4650-8bd2-1902fb4cda1a","order_by":6,"name":"Yufeng Xu","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAsElEQVRIiWNgGAWjYHAC5gcfKmx4+PkbiNfCZjjjTJqM5IwDJFgjzdt22MagIYFI5QbXzhgY8Lad5zFgOMD44WMOMVpupyU8kDh3m8ecuYFZcuY2IrSY3U4+YGBQdpvHsuEAGzMvcVoSGyQS2M7xGBxIIFpL8gGJA20HSNBifzstzbDhTDKP5IyDzcT5RXJ2jvHjPxV29vz8zQc/fCRGCxJgbCBN/SgYBaNgFIwC3AAAd8045LAYj3AAAAAASUVORK5CYII=","orcid":"","institution":"Nanjing Integrated Traditional Chinese and Western Medicine Hospital","correspondingAuthor":true,"prefix":"","firstName":"Yufeng","middleName":"","lastName":"Xu","suffix":""}],"badges":[],"createdAt":"2026-03-07 15:08:27","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-9059350/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-9059350/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":108181279,"identity":"dd3d72e3-7d0a-4151-9858-6123ee61a016","added_by":"auto","created_at":"2026-04-30 08:58:30","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":4713105,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eStemness-Related Clustering and Immune Landscape in TCGA-BLCA Cohort.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A) Violin plot illustrating the distribution of the Stemness Index in normal (n=51) and tumor (n=453) samples from the TCGA-BLCA cohort. The Stemness Index, calculated using ssGSEA and PCA, shows significantly higher stemness in tumor samples compared to normal tissues (**P \u0026lt; 0.01).(B) Circular heatmap depicting the consensus clustering results for the 35 stemness-related genes across all samples. The heatmap visualizes the expression patterns and clustering stability, highlighting distinct gene expression profiles in different clusters.(C)Bar plot showing the consensus index for different cluster numbers (k=2 to k=4). The optimal number of clusters (k=4) was determined based on the highest consensus index and stability of sample assignments.(D) Dendrogram representing the hierarchical clustering of samples based on the expression of stemness-related genes, forming two major clusters (Cluster 1: n=173, Cluster 2: n=230). The dendrogram indicates the similarity between samples within each cluster.(E) Heatmap showing the expression levels of the 35 stemness-related genes in the two identified clusters. The heatmap reveals significant differences in gene expression between Cluster 1 and Cluster 2, with Cluster 1 exhibiting higher expression of stemness-related genes, suggesting a stronger stemness phenotype.\u003c/p\u003e","description":"","filename":"FIG1.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/503487411da71f402e4a374a.png"},{"id":108071257,"identity":"98967387-8f37-4f8c-a7e2-8039d9cca30c","added_by":"auto","created_at":"2026-04-29 06:07:31","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":1968943,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eImmune Landscape and Stemness-Related Gene Expression in TCGA-BLCA Cohort.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A) A correlation heatmap showing the relationships between various immune cell types and stemness-related genes. The color gradient represents the strength and direction of the correlation, with darker colors indicating stronger correlations.(B) A violin plot depicting the distribution of the Tumor Immune Dysfunction and Exclusion (TIDE) scores between Cluster 1 and Cluster 2. The plot reveals significant differences in TIDE scores between the two clusters, with Cluster 2 showing higher TIDE scores, suggesting a more immunosuppressive microenvironment.(C) Violin plots illustrating the scores of immune cells including CD274, CD8, IFNG, and MSI.Expr.Sig across Cluster 1 and Cluster 2. The plots show significant differences in immune cell scores between the two clusters, with Cluster 2 generally exhibiting higher scores for these immune cells.(D) Violin plots depicting the scores of immunosuppressive cells including MDSC, CAF, and TAM.M2 across Cluster 1 and Cluster 2. The plots indicate that Cluster 1 has higher scores for these suppressor cells, suggesting a more immunoactive microenvironment in Cluster 1 compared to Cluster 2.\u003c/p\u003e","description":"","filename":"FIG2.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/79e1b3be986007c8967ab0e7.png"},{"id":108071258,"identity":"7f368b79-8968-45ff-9e01-a6ef7850d1d2","added_by":"auto","created_at":"2026-04-29 06:07:31","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":9647313,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDifferential Gene Expression and Functional Enrichment Analysis.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A) A volcano plot showing the distribution of differentially expressed genes (DEGs) between the two clusters. The plot highlights the significant upregulated (red dots) and downregulated (blue dots) genes, with the x-axis representing the log2 fold change and the y-axis representing the negative log10 of the adjusted p-value.(B) A bubble plot depicting the results of the pathway enrichment analysis for the differentially expressed genes. The plot shows the normalized enrichment score (NES) for each pathway, with the size of the bubbles representing the number of genes in the pathway and the color indicating the significance level.(C) A circular plot illustrating the GO enrichment analysis results, which categorizes genes into different biological processes, molecular functions, and cellular components. The plot shows the distribution of enriched GO terms across the three categories, with the color gradient representing the significance level.(D) A bubble plot showing the results of the KEGG pathway enrichment analysis for the differentially expressed genes. The plot highlights significant pathways with the size of the bubbles representing the number of genes in the pathway and the color indicating the significance level.\u003c/p\u003e","description":"","filename":"FIG3.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/68eaad777ff9fa981a0be49e.png"},{"id":108977028,"identity":"c9dca594-03db-4e4e-aa03-65ed7a611ba5","added_by":"auto","created_at":"2026-05-11 11:30:01","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":4612745,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003ePrognostic Modeling and Validation in TCGA-BLCA Cohort.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A)The plot shows the coefficients of various genes along the LASSO-Cox regularization path. Genes with non-zero coefficients at the optimal lambda value are selected for the prognostic model.(B) This plot illustrates the partial likelihood deviance for a range of lambda values, highlighting the lambda.min that provides the best model fit with cross-validation.(C) Forest plot displaying the hazard ratios (HR) and 95% confidence intervals for each gene in the prognostic model. Significant genes with HR \u0026gt; 1 or \u0026lt; 1 indicate their prognostic value.(D) The plot shows the distribution of the risk scores and their association with survival outcomes. The cutoff value of 0.07 is used to separate samples into high and low risk groups. The box plot below the scatter plot displays the expression levels of the genes in the high and low risk groups.(E) Kaplan-Meier curves demonstrating the overall survival of samples divided by risk score (left panel), and within each identified cluster (middle and right panels). The log-rank test p-values indicate the significance of the survival difference between the high and low risk groups.\u003c/p\u003e","description":"","filename":"FIG4.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/9b656f57a5ac2b69f16a0475.png"},{"id":108071260,"identity":"c3f9831a-41a6-41c7-9a03-4abb40522a49","added_by":"auto","created_at":"2026-04-29 06:07:31","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":2796266,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eValidation of Prognostic Model and Individualized Risk Prediction.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A) Kaplan-Meier survival curves for the GSE13507 dataset, comparing high and low risk groups based on the prognostic model. The log-rank test p-value indicates a significant difference in survival (p = 0.032).(B) Kaplan-Meier survival curves for the GSE32894 dataset, similarly comparing high and low risk groups. The survival difference is also statistically significant (p = 0.0017).(C) A nomogram integrating risk score, age, and stage to predict individual patient survival probability. The total points are calculated by summing the points assigned to each variable.(D) A calibration plot assessing the accuracy of the prognostic model, showing the observed probability against the predicted probability. The plot indicates good agreement between observed and predicted outcomes.(E) A decision curve analysis plot demonstrating the clinical utility of the prognostic model compared to other factors including age, stage, and the number of samples. The model shows a higher net benefit across a range of risk thresholds.(F) An ROC curve illustrating the model's ability to discriminate between high and low risk groups, with an area under the curve (AUC) of 0.751, indicating good predictive performance.\u003c/p\u003e","description":"","filename":"FIG5.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/79d9c91e6825d9116ff9391b.png"},{"id":108181991,"identity":"90b6e9f0-6a21-438d-bb51-eef33cda55a5","added_by":"auto","created_at":"2026-04-30 08:59:03","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":5062440,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eAnalysis of Immune Cell Infiltration and Correlation with Stemness-Related Genes.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A) A heatmap depicting the correlation between various immune cell types and stemness-related genes. The color gradient represents the strength and direction of the correlation, with darker colors indicating stronger correlations.(B) Box plots showing the composition of different immune cell types in the tumor microenvironment (TME) for samples with high and low stemness scores. Statistical significance is indicated by asterisks.(C) A matrix displaying the correlation between stemness-related genes (XPOT, UP6, RRAS2, POLE2, PLAA, DIAPH3) and various immune cell types. The color gradient represents the correlation coefficient, and the stars denote statistical significance (p \u0026lt; 0.05, *p \u0026lt; 0.01, **p \u0026lt; 0.001).\u003c/p\u003e","description":"","filename":"FIG6.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/9dd04e3fe368a3573a7d58cb.png"},{"id":108071261,"identity":"42cf90cc-9a90-4833-9175-e37abedb4b80","added_by":"auto","created_at":"2026-04-29 06:07:32","extension":"png","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":3371485,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDrug Sensitivity Prediction and Gene-Drug Correlation Analysis.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e(A)The volcano plot displays the differential expression of genes between two conditions, with red dots indicating upregulated genes and grey dots representing no significant change. The plot highlights genes such as SB-743921 and docetaxel that show significant differential expression.(B) Violin plots for paclitaxel, SB-743921, and docetaxel showing the expression levels of these drugs in Cluster 1 and Cluster 2. The plots indicate significant differences in expression between the two clusters for each drug, with asterisks denoting statistical significance.(C) catter plots depicting the relationship between the expression levels of PLAA, CENPH, and XPOT genes and the sensitivity to paclitaxel, SB-743921, and docetaxel, respectively. The plots show negative correlations, suggesting that higher expression levels of these genes are associated with lower drug sensitivity.\u003c/p\u003e","description":"","filename":"FIG7.png","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/b1068d1b24030d1ca1f47b54.png"},{"id":108979583,"identity":"5d76601c-05a6-4c95-b122-d3e27044c779","added_by":"auto","created_at":"2026-05-11 11:59:49","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":28723840,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-9059350/v1/3d6cb008-6005-4cad-ab8f-0005e946a01a.pdf"}],"financialInterests":"No competing interests reported.","formattedTitle":"Integrating Stemness Features and Immune Microenvironment in a Prognostic Model for Bladder Cancer Treatment","fulltext":[{"header":"Introduction","content":"\u003cp\u003eBladder cancer (BLCA) is one of the most common malignancies of the urinary system, with approximately 570,000 new cases and 210,000 deaths annually worldwide, and its occurrence rate continues to rise [\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e]. In China, the age-standardized incidence rate of BLCA has reached 6.8 per 100,000[\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e], with muscle-invasive bladder cancer (MIBC) accounting for ~\u0026thinsp;25% of new cases but contributing to over 50% of cancer-specific deaths [\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e]. Although cisplatin-based neoadjuvant chemotherapy combined with radical cystectomy remains the standard treatment for MIBC, ~\u0026thinsp;50% of patients develop distant metastases within two years post-surgery, with a 5-year overall survival rate below 50% [\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e, \u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e]. In recent years, immune checkpoint blockers (ICBs) have significantly transformed the therapeutic landscape for advanced BLCA, with programmed death-1/programmed death-ligand 1 (PD-1/PD-L1) monoclonal antibodies approved for second- and first-line treatment after chemotherapy failure [\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e, \u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e, \u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e]. However, real-world data indicate that only 20%\u0026ndash;30% of patients exhibit durable responses to ICBs, and no reliable biomarkers currently predict long-term benefits [9.10]. Thus, elucidating the molecular basis of BLCA heterogeneity, uncovering mechanisms driving recurrence and metastasis, and developing integrated prognostic-therapeutic models remain critical challenges [\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e, \u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eThe concept of \"cancer stemness\" provides a novel perspective for understanding tumor heterogeneity and treatment failure [\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e, \u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e]. Stemness does not refer to a single cell population but rather encompasses a phenotypic spectrum\u0026mdash;including self-renewal, unlimited proliferation, multipotent differentiation, and stress resistance\u0026mdash;regulated by specific genetic-epigenetic networks [\u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e, \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e]. Single-cell transcriptomic studies confirm the existence of rare stemness-associated subpopulations in BLCA, characterized by high expression of core transcription factors (e.g., SOX2, NANOG, OCT4) and enrichment in PI3K/AKT/mTOR, Hippo-YAP/TAZ, and epithelial-mesenchymal transition (EMT) pathways [\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e, \u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e]. Functionally, these cells exhibit marked resistance to conventional chemo- and radiotherapy and evade immune surveillance via PD-L1/CD47 upregulation and TGF-β secretion [\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e]. Clinically, large-scale TCGA-BLCA analyses demonstrate that stemness indices (e.g., mRNAsi, ssGSEA stemness score) correlate significantly with pathological grade, T-stage, lymph node metastasis, and poor prognosis, suggesting stemness as an independent risk factor for BLCA progression [\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e]. However, prior studies have focused on isolated stemness-related genes or pathways, lacking systematic integration of stemness features with tumor microenvironment (TME) interactions or validation in prospective cohorts or drug screening platforms [\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eBuilding on recent advances in BLCA and stemness research, we hypothesize that BLCA stemness features may serve as independent prognostic-therapeutic biomarkers beyond conventional clinical parameters. To test this, we integrated TCGA-BLCA, GEO, applying machine learning algorithms (consensus clustering, LASSO-Cox, Random Forest) to construct and validate a stemness-related prognostic index (SRPI) comprising nine core genes [\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e]. Mechanistically, the SRPI high-risk group exhibited an \"immune-desert\" phenotype with elevated tumor mutational burden (TMB) but reduced CD8\u003csup\u003e+\u003c/sup\u003eT-cell infiltration, indicative of ICB resistance [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e]. Collectively, this study translates BLCA stemness into a quantifiable prognostic-therapeutic decision tool, offering novel insights to overcome current treatment limitations .\u003c/p\u003e"},{"header":"Materials and Methods","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003e2.1. Data Collection and Integration\u003c/h2\u003e \u003cp\u003eIn this study, transcriptomic data from multiple public databases were comprehensively analyzed. The primary dataset was obtained from The Cancer Genome Atlas Program (TCGA) bladder cancer (BLCA) project, comprising 453 bladder cancer tumor tissue samples and 51 adjacent normal tissue samples. These data provided transcriptomic expression profiles and corresponding clinical information for bladder cancer tissues.\u003c/p\u003e \u003cp\u003eTo validate our findings, two independent datasets were used as external validation sets. The first validation set, GSE13507, was downloaded from the Gene Expression Omnibus (GEO) database and included gene expression profiles of bladder cancer tumor tissues and normal tissues. The second validation set, GSE32894, also sourced from GEO, contained gene expression data from bladder cancer and normal tissue samples.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec4\" class=\"Section2\"\u003e \u003ch2\u003e2.2. Stemness Feature Calculation and Gene Screening\u003c/h2\u003e \u003cp\u003eIn the TCGA-BLCA transcriptomic FPKM dataset, pathway activity scores were computed using the GSVA package via single-sample gene set enrichment analysis (ssGSEA) based on the stemness gene set defined by Miranda et al. [\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e]. Principal component analysis (PCA) was employed to integrate multi-pathway information, with the first principal component extracted as the comprehensive stemness index. Subsequently, gene screening was performed by calculating Spearman correlations between all gene expressions and the stemness index. Genes with significant correlations (|cor| \u0026gt; 0.5 and P\u0026thinsp;\u0026lt;\u0026thinsp;0.01) were defined as stemness-associated genes.\u003c/p\u003e \u003cp\u003eSamples were stratified into tumor (01) and normal (11) groups based on the 14th\u0026ndash;15th digits of TCGA barcodes. A two-sided Wilcoxon rank-sum test was applied to compare stemness index differences between groups. Statistical analyses were conducted using the ggpubr package, with significance thresholds set at P\u0026thinsp;\u0026lt;\u0026thinsp;0.001, P\u0026thinsp;\u0026lt;\u0026thinsp;0.01, and P\u0026thinsp;\u0026lt;\u0026thinsp;0.05 (uncorrected). All analyses were performed in R 4.3.2.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec5\" class=\"Section2\"\u003e \u003ch2\u003e2.3. Subtype Identification and Cluster Analysis\u003c/h2\u003e \u003cp\u003eTo elucidate heterogeneity among urothelial carcinoma (UC) patients, consensus clustering analysis was performed to stratify UC patients into distinct subtypes. In the development set (GSE92415 UC samples), the ConsensusClusterPlus R package and K-means algorithm were applied for consensus clustering based on key genes [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e]. The optimal number of clusters (ranging from 2 to 10) was determined, with two clusters generally exhibiting the best separation. Heatmaps generated using the ggplot2 package visualized expression differences of key genes across subtypes (P\u0026thinsp;\u0026lt;\u0026thinsp;0.05).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec6\" class=\"Section2\"\u003e \u003ch2\u003e2.4. Tumor Immune Microenvironment Assessment\u003c/h2\u003e \u003cp\u003eThe Tumor Immune Dysfunction and Exclusion (TIDE) platform (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://tide.dfci.harvard.edu\u003c/span\u003e\u003cspan address=\"http://tide.dfci.harvard.edu\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) was utilized to evaluate the immune status of TCGA-BLCA samples [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e]. The TIDE algorithm computed a composite TIDE score for each sample, along with infiltration abundances of immune cells (e.g., CD8\u0026thinsp;+\u0026thinsp;T cells, regulatory T cells), expression levels of immune checkpoint molecules (e.g., CD274/PD-L1, PDCD1/PD-1, CTLA4), and feature scores of immunosuppressive factors (e.g., myeloid-derived suppressor cells [MDSCs], cancer-associated fibroblasts [CAFs], M2 tumor-associated macrophages [TAMs]). Spearman or Pearson correlation analyses assessed associations between immune features (including TIDE scores) and stemness markers, visualized via heatmaps.\u003c/p\u003e \u003cp\u003eSubtype-specific comparisons employed the Mann-Whitney U test to evaluate differences in TIDE scores, immune effector cell infiltration (e.g., CD8\u0026thinsp;+\u0026thinsp;T cells), immune checkpoint expression, immunosuppressive factor scores, and immune activity markers (e.g., interferon-gamma [IFNG]-related signaling). Statistical significance was set at P\u0026thinsp;\u0026lt;\u0026thinsp;0.05. This analysis systematically untangled immune microenvironment composition, activation states, and suppression levels across BLCA subtypes, offering insights into immune evasion mechanisms.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec7\" class=\"Section2\"\u003e \u003ch2\u003e2.5. Differential Gene Expression Analysis\u003c/h2\u003e \u003cp\u003eThe DESeq2 R package (v3.5.2) was used to analyze differential expression in the TCGA colorectal cancer dataset (453 tumors vs. 51 normal tissues). Differentially expressed genes (DEGs) were identified with thresholds of |log2FoldChange| \u0026ge; 1 and P\u0026thinsp;\u0026lt;\u0026thinsp;0.05. Volcano plots (ggplot2) illustrated DEG distributions, while standardized (Z-score) expression profiles of significant DEGs were visualized via heatmaps (ComplexHeatmap).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003ch2\u003e2.6. Functional Pathway and Gene Set Enrichment Analysis\u003c/h2\u003e \u003cp\u003eGene set enrichment analysis (GSEA) was performed using the clusterProfiler R package with the MSigDB human Canonical Pathways dataset (3,917 gene sets) [\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e]. Results were clustered by similarity (aPEAR) and presented as network plots. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses annotated candidate genes across biological processes (BP), cellular components (CC), and molecular functions (MF). KEGG elucidated metabolic pathways and signal transduction mechanisms.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec9\" class=\"Section2\"\u003e \u003ch2\u003e2.7. Risk Score Model Construction and Validation\u003c/h2\u003e \u003cp\u003eIn the integrated TCGA-BLCA expression and clinical matrix (453 tumors, 51 normals; OS\u0026thinsp;\u0026gt;\u0026thinsp;0, missing values excluded), LASSO-Cox regression (glmnet v4.1-8) with 10-fold cross-validation (nfolds\u0026thinsp;=\u0026thinsp;10, α\u0026thinsp;=\u0026thinsp;1, nlambda\u0026thinsp;=\u0026thinsp;1000) selected seven non-zero-coefficient genes (lambda.min). Risk scores (Σ(βi \u0026times; gene expressioni)) stratified samples into high/low groups for prognostic nomograms, Kaplan-Meier (KM) curves, and subtype-specific survival validation.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec10\" class=\"Section2\"\u003e \u003ch2\u003e2.8. Survival Prognosis and Nomogram Development\u003c/h2\u003e \u003cp\u003ePrognostic performance was validated in GSE13507 and GSE32894 using KM survival curves. A nomogram (rms package) integrating key genes and clinical data predicted BLCA risk probabilities. Variables were mapped to \"Points,\" summed to \"Total Points,\" and translated to risk probabilities.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003e2.9. Immune Infiltration and Correlation Analysis\u003c/h2\u003e \u003cp\u003eThe CIBERSORT algorithm estimated 22 immune cell infiltration levels in BLCA subtypes (samples with P\u0026thinsp;\u0026gt;\u0026thinsp;0.05 excluded). Wilcoxon tests identified subtype-specific differences (P\u0026thinsp;\u0026lt;\u0026thinsp;0.05), visualized via ggplot2. Spearman correlations (psych package) assessed immune cell-gene relationships (|cor| \u0026gt; 0.3, P\u0026thinsp;\u0026lt;\u0026thinsp;0.05), plotted as heatmaps (pheatmap).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003e2.10. Drug Sensitivity Prediction and Gene Association Analysis\u003c/h2\u003e \u003cp\u003eThe oncoPredict package (v0.2) leveraged GDSC2022 data to predict IC50 values for 265 compounds in BLCA subtypes [\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e]. limma (v3.54) identified differentially sensitive drugs (adj. P\u0026thinsp;\u0026lt;\u0026thinsp;0.01, |log2FC| \u0026gt; 1). Spearman correlations (|cor| \u0026gt; 0.4, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001) linked key genes to drug IC50 values, highlighting candidate therapeutics.\u003c/p\u003e \u003c/div\u003e"},{"header":"Results","content":"\u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003e3.1 Identification of Stemness Features and Molecular Subtypes in Bladder Cancer\u003c/h2\u003e \u003cp\u003eBased on the TCGA-BLCA transcriptomic dataset, this study employed the stemness gene set defined by Miranda et al. to calculate pathway activity scores via ssGSEA and integrated principal component analysis (PCA) to derive a composite stemness index. The analysis revealed that the stemness index was significantly higher in the tumor group (TCGA barcode positions 14\u0026ndash;15: \"01\") compared to the normal group (\"11\") (Wilcoxon rank-sum test, P\u0026thinsp;\u0026lt;\u0026thinsp;0.01; Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eA). Further Spearman correlation analysis identified 35 stemness signature genes (|ρ| \u0026gt; 0.5 and P\u0026thinsp;\u0026lt;\u0026thinsp;0.01), which were enriched in multiple biological functions (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eB) and subsequently used for subtype classification.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eFunctional enrichment analysis of the stemness signature genes using ClusterProfiler revealed significant enrichment in biological processes such as cell movement, cell proliferation, and signal transduction, all of which are closely associated with the maintenance of tumor stemness and progression. In the molecular function category, enriched terms included protein binding, enzyme activity, and transcription factor activity, which are critical for gene expression regulation and cell fate determination. Cellular component analysis showed predominant enrichment in the nucleus, cytoplasm, and cell membrane, reflecting the broad functional roles of stemness genes in cellular architecture and function.\u003c/p\u003e \u003cp\u003eConsensus clustering analysis of TCGA-BLCA tumor samples based on these stemness signature genes (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eC) revealed two major subtypes (Cluster 1: n\u0026thinsp;=\u0026thinsp;173; Cluster 2: n\u0026thinsp;=\u0026thinsp;230; Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eD). Heatmap analysis demonstrated significant heterogeneity in stemness gene expression between the two subtypes, with Cluster 1 exhibiting markedly higher expression than Cluster 2, suggesting that Cluster 1 may possess stronger stemness characteristics and greater tumorigenic potential. This differential expression pattern may correlate with tumor aggressiveness, metastatic capacity, and treatment resistance, indicating that Cluster 1 represents a more aggressive and therapeutically challenging subtype. Conversely, Cluster 2 may correspond to a less stem-like subtype with potentially better treatment responsiveness. These findings provide novel insights into the role of tumor stemness in bladder cancer progression and may facilitate the development of personalized therapeutic strategies for distinct molecular subtypes.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section2\"\u003e \u003ch2\u003e3.2 Immune Microenvironment Profiling of Bladder Cancer Subtypes via TIDE\u003c/h2\u003e \u003cp\u003eIn the TCGA-BLCA cohort (n\u0026thinsp;=\u0026thinsp;453), TIDE platform analysis revealed distinct immune signatures between the two subtypes. Spearman correlation coefficients were calculated for eight core immune metrics (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eA), showing a strong positive correlation between CD8⁺ T-cell infiltration and IFN-γ signaling (ρ\u0026thinsp;=\u0026thinsp;0.813), both of which exhibited moderate positive correlations with PD-L1 (CD274) (ρ\u0026thinsp;=\u0026thinsp;0.611 and 0.755, respectively), suggesting IFN-γ-driven adaptive immune upregulation alongside PD-L1-mediated feedback inhibition. Cancer-associated fibroblasts (CAFs) positively correlated with TIDE score (ρ\u0026thinsp;=\u0026thinsp;0.405), whereas TAM-M2 showed a strong negative correlation with CD274 (ρ = \u0026minus;\u0026thinsp;0.763), indicating that pro-fibrotic CAFs and M2-polarized macrophages may cooperatively suppress antitumor immunity through distinct mechanisms. MSI.Expr.Sig exhibited weak correlations with all immune metrics (|ρ| \u0026le; 0.40), only showing a negative correlation with CAFs (ρ = \u0026minus;\u0026thinsp;0.395), implying that tumor stemness exerts limited direct regulatory effects on the immune microenvironment and requires further investigation in conjunction with other stemness dimensions.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eThe composite TIDE score was significantly higher in Cluster 2 (n\u0026thinsp;=\u0026thinsp;230) than in Cluster 1 (n\u0026thinsp;=\u0026thinsp;173) (Mann-Whitney U test, P\u0026thinsp;=\u0026thinsp;4.3 \u0026times; 10⁻⁷; Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eB), indicating stronger immune dysfunction and immune exclusion in Cluster 2, which may confer immunotherapy resistance.\u003c/p\u003e \u003cp\u003eAll immunosuppressive components exhibited significant differential distribution between the two subtypes (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eC). Cluster 2 displayed higher myeloid-derived suppressor cell (MDSC) infiltration (median: 0.62 vs. 0.33; P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), elevated CAF signature score (0.35 vs. 0.28, P\u0026thinsp;=\u0026thinsp;0.009), and increased TAM-M2 score (0.63 vs. 0.45, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001). Combined with the positive correlation network in Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eA, these findings suggest that Cluster 2 establishes a more robust immune barrier through enrichment of the MDSC-CAF-TAM-M2 inhibitory axis, likely driving its high TIDE score. Further comparison of key immune effector molecules (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eD) revealed significantly upregulated PD-L1 (CD274) in Cluster 2 (median: 0.61 vs. 0.46, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), consistent with its positive correlation with TIDE score, indicating an amplified immune-suppressive feedback loop in high-TIDE subtypes. CD8⁺ T-cell infiltration was significantly higher in Cluster 1 (median: 0.40 vs. 0.25, P\u0026thinsp;\u0026lt;\u0026thinsp;0.01) but exhibited a concomitant decline in IFN-γ signaling (Cluster 1 vs. Cluster 2: 0.81 vs. 0.75, P\u0026thinsp;=\u0026thinsp;0.028), suggesting that although Cluster 1 harbors more effector T cells, their functionality may be partially suppressed by the microenvironment. MSI.Expr.Sig (stemness index) showed no significant inter-subtype difference (P\u0026thinsp;=\u0026thinsp;0.17), further supporting that stemness per se is not the primary driver of immune divergence but may act indirectly via immunosuppressive pathways.\u003c/p\u003e \u003cp\u003eIn summary, our multi-dimensional immune profiling delineated two polarized BLCA subtypes:\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e- Cluster 1 (Immune-Responsive Subtype): Low TIDE score, high CD8⁺ T-cell infiltration, active IFN-γ signaling, and potential antitumor immunity.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- Cluster 2 (Immune-Desert Subtype): High TIDE score, enrichment of the MDSC-CAF-TAM-M2 inhibitory axis, PD-L1-mediated T-cell exhaustion, and predicted immunotherapy resistance. This stark contrast provides a molecular basis for subtype-specific therapeutic strategies.\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec16\" class=\"Section2\"\u003e \u003ch2\u003e3.3 Differential Gene Expression and Functional Characteristics of Tumor Subgroups\u003c/h2\u003e \u003cp\u003eTo elucidate molecular differences between the two stemness subtypes, we performed DESeq2 differential expression analysis on TCGA-BLCA samples (Cluster 2: n\u0026thinsp;=\u0026thinsp;230 vs. Cluster 1: n\u0026thinsp;=\u0026thinsp;173), applying thresholds of |log₂FC| \u0026ge; 1 and P\u0026thinsp;\u0026lt;\u0026thinsp;0.05. A total of 3,872 differentially expressed genes (DEGs) were identified, including 2,097 upregulated and 1,775 downregulated genes (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eA). Integrated GSEA, functional, and pathway enrichment analyses revealed that upregulated genes were significantly enriched in keratinization, extracellular matrix remodeling, and lipid metabolic pathways, whereas downregulated genes were associated with neuroactive ligand-receptor interactions and inhibition of cellular differentiation. Collectively, Cluster 2 exhibited a \"high-keratinization, high-metabolism, low-differentiation\" stemness-enhanced phenotype, whereas Cluster 1 retained more neuroendocrine and stromal homeostasis features, providing a gene expression foundation for subtype-specific therapeutic targeting.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003cb\u003e3.4 Construction and Cross-Subtype Validation of a 7-Gene Prognostic Model for Bladder Cancer Based on LASSO-Cox Regression\u003c/b\u003e \u003c/p\u003e \u003cp\u003eUsing the integrated TCGA-BLCA expression-clinical dataset (n\u0026thinsp;=\u0026thinsp;395 tumor samples, excluding cases with missing survival data), candidate genes were initially screened via univariate Cox regression, followed by 10-fold LASSO cross-validation (α\u0026thinsp;=\u0026thinsp;1, nλ\u0026thinsp;=\u0026thinsp;1,000) to construct the prognostic model. At λ\u0026thinsp;=\u0026thinsp;λ\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026thinsp;\u0026gt;\u0026thinsp;min\u0026lt;/sub\u0026gt;, seven genes were identified as key prognostic factors: POLE2, UTP6, and XPOT (risk coefficients\u0026thinsp;\u0026gt;\u0026thinsp;0) as well as RRAS2, PLAA, CENPH, and DIAPH3 (risk coefficients\u0026thinsp;\u0026lt;\u0026thinsp;0) (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eA\u0026ndash;B). The risk score formula was defined as:\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003erisk_score\u0026thinsp;=\u0026thinsp;exp(0.468\u0026times;POLE2\u0026thinsp;+\u0026thinsp;0.380\u0026times;UTP6\u0026ndash;0.183\u0026times;RRAS2\u0026ndash;0.439\u0026times;PLAA\u0026thinsp;+\u0026thinsp;0.237\u0026times;DIAPH3\u0026ndash;0.417\u0026times;CENPH\u0026thinsp;+\u0026thinsp;0.429\u0026times;XPOT).\u003c/p\u003e \u003cp\u003eMultivariate Cox regression confirmed the model\u0026rsquo;s significant independent prognostic value (global test P\u0026thinsp;=\u0026thinsp;1.4\u0026times;10\u0026thinsp;\u0026lt;\u0026thinsp;sup\u0026thinsp;\u0026gt;\u0026thinsp;\u0026minus;\u0026thinsp;8\u0026lt;/sup\u0026gt;, C-index\u0026thinsp;=\u0026thinsp;0.690), with POLE2 (HR\u0026thinsp;=\u0026thinsp;1.63, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), XPOT (HR\u0026thinsp;=\u0026thinsp;1.30, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), and PLAA (HR\u0026thinsp;=\u0026thinsp;0.64, P\u0026thinsp;=\u0026thinsp;0.002) exhibiting the strongest associations (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eC). A nomogram integrating risk scores and clinical variables further demonstrated robust predictive performance (C-index\u0026thinsp;=\u0026thinsp;0.69; Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eD).\u003c/p\u003e \u003cp\u003eUsing the median risk score (0.07) as the cutoff, Kaplan-Meier analysis revealed significantly shorter median overall survival in the high-risk group (720 vs. 2,880 days, log-rank P\u0026thinsp;\u0026lt;\u0026thinsp;0.001; Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eE). Subtype-stratified validation confirmed the model\u0026rsquo;s prognostic discrimination in both Cluster 1 (P\u0026thinsp;\u0026lt;\u0026thinsp;0.001) and Cluster 2 (P\u0026thinsp;\u0026lt;\u0026thinsp;0.001) (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eF), demonstrating its robustness across stemness subtypes. In summary, the 7-gene LASSO-Cox risk model provides a reliable prognostic tool for bladder cancer patients and lays the foundation for targeted intervention strategies.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003e3.5 Cross-Cohort Validation and Development of an Individualized Survival Prediction Tool\u003c/h2\u003e \u003cp\u003eThe 7-gene risk model (risk_score\u0026thinsp;=\u0026thinsp;Σβ\u003csub\u003ei\u003c/sub\u003e\u0026times;gene expression\u003csub\u003ei\u003c/sub\u003e) was validated in two independent external cohorts: GSE13507 (n\u0026thinsp;=\u0026thinsp;167) and GSE32894 (n\u0026thinsp;=\u0026thinsp;308). Survival analysis showed that in GSE13507, the high-risk group (n\u0026thinsp;=\u0026thinsp;83) had a median survival of 24 months, significantly lower than the low-risk group (n\u0026thinsp;=\u0026thinsp;84; 48 months, log-rank P\u0026thinsp;=\u0026thinsp;0.032; Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eA). Similarly, in GSE32894, the high-risk group (n\u0026thinsp;=\u0026thinsp;154) exhibited further reduced median survival (48 vs. 96 months, P\u0026thinsp;=\u0026thinsp;0.0017; Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eB), with both cohorts showing hazard ratios (HR)\u0026thinsp;\u0026gt;\u0026thinsp;1.5 and sustained survival curve separation, confirming the model\u0026rsquo;s cross-platform robustness.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eAn individualized nomogram integrating risk score, age, and pathological stage (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eC) demonstrated high concordance between predicted and observed 1\u0026ndash;5-year survival rates (mean absolute error\u0026thinsp;=\u0026thinsp;0.014; Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eD). Decision curve analysis revealed that the model significantly improved clinical net benefit over traditional TNM staging within the 0.2\u0026ndash;0.8 threshold probability range (AUC\u0026thinsp;=\u0026thinsp;0.751; Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eE). The \"Points\u0026rarr;Total Points\u0026rarr;Survival Probability\" mapping (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eC) enables rapid quantification of individual mortality risk (e.g., a 60-year-old Stage III patient with risk_score\u0026thinsp;=\u0026thinsp;0.5 has an estimated 5-year mortality probability of \u0026asymp;\u0026thinsp;60%), facilitating tailored follow-up or combination therapy. Thus, cross-cohort validation and nomogram implementation enable precise risk stratification and clinical translation of bladder cancer prognosis.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec18\" class=\"Section2\"\u003e \u003ch2\u003e3.6 Immune Microenvironment Characterization and Risk Gene\u0026ndash;Immune Cell Interaction Network in Bladder Cancer\u003c/h2\u003e \u003cp\u003eCIBERSORT-based immune infiltration analysis of the TCGA-BLCA development set (n\u0026thinsp;=\u0026thinsp;395) revealed significant immune polarization between Cluster 2 (n\u0026thinsp;=\u0026thinsp;230) and Cluster 1 (n\u0026thinsp;=\u0026thinsp;173). Cluster 2 exhibited effector cell exhaustion (CD8⁺ T cells: 1.5% vs. 2.8%, P\u0026thinsp;=\u0026thinsp;0.002; activated NK cells: 0.4% vs. 1.1%, P\u0026thinsp;=\u0026thinsp;0.007) and immunosuppressive cell enrichment (M2 macrophages: 15.6% vs. 8.3%, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001; Tregs: 3.4% vs. 2.0%, P\u0026thinsp;=\u0026thinsp;0.004), whereas Cluster 1 maintained antigen-presenting capacity (resting DCs: 1.1% vs. 0.5%, P\u0026thinsp;=\u0026thinsp;0.011), aligning with TIDE-classified \"immune-desert\" (Cluster 2) and \"immune-responsive\" (Cluster 1) phenotypes (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eA).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eSpearman correlation analysis (|ρ| \u0026gt; 0.3, P\u0026thinsp;\u0026lt;\u0026thinsp;0.05) further untangled the regulatory network between the seven risk genes and immune cells:\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e- Prognosis-promoting genes (XPOT, UTP6, POLE2) correlated positively with M2 macrophages (ρ\u0026thinsp;=\u0026thinsp;0.37\u0026ndash;0.42) and suppressed CD8⁺ T cells (ρ = \u0026minus;0.31 to \u0026minus;\u0026thinsp;0.29).\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- Prognosis-suppressing genes (RRAS2, PLAA, CENPH) enhanced activated NK cell activity (ρ\u0026thinsp;=\u0026thinsp;0.33\u0026ndash;0.36) and attenuated Treg suppression (ρ = \u0026minus;0.32 to \u0026minus;\u0026thinsp;0.28).\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- DIAPH3 promoted antigen presentation via resting DC activation (ρ\u0026thinsp;=\u0026thinsp;0.34) (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eB).\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003cp\u003eThe risk score positively correlated with M2 macrophage abundance (ρ\u0026thinsp;=\u0026thinsp;0.46, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001) and negatively with CD8⁺ T cells (ρ = \u0026minus;0.39, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001). After adjusting for immune variables, multivariate Cox regression retained its independent prognostic value (HR\u0026thinsp;=\u0026thinsp;2.27, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), indicating that risk genes drive poor outcomes by reshaping the \"M2/Treg immunosuppressive axis\" while suppressing the \"CD8⁺/NK effector axis\", concurrently capturing stemness-related regulatory signals beyond the immune microenvironment (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eC).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec19\" class=\"Section2\"\u003e \u003ch2\u003e3.7 Stemness Subtype-Specific Drug Resistance Mechanisms and Risk Gene\u0026ndash;Drug Interaction Profiles\u003c/h2\u003e \u003cp\u003eUsing the oncoPredict algorithm and GDSC2022 database, drug sensitivity to 265 compounds was predicted in the TCGA-BLCA cohort (n\u0026thinsp;=\u0026thinsp;453). Differential analysis (adj. P\u0026thinsp;\u0026lt;\u0026thinsp;0.01, |log\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;2\u0026lt;/sub\u0026thinsp;\u0026gt;\u0026thinsp;FC| \u0026gt; 1) revealed that the high-stemness Cluster 2 exhibited marked resistance to three chemotherapeutics: Aurora B inhibitor SB-743921 (log\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;2\u0026lt;/sub\u0026thinsp;\u0026gt;\u0026thinsp;FC\u0026thinsp;=\u0026thinsp;1.12, P\u0026thinsp;=\u0026thinsp;2.7\u0026times;10\u0026thinsp;\u0026lt;\u0026thinsp;sup\u0026thinsp;\u0026gt;\u0026thinsp;\u0026minus;\u0026thinsp;11\u0026lt;/sup\u0026gt;), microtubule stabilizer docetaxel (log\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;2\u0026lt;/sub\u0026thinsp;\u0026gt;\u0026thinsp;FC\u0026thinsp;=\u0026thinsp;0.98, P\u0026thinsp;=\u0026thinsp;4.8\u0026times;10\u0026thinsp;\u0026lt;\u0026thinsp;sup\u0026thinsp;\u0026gt;\u0026thinsp;\u0026minus;\u0026thinsp;14\u0026lt;/sup\u0026gt;), and paclitaxel (log\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;2\u0026lt;/sub\u0026thinsp;\u0026gt;\u0026thinsp;FC\u0026thinsp;=\u0026thinsp;0.85, P\u0026thinsp;=\u0026thinsp;3.2\u0026times;10\u0026thinsp;\u0026lt;\u0026thinsp;sup\u0026thinsp;\u0026gt;\u0026thinsp;\u0026minus;\u0026thinsp;13\u0026lt;/sup\u0026gt;), with systematic rightward shifts in IC\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;50\u0026lt;/sub\u0026gt; values (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eA\u0026ndash;B), suggesting cell cycle\u0026ndash;microtubule dynamics and Aurora kinase activation jointly mediate chemotherapy resistance.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eSpearman correlation (|ρ| \u0026gt; 0.4, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001) between risk genes (PLAA, CENPH, XPOT, RRAS2) and drug sensitivity revealed:\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e- PLAA overexpression reduced SB-743921 IC\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;50\u0026lt;/sub\u0026gt; (ρ = \u0026minus;0.33), indicating enhanced Aurora inhibitor sensitivity.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- CENPH overexpression inversely correlated with docetaxel resistance (ρ = \u0026minus;0.35), implicating cell cycle regulation in overcoming microtubule drug evasion.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- XPOT overexpression reversed paclitaxel resistance (ρ = \u0026minus;0.37), suggesting nuclear-cytoplasmic transport dysfunction exacerbates insensitivity (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eC).\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003cp\u003eThese findings demonstrate that Cluster 2 establishes cross-resistance via synergistic stemness pathway activation and risk gene dysregulation, whereas Cluster 1 may benefit from targeted interventions, providing a computational framework for dual stemness-risk gene biomarker-guided precision therapy.\u003c/p\u003e \u003c/div\u003e"},{"header":"Discussion","content":"\u003cp\u003eThis study systematically elucidates the core drivers of bladder cancer heterogeneity from a \"stemness\u0026ndash;immune\u0026ndash;prognosis\u0026ndash;drug\" integrative perspective, establishing a clinically actionable risk stratification and targeting system to address ICIs resistance and chemotherapy failure.\u003c/p\u003e \u003cp\u003eConsensus clustering of 35 stemness signature genes divided TCGA-BLCA into Cluster 1 (immune-responsive) and Cluster 2 (immune-desert). Cluster 2 displayed a \"high-keratinization\u0026ndash;high-metabolism\u0026ndash;low-differentiation\" transcriptional program with upregulated EMT, lipid metabolic reprogramming, and keratinization pathways, whereas Cluster 1 retained neuroendocrine and stromal homeostasis features. This classification was validated not only transcriptionally but also immunologically: Cluster 2 was enriched with M2 macrophages, CAFs, and Tregs but depleted of CD8⁺ T and activated NK cells, forming an immunosuppressive loop. Notably, limited correlation between stemness indices and immune metrics suggests stemness indirectly shapes the immune landscape via extracellular matrix remodeling and inhibitory cell recruitment, expanding the traditional \"stemness\u0026ndash;immune\" unidirectional model to emphasize microenvironment-mediated mechanisms.\u003c/p\u003e \u003cp\u003eThe 7-gene LASSO-Cox model (POLE2, UTP6, XPOT [prognosis-promoting]; RRAS2, PLAA, CENPH, DIAPH3 [prognosis-suppressing]) robustly discriminated survival in the TCGA training set (C-index\u0026thinsp;=\u0026thinsp;0.69) and external cohorts (HR\u0026thinsp;\u0026gt;\u0026thinsp;1.5, P\u0026thinsp;\u0026lt;\u0026thinsp;0.05), confirming cross-platform generalizability. Its independence from stemness subtypes, age, or stage implies capture of biological features beyond conventional clinical variables. Mechanistically, XPOT/POLE2 activate the \"M2/Treg axis\" to impair antitumor immunity, whereas PLAA/CENPH sustain the \"CD8⁺/NK axis\", forming a dual-axis regulatory network. Correlation with immune cells further reveals the model\u0026rsquo;s integration of stemness-cycling and microenvironmental signals, ensuring prognostic precision.\u003c/p\u003e \u003cp\u003eCluster 2\u0026rsquo;s resistance to microtubule inhibitors (paclitaxel, docetaxel) and SB-743921 (elevated IC\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;50\u0026lt;/sub\u0026gt;) implicates cell cycle\u0026ndash;microtubule dynamics and Aurora kinase activation in therapeutic escape. Key risk genes negatively correlated with drug resistance: PLAA enhanced SB-743921 sensitivity (ρ = \u0026minus;0.33), CENPH reversed docetaxel resistance (ρ = \u0026minus;0.35), and XPOT reduced paclitaxel IC\u0026thinsp;\u0026lt;\u0026thinsp;sub\u0026gt;50\u0026lt;/sub\u0026gt; (ρ = \u0026minus;0.37). These findings validate risk gene functionality and propose a dual \"stemness-risk gene\" biomarker framework:\u003c/p\u003e \u003cp\u003e \u003cul\u003e \u003cli\u003e \u003cp\u003e- Cluster 2: Aurora inhibitor\u0026thinsp;+\u0026thinsp;microtubule stabilizer combinations to overcome resistance.\u003c/p\u003e \u003c/li\u003e \u003cli\u003e \u003cp\u003e- Low-risk patients: Mono-/immunotherapy strategies.\u003c/p\u003e \u003c/li\u003e \u003c/ul\u003e \u003c/p\u003e \u003cp\u003eAdditionally, risk gene overexpression in immune-desert subtypes suggests their targeting may synergize with ICIs to disrupt immunosuppression, enabling \"stemness\u0026ndash;immune\" dual intervention.\u003c/p\u003e\n\u003ch3\u003eLimitations and Future Directions\u003c/h3\u003e\n\u003cp\u003eThe retrospective nature of public data warrants prospective trials to validate the model\u0026rsquo;s predictive power and targeting safety. Single-cell multi-omics could delineate spatiotemporal dynamics of risk gene\u0026ndash;immune cell crosstalk and optimize combination therapy timing/dosing.\u003c/p\u003e"},{"header":"Conclusion","content":"\u003cp\u003eThis study establishes the first \"stemness subtyping\u0026ndash;prognostic modeling\u0026ndash;targeted therapy\" decision framework for bladder cancer. The 7-gene risk score overcomes ICIs response prediction barriers, while Cluster 2\u0026rsquo;s drug sensitivity profile unveils novel targets to reverse chemotherapy resistance, advancing precision therapy from \"ICI monotherapy\" to \"stemness\u0026ndash;immune\u0026ndash;targeted\" tripartite strategies.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e \u003ch2\u003eEthics approval and consent to participate\u003c/h2\u003e \u003cp\u003eNot applicable.\u003c/p\u003e \u003c/p\u003e \u003cp\u003e \u003cstrong\u003eCompeting interests\u003c/strong\u003e \u003cp\u003eThe authors declare that they have no competing interests.\u003c/p\u003e \u003c/p\u003e\u003ch2\u003eFunding\u003c/h2\u003e \u003cp\u003eNot applicable.\u003c/p\u003e\u003ch2\u003eAuthor Contribution\u003c/h2\u003e\u003cp\u003eZ L: Writing\u0026ndash;original draft, Writing\u0026ndash;review \u0026amp; editing, Software, Methodology, Conceptualization, Data curation, Formal analysis. Q Z: Validation, Writing\u0026ndash;original draft, Writing\u0026ndash;review \u0026amp;\u0026nbsp;editing, Methodology, Visualization. D C: Resources, Writing\u0026ndash;review \u0026amp; editing, Software, Visualization, Formal analysis. J L: Resources, Writing\u0026ndash;review \u0026amp; editing, Software, Visualization, Formal analysis. Z X:Writing\u0026ndash;review \u0026amp; editing, Resources. D L: Writing\u0026ndash;review \u0026amp; editing, Resources. Y X :Supervision, Conceptualization, Writing\u0026ndash;review \u0026amp; editing.\u003c/p\u003e\u003ch2\u003eAcknowledgements\u003c/h2\u003e \u003cp\u003eNot applicable.\u003c/p\u003e\u003ch2\u003eData Availability\u003c/h2\u003e\u003cp\u003eAll data were sourced from public databases: TCGA-BLCA (via GDC Portal), GSE13507, and GSE32894 (via GEO). Drug sensitivity data were derived from GDSC2022 (oncoPredict R package). TIDE scores were computed via http://tide.dfci.harvard.edu. Data access complies with original database policies, ensuring transparency and reproducibility. The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eSiegel RL, Miller KD, Fuchs HE, Jemal A, Cancer statistics. 2023. CA: A Cancer Journal for Clinicians. 2023;73(1):17\u0026ndash;48.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSun L, Zhao K, Liu X, Meng X. Global, regional, and national burden of esophageal cancer using the 2019 global burden of disease study. Sci Rep. 2025;15(1):3284.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBoegemann M, Krabbe LM. Prognostic Implications of Immunohistochemical Biomarkers in Non-muscle-invasive Blad Cancer and Muscle-invasive Bladder Cancer. Mini Rev Med Chem. 2020;20(12):1133\u0026ndash;52.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFunt SA, Rosenberg JE. Systemic, perioperative management of muscle-invasive bladder cancer and future horizons. Nat Rev Clin Oncol. 2017;14(4):221\u0026ndash;34.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKim KH, Lee HW, Ha HK, Seo HK. Perioperative systemic therapy in muscle invasive bladder cancer: Current standard method, biomarkers and emerging strategies. Investig Clin Urol. 2023;64(3):202\u0026ndash;18.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFunt SA, Lattanzi M, Whiting K, et al. Neoadjuvant atezolizumab with gemcitabine and cisplatin in patients with muscle-invasive bladder cancer: a multicenter, single-arm, phase II trial. J Clin Oncol. 2022;40(12):1312\u0026ndash;22.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen M, Jiang M, Wang X, Shen L, Li J. Envafolimab - first PD-1/PD-L1 antibody to be administered by subcutaneous injection for microsatellite instability-high or deficient mismatch repair advanced solid tumors. Expert Opin Biol Ther. 2022;22(10):1227\u0026ndash;32.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWu Y, Mou H. Complete remission in a patient with advanced renal pelvis carcinoma with lung metastasis treated with durvalumab immunotherapy: a case report. Ann Palliat Med. 2020;9(4):2408\u0026ndash;13.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePowles T, Rosenberg JE, Sonpavde CP, et al. Enfortumab vedotin in previously treated advanced urothelial carcinoma. N Engl J Med. 2021;384(12):1125\u0026ndash;35.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFan B, Zheng X, Wang Y, Hu X. Predicting prognosis and clinical efficacy of immune checkpoint blockade therapy via interferon-alpha response in muscle-invasive bladder cancer. Pathol Oncol Res. 2023;29:1611117.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMariathasan S, Turley SJ, Nickles D, et al. TGF-β attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature. 2018;554(7693):544\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen Y, Yan Z, Li L, Liang Y, Wei X, Zhao Y, Cao Y, Zhang H, Tang L. Identification and validation of a 9-RBPs-related gene signature associated with prognosis and immune infiltration in bladder cancer based on bioinformatics analysis and machine learning. Transl Androl Urol. 2025;14(4):1066\u0026ndash;81.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLoh JJ, Ma S. Hallmarks of cancer stemness. Cell Stem Cell. 2024;31(5):617\u0026ndash;39.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSato R, Semba T, Saya H, Arima Y. Concise Review: Stem Cells and Epithelial-Mesenchymal Transition in Cancer: Biological Implications and Therapeutic Targets. Stem Cells. 2016;34(8):1997\u0026ndash;2007.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLin G, Lin J, Wang X, et al. Characterization of the stem cell landscape and identification of a stemness-associated prognostic signature in bladder cancer. Cancer Cell Int. 2024;24:45.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHass R, von der Ohe J, Ungefroren H. Impact of the Tumor Microenvironment on Tumor Heterogeneity and Consequences for Cancer Cell Plasticity and Stemness. Cancers (Basel). 2020;12(12):3716.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang Y, Zhang X, Li H, et al. Tumor stemness score to estimate epithelial-to-mesenchymal transition and cancer stem cells characterization and to predict prognosis and immunotherapy response in bladder urothelial carcinoma. Stem Cell Res Ther. 2023;14(1):15.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSathe A, Guerth F, Cronauer MV, Heck MM, Thalgott M, Gschwend JE, Retz M, Nawroth R. Mutant PIK3CA controls DUSP1-dependent ERK 1/2 activity to confer response to AKT target therapy. Br J Cancer. 2014;111(11):2103\u0026ndash;13.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu S, Tan Z, Shi H, et al. Development of a stemness-related prognostic index to provide therapeutic strategies for bladder cancer. NPJ Precision Oncol. 2024;8(1):14.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMiranda A, Hamilton PT, Zhang AW, et al. Cancer stemness, intratumoral heterogeneity, and immune response across cancers. Proc Natl Acad Sci U S A. 2019;116(18):9020\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLai Y, Lin P, Lin F, et al. Identification of immune microenvironment subtypes and signature genes for Alzheimer's disease diagnosis and risk prediction based on explainable machine learning. Front Immunol. 2022;13:1046410.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu S, Tan Z, Shi H, Chen J, Zhang Y, Guo C, Feng W, Xu H, Wang J, Wang H. Development of a stemness-related prognostic index to provide therapeutic strategies for bladder cancer. NPJ Precis Oncol. 2024;8(1):14. \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003e10.1038/s41698-024-00510-3\u003c/span\u003e\u003cspan address=\"10.1038/s41698-024-00510-3\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMellman I, Chen DS, Powles T, Turley SJ. The cancer-immunity cycle: Indication, genotype, and immunotype. Immunity. 2023;56(10):2188\u0026ndash;205.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu J, Li K, Zhang W, et al. Large-scale public data reuse to model immunotherapy response and resistance. Genome Med. 2020;12(1):21.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiberzon A, Birger C, Thorvaldsd\u0026oacute;ttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst. 2015;1(6):417\u0026ndash;25.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMaeser D, Gruener RF, Huang RS. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform. 2021;22(6):bbab260.\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"Bladder cancer (BLCA), Cancer stemness, Molecular subtyping, Prognostic model, Tumor immune microenvironment (TME), Targeted therapy","lastPublishedDoi":"10.21203/rs.3.rs-9059350/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-9059350/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003ch2\u003eBackground\u003c/h2\u003e \u003cp\u003eThe high heterogeneity and treatment resistance of bladder cancer (BLCA) are closely associated with cancer stemness, yet an integrated prognostic-therapeutic model incorporating stemness features and the immune microenvironment remains lacking.\u003c/p\u003e\u003ch2\u003eMethods\u003c/h2\u003e \u003cp\u003eWe integrated data from the TCGA-BLCA (n\u0026thinsp;=\u0026thinsp;453), GSE13507 (n\u0026thinsp;=\u0026thinsp;167), and GSE32894 (n\u0026thinsp;=\u0026thinsp;308) cohorts. Stemness indices were calculated via ssGSEA/PCA, and consensus clustering was employed for subtyping. A prognostic model was constructed using LASSO-Cox regression, while immune microenvironment analysis was performed via CIBERSORT/TIDE. Drug sensitivity was predicted using oncoPredict.\u003c/p\u003e\u003ch2\u003eResults\u003c/h2\u003e \u003cp\u003eBased on 35 stemness-related genes, BLCA was classified into Cluster 1 (immune-responsive, n\u0026thinsp;=\u0026thinsp;173) and Cluster 2 (immune-desert, n\u0026thinsp;=\u0026thinsp;230). The latter exhibited significantly elevated stemness indices, TIDE scores, and immunosuppressive cell infiltration (MDSC/CAF/TAM-M2, P\u0026thinsp;\u0026lt;\u0026thinsp;0.01). A 7-gene prognostic model (POLE2/UTP6/XPOT/RRAS2/PLAA/CENPH/DIAPH3) was established, demonstrating that high-risk patients had a threefold shorter median survival (720 vs. 2880 days, P\u0026thinsp;\u0026lt;\u0026thinsp;0.001), with consistent validation across cohorts (HR\u0026thinsp;\u0026gt;\u0026thinsp;1.5, P\u0026thinsp;\u0026lt;\u0026thinsp;0.05). Cluster 2 exhibited resistance to microtubule inhibitors (paclitaxel/docetaxel) and Aurora B inhibitors (log₂FC\u0026thinsp;\u0026gt;\u0026thinsp;0.85, P\u0026thinsp;\u0026lt;\u0026thinsp;10⁻\u0026sup1;\u0026sup1;), whereas low PLAA expression enhanced sensitivity to SB-743921 (ρ=-0.33). Risk genes were found to drive poor prognosis by modulating the \"M2/Treg suppression axis\" (XPOT/POLE2 positively correlated, ρ\u0026thinsp;\u0026gt;\u0026thinsp;0.37) and the \"CD8⁺/NK effector axis\" (PLAA/CENPH negatively correlated, ρ\u0026lt;-0.32).\u003c/p\u003e\u003ch2\u003eConclusion\u003c/h2\u003e \u003cp\u003eThis study establishes the first stemness-immune classification system for BLCA and a cross-subtype prognostic model, offering novel targeted therapeutic strategies (e.g.Aurora inhibitor combination therapy) for immune checkpoint inhibitor (ICI)-resistant patients.\u003c/p\u003e","manuscriptTitle":"Integrating Stemness Features and Immune Microenvironment in a Prognostic Model for Bladder Cancer Treatment","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2026-04-29 06:07:22","doi":"10.21203/rs.3.rs-9059350/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"43255113-8ac8-43bc-8bfb-bbc3af39dc9e","owner":[],"postedDate":"April 29th, 2026","published":true,"recentEditorialEvents":[{"type":"decision","content":"Rejected","date":"2026-05-02T18:20:46+00:00","index":"","fulltext":""}],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[],"tags":[],"updatedAt":"2026-05-02T18:24:41+00:00","versionOfRecord":[],"versionCreatedAt":"2026-04-29 06:07:22","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-9059350","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-9059350","identity":"rs-9059350","version":["v1"]},"buildId":"XKTyCvWXoU3ODBz1xrDgd","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}