Zhaokai Zhou1
Yajun Chen2
Zhan Wang3
Shuai Yang4
Zhengrui Li5Run Shi6Ruizhi Wang3Kui Liu7Xiaojuan Tang8
Qi Li9Ran Xu1*- 1Department of Urology, The Second Xiangya Hospital of Central South University, Changsha, China
- 2College of Traditional Chinese Medicine, Guangzhou University of Chinese Medicine, Guangzhou, China
- 3Department of Urology, The First Affiliated Hospital of Zhengzhou University, Zhengzhou, China
- 4Department of Pediatric Urology, Guangzhou Women and Children’s Medical Center, Guangzhou Medical University, Guangzhou, China
- 5Department of Oral and Maxillofacial-Head and Neck Oncology, Shanghai Ninth People’s Hospital, Shanghai Jiao Tong University School of Medicine, Shanghai, China
- 6Department of Oncology, The First Affiliated Hospital of Nanjing Medical University, Nanjing, China
- 7Department of Pediatric Surgery, The First Affiliated Hospital of Henan University of Science and Technology, Luoyang, China
- 8Department of Plastic and Reconstructive Surgery, The First Affiliated Hospital of Zhengzhou University, Zhengzhou, Henan, China
- 9Department of Pediatric Surgery, The First Affiliated Hospital of Zhengzhou University, Zhengzhou, China
Background: Bladder cancer (BLCA) continues to be a significant cause of cancer mortality in the urinary tract, with therapeutic resistance representing a major barrier to improving patient outcomes. Within the tumor microenvironment (TME), cancer-associated fibroblasts (CAFs) are pivotal drivers of BLCA progression, contributing to immune evasion and therapy resistance. This study leverages single-cell analysis to delineate CAF subclusters and explore the immune characteristics of CAFs-based BLCA classification.
Materials and methods: Signal-cell RNA sequencing (scRNA-seq) datasets were used to identify CAF subpopulations in BLCA, and bulk RNA-seq datasets were used to construct CAFs-based BLCA classification. Next, we comprehensively explored the distinct heterogeneity and characteristics for four CAFs-based BLCA subtypes. Moreover, machine learning algorithms were applied to identify novel potential targets for each subtype, and experimentally validate their effects.
Results: This study identified CAFs closely associated with BLCA development based on scRNA-seq datasets. Through further systematic clustering and functional analysis of CAFs, we successfully identified 10 distinct CAF sub-clusters, including PSCA+ Pericyte, ISG15+ Pericyte, ACTA2+ Smooth muscle cell (SMC), ACTG2+ SMC, CCL21+ inflammatory Pericyte, CD74+ apCAF, STMN1+ pCAF, CXCL14+ mCAF, APOD+ iCAF, CFD+ iCAF. The study identified four pCAFs-based BLCA distinct subtypes with different molecular, functional, and immunologic characteristics. C3 exhibited an immune-rich subtype accompanied by poor clinical prognosis, cell death pathway enrichment, higher expression of MHC molecules and co-stimulatory/co-inhibitory molecules. Conversely, C4 subtype has a smaller number of patients and an optimal prognosis, associated with lower levels of cell death pathway enrichment, lower frequency of tumor mutations, and an “immune desert” TME. C1 is mainly enriched in metabolism-related pathways, and C2 is mainly enriched in the activation of genome instability pathways, accompanied by more frequent mutations and higher Atezolizumab response. Furthermore, this study identified potential target genes or prognostic markers for each subtype.
Conclusion: Various heterogeneous CAF subgroups exist in BLCA, which is closely associated with the development of BLCA. This study identified a promising platform for understanding heterogeneity of CAFs-based BLCA subtypes, providing novel insights into the intricate molecular mechanisms of BLCA. Potential target genes for each subtype provide a basis for diagnosis and screening of BLCA patients.
1 Introduction
Bladder cancer (BLCA) represents a globally significant health burden with pronounced epidemiological impact, particularly among males (1, 2). BLCA not only affects the quality of life of patients, but also places a heavy financial burden on the healthcare system. The molecular mechanisms underlying the pathogenesis of BLCA have not been fully clarified to date. Emerging studies underscore the critical role of the tumor microenvironment (TME) in BLCA. The TME is mainly composed of cancer-associated fibroblasts (CAFs), immune cells, extracellular matrix, and numerous signaling molecules, interacting with each other during tumorigenesis to promote tumor progression and drug resistance. Of note, CAFs play the multi-faceted roles in driving BLCA progression, facilitating immune evasion, and conferring resistance to therapies (3–6). CAFs constitute a remarkably heterogeneous and functionally diverse cell population within the TME, originating from various sources and characterized by distinct molecular markers and biological activities (3, 7–9). Through single-cell RNA sequencing (scRNA-seq), Xu et al. identified a new CAF subgroup, MMP11+ myofibroblast (mCAF), which gradually accumulates during BLCA progression and is significantly associated with poor prognosis (10). Utilizing scRNA-seq and spatial transcriptomics, Zheng et al. also identify a previously unknown CAF subset in early-stage BLCA, which is associated with lymphovascular invasion, lymph node metastasis, and poor prognosis in early-stage BLCA (11). In short, CAFs play a crucial role in tumor progression, and their heterogeneity influences treatment response and prognosis.
The increasing importance of immunotherapy and targeted therapies has created an urgent need for novel molecular subtyping of BLCA to develop personalized treatment strategies (12, 13). To elucidate the molecular mechanisms of BLCA and guide individualized patient treatment, various gene expression-based classification systems have classified BLCA into molecular and clinical subtypes (13–15). These findings emphasize the critical role of molecular classifications in prognostic prediction for BLCA patients. Clinical decision-making and tailored therapeutics will benefit from these results, providing a valuable reference for further treatment. However, their specific roles, molecular diversity, and CAF-based classifications in BLCA have remained underexplored.
The bulk RNA-seq analysis intrinsically limits the resolution needed to accurately delineate and characterize these functionally distinct CAF subclusters in TME. Recent advances in high-resolution scRNA-seq and transcriptomic analysis could systematically unravel the heterogeneity of CAFs (16). Therefore, this study leverages scRNA-seq data to comprehensively define molecularly distinct CAF subpopulations within BLCA. The novel pCAF-based classification system was subsequently developed and rigorously validated using different datasets of bulk RNA-seq. Finally, machine learning algorithms were applied to identify novel potential targets for each subtype and experimentally validate their effects.
2 Materials and methods
2.1 Data sources and single-cell RNA-seq processing
The single-cell RNA-seq datasets, including GSE130001, GSE129845, GSE135337, GSE190888, and GSE222315 were downloaded from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). The R package Seurat (version 4.4.0) is used for single-cell data preprocessing. Potential doublets were predicted and removed using DoubletFinder (17). Low-quality cells (datasets with fewer than 1000 cells, cells with fewer than 500 genes, or more than 20% mitochondrial content) were removed. Normalization was performed using the LogNormalize method with a scale factor set to 1e5. For downstream analysis, the top 2000 highly variable genes were identified using the FindVariableFeatures function. The unwanted sources of variation were regressed out using the ScaleData function. Dimensionality reduction was performed using principal component analysis (PCA) with the first 30 components. The Harmony algorithm was used for batch correction before clustering analysis (18). The FindNeighbors function was used to construct a shared nearest neighbor (SNN) graph based on edge weights between cells. Subsequently, cell clusters were identified using the FindClusters function and visualized using the uniform manifold approximation and projection (UMAP) algorithm (19). In the initial round of annotation, each cluster was annotated based on known markers collected from established literature (11, 20–22), including Epithelial cells (EPCAM, CD24, KRT18, KRT8, KRT19, CLDN4), Endothelial cells (VWF, PECAM1, CDH5, ENG, CLDN5, ACKR1, PLVAP), Stromal cell (COL1A1, COL1A2, COL3A1, MYH11, ACTA2, DCN), T cells (CD2, CD3D, CD3E, TRAC, TRBC1, CD4, CD8A, CD8B, IL7R), NK cells (PRF1, KLRF1, KLRD1, FGFBP2, NKG7, XCL2), B cells (CD19, CD79A, CD79B, MS4A1), Plasma cells (TNFRSF17, MZB1, IGHG1, IGHA1), Mast cells (TPSAB1, TPSB2, MS4A2), and Myeloid cells (CD14, CD68, CD163, LYZ, S100A8, FCGR3A).
2.2 Cell type distribution
To explore the heterogeneity of the TME in each sample, this study demonstrated the composition of the broad cell types in each sample. Next, we further investigated the association between the proportion of cell types and tumor grade to explore potential clinical associations. In addition, we further compared the cell proportion profiles in non-muscle invasive bladder cancer (NMIBC) samples and normal samples. The cell proportion profile further confirmed the importance of CAFs (containing mural cells) in tumorigenesis and development. Subsequently, we clustered all the fibroblasts (containing mural cells) into sub-clusters. Likewise, the top 2000 highly variable genes were identified using the FindVariableFeatures function. The unwanted sources of variation were regressed out using the ScaleData function. Dimensionality reduction was performed using PCA with the first 30 components. The Harmony algorithm was used for batch correction before clustering analysis (18). The FindNeighbors function was used to construct an SNN graph based on edge weights between cells. Subsequently, cell clusters were identified using the FindClusters function and visualized using the UMAP algorithm (19). Based on unsupervised clustering analysis, we chose the medium clustering resolution parameter (RNA_snn_res.0.2) for cellular clustering of fibroblasts (containing mural cells). In the second round of annotation, each cluster was annotated based on known markers collected from established literature (11, 22, 23), including Mural cell (CALD1, NDUFA4L2, COX4I2, RGS5, NOTCH3, DES, MYH11), inflammatory CAF (iCAF) (CCL2, COL14A1, CXCL14, CXCL1, CXCL12, CXCL2, CXCL8, CCL13, IL6, IL11, IL8, CCL5, CCL22, CLEC3B, COL14A1, LY6C, CCL17, CXCL1, CXCL2, CXCL12, CXCL14, LIF, HGF, APOD, IGF1, C3, C7, ITM2A, MGP, CCL11), mCAF (AOC3, COL12A1, COL15A1, COL1A2, COL1A1, COL8A1, FAP, MYL9, TAGLN, MYLK, TPM1, TPM2, POSTN, HOPX, PDGFA, COL5A1, VIM, COL6A1, αSMA, TGFβ, TAGN, CTHRC1,THY1, CTGF, CTA2, FAP, COL12A1, FBLN1, SERPINF1, VCAN, COL11A1, COL10A1, ACTC2), vascular CAF (vCAF) (VEGFA, FGF, PDGF, PDGFα, PDFGRβ, ANG2, CLIC3, ACTA2, ANG1, MMP1, MMP2, MMP3, MMP9, MMP11, MMP14, PDPN, IMF1), antigen-presenting CAF (apCAF) (CD74, HLA-DRA, HLA-DRA1, HLA-DPA1, HLA-DQA1, CD52, IGFBP3, SLP1, SAA3, FSP1), proliferating CAF (pCAF) (NUSAP1, TOP2A, CENPF, PTTG1, SCG2, STMN1, PLAU, TUBA1B, IFI27, IGFBP2). To identify the subpopulations most relevant to tumor grading, we again explored the proportional profile of the CAF sub-clusters.
2.3 Enrichment analysis and deconvolution algorithm
Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) were used to analyze the biological process and pathway enrichment of each CAF sub-clusters through the “clusterProfiler” R package. Gene set variation analysis (GSVA) was performed on hallmark pathways obtained in the Molecular Signatures Database (MSigDB, version 7.0).
Bulk RNA-seq data and clinical information about BLCA were sourced from the TCGA database (https://xenabrowser.net/datapages/). The deconvolution algorithm serves as a tool to infer TME composition (different cell populations) and gene expression in bulk RNA data based on scRNA-seq data. In this study, we analyzed the relationship between CAF sub-clusters and survival by deconvolution algorithms such as Multi-subject single cell deconvolution (MuSiC), Convolutional pose machine (CPM), Dampened weighted least squares (DWLS), Bayesprism, and InstaPrism. Survival analysis was performed for each CAF sub-cluster using the optimal cutoff value.
2.4 Construction and verification of STMN1+ pCAF-based subtypes
This study performed unsupervised consensus clustering (k-means algorithm, 1,000 iterations) using 1,215 marker genes in the STMN1+ pCAF (adjusted p-value < 0.05, mean log2 fold-change > 1) via the “FindMarkers” function. Subsequently, we applied the Nearest Template Prediction (NTP) algorithm to classify TCGA-BLCA tumors into four subtypes. NTP computes subtype probabilities by projecting bulk transcriptomes onto predefined gene signatures through a cosine similarity metric, inherently mitigating technical batch effects. The reproducibility of this subtyping was rigorously validated in four independent cohorts (GSE13507, GSE48075, E−MTAB−4321, and IMvigor210 cohort), and a Meta-data cohort (including GSE13507, GSE19423, and GSE37815). IMvigor210 cohort that was additionally downloaded as external datasets using the “IMvigor210CoreBiologies” R package and the GEO cohorts were used to confirm the stability of STMN1+ pCAF-derived subtypes.
2.5 Association of four subtypes with clinicopathological characteristics
Four distinct molecular subtypes (C1, C2, C3, C4) were delineated through consensus clustering. To characterize their clinical relevance, we evaluated each subtype association with clinicopathological parameters: age, sex, TNM stage, histological grade, and growth pattern.
2.6 Cell death-related modalities and functional enrichment analysis for four subtypes
To investigate the cell death patterns of each subtype, we downloaded 56 pathways of seven cell death types including apoptosis (n = 12), autophagy (n = 19), necrosis (n = 10), lysosome-dependent cell death (n = 11), pyroptosis (n = 2), necroptosis (n = 1), and ferroptosis (n = 1) from MSigDB (https://www.gsea-msigdb.org/gsea/msigdb/human/search.jsp) and then assessed their enrichment in the four subtypes (24). What’s more, we used the online tool Metascape to construct a network of statistically significant enrichment profiles in the four subtypes to explore the functional differences and the associated pathways. The input genes were screened specifically for each subgroup and have been pre-screened by the COSG algorithm. Meanwhile, KEGG pathway enrichment analyses were performed for the four pCAF-based BLCA subtypes.
2.7 Immune microenvironment and genomic alteration analysis
This study quantified the relative abundance of 22 immune cell types across samples from the TCGA-BLCA cohort using CIBERSORT analysis, an algorithm widely used to assess immune cell infiltration. Moreover, unique scores for gene sets in specific samples were calculated by single-sample gene set enrichment analysis (ssGSEA), a method commonly used to assess immune cell infiltration abundance. Based on the R language package “GSVA”, we used ssGSEA to obtain the infiltration status of 28 immune cells within each sample from the TCGA-BLCA data, which enabled the characterization of immune cell properties among different subtypes. Meanwhile, we analyzed the expression levels of 27 immune checkpoints in the BLCA cohort and compared the differential expression among the four groups. To assess the predictive power of cancer immunotherapy risk scores, we included inflammatory genes associated with immune checkpoint blockade (ICB) response. Differences in mutation patterns among the four groups were examined using the “maftools” software package. Additionally, tumor neoantigen load, an indicator of genomic instability, was evaluated.
2.10 Machine learning-driven identification of hub genes in four subtypes
To further identify pivotal genes associated with four subtypes, we used two machine learning algorithms, namely Random Forest (RF) and Least Absolute Shrinkage and Selection Operator (LASSO) regression analyses. The RF models were analyzed using the ‘randomForest’ R package (version 4.7-1.1), and genes with MeanDecreaseGini values greater than 0.6 were selected as candidate genes. LASSO regression using the “glmnet” R package (version 4.1-8). The genes cross-identified by two algorithms were considered as potential hub genes. Thereafter, the hub genes were analyzed for single-gene GSEA enrichment using the “GSEABase” R package (version 1.60.0). The relationship between the three hub genes and 22 immune cells was analyzed using the “corrplot” R package (version 0.92).
2.11 Survival analysis of hub genes
The prognostic role of hub genes was evaluated in the TCGA-BLCA cohort and patients were categorized into high and low expression groups using the median expression threshold for each gene to plot the Kaplan-Meier survival curve. Furthermore, differences in stage and grade between the two groups were also compared. Staging and grading differences between the two groups were also compared.
2.12 Protein extraction and western blot analysis
Protein extraction and Western blot analysis refer to the steps of our previous studies (14). Among them, primary antibodies included: FAT4 (1:1000; PA5-72970, Invitrogen, China), RPL37P1 (1:1000; ab228542, Abcam, China), FGFR1 (1:1000; R381166, Zenbio, China), RNASEH1 (1:1000; 82771-1-RR, Proteintech, China), AHNAK2 (1:1000; 680021, Zenbio, China), SLC9B2 (1:1000; 24065-1-AP, Proteintech, China), MN1 (1:1000; 24697-1-AP, Proteintech, China), TTLL3 (1:1000; PA5-70598, Invitrogen, China), FABP6 (1:1000; 126828, Zenbio, China), TBC1D3 (1:1000; DF3346, Affinity, China).
2.14 Statistical analysis
All data processing was done using R 4.2.2 software. The “survival” R package was utilized to assess disease-specific survival (DSS), progression-free survival (PFS), and overall survival (OS) among subtypes using the Kaplan-Meier method and log-rank test. Specifically, statistical analyses were performed as follows: 1) Continuous variables: Kruskal-Wallis test for multi-group comparisons; Wilcoxon rank-sum or Student’s t-test for pairwise comparisons; 2) Categorical variables: Fisher’s exact or Chi-square test; 3) Pathway enrichment: Mann-Whitney-Wilcoxon Gene Set Test. Statistical significance was defined as two-sided P< 0.05.
3 Results
3.1 Single-cell atlas reveals CAFs heterogeneity in BLCA
In this study, T-SNE and UMAP dimensionality reduction resolved seven major cell clusters within BLCA. The cell types in the clusters were identified by characterized genes and cell type markers in the established literature (11, 20–22), including epithelial cells, fibroblasts, endothelial cells, myeloid cells, TNK cells, plasma cells, and mural cells. Figure 1A presents the expression of significant markers in all cell types. Analysis of cell type distribution revealed pronounced inter-patient heterogeneity (Figure 1B), corroborating significant heterogeneity of TME among individuals. Of note, the proportion of the seven cell subgroups varied among staging (Figure 1C), implicating that changes in cell type may drive tumor progression. Since there were only two cases of muscle-invasive bladder cancer (MIBC) samples, in order to reduce the error caused by too small amount of samples, we compared the cellular composition between NMIBC and paracancerous tissues. Numerous studies have pointed out that CAFs significantly impact tumor progression (20). Likewise, in our study, the proportion of CAFs was significantly higher in tumor tissues compared with paracancerous tissues (Figure 1D), further demonstrating that CAFs were intimately associated with tumor progression. The percentage of endothelial cells was also higher in the tumor tissues compared to the paracancerous tissues. Interestingly, T/NK and B cells were higher in the paracancerous tissues, which may suggest that the tumor margin region is the primary place for adaptive immune responses.
Figure 1. Cellular atlas and identification of CAFs in BLCA. (A) The tSNE (left) and UMAP (right) graphic in BLCA and paracancer tissues, color-coded by cell types. Bubble diagram for expression of marker genes in each cell type (bottom). (B) The distribution plots of cell types in every samples. (C) The distribution plots of cell types in different tumor stages (left) and tissues (right). (D) The relative proportions of 7 cell types in NMIBC and paracancer tissues. (E) The CAFs’ tSNE (left) and UMAP (right) graphic, color-coded by cell types. Bubble diagram for expression of marker genes in each cell type (bottom). (F) The tSNE (left) and UMAP (right) graphic in BLCA and paracancer tissues, color-coded by CAFs sub-cluster. The relative proportions of 5 CAFs in NMIBC and paracancer tissues. (G) The survival curve based on pCAF in TCGA BLCA cohort (Bayesprism). (H) Pathway enrichment analysis based on pCAF marker genes (HAllMARK). Adj, adjacent; N, normal; RT, recurrent tumor; T, tumor.
Unsupervised clustering of CAFs identified five general clusters: comprising apCAF, iCAF, mCAF, pCAF, and Mural cells (Figure 1E). Specifically, this study annotated each of the sub-clusters based on their characteristic genes and functional enrichment, including PSCA+ Pericyte, ISG15+ Pericyte, ACTA2+ Smooth muscle cell (SMC), ACTG2+ SMC, CCL21+ inflammatory Pericyte (iPericyte), CD74+ apCAF, STMN1+ pCAF, CXCL14+ mCAF, APOD+ iCAF, CFD+ iCAF (Figure 1F). Heatmap showed each cluster-specific gene expression pattern. Among sub-clusters, STMN1+ pCAFs were significantly enriched in NMIBC samples compared to paracancerous tissues (Figure 1F and Supplementary Figure S1A). Deconvolution analysis showed that patients with high STMN1+ pCAFs abundance exhibited significantly reduced overall survival (Figure 1G and Supplementary Figure S1B), suggesting a strong association between STMN1+ pCAFs and aggressive progression of BLCA. GSEA analysis revealed the enrichment of Hallmark_e2f_targets, Hallmark_g2m_checkpoint, Cytochrome-c oxidase activity, Cytochrome complex, Glutathione metabolism, Oxidative phosphorylation, P53 signaling pathway (Figure 1H and Supplementary Figures S1C, D). These pathways collectively explain pCAF-driven proliferative characteristics, which may potentially contribute to BLCA progression.
3.2 Identification and validation of STMN1+ pCAFs-based BLCA heterogeneity subtypes
Firstly, we obtained 1215 STMN1+ pCAF-specific expressed genes with adj_pval < 0.05 and log2FC > 1 (Supplementary Table S1), further screened 85 STMN1+ pCAF-specific prognosis-related genes (Supplementary Table S2), and finally obtained the expression of specific prognosis-related genes in TCGA-BLCA for clustering analysis. Consensus clustering was performed to define tumor subtypes, with sample similarity quantified by proximity metrics. The optimal cluster number (k=4) was determined by convergence of the delta area plot and consensus cumulative distribution function (CDF) curve (Figure 2A), resulting in the identification of four pCAF-based BLCA subtypes in the TCGA cohort. A total of 98 patients were assigned to C1 (31.2%), 88 to C2 (28.0%), 91 to C3 (29.0%), and 37 to C4 (11.8%). Figure 2B showed that the survival analysis was statistically significant (n=314, OS, P = 0.007; DSS, P = 0.003; PFS, P = 0.012). C3 subtype demonstrated poor survival outcomes. Conversely, C4 was characterized by superior survival outcomes and a lower recurrence rate. Next, we extracted the characteristic genes for each subtype, namely 90 genes for C1, 86 for C2, 54 for C3, and 51 for C4 (Supplementary Table S3). Subsequently, we verified the classification capabilities of these characteristic genes in multiple external cohorts using NTP. Similarly, in the IMvigor210 cohort, four subtypes can be clearly distinguished. C3 displayed the worst prognosis and C4 showed superior prognosis (n=213, P = 0.032, Figure 2C), in agreement with the aforementioned results, validating the accuracy and consistency of the classification. C3 with the worst prognosis had the highest proportion of non-responders to atezolizumab in the IMvigor210 cohort. Interestingly, C2 exhibited the most favorable response to atezolizumab in the IMvigor210 cohort, suggesting that C2 may be more sensitive to chemotherapy (Figure 2D). Likewise, in the meta-cohort (n=106, P = 0.043, Figure 2E), GSE13507 (n=81, P = 0.019, Supplementary Figure S2A), GSE48075 (n=39, P = 0.016, Supplementary Figure S2B), and E−MTAB−4321 (n=219, P < 0.01, Supplementary Figure S2C), the survival analysis showed a comparable trend with the training cohort, confirming the rationality of the subtyping. Collectively, the robustness and generalization of pCAF-based BLCA classification were validated by large-scale RNA-seq data.
Figure 2. The identification and verification of pCAFs-based subtypes in bulk cohorts and clinical characteristics for four subtypes. (A) Identifying four heterogeneity subtypes through consensus clustering. (B) OS, DSS, and PFS among four subtypes in TCGA-BLCA cohort. (C) OS for four subtypes and NTP results in IMvigor210-cohort. (D) The distribution plots of immunotherapy sensitivity among four subtypes in IMvigor210-cohort. (E) OS for four subtypes and NTP results in GEO Meta cohort. (F) The distinct clinical characteristics among four BLCA subtypes. (G) The distribution plots of age and tumor growth pattern among four subtypes in TCGA cohort. (H) Sankey diagram of pCAF-based BLCA classifications and previously described consensus molecular classifications, comprising NMIBC_class, MIBC_class, CIT.subtype, Lund.subtype, MDA.subtype, TCGA.subtype, and UNC.subtype. nsP > 0.05, *P < 0.05, **P < 0.01, ***P < 0.001.
3.3 Comparisons of clinicopathological characteristics across four pCAF subtypes
Clinical characteristics from the TCGA-BLCA cohort were analyzed to assess the clinical significance of the four subtypes. Age, Ts, and growth pattern varied significantly across the four subtypes (P < 0.05, Figures 2F, G and Supplementary Figure S2D). The C3 subtype with the worst prognosis had the highest proportion of patients greater than 65 years of age and was associated with a higher proportion of non-papillary growth patterns (Figure 2F). Moreover, the proportion of C3 subtype at stage T1/T2 was the lowest (Figure 2G). These findings may be related to the worst prognosis of C3 subtype. To evaluate the correlation between our classification and other BLCA consensus classifiers, the pCAF-based BLCA classification of 314 samples was compared with various classifiers (Figure 2H). As expected, CAFs-based BLCA subtypes have limited overlap with other classifications, indicating the novelty of our subtypes and their complementary role to previous BLCA classifications.
3.4 Identification of cell death-related modalities among four BLCA subtypes
As is well documented, cell death significantly affects tumor cell proliferation, invasion, metastasis, and the efficacy of chemotherapy, thereby influencing the prognosis of BLCA (25–27). The results showed that most of cell death-related pathways were significantly up-regulated in the C3 subtype, while most of these pathways were only weakly enriched in C4 (Figure 3A), suggesting that these cell death-related pathways are closely associated with prognosis of BLCA. The expression levels of C1 and C2 in the cell death-related pathways were located between C3 and C4. The hypothesis can be drawn that the enrichment of stronger cell death pathways may promote tumor malignancy (28).
Figure 3. Distinctive cell death-related modalities and biological features among four BLCA subtypes. (A) Heatmap demonstrated differences in cell death signaling pathways among the four subtypes. Enrichment map network and KEGG enrichment analysis in C1 subtype (B), C2 subtype (C), C3 subtype (D), C4 subtype (E).
3.5 Distinctive biological features of four subtypes
Functional enrichment analyses revealed that C1 is mainly enriched in molecular metabolism-related pathways such as gaba_synthesis_release_reuptake_and_degradation, alanine_aspartate_and_glutamate_metabolism, and butanoate_metabolism (Figure 3B). C2 are mainly enriched in cell cycle and DNA repair pathways, exhibiting concerted activation of genome instability pathways, including mitotic_spindle_checkpoint, inhibition_of_replication_initiation_of_damaged_DNA_by_rb1_e2f1, mitotic_prometaphase, trna_processing_in_the_nucleus, and resolution_of_sister_chromatid_cohesion (Figure 3C). C3 subtype identified unique enrichment in keratinization and immune-modulatory pathways such as formation_of_the_cornified_envelope and antiviral_mechanism_by_ifn_stimulated_genes (Figure 3D). C4 subtype are mainly enriched in substance transport and metabolism pathways, including transport_of_organic_anions, eicosanoid_ligand_binding_receptors, recycling_of_bile_acids_and_salts, and triglyceride_catabolism (Figure 3E). Overall, the C1-C4 subtypes exhibit different bio-functional states.
3.6 Immune landscape for the four pCAF-based BLCA heterogeneity subtypes
Considering that the tumor immune microenvironment plays a key role in BLCA progression, we comprehensively explored the infiltration levels of immune cell subsets and the expression of various immune-related markers. Cibersort and ssGSEA was performed to assess the infiltration abundance of different immune cells across four subtypes (Figure 4A, Supplementary Figures S2E, F). The degree of infiltration of immune cells, including natural killer cells, CD8+ cells, CD4+ cells, neutrophil, regulatory T cells, and myeloid-derived suppressor cells (MDSCs), is the highest in the C3 subtype. The C4 subtype had the lowest degree of immune infiltration among the four subtypes. Immuno-radargrams further demonstrated that C3 is an immune-infiltrating tumor and C4 is an “immune-desert” tumor (Figure 4B). We further compared the immunogenicity-associated indicators of the four subtypes to further explore the characteristics of the immune microenvironment (Figure 4C). C2 subtype had higher immunogenicity and higher genomic instability, such as single-nucleotide variant (SNV) neoantigens, indel neoantigens, loss of heterozygosity (LOH) segments, LOH fraction alterations, and homologous recombination defects, accompanied by a higher T-cell receptor (TCR) richness (TCR richness and Shannon). The C1 and C4 subtypes have lower levels of these indicators, which may reflect poor survivor outcomes. Nevertheless, C4 had higher SNV neoantigens and Indel neoantigens, suggesting that C4, as an immune-desert tumor, has higher neoantigens at the same time, but fails to elicit an effective immune response.
Figure 4. The immune landscape across distinct BLCA subtypes. (A) Heatmap of immune cell infiltration by SsGSEA in four subtypes. (B) Radar plot of immune gene expression levels. (C) Comparison of immune indicators in the four subtypes. (D) The relative expression of co-stimulatory and co-inhibitory immune check points and specific MHC molecules. (E) Heatmap of immune checkpoints of four BLCA subtypes. (F) The relative expression of immune checkpoints of the four subtypes. (G) The tumor inflammation signature (TIS) score, antigen presentation score (APS), and neoantigen in BLCA subtypes. (H) Mutation landscape in TCGA BLCA cohort. (I) Expression heatmap of the top 20 highly mutated genes in the four subtypes. (J) Tumor mutational burden, single nucleotide polymorphism, and insertion-deletion in BLCA subtypes. ns P > 0.05, *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001.
Next, we explore the characteristics of the immune molecules of the four subtypes. C2 and C3 showed increased activation of co-stimulatory immune checkpoint inhibitors (ICIs) such as CD226, CD27, CD28, and TNFRSF14, and co-inhibitory receptors such as CD160, CD244, CD48, CD274, CD276, LAG3, LGALS9, and CTLA4, MHC molecules like HLA-A/B/C and HLA-DQA1/2 (Figure 4D). However, expressions of immune checkpoints were scarce in the C1 and C4 subtypes, suggesting that ICI therapy might not be effective (Figures 4E, F). Furthermore, antigen presentation efficiency was further assessed using the antigen processing and presenting machinery score (APS) (Figure 4G). It has been demonstrated in numerous studies that tumor inflammation signatures (TIS), APS, and neoantigen score are all predictive biomarkers for tumor immunotherapy, where higher values generally correlate with better immunotherapy outcomes (24). Similarly, higher TIS, APS score, score, and neoantigen score further demonstrate that the C2 and C3 subtypes could achieve greater clinical efficacy with immunotherapy (Figure 4G).
Moreover, this study further explored the gene mutations in the four subtypes (Figures 4H–J). Figure 4I shows that C2 and C3 have significantly higher TP53 and RB1 mutation frequencies than C1 and C4. TP53 mutations could cause DNA damage repair defects and RB1 mutations could cause chromosome segregation errors, thereby promoting the development of multiple tumors. It was demonstrated that RB1 and TP53 co-mutations are strongly associated with genomic biomarkers of ICI response in urothelial BLCA (29). The higher mutation frequencies tend to portend a poor prognosis. However, studies have shown that high tumor mutational load (TMB) is associated with longer survival after ICI therapy (30). C2 and C3 have higher single-nucleotide polymorphism (SNP), insertion-deletion (INDEL), and TMB compared to C1 and C4, which further demonstrates the suitability of C2 and C3 for immunotherapy.
3.7 Identification and validation of hub genes for C1 subtypes
To explore the molecular mechanisms and prognosis-related markers of the four subtypes, machine-learning algorithms were used to identify hub genes. For C1 subtype, the randomForest algorithm analysis screened two potential hub genes using a node size threshold of 3. LASSO regression analysis finally screened five target genes. By taking the intersection of the genes identified by the two algorithms, we identified two key genes, FAT4 and RPL37P1, that are closely associated with the malignant progression of BLCA (Figure 5A). The survival curves illustrated that the high-expression FAT4 had a poor prognosis and the high-expression RPL37P1 had a favorable prognosis (Figure 5B). Moreover, significant differences were observed in the expression levels of hub genes across different stages, with higher expression levels of FAT4 in patients at higher stages and higher expression levels of RPL37P1 in patients at lower stages (Figure 5B). Likewise, high-grade BLCA patients had higher gene expression levels of FAT4 and lower gene expression levels of RPL37P1, further corroborating that FAT4 gene was positively correlated with the degree of BLCA malignancy and RPL37P1 gene was protective (Figure 5B). Next, their functions and correlations with immune cells were examined. FAT4 is predominantly enriched in pathways associated with driving malignant phenotypic transformation and microenvironmental remodeling such as HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION, HALLMARK_ANGIOGENESIS, HALLMARK_TGF_BETA_SIGNALING, and HALLMARK_KRAS_SIGNALING_UP, while RPL37P1 is predominantly enriched in pathways related to metabolic reprogramming and genomic destabilization, such as HALLMARK_DNA_REPAIR, HALLMARK_P53_PATHWAY, and HALLMARK_GLYCOLYSIS (Figure 5C). Moreover, FAT4 was significantly positively correlated with mast cells resting, T cells CD4 memory resting, Macrophages M2, and monocytes. RPL37P1 was significantly positively correlated with T cells CD4 naive, T cells gamma delta, and B cells memory, closely related to immune response (Figure 5D). The results of Western blot analysis showed that higher stages were correlated with higher gene expression levels of FAT4. However, the relationship between RPL37P1 expression level and staging was not significant (Figure 5E).
Figure 5. Machine learning algorithms to obtain hub targets of C1 subtype. (A) The identification of potential targets using the randomForest (left) and Lasso (right) algorithm. Venn diagram for identifying hub genes of C1 subtype. (B) Survival analysis in TCGA BLCA cohort, respectively, for the grouping of the median expression of FAT4 and RPL37P1. The relative expression of different stages (left) and grades (right) for FAT4 and RPL37P1, respectively. (C) GSEA enrichment analysis for FAT4 and RPL37P1 based on hallmark gene set. (D) Correlation between 22 infiltrating immune cells and FAT4 and RPL37P1, respectively. (E) WB and the relative expression of FAT4 and RPL37P1 proteins in different tumor stages. *P < 0.05, **P < 0.01, ***P < 0.001.
3.8 Identification and validation of hub genes for C2 subtypes
For C2 subtype, randomForest algorithm analysis screened three potential hub genes and LASSO regression analysis screened 14 target genes, finally identifying three key genes, namely FGFR1, RNASEH1, DDX1 (Figure 6A). Among them, the two hub genes FGFR1 and RNASEH1, are associated with survival in BLCA patients. The KM curves illustrated that the high-expression FAT4 and RNASEH1 had poor prognosis (Figure 6B) and patients with advanced stages and grades have higher levels of gene expression of FGFR1 and RNASEH1 (Figure 6C), indicating that these genes were malignant and protumorigenic. Likewise, FGFR1 was also predominantly enriched in pathways associated with malignant phenotypic transformation and microenvironmental remodeling, while RNASEH1 was predominantly enriched in pathways related to genome stability (Figure 6C). Moreover, FAT4 and RNASEH1 were significantly positively correlated with macrophages (Figure 6D). The results of western blot analysis showed that higher stages were correlated with higher gene expression levels (Figure 6E).
Figure 6. Machine learning algorithms to obtain hub targets of C2 subtype. (A) The identification of potential targets using the randomForest and Lasso algorithm. Venn diagram for identifying hub genes of C2 subtype. (B) Survival analysis in TCGA BLCA cohort, respectively, for the grouping of the median expression of FGFR1 and RNASEH1. The relative expression of different stages and grades for FGFR1 and RNASEH1. (C) GSEA enrichment analysis for FGFR1 and RNASEH1 based on hallmark gene set. (D) Correlation between 22 infiltrating immune cells and FGFR1 and RNASEH1, respectively. (E) WB and the relative expression of FGFR1 and RNASEH1 proteins in different tumor stages. *P < 0.05, **P < 0.01, ***P < 0.001.
3.9 Identification and validation of hub genes for C3 subtypes
Similarly, three key genes (AHNAK2, SLC9B2, and MN1) in C3 subtype were identified (Figure 7A). The high-expression of these genes is strongly associated with a poor prognosis (Figure 7B). What’s more, these genes are expressed at higher levels in advanced and highly graded patients (Figure 7B). Interestingly, They are all predominantly enriched in HALLMARK_EPITHELIAL_MESENCHYMAL_TRANSITION and are closely related to Macrophages (Figures 7C, D). Western blot analysis demonstrated a progressive increase in the expression of these genes with staging (Figure 7E).
Figure 7. Machine learning algorithms to obtain hub targets of C3 subtype. (A) The identification of potential targets using the randomForest and Lasso algorithm in C3 subtype. (B) Survival analysis for AHANK2, SLC9B2, and MN1. The relative expression of different stages and grades for AHANK2, SLC9B2, and MN1, respectively. (C) GSEA enrichment analysis for AHANK2, SLC9B2, and MN1 based on hallmark gene set. (D) Correlation between 22 infiltrating immune cells and AHANK2, SLC9B2, and MN1, respectively. (E) WB and the relative expression of AHANK2, SLC9B2, and MN1 proteins in different tumor stages. ns P > 0.05, *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001.
3.10 Identification and validation of hub genes for C4 subtypes
In C4 subtype, this study illustrated that the high expression of four key genes (TTLL3, FABP6, TBC1D3, and CYP4F35P) are associated with the favorable prognosis in the TCGA cohort (Figures 8A, B). These genes are more lowly expressed in patients with high staging and advanced grade (Figure 8B). Enrichment analysis indicated that Each gene has a different enrichment function and are all closely related to T cells CD4 naive and T cells gamma delta (Figures 8C, D). To further validate the aforementioned results, the expression levels of three hub genes (No CYP4F35P antibody) were quantified across the four different stages of BLCA tissues (Figure 8E). Taken together, we identified the corresponding hub genes for each subtype. The specific role of these genes remains to be validated in further clinical practice.
Figure 8. Machine learning algorithms to obtain hub targets of C4 subtype. (A) The identification of potential targets. (B) Survival analysis in TCGA-BLCA cohort for the TTLL3, FABP6, TBC1D3, and CYP4F35P. The relative expression of different stages and grades for TTLL3, FABP6, TBC1D3, and CYP4F35P. (C) GSEA enrichment analysis for TTLL3, FABP6, TBC1D3, and CYP4F35P. (D) Correlation between 22 infiltrating immune cells and TTLL3, FABP6, TBC1D3, and CYP4F35P, respectively. (E) WB and the relative expression of TTLL3, FABP6, TBC1D3, and CYP4F35P proteins in different tumor stages. ns P > 0.05, *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001.
4 Discussion
CAFs have been shown to have various cellular communications with tumor cells and immune cells, which not only secrete cytokines and mediate extracellular matrix remodeling to indirectly affect tumor cells and immune cells, but also promote tumor cell proliferation, invasion, and metastasis through direct cell-to-cell interactions (3, 7, 8, 31), further promoting immune escape, chemotherapy resistance, and tumor recurrence. This study identified CAFs closely associated with BLCA development based on scRNA-seq datasets. Through further systematic clustering and functional analysis of CAFs, we successfully identified 10 distinct CAF sub-clusters. The novel pCAF-based BLCA classification system was subsequently developed and rigorously validated using different datasets of tumor bulk RNA-seq data with excellent robustness and generalizability.
The four CAF-based BLCA subtypes have different molecular, functional and immunological characteristics. The C3 subtype exhibits a malignant phenotype with poor clinical prognosis, with a higher enrichment of the cell death pathway, a significant immune infiltration, and higher expression of MHC molecules and co-stimulatory/co-inhibitory molecules, suggesting that C3 is an immune-enriched phenotype and may be sensitive to immunotherapy. The higher TIS, APS, neoantigen score, and TMB further implied that C3 was potentially responsive to immunotherapy. Moreover, some studies indicate that cell death signals play an important role in the immune environment, and targeting cell death signaling pathways may enhance the efficacy of tumor immunotherapy (28, 32). Therefore, it is possible that cell death inducers and immunotherapy used together will achieve better results for C3. In contrast, the C4 subtype has a smaller population and a better prognosis, as characterized by a lower enrichment of cell death pathways, a lower frequency of tumor mutations, and an ‘immune desert’ TME. C4 had lower expression levels of immune molecules, APS, and TMB, further suggesting resistance to immunotherapy. C1 subtype was mainly enriched in metabolism-related pathways, whereas C2 subtype was mainly enriched in genome instability pathway activation, with a higher frequency of mutations and better response to atezolizumab treatment. The expression levels of cell death-related pathways and the degree of immune infiltration of C1 and C2 were intermediate between C3 and C4. Collectively, CAF-based BLCA classification effectively distinguished independent BLCA groups with low overlap with existing transcriptional classification methods, demonstrating unprecedented predictive value. Besides, this study identified two potential key genes (FAT4 and RPL37P1) as significantly associated with C1, two potential target genes (FGFR1 and RNASEH1) for C2, three potential target genes (AHNAK2, SLC9B2, and MN1) for C3, and four potential target genes (TTLL3, FABP6, TBC1D3 and CYP4F35P) for C4. These hub genes play critical roles in BLCA patients and provide potential molecular targets and clinical strategies for designing targeted therapeutic regimens for BLCA patients.
FAT4 (FAT atypical cadherin 4) is a member of the cadherin family and is involved in cell adhesion, epithelial-mesenchymal transition (EMT), and the regulation of the TME (33–36). It has been verified that FAT4 could affect tumor metastasis and regulate TMB and microsatellite instability to influence the immunotherapy response through calcium signaling pathways and chemokine signaling pathways (34). Additionally, the interaction between FAT4 and CAFs could reshape the TME and promote BLCA progression (34, 35). Likewise, in our study, FAT4 shows that high expression is associated with poor prognosis. RPL37P1 (Ribosomal Protein L37 pseudogene 1) is a pseudogene related to ribosome biogenesis (37). Its parent gene, RPL37, plays a crucial role as a ribosomal protein in protein translation and participates in regulating the tumor suppression (38). RPL37 can bind and inhibit the E3 ubiquitin ligase activity of MDM2, thereby stabilizing p53 and inducing cell cycle arrest or apoptosis (38), which is closely related to the prognosis of patients in various cancers. Moreover, DNA damage (such as exposure to ultraviolet rays or treatment with cisplatin) can trigger the proteasomal degradation of RPL37, and then activate p53 pathway through a L11-dependent manner, indicating that RPL37 and its pseudogene may be involved in the DNA damage stress response (39). In summary, FAT4 and RPL37P1 play important roles in reprogramming the TME and have the potential to become a cancer prognostic marker and therapeutic target for C1 subtype.
FGFR1 (fibroblast growth factor receptor 1) is a transmembrane receptor with tyrosine kinase activity, belonging to the FGFR family (FGFR1-4) (40). FGFR1 regulates cell proliferation, differentiation, and migration by activating downstream signaling pathways (such as MAPK/PLCγ and PI3K/AKT) (41, 42). Overexpression or abnormal splicing of FGFR1 can promote EMT, enhance invasion and metastasis ability, and drive cell proliferation through MAPK/ERK or COX-2 pathways, which is closely related to tumor malignancy progression (43–45). It has been found that FGFR1 overexpression is an independent factor for poor prognosis in patients with MIBC, and is associated with shorter recurrence time and overall survival (43). Moreover, the abnormal activation of FGFR1 can be further strengthened by collaborating with molecules like TTYH3, thereby reinforcing downstream oncogenic signals (46). Although pan-FGFR inhibitors like erdafitinib are effective for BLCA with FGFR3 mutations/fusions, their efficacy is limited for FGFR1-dependent tumors (44, 47). RNASEH1 (Ribonuclease H1) is a nucleolytic endonuclease that specifically degrades RNA-DNA hybrids (R-loops), which plays a crucial role in maintaining genomic stability, regulating DNA replication and repair (48–50). RNASEH1 could prevent replication stress and DNA damage by removing R-loops, and dysfunction of RNASEH1 is closely related to the progression of various cancers, including gastric cancer, B-cell lymphoma, and melanoma (51–54). Some cancer cells may rely on the telomerase-independent ‘Alternative Lengthening of Telomeres’ pathway to maintain telomere length (49). RNASEH1 regulates the balance of hybrid structures between telomere DNA and long non-coding RNA TERRA, affecting homologous recombination efficiency. Abnormal expression of RNASEH1 may lead to telomere dysfunction or excessive recombination, promoting the immortalization and invasiveness of cancer cells. A pan-cancer analysis revealed that overexpression of RNASEH1 is associated with the regulation of the TME and poor prognosis, suggesting that RNASEH1 may promote BLCA immune escape by inhibiting anti-tumor immunity (48). Our study also showed that high expressions of FGFR1 and RNASEH1 related to poor prognosis. In conclusion, FGFR1 and RNASEH1 may affect BLCA progression through multiple mechanisms and may become prognostic markers or therapeutic targets for BLCA.
AHNAK2 is a large molecular nuclear protein with a molecular weight exceeding 600 kDa, belonging to the AHNAK protein family, regulating tumor progression by activating signaling pathways such as ERK, MAPK, Wnt, and MEK, as well as promoting EMT (55). AHNAK2 is abnormally highly expressed in numerous cancers and is closely associated with poor prognosis (55, 56). Several researchers have revealed that high expression of AHNAK2 is significantly correlated with the malignancy of BLCA (56–58). For patients who underwent radical cystectomy and have high AHNAK2 expression in tumor tissues, their disease-free survival and cancer-specific survival were significantly shortened, and multivariate analysis confirmed that AHNAK2 was an independent poor prognostic factor (56). From the molecular mechanism perspective, the expression of AHNAK2 in the tumor hypoxic microenvironment is regulated by HIF-1α, which may drive tumor progression by inducing EMT and enhancing tumor stem cell characteristics (59). It is noteworthy that the level of AHNAK2 protein in the urine of BLCA patients is significantly elevated, especially in MIBC, making it a potential non-invasive diagnostic marker (58). SLC9B2 (also known as NHA2) is a unique Na/H antiporter, belonging to the SLC9B subfamily, with 14 transmembrane segments and forming a unique lipid-sensitive dimer structure (60, 61). Nevertheless, there are relatively few direct studies on SLC9B2 in BLCA. In terms of physiological functions, SLC9B2 is highly expressed in the distal convoluted tubules of the kidney, participating in blood pressure homeostasis and electrolyte balance by regulating the WNK4-NCC signaling pathway (62). Interestingly, SLC9B2 exhibits abnormal expression in various pathological conditions. In polycystic kidney disease, SLC9B2 is significantly upregulated and positively correlated with cyst size, regulated by the polycystin-1/Ca2+/NFAT signaling axis, which can be induced by vasopressin and methylxanthine drugs (63). What’s more, SLC9B2 is related to the pathogenesis of diabetes and hypertension, and is crucial for sperm motility and fertility (64, 65). The MN1 (Meningioma 1) gene is an important transcriptional co-regulatory factor that is initially discovered in meningiomas and has been confirmed to be involved in the occurrence and development of various malignant tumors in recent years (66, 67). MN1 enhances mRNA stability through m6A methylation modification mediated by METTL14, promoting tumor progression and chemotherapy resistance in osteosarcoma (68). In astrocytoma and other nervous system tumors, MN1 gene rearrangements (such as MN1-BEND2 and MN1-CXXC5) drive tumor occurrence via activating the PDGFRα signaling pathway (69). However, the high expression of MN1 in low-grade gliomas actually predicts a better prognosis, and this tissue-specific difference suggests that the regulatory network of MN1 is highly complex (66). Of note, MN1 exerts a pro-tumor effect through the XIST/miR-15a-5p/MN1/FZD2 signaling axis, and its significantly high expression is observed in female patients with poor prognosis, which may be an important molecular basis for the gender difference in BLCA (70). The expression regulation of MN1 may involve epigenetic regulation of long non-coding RNAs, and the activation of its downstream effector molecule FZD2 may affect the remodeling of TME (70). In C3, high expression of AHNAK2, SLC9B2, and MN1 is correlated with poor survival and may be potential therapeutic targets. These findings provide a novel perspective for understanding the molecular mechanism underlying the poor prognosis of BLCA.
TTLL3 is a tubulin glycine ligase that regulates the assembly of cilia, which plays a crucial role in maintaining the structure and function of cilia (71). Studies have shown that TTLL3 is expressed in colon tissues, and its absence or decreased levels of its expression lead to the deficiency of tubulin glycosylation and a reduction in primary cilia, thereby promoting the development of colorectal cancer (72). FABP6 (Fatty Acid Binding Protein 6) is an intracellular lipid transporter protein, which precisely regulates metabolic pathways, signal transduction, and gene expression (73). It is abnormally expressed in various cancers and is closely related to the occurrence and progression of cancer. In colorectal cancer, FABP6 is negatively correlated with immune infiltration, and its downregulation can enhance tumor immunogenicity and promote the recruitment of CD8+ T cells (74). FABP6 is regulated by the transcription factor REST and participates in gastric cancer progression by influencing autophagy and the Akt/mTOR pathways (75). Knockdown of FABP6 can significantly inhibit the proliferation and motility of BLCA cells, while downregulating the expression of cell cycle proteins such as CDK2 and CDK4 and blocking the AKT-mTOR signaling pathway, enhancing the therapeutic effect of cisplatin (76). Interestingly, some researchers indicate that FABP6 has been identified as a protective gene in the prognostic model of BLCA, and its expression level is related to the survival period of patients (77). TBC1D3, a hominoid-specific oncogene belonging to the TBC1 domain protein family, is encoded by a cluster of paralogues located on chromosome 17q12 (78). TBC1D3 enhances cell migration ability in breast cancer by activating the TNFα/NF-κB signaling pathway and upregulating the expression of OLR1 (79). The oncogenic mechanism of TBC1D3 involves the regulation of multiple key signaling pathways, including maintaining the stability of EGFR through interaction with the microtubule network, delaying EGFR degradation, and enhancing Ras activity to promote cell proliferation, as well as affecting gene expression through epigenetic regulation (80, 81). However, there are relatively few direct studies on TBC1D3 in BLCA. CYP4F35P is an lncRNA derived from a pseudogene and belongs to the cytochrome P450 family (82, 83). In recent years, it has been found to have regulatory effects in various cancers. In laryngeal squamous cell carcinoma and tongue squamous cell carcinoma, CYP4F35P shows differential expression patterns and forms a co-expression network with genes such as MUC21 and CEACAM1, suggesting that it may affect tumor progression by regulating cell adhesion and signal transduction processes (83, 84). In BLCA, CYP4F35P has been identified as one of the seven key prognostic-related lncRNAs, and its expression characteristics are significantly correlated with the survival rate of patients (82). The risk scoring model based on it can effectively predict the prognosis of BLCA patients (82). These findings suggest that TTLL3, FABP6, TBC1D3, and CYP4F35P may become potential biomarkers for the diagnosis and prognosis assessment of BLCA. Future research needs to further verify the specific functional mechanism of these genes in BLCA and their potential value as a diagnostic marker and therapeutic target.
5 Limitations and conclusions
Several limitations need to be noted. Firstly, CAF clustering and CAF-based BLCA classification were constructed from retrospective data from public databases. And the inconsistent number of patient stage distribution and the lack of prior treatment history in these public databases further limited our study. Thus, future validation in more prospective, multi-center BLCA cohorts is needed. Secondly, we only initially explored the key genes for each subtype, and subsequent studies will need to delve deeper into their potential mechanisms of action in the development of BLCA. In conclusion, this study identified a promising platform for understanding CAF heterogeneity and BLCA classification, which could provide novel insights into the complex molecular mechanisms of BLCA. Each heterogeneous subtype possesses different unique molecular, functional and immunological characteristics, implying a different therapeutic strategy, further facilitating personalized medicine for BLCA patients.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.
Ethics statement
The studies involving humans were approved by The Ethics Committee of Zhengzhou University’s First Affiliated Hospital (2024-KY-1687). The studies were conducted in accordance with the local legislation and institutional requirements. The human samples used in this study were acquired from primarily isolated as part of your previous study for which ethical approval was obtained. Written informed consent for participation was not required from the participants or the participants’ legal guardians/next of kin in accordance with the national legislation and institutional requirements.
Author contributions
ZW: Conceptualization, Data curation, Investigation, Methodology, Software, Validation, Visualization, Writing – original draft. YC: Investigation, Methodology, Project administration, Software, Supervision, Writing – review & editing. ZZ: Formal Analysis, Investigation, Methodology, Project administration, Resources, Software, Supervision, Visualization, Writing – review & editing. SY: Investigation, Methodology, Project administration, Software, Supervision, Writing – review & editing. ZL: Investigation, Project administration, Resources, Software, Validation, Visualization, Writing – review & editing. RS: Investigation, Methodology, Project administration, Software, Validation, Writing – review & editing. RW: Investigation, Methodology, Project administration, Software, Visualization, Writing – review & editing. KL: Data curation, Formal Analysis, Methodology, Project administration, Resources, Writing – review & editing. XT: Investigation, Methodology, Project administration, Resources, Software, Visualization, Writing – review & editing. QL: Conceptualization, Investigation, Methodology, Project administration, Resources, Supervision, Validation, Visualization, Writing – review & editing. RX: Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Validation, Visualization, Writing – review & editing.
Funding
The author(s) declare that no financial support was received for the research and/or publication of this article.
Acknowledgments
Bladder and experiment cartoon images used in figures were obtained from https://smart.servier.com/.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that no Generative AI was used in the creation of this manuscript.
Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2025.1580986/full#supplementary-material
Supplementary Figure 1 | Identification of CAFs in BLCA.
Supplementary Figure 2 | The verification of CAFs-based BLCA subtypes in bulk cohorts and the immune landscape across distinct BLCA subtypes.
Supplementary Table 1 | STMN1+ pCAF-specific expressed genes.
Supplementary Table 2 | STMN1+ pCAF-specific prognosis-related genes.
Supplementary Table 3 | Marker_genes of four subtype.
Abbreviations
BLCA, Bladder cancer; TME, tumor microenvironment; CAFs, cancer-associated fibroblasts; scRNA-seq, single-cell RNA sequencing; mCAF, myofibroblast; GEO, Gene Expression Omnibus; PCA, principal component analysis; SNN, shared nearest neighbor; UMAP, uniform manifold approximation and projection; NMIBC, non-muscle invasive bladder cancer; iCAF, inflammatory CAF; vCAF, vascular CAF; apCAF, antigen-presenting CAF; pCAF, proliferating CAF; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; GSVA, Gene set variation analysis; MuSiC, Multi-subject single cell deconvolution; CPM, Convolutional pose machine; DWLS, Dampened weighted least squares; NTP, Nearest Template Prediction; ssGSEA, single-sample gene set enrichment analysis; ICB, immune checkpoint blockade; RF, Random Forest; LASSO, Least Absolute Shrinkage and Selection Operator; DSS, disease-specific survival; PFS, progression-free survival; OS, overall survival; MIBC, muscle invasive bladder cancer; SMC, Smooth muscle cell; CDF, cumulative distribution function; MDSC, myeloid-derived suppressor cell; SNV, single-nucleotide variant; LOH, loss of heterozygosity; TCR, T-cell receptor; ICI, immune checkpoint inhibitor; APS, antigen processing and presenting machinery score; TMB, tumor mutational load; SNP, single nucleotide polymorphism; INDEL, insertion-deletion; EMT, epithelial-mesenchymal transition.
References
1. Huang X, Du G, Yang Y, Su P, Chen S, Cai C, et al. Advancing bladder cancer management: development of a prognostic model and personalized therapy. Front Immunol. (2024) 15:1430792. doi: 10.3389/fimmu.2024.1430792
2. Su Y, Wang B, Huang J, Huang M, and Lin T. YTHDC1 positively regulates PTEN expression and plays a critical role in cisplatin resistance of bladder cancer. Cell Prolif. (2023) 56:e13404. doi: 10.1111/cpr.13404
3. Caramelo B, Zagorac S, Corral S, Marqués M, and Real FX. Cancer-associated fibroblasts in bladder cancer: origin, biology, and therapeutic opportunities. Eur Urol Oncol. (2023) 6:366–75. doi: 10.1016/j.euo.2023.02.011
4. Zheng Y, Wang X, Yang X, and Xing N. Single-cell RNA sequencing reveals the cellular and molecular characteristics of high-grade and metastatic bladder cancer. Cell Oncol (Dordr). (2023) 46:1415–27. doi: 10.1007/s13402-023-00820-x
5. Qin Y, Zu X, Li Y, Han Y, Tan J, Cai C, et al. A cancer-associated fibroblast subtypes-based signature enables the evaluation of immunotherapy response and prognosis in bladder cancer. iScience. (2023) 26:107722. doi: 10.1016/j.isci.2023.107722
6. Shi ZD, Liu Y, Tao ZQ, Chao L, Zhang ZG, Sun F, et al. Apigenin inhibits recurrent bladder cancer progression by targeting VEGF-β. Cancer Lett. (2025) 620:217676. doi: 10.1016/j.canlet.2025.217676
7. Ma Z, Li X, Mao Y, Wei C, Huang Z, Li G, et al. Interferon-dependent SLC14A1(+) cancer-associated fibroblasts promote cancer stemness via WNT5A in bladder cancer. Cancer Cell. (2022) 40:1550–1565.e7. doi: 10.1016/j.ccell.2022.11.005
8. Burley A, Rullan A, and Wilkins A. A review of the biology and therapeutic implications of cancer-associated fibroblasts (CAFs) in muscle-invasive bladder cancer. Front Oncol. (2022) 12:1000888. doi: 10.3389/fonc.2022.1000888
9. Deng S, Chen YZ, and Du J. Decoding the plasticity of cancer-associated fibroblasts: mechanistic insights and precision targeting strategies in gastric cancer progression and therapeutic resistance. Cell Prolif. (2025) p:e70094. doi: 10.1111/cpr.70094
10. Xu W, Liang T, Fang H, Fu L, Deng D, Tan X, et al. Single-cell RNA sequencing identifies MMP11(+) cancer-associated fibroblasts as drivers of angiogenesis and bladder cancer progression. Adv Sci (Weinh). (2025) p:e02774. doi: 10.1002/advs.202502774
11. Zheng H, An M, Luo Y, Diao X, Zhong W, Pang M, et al. PDGFRα(+)ITGA11(+) fibroblasts foster early-stage cancer lymphovascular invasion and lymphatic metastasis via ITGA11-SELE interplay. Cancer Cell. (2024) 42:682–700.e12. doi: 10.1016/j.ccell.2024.02.002
12. Sun L, Liu XP, Yan X, Wu S, Tang X, Chen C, et al. Identification of molecular subtypes based on liquid-liquid phase separation and cross-talk with immunological phenotype in bladder cancer. Front Immunol. (2022) 13:1059568. doi: 10.3389/fimmu.2022.1059568
13. Zhang Y, Chen Y, and Wen W. Four types of adenine-related RNA modification writers -mediated molecular subtypes contribute to predicting clinical outcomes and treatment options in bladder cancer. Front Immunol. (2023) 14:1152806. doi: 10.3389/fimmu.2023.1152806
14. Wang Z, Zhou Z, Yang S, Li Z, Shi R, Wang R, et al. Pathway-based cancer transcriptome deciphers a high-resolution intrinsic heterogeneity within bladder cancer classification. J Transl Med. (2025) 23:666. doi: 10.1186/s12967-025-06682-1
15. Dong Y, Xu C, Su G, Li Y, Yan B, Liu Y, et al. Clinical value of anoikis-related genes and molecular subtypes identification in bladder urothelial carcinoma and in vitro validation. Front Immunol. (2023) 14:1122570. doi: 10.3389/fimmu.2023.1122570
16. Zhang X, Zhu R, Yu D, Wang J, Yan Y, and Xu K. Single-cell RNA sequencing to explore cancer-associated fibroblasts heterogeneity: “Single” vision for “heterogeneous” environment. Cell Prolif. (2024) 57:e13592. doi: 10.1111/cpr.13592
17. McGinnis CS, Murrow LM, and Gartner ZJ. DoubletFinder: doublet detection in single-cell RNA sequencing data using artificial nearest neighbors. Cell Syst. (2019) 8:329–337.e4. doi: 10.1016/j.cels.2019.03.003
18. Korsunsky I, Millard N, Fan J, Slowikowski K, Zhang F, Wei K, et al. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat Methods. (2019) 16:1289–96. doi: 10.1038/s41592-019-0619-0
19. Becht E, McInnes L, Healy J, Dutertre CA, Kwok IWH, Ng LG, et al. Dimensionality reduction for visualizing single-cell data using UMAP. Nat Biotechnol. (2018). doi: 10.1038/nbt.4314
20. Chen Z, Zhou L, Liu L, Hou Y, Xiong M, Yang Y, et al. Single-cell RNA sequencing highlights the role of inflammatory cancer-associated fibroblasts in bladder urothelial carcinoma. Nat Commun. (2020) 11:5077. doi: 10.1038/s41467-020-18916-5
21. Lai H, Cheng X, Liu Q, Luo W, Liu M, Zhang M, et al. Single-cell RNA sequencing reveals the epithelial cell heterogeneity and invasive subpopulation in human bladder cancer. Int J Cancer. (2021) 149:2099–115. doi: 10.1002/ijc.33794
22. Ma C, Yang C, Peng A, Sun T, Ji X, Mi J, et al. Pan-cancer spatially resolved single-cell analysis reveals the crosstalk between cancer-associated fibroblasts and tumor microenvironment. Mol Cancer. (2023) 22:170. doi: 10.1186/s12943-023-01876-x
23. Gao Y, Li J, Cheng W, Diao T, Liu H, Bo Y, et al. Cross-tissue human fibroblast atlas reveals myofibroblast subtypes with distinct roles in immune modulation. Cancer Cell. (2024) 42:1764–1783.e10. doi: 10.1016/j.ccell.2024.08.020
24. Zhou Z, Zhang Y, Li J, Weng S, Li J, Chen S, et al. Crosstalk between regulated cell death and immunity in redox dyshomeostasis for pancreatic cancer. Cell Signal. (2023) 109:110774. doi: 10.1016/j.cellsig.2023.110774
25. Zhang Q, Tan Y, Zhang J, Shi Y, Qi J, Zou D, et al. Pyroptosis-related signature predicts prognosis and immunotherapy efficacy in muscle-invasive bladder cancer. Front Immunol. (2022) 13:782982. doi: 10.3389/fimmu.2022.782982
26. Liu K, Chen H, Li Y, Wang B, Li Q, Zhang L, et al. Autophagy flux in bladder cancer: Cell death crosstalk, drug and nanotherapeutics. Cancer Lett. (2024) 591:216867. doi: 10.1016/j.canlet.2024.216867
27. Xu W, Tang HJ, Anwaier A, Liu W, Tian X, Su J, et al. Immunogenomic characteristics of cell-death-associated genes with prognostic implications in bladder cancer. Front Immunol. (2022) 13:909324. doi: 10.3389/fimmu.2022.909324
28. Tong X, Tang R, Xiao M, Xu J, Wang W, Zhang B, et al. Targeting cell death pathways for cancer therapy: recent developments in necroptosis, pyroptosis, ferroptosis, and cuproptosis research. J Hematol Oncol. (2022) 15:174. doi: 10.1186/s13045-022-01392-3
29. Manzano RG, Catalan-Latorre A, and Brugarolas A. RB1 and TP53 co-mutations correlate strongly with genomic biomarkers of response to immunity checkpoint inhibitors in urothelial bladder cancer. BMC Cancer. (2021) 21:432. doi: 10.1186/s12885-021-08078-y
30. Samstein RM, Lee CH, Shoushtari AN, Hellmann MD, Shen R, Janjigian YY, et al. Tumor mutational load predicts survival after immunotherapy across multiple cancer types. Nat Genet. (2019) 51:202–6. doi: 10.1038/s41588-018-0312-8
31. Li S, Zhou X, Feng H, Huang K, Chen M, Lin M, et al. Deciphering the immunomodulatory function of GSN(+) inflammatory cancer-associated fibroblasts in renal cell carcinoma immunotherapy: insights from pan-cancer single-cell landscape and spatial transcriptomics analysis. Cell Prolif. (2025) p:e70062. doi: 10.1111/cpr.70062
32. Yu L, Huang K, Liao Y, Wang L, Sethi G, and Ma Z. Targeting novel regulated cell death: Ferroptosis, pyroptosis and necroptosis in anti-PD-1/PD-L1 cancer immunotherapy. Cell Prolif. (2024) 57:e13644. doi: 10.1111/cpr.13644
33. Wang D, Wu S, He J, Sun L, Zhu H, Zhang Y, et al. FAT4 overexpression promotes antitumor immunity by regulating the β-catenin/STT3/PD-L1 axis in cervical cancer. J Exp Clin Cancer Res. (2023) 42:222. doi: 10.1186/s13046-023-02758-2
34. Mao W, Zhou J, Hu J, Zhao K, Fu Z, Wang J, et al. A pan-cancer analysis of FAT atypical cadherin 4 (FAT4) in human tumors. Front Public Health. (2022) 10:969070. doi: 10.3389/fpubh.2022.969070
35. Yang Y, Li Y, Yang Q, Liu Z, Chang X, Yang H, et al. FAT4 activation inhibits epithelial-mesenchymal transition (EMT) by promoting autophagy in H2228/Cer cells. Med Oncol. (2022) 40:64. doi: 10.1007/s12032-022-01934-2
36. Huang FY, Wong DK, Mak LY, Cheung TT, Zhang SS, Chau HT, et al. FAT4 loss initiates hepatocarcinogenesis through the switching of canonical to noncanonical WNT signaling pathways. Hepatol Commun. (2023) 7(12):e0338. doi: 10.1097/HC9.0000000000000338
37. Daftuar L, Zhu Y, and Prives C. Ribosomal protein L37 and the p53 network. Cell Cycle. (2010) 9:4259. doi: 10.4161/cc.9.21.13527
38. Daftuar L, Zhu Y, Jacq X, and Prives C. Ribosomal proteins RPL37, RPS15 and RPS20 regulate the Mdm2-p53-MdmX network. PloS One. (2013) 8:e68667. doi: 10.1371/journal.pone.0068667
39. Llanos S and Serrano M. Depletion of ribosomal protein L37 occurs in response to DNA damage and activates p53 through the L11/MDM2 pathway. Cell Cycle. (2010) 9:4005–12. doi: 10.4161/cc.9.19.13299
40. Katoh M and Nakagama H. FGF receptors: cancer biology and therapeutics. Med Res Rev. (2014) 34:280–300. doi: 10.1002/med.21288
41. Tomlinson DC, Baxter EW, Loadman PM, Hull MA, and Knowles MA. FGFR1-induced epithelial to mesenchymal transition through MAPK/PLCγ/COX-2-mediated mechanisms. PloS One. (2012) 7:e38972. doi: 10.1371/journal.pone.0038972
42. Carneiro BA, Meeks JJ, Kuzel TM, Scaranti M, Abdulkadir SA, and Giles FJ. Emerging therapeutic targets in bladder cancer. Cancer Treat Rev. (2015) 41:170–8. doi: 10.1016/j.ctrv.2014.11.003
43. Lim S, Koh MJ, Jeong HJ, Cho NH, Choi YD, Cho do Y, et al. Fibroblast growth factor receptor 1 overexpression is associated with poor survival in patients with resected muscle invasive urothelial carcinoma. Yonsei Med J. (2016) 57:831–9. doi: 10.3349/ymj.2016.57.4.831
44. Cheng T, Roth B, Choi W, Black PC, Dinney C, and McConkey DJ. Fibroblast growth factor receptors-1 and -3 play distinct roles in the regulation of bladder cancer growth and metastasis: implications for therapeutic targeting. PloS One. (2013) 8:e57284. doi: 10.1371/journal.pone.0057284
45. Tomlinson DC and Knowles MA. Altered splicing of FGFR1 is associated with high tumor grade and stage and leads to increased sensitivity to FGF1 in bladder cancer. Am J Pathol. (2010) 177:2379–86. doi: 10.2353/ajpath.2010.100354
46. Biswas PK, Kwak Y, Kim A, Seok J, Kwak HJ, Lee M, et al. TTYH3 modulates bladder cancer proliferation and metastasis via FGFR1/H-ras/A-raf/MEK/ERK pathway. Int J Mol Sci. (2022) 23(18):10496. doi: 10.3390/ijms231810496
47. Katoh M, Loriot Y, Brandi G, Tavolari S, Wainberg ZA, and Katoh M. FGFR-targeted therapeutics: clinical activity, mechanisms of resistance and new directions. Nat Rev Clin Oncol. (2024) 21:312–29. doi: 10.1038/s41571-024-00869-z
48. Yi C, Yang J, Zhang T, Xie S, Li W, Qin L, et al. A pan-cancer analysis of RNASEH1, a potential regulator of the tumor microenvironment. Clin Transl Oncol. (2023) 25:2569–86. doi: 10.1007/s12094-023-03142-4
49. Arora R, Lee Y, Wischnewski H, Brun CM, Schwarz T, and Azzalin CM. RNaseH1 regulates TERRA-telomeric DNA hybrids and telomere maintenance in ALT tumour cells. Nat Commun. (2014) 5:5220. doi: 10.1038/ncomms6220
50. Sun J, Li Y, Tian H, Chen H, Li J, and Li Z. Comprehensive analysis identifies long non-coding RNA RNASEH1-AS1 as a potential prognostic biomarker and oncogenic target in hepatocellular carcinoma. Am J Cancer Res. (2024) 14:996–1014. doi: 10.62347/JPHF4071
51. Silva B, Pentz R, Figueira AM, Arora R, Lee YW, Hodson C, et al. FANCM limits ALT activity by restricting telomeric replication stress induced by deregulated BLM and R-loops. Nat Commun. (2019) 10:2253. doi: 10.1038/s41467-019-10179-z
52. Nie X, Zhao C, Zhang Y, Huang W, Zhou Y, Liu H, et al. Nucleolar protein treacle ribosome biogenesis factor 1 maintains gastric cancer cell proliferation by regulating R-loop associated DNA replication stress. J Gastroenterol Hepatol. (2023) 38:1170–80. doi: 10.1111/jgh.16181
53. Zong C, Zhang Z, Gao L, He J, Wang Y, Li Q, et al. APOBEC3B coordinates R-loop to promote replication stress and sensitize cancer cells to ATR/Chk1 inhibitors. Cell Death Dis. (2023) 14:348. doi: 10.1038/s41419-023-05867-0
54. Shen YJ, Ho SSW, Tan NY, Koo CX, Khatoo M, Cheung FSG, et al. Genome-derived cytosolic DNA contributes to type I interferon expression and immunogenicity of B-cell lymphoma cells. Cytokine. (2015) 76:581–2. doi: 10.1016/j.cyto.2015.05.024
55. Zhang S, Cai Z, and Li H. AHNAKs roles in physiology and Malignant tumors. Front Oncol. (2023) 13:1258951. doi: 10.3389/fonc.2023.1258951
56. Koguchi D, Matsumoto K, Shimizu Y, Kobayashi M, Hirano S, Ikeda M, et al. Prognostic impact of AHNAK2 expression in patients treated with radical cystectomy. Cancers (Basel). (2021) 13(8):1748. doi: 10.3390/cancers13081748
57. Yan Y, Cao X, Wang Z, Huang Z, Cai J, Tang P, et al. Pyroptosis-related patterns predict tumor immune landscape and immunotherapy response in bladder cancer. Front Mol Biosci. (2022) 9:815290. doi: 10.3389/fmolb.2022.815290
58. Komina S, Petrusevska G, Jovanovic R, Kunovska SK, Stavridis S, Dohcev S, et al. AHNAK2 urinary protein expression as potential biomarker for bladder cancer detection: A pilot study. Turk J Urol. (2022) 48:423–30. doi: 10.5152/tud.2022.22132
59. Wang M, Li X, Zhang J, Yang Q, Chen W, Jin W, et al. AHNAK2 is a novel prognostic marker and oncogenic protein for clear cell renal cell carcinoma. Theranostics. (2017) 7:1100–13. doi: 10.7150/thno.18198
60. Ho TM, Berger S, Müller P, Simonin C, Reymond JL, Von Ballmoos C, et al. Physiological and molecular function of the sodium/hydrogen exchanger NHA2 (SLC9B2). Chimia (Aarau). (2022) 76:1019–24. doi: 10.2533/chimia.2022.1019
61. Anderegg MA, Gyimesi G, Ho TM, Hediger MA, and Fuster DG. The less well-known little brothers: the SLC9B/NHA sodium proton exchanger subfamily-structure, function, regulation and potential drug-target approaches. Front Physiol. (2022) 13:898508. doi: 10.3389/fphys.2022.898508
62. Anderegg MA, Albano G, Hanke D, Deisl C, Uehlinger DE, Brandt S, et al. The sodium/proton exchanger NHA2 regulates blood pressure through a WNK4-NCC dependent pathway in the kidney. Kidney Int. (2021) 99:350–63. doi: 10.1016/j.kint.2020.08.023
63. Prasad H, Dang DK, Kondapalli KC, Natarajan N, Cebotaru V, and Rao R. NHA2 promotes cyst development in an in vitro model of polycystic kidney disease. J Physiol. (2019) 597:499–519. doi: 10.1113/JP276796
64. Matsuoka R, Fudim R, Jung S, Zhang C, Bazzone A, Chatzikyriakidou Y, et al. Structure, mechanism and lipid-mediated remodeling of the mammalian Na(+)/H(+) exchanger NHA2. Nat Struct Mol Biol. (2022) 29:108–20. doi: 10.1038/s41594-022-00738-2
65. Gardner CC and James PF. Na(+)/H(+) exchangers (NHEs) in mammalian sperm: essential contributors to male fertility. Int J Mol Sci. (2023) 24(19):14981. doi: 10.3390/ijms241914981
66. Saini M, Jha AN, Tangri R, Qudratullah M, and Ali S. MN1 overexpression with varying tumor grade is a promising predictor of survival of glioma patients. Hum Mol Genet. (2021) 29:3532–45. doi: 10.1093/hmg/ddaa231
67. Chen W, Soon YY, Pratiseyo PD, Sutanto R, Hendriansyah L, Kuick CH, et al. Central nervous system neuroepithelial tumors with MN1-alteration: an individual patient data meta-analysis of 73 cases. Brain Tumor Pathol. (2020) 37:145–53. doi: 10.1007/s10014-020-00372-0
68. Li HB, Huang G, Tu J, Lv DM, Jin QL, Chen JK, et al. METTL14-mediated epitranscriptome modification of MN1 mRNA promote tumorigenicity and all-trans-retinoic acid resistance in osteosarcoma. EBioMedicine. (2022) 82:104142. doi: 10.1016/j.ebiom.2022.104142
69. Shi Y, Sun Q, Jia F, Xie X, Zhou X, Guo R, et al. Oncogenic fusions converge on shared mechanisms in initiating astroblastoma. Nature. (2025) 643:551–61. doi: 10.1038/s41586-025-08981-5
70. Cai F, Xu S, Li Y, He Q, Su Q, Chen H, et al. The role of the LncRNA XIST/miR-15a-5p/MN1 signaling axis in gender disparities in bladder cancer prognosis. Front Immunol. (2025) 16:1554829. doi: 10.3389/fimmu.2025.1554829
71. Wloga D, Webster DM, Rogowski K, Bré MH, Levilliers N, Jerka-Dziadosz M, et al. TTLL3 Is a tubulin glycine ligase that regulates the assembly of cilia. Dev Cell. (2009) 16:867–76. doi: 10.1016/j.devcel.2009.04.008
72. Rocha C, Papon L, Cacheux W, Marques Sousa P, Lascano V, Tort O, et al. Tubulin glycylases are required for primary cilia, control of cell proliferation and tumor development in colon. EMBO J. (2014) 33:2247–60. doi: 10.15252/embj.201488466
73. Wu L, Ou GL, Zhang W, Ma HX, Li XY, Zhen YH, et al. Fatty acid-binding proteins in cancers. Int J Surg. (2025). doi: 10.1097/JS9.0000000000003049
74. Lian W, Wang Z, Ma Y, Tong Y, Zhang X, Jin H, et al. FABP6 expression correlates with immune infiltration and immunogenicity in colorectal cancer cells. J Immunol Res. (2022) 2022:3129765. doi: 10.1155/2022/3129765
75. Luo J, Yu H, Yuan Z, Tang W, and Wang C. REST promotes autophagy in gastric cancer by transcriptionally activating FABP6 to inhibit the akt/mTOR signaling pathway. Front Biosci (Landmark Ed). (2024) 29:212. doi: 10.31083/j.fbl2906212
76. Lin CH, Chang HH, Lai CR, Wang HH, Tsai WC, Tsai YL, et al. Fatty acid binding protein 6 inhibition decreases cell cycle progression, migration and autophagy in bladder cancers. Int J Mol Sci. (2022) 23(4):2154. doi: 10.3390/ijms23042154
77. Lin J, Yang J, Xu X, Wang Y, Yu M, and Zhu Y. A robust 11-genes prognostic model can predict overall survival in bladder cancer patients based on five cohorts. Cancer Cell Int. (2020) 20:402. doi: 10.1186/s12935-020-01491-6
78. Hodzic D, Kong C, Wainszelbaum MJ, Charron AJ, Su X, and Stahl PD. TBC1D3, a hominoid oncoprotein, is encoded by a cluster of paralogues located on chromosome 17q12. Genomics. (2006) 88:731–6. doi: 10.1016/j.ygeno.2006.05.009
79. Wang B, Zhao H, Zhao L, Zhang Y, Wan Q, Shen Y, et al. Up-regulation of OLR1 expression by TBC1D3 through activation of TNFα/NF-κB pathway promotes the migration of human breast cancer cells. Cancer Lett. (2017) 408:60–70. doi: 10.1016/j.canlet.2017.08.021
80. He Z, Tian T, Guo D, Wu H, Chen Y, Zhang Y, et al. Cytoplasmic retention of a nucleocytoplasmic protein TBC1D3 by microtubule network is required for enhanced EGFR signaling. PloS One. (2014) 9:e94134. doi: 10.1371/journal.pone.0094134
81. Wainszelbaum MJ, Charron AJ, Kong C, Kirkpatrick DS, Srikanth P, Barbieri MA, et al. The hominoid-specific oncogene TBC1D3 activates Ras and modulates epidermal growth factor receptor signaling and trafficking. J Biol Chem. (2008) 283:13233–42. doi: 10.1074/jbc.M800234200
82. Li WX and Zhang YL. Novel long non-coding RNA markers for prognostic prediction of patients with bladder cancer. Chin Med Sci J. (2020) 35:239–47. doi: 10.24920/003662
83. Yang M, Zhang Y, Zhou Y, Zhao T, Li Z, Yue H, et al. Analysis of the expression profiles of long noncoding RNAs and messenger RNAs in tongue squamous cell carcinoma. Oral Surg Oral Med Oral Pathol Oral Radiol. (2022) 134:230–9. doi: 10.1016/j.oooo.2022.01.006
Keywords: bladder cancer, single-cell RNA-seq, cancer-associated fibroblasts, molecular subtypes, immune microenvironment
Citation: Zhou Z, Chen Y, Wang Z, Yang S, Li Z, Shi R, Wang R, Liu K, Tang X, Li Q and Xu R (2025) Characterization of cancer-related fibroblasts in bladder cancer and construction of CAFs-based bladder cancer classification: insights from single-cell and multi-omics analysis. Front. Immunol. 16:1580986. doi: 10.3389/fimmu.2025.1580986
Received: 21 February 2025; Accepted: 14 August 2025;
Published: 11 September 2025.
Edited by:
Nan-Shan Chang, China Medical University, TaiwanReviewed by:
Debbie Robbrecht, Erasmus Medical Center, NetherlandsHieu-Huy Nguyen-Tran, China Medical University, Taiwan
Copyright © 2025 Zhou, Chen, Wang, Yang, Li, Shi, Wang, Liu, Tang, Li and Xu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Ran Xu, eHVyYW5AY3N1LmVkdS5jbg==