Wenqian Zhou1†
Zhiheng Lin- 1Tongji Hospital, School of Medicine, Tongji University, Shanghai, China
- 2Shandong University of Traditional Chinese Medicine, Jinan, Shandong, China
- 3Xiangya Boai Rehabilitation Hospital, Changsha, Hunan, China
Background: Clear cell renal cell carcinomas (ccRCCs) epitomize the most formidable clinical subtype among renal neoplasms. While the impact of tumor-associated fibroblasts on ccRCC progression is duly acknowledged, a paucity of literature exists elucidating the intricate mechanisms and signaling pathways operative at the individual cellular level.
Methods: Employing single-cell transcriptomic analysis, we meticulously curated UMAP profiles spanning substantial ccRCC populations, delving into the composition and intrinsic signaling pathways of these cohorts. Additionally, Myofibroblasts were fastidiously categorized into discrete subpopulations, with a thorough elucidation of the temporal trajectory relationships between these subpopulations. We further probed the cellular interaction pathways connecting pivotal subpopulations with tumors. Our endeavor also encompassed the identification of prognostic genes associated with these subpopulations through Bulk RNA-seq, subsequently validated through empirical experimentation.
Results: A notable escalation in the nFeature and nCount of Myofibroblasts and EPCs within ccRCCs was observed, notably enriched in oxidation-related pathways. This phenomenon is postulated to be closely associated with the heightened metabolic activities of Myofibroblasts and EPCs. The Myofibroblasts subpopulation, denoted as C3 HMGA1+ Myofibroblasts, emerges as a pivotal subset, displaying low differentiation and positioning itself at the terminal point of the temporal trajectory. Intriguingly, these cells exhibit a high degree of interaction with tumor cells through the MPZ signaling pathway network, suggesting that Myofibroblasts may facilitate tumor progression via this pathway. Prognostic genes associated with C3 were identified, among which TUBB3 is implicated in potential resistance to tumor recurrence. Finally, experimental validation revealed that the knockout of the key gene within the MPZ pathway, MPZL1, can inhibit tumor activity, proliferation, invasion, and migration capabilities.
Conclusion: This investigation delves into the intricate mechanisms and interaction pathways between Myofibroblasts and ccRCCs at the single-cell level. We propose that targeting MPZL1 and the oxidative phosphorylation pathway could serve as potential key targets for treating the progression and recurrence of ccRCC. This discovery paves the way for new directions in the treatment and prognosis diagnosis of ccRCC in the future.
Introduction
Renal neoplasia stands as one of the most prevalent genitourinary malignancies globally (1). According to data disseminated by the World Health Organization, the annual mortality toll attributed to renal cancer reaches a staggering 140,000-170,000 patients (2, 3). Among these cases, clear cell renal cell carcinomas (ccRCCs) emerge as the predominant histological subtype, contributing to approximately 80% of renal cancers and representing the clinical variant with the highest mortality rates (3–5). Presently, localized instances of ccRCCs are conventionally addressed through surgical resection (6, 7). However, due to the latent nature of symptoms during the early stages and the proclivity for tumor metastasis by the time of detection, reliance solely on surgical intervention becomes less efficacious (8). Distinguishing ccRCCs from other genitourinary tumors lies in their refractoriness to chemotherapy and radiotherapy (9). The prevailing use of chemotherapy and cytokine therapy involving interferon-alpha (INFα) and IL-2 not only presents challenges in ensuring efficacy but also manifests susceptibility to drug resistance and adverse effects (10). Furthermore, the overall survival rate for advanced ccRCCs languishes below 30% (11, 12). Consequently, an imperative exists to investigate the fundamental mechanisms and signaling pathways governing the progression and recurrence of ccRCCs, identify novel therapeutic targets, and formulate innovative prognostic models.
Cancer-associated fibroblasts (CAFs) are activated fibroblasts found within or in the vicinity of cancerous tissues, closely associated with cancer progression. Currently, it is speculated that they primarily originate from resident fibroblasts and can modulate certain cellular mechanisms of tumors, including functions related to tumor proliferation, invasion, and metastasis, leading to adverse prognoses in affected patients (13–17). Existing research has revealed that fibroblasts can promote the tumor drug resistance of ccRCCs through the TDO/Kyn/AhR signaling pathway (18), and they are positively correlated with tumor initiation and adverse prognosis (19). There have been documented interactions between fibroblasts and tumors, enabling them to regulate tumor cell metabolism and modulate associated immune functions, making them conducive to tumor growth (20, 21). Additionally, numerous studies have found that tumor-associated fibroblasts can participate in regulating the stemness and epithelial-mesenchymal transition of tumor cells, thereby promoting tumor progression (22–24). However, there is relatively limited literature available on the interaction between fibroblasts and ccRCCs tumor cells at the single-cell level. This study delves into elucidating the functional roles and underlying mechanisms of Myofibroblasts and EPCs in ccRCCs at the single-cell level. Through meticulous subpopulation stratification of Myofibroblasts and EPCs, we explore intricate cell-cell interactions, specifically focusing on those occurring between the crucial subpopulation C3 and tumor cells within EPCs. Simultaneously, we construct prognostic models linked to these interactions, subsequently subjecting them to rigorous experimental validation. Based on these research findings, we hypothesize a strong correlation between the progression of ccRCCs and the interplay pathways of oxidative phosphorylation and Myofibroblasts, with the pathway gene MPZL1 playing a pivotal role. These studies contribute to a deeper understanding of ccRCCs’ proliferation and recurrence at the single-cell level, aiming to provide insights for the design of better therapeutic targets and prognostic models for ccRCCs treatment.
Methods
Data processing and download
We acquired the neoplastic specimens and adjacent tissues from six patients diagnosed with Clear Cell Renal Cell Carcinomas (ccRCCs) through the European Nucleotide Archive (ENA) database (https://www.ebi.ac.uk/ena/browser/home). The Project number is PRJNA705464. Single-cell samples were extracted. Gene expression quantification using RNA-Seq (HTSeqFPKM) and clinical sample data for ccRCCs were retrieved from The Cancer Genome Atlas (TCGA) website (https://portal.gdc.cancer.gov) and subsequently processed and normalized utilizing R software (R 4.3.0).
Quality control
The R package DoubletFinder was employed for the filtration and removal of low-quality cellular data, including instances of doublet cells and multi-cellularity in the samples (25, 26). The filtering criteria encompassed the following parameters: 1) Total cellular gene transcripts (nCount) ranging from 500 to 80,000; 2) Total cellular gene counts (nFeature) between 200 and 6,000; 3) The proportion of mitochondrial genes less than 25%; 4) The proportion of erythrocyte gene counts below 5%, with a count proportion also below 5%.
Clustering and annotation
We utilized the NormalizeData function within the R package Seurat to standardize the curated single-cell data (27). Subsequently, we computed the variance and standard deviation for each gene through the FindVariableFeatures function, selecting the top 2000 genes as highly variable based on the genes’ average expression and the extent of their dispersion (28, 29). All genes underwent centering via the ScaleData function, and the cyclical effects of distinct cells were assessed using the CellCycleFeatures function. The RunPCA function facilitated dimensionality reduction based on the expression of the top 2000 high-margin genes, followed by the mitigation of batch effects from the samples using the RunHarmony function from the R package harmony.
For clustering the data after dimensionality reduction, the FindNeighbors and FindClusters functions from the Seurat package were employed. Integration and annotation of distinct cell clusters post-clustering were conducted with reference to the CellMarker database (http://bio-bigdata.hrbmu.edu.cn/CellMarker/) and the SingleR function. Subsequent corrections were applied through a meticulous review of pertinent literature. Sample source analysis for diverse large cell clusters was executed, and the identification and annotation of differentially expressed MARKER genes among these clusters were performed using the FindAllMarkers function (30, 31).
Analysis of cell sources and enrichment of large populations
We scrutinized the cellular origins and relative sizes of the principal clusters within ccRCCs. A comparative analysis was conducted on the total number of transcribed genes and the transcriptional profile across different major clusters, visually represented through column line graphs and UMAP (32, 33). Furthermore, we identified differential genes (DEGs) among the cells of these major clusters, necessitating DEGs to be detected in a minimum of 25% of the cells, with a false discovery rate (FDR) below 0.05, | logFCfilter | exceeding 1. The DEGs identified in the major clusters underwent Gene Ontology (GO) enrichment analysis.
Additionally, we independently computed DEGs between the two cellular macropopulations of Myofibroblasts and Epithelial Cells (EPCs) and other macropopulations. This analysis was carried out using the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (c2.cp.kegg.v7.5.1.symbols.gmt). The marker genes retrieved from the database underwent filtration and analysis for Gene Set Enrichment Analysis (GSEA), with statistical significance set at FDR < 0.05.
Analysis of cellular metabolism in large populations
We computed the expression levels of stemness genes across various major clusters within ccRCCs. Additionally, we evaluated cellular metabolic pathways within the major clusters of ccRCCs and distinct cellular sources using the R package scMetabolism. The outcomes were presented in the form of heatmaps, showcasing the top 20 pathways. Further exploration involved the selection of the top three pathways with the highest scores, visualizing their expression across major populations through UMAP. A comparative analysis was performed on the top three metabolic pathways within major populations and diverse sources of ccRCCs, illustrated in the form of violin plots.
Myofibroblasts subgroup correlation analysis
We employed the FindNeighbors function and FindClusters within the Seurat package to reclustering the downscaled Myofibroblasts cell clusters. Subsequently, we identified marker genes for each subgroup based on the percentage of expression and the frequency of expression, and conducted chromosome copy variation analysis (inferCNV) in different subpopulations of Myofibroblasts. We delved into the tissue origin and cellular staging of cells within each subpopulation, visualizing the findings through UMAP. Furthermore, we conducted a comparative analysis of the tissue origin and cellular staging ratios among different subpopulations. Staging scores, including G2M.Score and S.Score, were computed and compared across subpopulations. Additionally, we explored and compared the nFeature and nCount of genes within each subgroup. Disparities in G2M.Score, S.Score, nFeature, and nCount were calculated and compared across various subpopulations.
Trajectory and correlation analysis of myofibroblasts subgroups
To investigate the differentiation and evolutionary processes among cells within each subpopulation of Myofibroblasts, we conducted a proposed time-series analysis. Initially, we predicted the relative differentiation status of each cell subpopulation using the R package Cytotrace, assigning scores within the interval of 0-1, where higher scores indicate heightened cell stemness. Subsequently, the developmental time of each cell subpopulation was calculated and organized by the R package monocle, unveiling the temporal transcriptional dynamics of marker genes within each cell subpopulation.
Following chronological ordering, cells were categorized into three periods, and the sequential changes of subpopulation cells, along with cells originating from different sources within the subpopulation, were computed for each period. Ultimately, leveraging slingshot, cell clustering clusters, and spatial dimensionality reduction information, we executed cell differentiation genealogy construction and pseudo-temporal extrapolation. This approach validated the temporal order of different subpopulations at an elevated level and calculated the proportionality between cell staging and diverse subpopulations according to slingshot.
We further computed and compared the stemness genes of cells within each subpopulation, visualizing the results. Through an amalgamation of tissue origin and cell staging information of subpopulation cells, coupled with the analysis of subpopulation cell stemness and associated genes, we identified the C3 subpopulation as a potential key subpopulation. Subsequently, using the pyscenic algorithm, we calculated the top 5 transcription factors for each subpopulation, delving into the overall distribution and regulatory regions of the pivotal C3 subpopulation.
Subpopulation cell interaction analysis
To elucidate the cell-cell interactions among subpopulations and larger populations, and to delve further into the interaction dynamics between the pivotal C3 subpopulation of Myofibroblasts and tumor cells, we conducted Copy Number Variation (CNV) analysis on EPCs cells, distinguishing tumor cells within EPCs based on the count of copy variants. Subsequently, utilizing the R package cellchat, we analyzed the overall count of cell-cell interactions and their respective strengths. These interactions were further categorized into those involving tumor cells and those with other cell types. Based on the interaction results, we scrutinized the crucial interaction pathway MPZ between C3 and tumor subpopulations.
Clinical correlation and independent prognosis analysis of C3 subpopulations
The R software was employed to determine the intersection of the identified marker genes within the C3 subgroup with both ccRCC tumor tissues and normal tissues. Samples from ccRCCs with incomplete clinical data were excluded. Subsequently, the intersected genes were merged with the standardized clinical data of ccRCCs. Univariate COX risk regression analysis was conducted using the coxph function within the R package survival. This analysis was validated through the application of the Least Absolute Shrinkage and Selection Operator (LASSO)-penalized Cox regression (34–36). Finally, multivariate COX risk regression analysis was performed to identify prognostic differential genes.
The risk score for each sample was calculated (risk score = Xλ, where X represents the relative expression level of prognostic genes and coefλ is the coefficient) (37, 38). Samples were categorized into high and low-risk groups based on their risk scores. Coefficient values of prognostic genes were computed, and the distribution of these genes in the high and low-risk groups was visualized through a heatmap. Kaplan-Meier curves illustrated the survival at different time points in the high- and low-risk groups (39, 40), while the specificity and sensitivity of prognostic genes were further validated using time-dependent receiver operating characteristic (ROC) curves (41, 42).
Additionally, the correlation between specific prognostic genes and risk scores was calculated and demonstrated. To construct the prognostic model, information on race, gender, age, and tumor stage in the samples was incorporated. A nomogram was created using the R package rms, combining the risk scores to predict the prognosis of ccRCC patients (43). This model was finally validated through ROC curves and decision curve analysis (DCA).
Immune correlation analysis and enrichment analysis
We computed and authenticated immune infiltration in the high and low-risk groups utilizing the xCell and CIBERSORT deconvolution algorithms to compare and validate immune infiltration, respectively. The correlation of immune cells with the risk score and prognostic genes was explored. TIDE scores were calculated for high and low-risk groups, and the correlation of prognostic genes and risk scores with immune checkpoint genes was investigated. Tumor immune microenvironment-specific scores of the high and low-risk groups were also calculated and compared.
To delve into the correlation mechanism of the high and low-risk groups, we further filtered the differential genes using the R package limma, employing filtering conditions of |log FC| > 1 and a threshold value of FDR (BH) corrected P.adj <0.05. GO and KEGG enrichment analyses were conducted using the R package clusterProfiler, encompassing biological processes, cellular components, and molecular functions in GO enrichment analyses (44, 45).
Cell culture
The 786-O and CAKI-1 cell lines were procured from the American Type Culture Collection (ATCC). Both cell lines underwent cultivation in F12K medium supplemented with 10% fetal bovine serum (Gibco BRL, USA) and 1% streptomycin/penicillin. PRMI1640 medium (Gibco BRL, USA) was utilized for the CAKI-1 cell line. Standard conditions of incubation were maintained at 37°C, 5% CO2, and 95% humidity.
Cell transfection
MPZL1 knockdown was achieved through the utilization of small interfering RNA (siRNA) constructs obtained from GenePharma, Suzhou, China. The transfection protocol closely adhered to the steps outlined for Lipofectamine 3000RNAiMAX (Invitrogen, USA). Cells were seeded in 6-well plates at 50% confluence and subsequently transfected with negative control (si-NC) and knockdown (Si-MPZL1-1 and Si-MPZL1-2). Each transfection process was carried out utilizing Lipofectamine 3000RNAiMAX (Invitrogen, USA).
Cell viability assay
Cell viability of 786-O and CAKI-1 cells post-transfection was assessed using the CCK-8 assay. Cell suspensions were seeded at a density of 5*103 cells per well in 96-well plates and cultured for 24 hours. Subsequently, 10 μL of CCK-8 reagent (A311-01, Vazyme) was added to each well and incubated in the dark at 37°C for 2 hours. Cell viability was determined by measuring absorbance at 450 nm using an enzyme marker (A33978, Thermo) on days 1, 2, 3, and 4. The mean OD values were computed and presented on a line graph.
5-Ethyl-2 ' -deoxyuridine proliferation assay
Transfected 786-O and CAKI-1 cells were seeded in 6-well plates at a density of 5 × 103 cells per well and allowed to incubate overnight. A 2x EdU working solution was prepared by incorporating a 10 mM EdU solution in serum-free medium, which was subsequently added to the cell culture and incubated at 37°C for 2 hours. Following incubation, the medium was aspirated, and cells were washed with PBS before fixation with 4% paraformaldehyde for 30 minutes. Subsequent treatment involved exposure to glycine (2 mg/ml) and 0.5% Triton X-100 for 15 minutes. Cells were then incubated with 1 ml 1X Apollo and 1 ml 1X Hoechst 33342 for 30 minutes at room temperature. Cell proliferation was quantified through fluorescence microscopy.
Wound healing
The transfected cells were seeded in 6-well plates and cultured until reaching a cell density of 95%. A sterile 200 μL pipette tip was initially employed to create a straight-line scratch across the cell layer in the culture wells, followed by gentle rinsing with PBS. Subsequently, the medium was replaced to facilitate ongoing cell culture. Photographs of the scratch were captured at the identical location at 0 hours and 48 hours, and the width of the scratch was measured.
Transwell experiment
Cells underwent serum deprivation in medium lacking serum for a duration of 24 hours preceding the experiment. Following treatment involving the addition of matrix gel (BD Biosciences, USA), the cell suspension was introduced into the upper chamber, which contained Costar, while serum medium was introduced into the lower chamber. Subsequent to a 48-hour incubation period in the incubator, cells were fixed using 4% paraformaldehyde and stained with crystal violet to assess their invasive capability.
Results
ccRCCs large population cell categorization
We acquired the single-cell data of ccRCCs from the ENA database (PRJNA705464), encompassing Tumor, Normal, LymphNode, and Peripheral Blood (PBMC) samples. Utilizing the R package DoubletFinder, we eliminated duplicated cells and subsequently filtered out subpar cells in the samples (Supplementary Figure 1A). Cell staging tests revealed a concentrated distribution in the PCA plot, indicating minimal impact on the study results (Supplementary Figure 1B). Selecting the top 2000 highly variable genes based on expression and dispersion (Supplementary Figure 1C), we downscaled them using RunPCA and retained the top 30 dimensions for subsequent analysis (Supplementary Figures 1D, E). Exploration of the top 10 highly variable genes for each of the first 9 dimensions was visualized (Supplementary Figure 1F). We categorized ccRCCs into 37 clusters (Figure 1A) and annotated them using the R package singleR and CellMarker databases. After further correction based on relevant literature, we categorized them into 8 major cell clusters (Figure 1B), namely T_NK(81753), Myeloid_cells(34761), ECs(7652), Myofibroblasts(7959), Pericytes(3725), and B_Plasma(3781).
Figure 1 Classification of ccRCCs cells into large clusters. (A) 6 patients with ccRCCs ccRCCs tumor and surrounding tissues were classified into 38 clusters after single cell clustering; (B) According to different MARKER genes, the ovaries were annotated into T_NK, Myeloid_cells, ECs, Myofibroblasts, Pericytes, B_Plasma, EPCs, MCs respectively, a total of 8 (C) Distribution of each sample source. (D) Each population of cells exhibits the top five marker genes in ccRCCs.
During the study, Myofibroblasts, Epithelial Cells (EPCs)(4265), and MCs(1173) emerged as crucial cell groups in the tumor. Reviewing literature on the annotated cell populations, we identified EPCs as pivotal cells in the tumor, with a close association with Myofibroblasts. Disparities in sample distribution among different cells were observed, possibly linked to the individual source of the samples and the location from which tumor tissues were selected (Figure 1C). Examination of the top 5 MARKER genes in each subpopulation (Figure 1D) revealed that MARKER genes of EPCs also exhibited some expression in LymphNode cells.
Analysis of cell origin and gene transcription correlation in large populations
We investigated the spatial distribution of the cell population in ccRCCs (Figures 2A–C) and noted a relatively high proportion of tumor cells within Myofibroblasts. In contrast, EPCs displayed a mixed composition, including a significant presence of normal tissue cells. Acknowledging the necessity to extract tumor cells from EPCs in subsequent analyses, we proceeded accordingly. Examination of nCount and nFeature (Figures 2D–G) revealed higher total transcript numbers and gene quantities in EPCs and Myofibroblasts cells, indicating elevated cellular activities and functions in these two locations.
Figure 2 Source site and transcriptional analysis of the major clusters of cells. (A–C) ccRCCs major clusters of cells source site and distribution plot; (D, E) ccRCCs major clusters of cells nFeature and nCount UNMP distribution plot; (F, G) ccRCCs major clusters of cells nFeature and nCount expression level box line plot.
Enrichment function analysis
To explore the functional attributes of the cell populations in the large cohort of ccRCCs, we scrutinized the differential genes between these cells and others (Figures 3A–H). Subsequently, we conducted GO enrichment function analysis based on the identified differential genes (Figure 3I). Notably, both EPCs and Myofibroblasts exhibited enrichment in oxidative respiration, with pathways such as Aerobic respiration and Oxidative phosphorylation significantly enriched in EPCs, while ATP metabolic process and Purine ribonucleoside triphosphate metabolic process were notably enriched in Myofibroblasts.
Figure 3 Enrichment analysis of major clusters of cells. (A–H) Differential gene distribution plot of ccRCCs major population cells versus other cells, in which the top 5 selected genes from each of the high and low expression genes are shown; (I) GO enrichment analysis plot of ccRCCs major population cells; (J–O) GSEA enrichment analysis plot of Myofibroblasts, with 5 pathways selected from each of the low expression group and the high expression group. (P–U) GSEA of EPCs Enrichment analysis plot with 5 pathways selected from each of the low expression group vs. high expression group.
To delve deeper, we performed GSEA enrichment analysis for EPCs and Myofibroblasts. In Myofibroblasts (Figures 3J–O), the low-expression group showed significant enrichment in pathways related to immune response-activating signal transduction, immune response-activating cell surface receptor signaling pathway, and leukocyte-mediated immunity. Conversely, the high-expression group exhibited enrichment in pathways like Glycolytic process, Nucleotide phosphorylation, and ATP generation from ADP. In EPCs (Figures 3P–U), the low expression group displayed enrichment in pathways related to the regulation of immune response, immune response-activating cell surface receptor signaling pathway, and regulation of leukocyte activation. In contrast, the high expression group showed enrichment in pathways associated with cellular respiration, aerobic respiration, and oxidative phosphorylation.
These findings suggest that Myofibroblasts and EPCs share similar GSEA-enriched pathways, hinting at a potential connection between the high- and low-expression groups. We posit that the oxidative response of Myofibroblasts may drive the suppression of immune function.
Cell stemness and metabolic analysis
We conducted an analysis of cell stemness genes in the large cell populations, identifying a high score for CD44 in Myofibroblasts (Figure 4A). Further exploration of metabolic function differences between large cell populations and cells of different origins (Figures 4B, C) revealed significantly higher scores for cellular metabolic pathways, including Oxidative phosphorylation, Glycolysis/Gluconeogenesis, and the Pentose phosphate pathway—all associated with oxidative respiration. To provide a comprehensive view, we compared the distribution and scores of these three pathways in each cell population (Figures 4D–L), observing their pronounced significance in both EPCs and Myofibroblasts. This aligns with our earlier findings from GO and GSEA enrichment analyses.
Figure 4 ccRCCs large population of cell stemness genes and metabolic pathways map. (A)Analysis of stemness genes in large populations of ccRCCs; (B) Top 20 metabolic pathways in large populations of ccRCCs; (C) Top 20 metabolic pathways in cells of different positional origins of Tumor, Normal, LymphNode, and PBMC; (D–L) Distribution of UNMPs for the three pathways “Oxidative phosphorylation, Glycolysis/Gluconeogenesis, Pentose phosphate pathway” and violin maps of the scores of the three pathways in the large population of ccRCCs with different locations of origin. **** means p <0.0001; ns means no statistical difference.
Subclustering of myofibroblasts
Distinct Myofibroblast subpopulations were identified as C0, C2 (FXYD2+), C3 (HMGA1+), C4 (ITGA1+), and C5 (PTTG1+) (Figure 5A). We calculated copy number variations (CNVs) of these subpopulations and the respective marker genes (Figures 5B, C). The cellular origin investigation revealed higher tumor tissue content in C0, C3, and C4 (Figures 5D, E). Assessing the proportionality of subpopulations with respect to cell origin (Figures 5J, K) indicated the highest amount and purity in C3 tumor tissues. Examining cell stage distribution (Figures 5F, G) showed a relatively high proportion of G2M and S stages in C3 and C5. The cell stage ratio graph (Figures 5H, I) highlighted significantly higher S stage in C3 compared to G2M and G1 stages, indicating vigorous cellular DNA synthesis, replication, and proliferative ability. G2M.Score and S.Score distributions (Figures 5L–O) revealed higher scores for both C3 and C5 subpopulations. Additionally, comparing nFeature and nCount values of different subpopulation cells (Figures 5P, Q) demonstrated significantly higher values in C3 compared to other subpopulations.
Figure 5 Subpopulation of ccRCCs cells and correlation analysis. (A) UNMP distribution of subpopulation cells; (B) Heatmap of CNV situation of subpopulation cells; (C) Expression of the first 5 marker genes in each subpopulation cells; (D, E) UMAP map of the distribution of source cells in different locations of Tumor, Normal, LymphNode; (F, G) UNMP map of the distribution of subpopulation of cells in each staging cells of G2M, S, G1; (H, I) G2M, S, G1 Proportional distribution between cells and subpopulations of cells in each stage; (J, K) Proportional distribution between subpopulations of cells and cells originating from different locations of Tumor, Normal, and LymphNode; (L) UNMP distribution of G2M.Score; (M) S.Score UNMP distribution; (N–Q) G2M.Score, S.core, nFeature and nCount Violin graph representing score.
Proposed time series analysis of subpopulation cells
We predicted the differentiation of individual subpopulations by the R package Cytotrace and found that the C3 group had the highest degree of differentiation, followed by C0 and C4 (Figures 6A, B). To understand the temporal trajectory relationship between cell subpopulations, we calculated the relative temporal order of the individual subpopulations by the R package monocle (Figure 6C), and found that the C1 subpopulation was at the starting point, whereas the C3 and C0 subpopulations were near the end point. We also analyzed the temporal trajectory relationship of marker genes of each subpopulation (Figure 6D) and explored the temporal trajectories of different cellular origins of each subpopulation and individual subpopulations (Figures 6E–G), and we found that subpopulations C1 and C2 were more at the beginning of the trajectory, while C0 and C3 were mostly located at the end of the trajectory. We also validated the intercellular trajectories by slingshot and found that they could be divided into three periods and there were two time trajectories (Figures 6H, I), while C3 was at the middle and end of trajectory 1 and at the tail end of trajectory 2. We also calculated the proportional relationship between different subpopulations and different periods with each other (Figures 6J, K), and found that the C3 subpopulation was in stages 2 and 3, and mainly existed in stage 3, indicating that the C3 subpopulation was at the end of the time trajectory compared with the C0 subpopulation, which mostly represented that the tumors had already completed the transition from normal cells to malignant cells, and might have a high degree of invasion and value-added ability, as well as drug resistance. With the above analysis, we consider the C3 subgroup as the key subgroup.
Figure 6 Proposed time-series analysis of subpopulation cells. (A, B) Cytotrace scores of subpopulation cells; (C, D) monocle-based time trajectory maps of each subpopulation and marker genes; (E, F) monocle-based time trajectory staging, which can be categorized into three periods; (G) distribution of subpopulation cells with time trajectories based on monocle; (H, I) slingshot-based time staging maps of subpopulation cells, which can be divided into three periods with two trajectories of merit; (J, K) distribution of the ratio of subpopulation cells to time staging period.
Analysis of cell stemness genes and transcription factors in subpopulations
In the analysis of cell stemness genes (Figures 7A–C), it was discerned that the C3 subpopulation exhibited the most robust scores for cell stemness genes, with CD44 showing heightened expression in both the C3 subpopulation and tumors. Furthermore, an examination of transcription factors in each subpopulation revealed that the top 5 transcription factors in the C3 subpopulation were HMGA1, PBX1, NFIC, BCLAF1, and RFX3 (Figure 7D). A detailed exploration of the expression of these 5 transcription factors across cellular subpopulations (Figures 7E–K) highlighted the prominent expression of the transcription factor HMGA1, particularly in the predominant C3 subpopulation.
Figure 7 Subpopulation cell stemness gene and transcription factor analysis. (A) Cellular stemness gene expression level violin plot for each subpopulation of cells in Myofibroblasts; (B) UMAP distribution of stemness gene expression level in subpopulation of cells; (C) Stemness gene expression in each subpopulation of cells in Myofibroblasts; (D) Top 5 transcription factors of each subpopulation of cells in Myofibroblasts; (E–K) Top 5 transcription factors UMAP distribution graph.
Analysis of cellular interactions
Through the computation of CNV values in EPCs, we pinpointed the tumor cell populations within EPCs (Supplementary Figure 2). Subsequent exploration into the number of cellular interactions within subpopulations of cells and the broader cell populations (Figures 8A–C) unveiled that tumor cells exhibited the highest degree of contact with the C3 subpopulation. Additionally, in the assessment of interaction strength (Figures 8D–F), tumor cells demonstrated more robust connections with the C3 subpopulation compared to other cells. Our deduction is that the C3, identified as a key subpopulation of Myofibroblasts, likely harbors the most potent interactions with tumor cells. Further scrutiny into the interactions involving the C3 subpopulation and other subpopulations (Figures 8G–I) brought to light that the MPZ signaling pathway network held particular significance among the subpopulations. Notably, MPZL1, a pivotal gene in this pathway, scored the highest within the C3 subpopulation.
Figure 8 Interaction analysis between subpopulation cells and large population cells. (A-C) Number of interactions between Myofibroblasts subpopulation cells and ccRCCs large population cells graph, where A is the number of interactions between all cells, B figure tumor cells as source, C figure tumor cells as target, the thicker line represents the more number; (D–F) Interaction between Myofibroblasts subpopulation cells and ccRCCs large population cells Intensity map, where D is the number of interactions among all cells, (E) means that Tumor cells are the source, and (F) means that Tumor cells are the target; thicker lines represent higher intensity; (G) MPZ signaling pathway network pathway interactions between subpopulations and large populations of cells; (H) MPZL1 gene expression in subpopulations and large populations; (I) Scores of action sites in MPZ signaling pathway network pathway in each subpopulation and large populations.
Cellular prognostic analysis
We conducted a comprehensive analysis of ccRCCs by intersecting differentially expressed genes from the TCGA database with C3 subpopulation marker genes. Through univariate COX risk regression analysis, we identified 46 prognostically relevant differential genes (Figure 9A). Validation with LASSO COX risk regression analysis confirmed the stability and reliability of these genes (Figures 9B, C). Subsequently, multivariate COX risk regression analysis yielded 22 final prognostic genes, including high-risk genes such as S100A16, FABP5, COL6A2, DCBLD2, KCNN4, COL6A1, TPM4, LINC00472, PYGB, ARHGAP29, TUBB3, and STEAP3, and low-risk genes such as ANXA2, MACF1, ACTN4, ZBTB38, LBH, NCKAP5, VGF, SPARC, AP1S1, and COL21A1. Patients were stratified into high- and low-risk groups based on risk scores (Figures 9D, E). The survival analysis revealed a progressive increase in the number of patient deaths over time, with the high-risk group exhibiting lower survival rates (Figures 9D, E). Examination of coefficient values (Figure 9F) showed TUBB3 with the highest score and ANXA2 with the lowest. The expression patterns of prognostic genes in high- and low-risk groups were explored (Figure 9G). Kaplan-Meier curve analysis demonstrated significantly lower survival rates in the high-risk group (Figure 9H), with a meaningful result (P < 0.0001). ROC curve analysis indicated stable and good predictive performance at 1, 3, and 5 years (Figure 9I). We further examined the correlation of different genes with risk score and prognosis (Figure 9J) and visualized six genes significantly associated with risk score (Figure 9K). Exploration of the correlation of clinical factors with the risk score and construction of a nomogram (Figures 10A–G) allowed for predicting the survival rate of ccRCC patients at 1, 3, and 5 years. Validation through ROC and DCA curves demonstrated the model’s high stability, sensitivity, and clinical utility (Figures 10H, I).
Figure 9 Independent prognostic analysis. (A) 46 prognosis-related differential genes; (B) Distribution of LASSO analysis coefficient spectra of 22 prognosis genes; (C) Optimal cross-validation of parameter selection in LASSO regression; (D) Survival time and survival status plots of patients in different risk groups over time; (E) Classification of patients into high- and low-risk groups based on risk scores; (F) Coef values of prognosis-related genes; (G) Heat map of prognosis-related gene distribution; (H) Kaplan-Meier prognostic analysis curves for high and low risk groups; (I) Time-dependent ROC curves with area under the curve (AUC) of 0.79, 0.74, 0.79 at 1, 3, and 5 years; (J). Correlation analyses between genes and risk scores and OS; (K) Dot plots of the top 6 prognostic genes that had a strong correlation with the risk score.
Figure 10 Clinical correlation analysis. (A–F) analysis of the correlation of risk scores with racial factors, gender, age, and tumor stage T, N and M (G) Nomogram plots of patients with ccRCCs at 1, 3, and 5 years; (H) time-dependent ROC plots, with AUCs of 0.86, 0.80, and 0.80 at 1, 3, and 5 years, respectively; (I) DCA analysis of prognostic models Figure.
Immune correlation analysis
We computed and compared immune cell infiltration in the high- and low-risk cohorts (Figures 11A, B) and examined the correlation of immune cells with risk scores and prognostic genes (Figures 11C, D). T cells CD4 memory activated showed a positive correlation with RISK, whereas T cells CD4 memory resting exhibited a negative correlation with RISK. Immune cell disparities between high and low-risk groups were also assessed (Figure 11E), revealing relatively high scores of Macrophages M2, particularly notable in the low-risk group. TIDE score analysis (Figure 11F) indicated a high score in the high-risk group, suggesting a potential for immune escape. Additionally, immune checkpoint score calculation (Figure 11G) revealed a predominantly negative correlation with RISK score. Comparing tumor microenvironment-related scores in the high- and low-risk groups (Figures 11H–J) showed higher scores in both high-risk cohorts.
Figure 11 Immune correlation analysis. (A) Heatmap of prognostic genes, tumor microenvironment and immune cell profile in high and low risk groups; (B) Difference in the proportion of immune cells in high and low risk groups; (C) Immune cell risk score relationship score; (D) Immune cell correlation with prognostic genes, risk score and OS; (E) Difference in immune cell expression between high and low risk groups; (F) TIDE score violin plot between high and low risk groups; (G) Situation point plot of the relationship between immune checkpoints and prognostic genes, risk scores; (H–J) Difference in immune scores, stromal scores, and total scores of tumor microenvironment between high and low risk groups. * means p <0.05; ** means p <0.01; *** means p <0.001; **** means p <0.0001; ns means no statistical difference.
Experimental result
We selected two distinct types of ccRCC cells for the experiment, conducting a comparative analysis between the control group and the knockdown infection group. Remarkably, both cell types exhibited significantly diminished cell viability levels upon MPZL1 knockdown (Figures 12A, B). In the cluster formation assay, the number of colonies showed a substantial reduction in both cell types following MPZL1 knockdown (Figure 12C). The Transwell assay results (Figure 12D) revealed markedly reduced staining regions in cells with MPZL1 knockdown, indicating a diminished invasive capability. In cell scratch assay experiments (Figure 12E), the 48-hour scratch width in both cell lines subjected to MPZL1 knockdown was significantly wider compared to the negative control group. The cell proliferation assay (Figure 12F) further demonstrated a lower cell count in both cell lines with MPZL1 knockdown compared to the negative control group. Thus, our experimental findings affirm that MPZL1 plays a promotive role in the proliferation, migration, and invasion of ccRCCs, and knockdown of MPZL1 effectively inhibits these processes, restraining the recurrence of ccRCCs.
Figure 12 In vitro experimental validation of MPZL1. (A, B) CCK-8 assay showed that cell viability was significantly reduced after MPZL1 knockdown; (C) Colony formation assay showed that the number of colonies in MPZL1 knockdown cells was significantly lower than that in si-NC group; (D) Transwell assay showed that MPZL1 knockdown significantly slowed down the growth of 786 -O cells and CAKI-1 cells invasion; (E) Scratch assay showed that MPZL1 knockdown significantly slowed down the migration of 786-O cells and CAKI-1 cells; (F) EdU staining results showed that MPZL1 knockdown inhibited the proliferation of 786-O cells and CAKI-1 cells. * means p <0.05; ** means p <0.01; *** means p <0.001.
Discussion
Renal cell carcinoma stands as the seventh most frequently diagnosed malignant tumor, exhibiting a gradual rise in incidence within developed countries in recent years (46). Among its common clinical subtypes, clear cell renal cell carcinomas (ccRCCs) are associated with diverse clinical prognostic outcomes, revealing metastasis in up to 30% of patients at the onset (47). Surgical intervention, though a common treatment approach, proves ineffective, and there remains a notable risk of recurrence post-surgery (48). Given the limited sensitivity of ccRCCs to radiotherapy drugs (9), the identification of mechanisms and therapeutic targets related to ccRCC progression becomes crucial.
Conducting enrichment analysis and metabolic pathway examination of Myofibroblasts and Endothelial Progenitor Cells (EPCs), we observed elevated scores in oxidation-related pathways and diminished immune function. Notably, the oxidative phosphorylation metabolic pathway exhibited the highest score in both EPCs and Myofibroblasts. A literature review revealed that redox homeostasis underpins normal cellular physiological activities and survival (49). It has been recognized that tumor cells display heightened oxidative processes compared to normal cells (50). Oxidative phosphorylation, a process associated with the generation of reactive oxygen species (ROS), poses a threat to cellular DNA and heightens cancer risk (51–54). Studies by Kuo (55) et al. have indicated that ROS activation correlates with increased blood vessel proliferation, suppression of immune microenvironmental functions, and intensified production of macrophage IFN and IL-6, ultimately hampering the immune response in the tumor microenvironment (56–58). Albiñana (59) et al. identified that reducing oxidative stress and ROS levels through targeting ADRB2 proves beneficial in controlling ccRCC progression. Additionally, Costa (60) et al. established that ROS can induce inflammation, fostering the transformation of fibroblasts into Myofibroblasts, exhibiting characteristics akin to invasive adenocarcinoma-associated fibroblasts.
Thus, we posit that the oxidative phosphorylation pathway and its affiliated pathways play a pivotal role in the progression of clear cell renal cell carcinomas (ccRCCs), impeding the immune microenvironment of these tumors and facilitating their transformation into Myofibroblasts. Consequently, this phenomenon propels the progression and augmentation of the tumor.
Through an exploration of cellular interactions, we discovered that the MPZ signaling pathway network exhibits substantial cellular interactions in crucial subpopulations of Myofibroblasts and tumor cells, with the gene MPZL1 identified as a key player. A comprehensive review of existing literature elucidates that MPZL1 (Myelin protein zero like 1) serves as a binding protein and substrate of tyrosine phosphatase SHP-2 (61), demonstrating widespread expression across various cell types and involvement in essential biological processes and molecular functions (62, 63). Furthermore, the biological impact of MPZL1 in diverse tumors, including ovarian, colorectal, and bladder cancers, has been investigated and corroborated (64–66), highlighting its promotional influence on tumor progression. For instance, Jia et al. (67) observed that MPZL1 promotes the migration of hepatocellular carcinoma cells by inducing the phosphorylation and activation of pre-metastatic proteins. Similarly, Chen (64) reported that elevated expression of MPZL1 stimulates the phosphorylation of Src kinase, facilitating the proliferation and migration of ovarian cancer cells. Additionally, Wang et al. (68) found that overexpression of MPZL1 correlates with the suppression of immune function in lung cancer studies. We also postulate that MPZL1 is intricately linked to the proliferation and migration of ccRCCs, exerting inhibitory effects on the immune microenvironment of tumors, thereby contributing to tumor progression or recurrence. Our experimental results substantiate this hypothesis.
TUBB3 (tubulin beta 3 class III), belonging to the β-tubulin protein family, intricately engages in neurogenesis and microtubule assembly within neuronal cells. Its dynamic interplay involves elongation or depolymerization throughout the cell cycle, exerting substantial influence on cellular morphology, division, and cytoskeletal architecture (69–71). A plethora of studies has consistently emphasized the notable correlation between elevated TUBB3 expression and unfavorable prognoses across various malignancies, including lung, ovarian, and gastric cancers, establishing it as a prognostic indicator for cancer recurrence (72–74). Moreover, certain studies have indicated that high expression of TUBB3 leads to increased drug resistance in tumor cells, which adversely affects patient prognosis (75, 76). Additionally, adenocarcinoma patients with negative TUBB3 expression tend to have longer survival times compared to those with positive expression (77). Furthermore, it has been observed that the TUBB3-associated protein network is involved in certain oxidative stress processes, which may enhance the vitality and drug resistance of tumor cells (78). Experimental research on anti-microtubule chemotherapy drugs has revealed that tumor cells with high TUBB3 expression exhibit stronger resistance to drugs like paclitaxel compared to those with lower expression (79, 80). Our research similarly hypothesizes a close relationship between TUBB3 and oxidative stress-related pathways, potentially promoting the proliferation and recurrence of ccRCCs through this process. This interaction may lead to ccRCCs developing resistance to relevant chemotherapy drugs, ultimately facilitating the progression of ccRCCs.
Conclusion
In this comprehensive study, a thorough examination of single-cell data related to clear cell renal cell carcinomas (ccRCCs) systematically elucidates the pivotal role of oxidative-related pathways within the diverse cellular milieu of ccRCCs. Simultaneously, we meticulously investigate the critical involvement of these pathways in myofibroblasts and endothelial progenitor cells (EPCs). Additionally, a detailed exploration is conducted into the interactive pathways between key subgroups, such as the C3 HMGA1+ Myofibroblasts in myofibroblasts, and tumor cells within EPCs. This comprehensive analysis extends to the investigation of the essential gene MPZL1, culminating in the construction of a ccRCCs prognostic model intricately linked to these pathways. We posit that both MPZL1 and the oxidative stress pathway harbor the potential to emerge as critical therapeutic targets for the treatment and resistance against recurrence in ccRCCs. This assertion is substantiated through empirical validation, specifically demonstrating the concrete impact of MPZL1 on ccRCCs. Our research has confirmed the significant roles played by Myofibroblasts and the key gene MPZL1 in the progression of ccRCCs, providing novel insights for future studies on ccRCCs. We have discovered that targeting MPZL1 and the oxidative phosphorylation pathway could serve as potential key targets for the treatment of ccRCCs progression and recurrence. This finding opens up new directions for the treatment and prognosis diagnosis of ccRCCs in the future. In conclusion, our discoveries explore the mechanisms underlying the proliferation and recurrence of ccRCCs and shed light on prognostic markers and therapeutic targets for ccRCCs.
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.
Author contributions
WZ: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Project administration, Resources, Software, Validation, Visualization, Writing – original draft, Writing – review & editing. ZL: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Project administration, Resources, Software, Validation, Visualization, Writing – original draft, Writing – review & editing. WT: Funding acquisition, Supervision, Writing – review & editing.
Funding
The author(s) declare that no financial support was received for the research, authorship, and/or publication of this article.
Acknowledgments
Thanks to everyone who have contributed to this article.
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.
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.2024.1374931/full#supplementary-material
Supplementary Figure 1 | Quality control and dimensionality reduction clustering of data. (A) Quality-controlled and filtered single-cell data of ccRCCs. (B) Filtered ccRCCs cell stage examination showing G1, G2M, and S stages. (C) Selection of the top 2000 highly variable genes based on gene expression and dispersion, with the top 10 highly variable genes indicated. (D) RunPCA dimensionality reduction plot showing the first 30 dimensions of the selected 50 dimensions. (E) PCA plot of different sample sources after dimensionality reduction. (F) Heatmap of the top ten highly variable genes in the first nine dimensions.
Supplementary Figure 2 | Heatmap of EPCs cell CNV (Copy Number Variation) status.
References
1. Siegel RL, Miller KD, Jemal A. Cancer statistics, 2020. CA Cancer J Clin. (2020) 70:7–30. doi: 10.3322/caac.21590
2. Delman KA. Introducing the "Virtual tumor board" series in CA: A cancer journal for clinicians. CA Cancer J Clin. (2020) 70:77. doi: 10.3322/caac.21598
3. Ferlay J, Soerjomataram I, Dikshit R, Eser S, Mathers C, Rebelo M, et al. Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012. Int J Cancer. (2015) 136:E359–86. doi: 10.1002/ijc.29210
4. Bray F, Ferlay J, Soerjomataram I, Siegel RL, Torre LA, Jemal A. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. (2018) 68:394–424. doi: 10.3322/caac.21492
5. Zhang Z, Lin E, Zhuang H, Xie L, Feng X, Liu J, et al. Construction of a novel gene-based model for prognosis prediction of clear cell renal cell carcinoma. Cancer Cell Int. (2020) 20:27. doi: 10.1186/s12935-020-1113-6
6. Atkins MB, Tannir NM. Current and emerging therapies for first-line treatment of metastatic clear cell renal cell carcinoma. Cancer Treat Rev. (2018) 70:127–37. doi: 10.1016/j.ctrv.2018.07.009
7. Motzer RJ, Jonasch E, Boyle S, Carlo MI, Manley B, Agarwal N, et al. NCCN guidelines insights: kidney cancer, version 1.2021. J Natl Compr Canc Netw. (2020) 18:1160–70. doi: 10.6004/jnccn.2020.0043
8. Xu S, Liu D, Chang T, Wen X, Ma S, Sun G, et al. Cuproptosis-associated lncRNA establishes new prognostic profile and predicts immunotherapy response in clear cell renal cell carcinoma. Front Genet. (2022) 13:938259. doi: 10.3389/fgene.2022.938259
9. Ljungberg B, Bensalah K, Canfield S, Dabestani S, Hofmann F, Hora M, et al. EAU guidelines on renal cell carcinoma: 2014 update. Eur Urol. (2015) 67:913–24. doi: 10.1016/j.eururo.2015.01.005
10. Hammers H. Immunotherapy in kidney cancer: the past, present, and future. Curr Opin Urol. (2016) 26:543–7. doi: 10.1097/MOU.0000000000000338
11. Du B, Zhou Y, Yi X, Zhao T, Tang C, Shen T, et al. Identification of immune-related cells and genes in tumor microenvironment of clear cell renal cell carcinoma. Front Oncol. (2020), 10:1770. doi: 10.3389/fonc.2020.01770
12. Motzer RJ, Escudier B, McDermott DF, George S, Hammers HJ, Srinivas S, et al. CheckMate 025 investigators. Nivolumab versus everolimus in advanced renal-cell carcinoma. N Engl J Med. (2015) 373:1803–13. doi: 10.1056/NEJMoa1510665
13. Errarte P, Guarch R, Pulido R, Blanco L, Nunes-Xavier CE, Beitia M, et al. The expression of fibroblast activation protein in clear cell renal cell carcinomas is associated with synchronous lymph node metastases. PloS One. (2016) 11:e0169105. doi: 10.1371/journal.pone.0169105
14. López JI, Errarte P, Erramuzpe A, Guarch R, Cortés JM, Angulo JC, et al. Fibroblast activation protein predicts prognosis in clear cell renal cell carcinoma. Hum Pathol. (2016) 54:100–5. doi: 10.1016/j.humpath.2016.03.009
15. López-Lago MA, Thodima VJ, Guttapalli A, Chan T, Heguy A, Molina AM, et al. Genomic deregulation during metastasis of renal cell carcinoma implements a myofibroblast-like program of gene expression. Cancer Res. (2010) 70:9682–92. doi: 10.1158/0008-5472.CAN-10-2279
16. Massari F, Ciccarese C, Santoni M, Lopez-Beltran A, Scarpelli M, Montironi R, et al. Targeting fibroblast growth factor receptor (FGFR) pathway in renal cell carcinoma. Expert Rev Anticancer Ther. (2015) 15:1367–9. doi: 10.1586/14737140.2015.1110488
17. Chen Y, McAndrews KM, Kalluri R. Clinical and therapeutic relevance of cancer-associated fibroblasts. Nat Rev Clin Oncol. (2021) 18:792–804. doi: 10.1038/s41571-021-00546-5
18. Chen LB, Zhu SP, Liu TP, Zhao H, Chen PF, Duan YJ, et al. Cancer associated fibroblasts promote renal cancer progression through a TDO/kyn/ahR dependent signaling pathway. Front Oncol. (2021) 11:628821. doi: 10.3389/fonc.2021.628821
19. Ambrosetti D, Coutts M, Paoli C, Durand M, Borchiellini D, Montemagno C, et al. CAFsin renal cell carcinoma: implication in prognosis and resistance to anti-angiogenic therapy. BJU Int. (2022) 129:80–92. doi: 10.1111/bju.15506
20. Errarte P, Larrinaga G, López JI. The role of CAFss in renal cell carcinoma. An example of tumor modulation through tumor/non-tumor cell interactions. J Adv Res. (2019), 21103–108. doi: 10.1016/j.jare.2019.09.004
21. Tao L, Huang G, Song H, Chen Y, Chen L. Cancer associated fibroblasts: An essential role in the tumor microenvironment. Oncol Lett. (2017) 14:2611–20. doi: 10.3892/ol.2017.6497
22. Erin N, Grahovac J, Brozovic A, Efferth T. Tumor microenvironment and epithelial mesenchymal transition as targets to overcome tumor multidrug resistance. Drug Resist Updat. (2020) 53:100715. doi: 10.1016/j.drup.2020.100715
23. Karagiannis GS, Poutahidis T, Erdman SE, Kirsch R, Riddell RH, Diamandis EP. Cancer-associated fibroblasts drive the progression of metastasis through both paracrine and mechanical pressure on cancer tissue. Mol Cancer Res. (2012) 10:1403–18. doi: 10.1158/1541-7786.MCR-12-0307
24. Kobayashi H, Gieniec KA, Wright JA, Wang T, Asai N, Mizutani Y, et al. The balance of stromal BMP signaling mediated by GREM1 and ISLR drives colorectal carcinogenesis. Gastroenterol. (2021) 160:1224–1239.e30. doi: 10.1053/j.gastro.2020.11.011
25. Xing J, Cai H, Lin Z, Zhao L, Xu H, Song Y, et al. Examining the function of macrophage oxidative stress response and immune system in glioblastoma multiforme through analysis of single-cell transcriptomics. Front Immunol. (2024) 14:1288137. doi: 10.3389/fimmu.2023.1288137
26. Zhao Z, Ding Y, Tran LJ, Chai G, Lin L. Innovative breakthroughs facilitated by single-cell multi-omics: manipulating natural killer cell functionality correlates with a novel subcategory of melanoma cells. Front Immunol. (2023) 14:1196892. doi: 10.3389/fimmu.2023.1196892
27. Lin Z, Sui X, Jiao W, Chen C, Zhang X, Zhao J. Mechanism investigation and experiment validation of capsaicin on uterine corpus endometrial carcinoma. Front Pharmacol. (2022) 13:953874. doi: 10.3389/fphar.2022.953874
28. Liu P, Deng X, Zhou H, Xie J, Kong Y, Zou Y, et al. Multi-omics analyses unravel DNA damage repair-related clusters in breast cancer with experimental validation. Front Immunol. (2023) 14:1297180. doi: 10.3389/fimmu.2023.1297180
29. Xie J, Deng W, Deng X, Liang JY, Tang Y, Huang J, et al. Single-cell histone chaperones patterns guide intercellular communication of tumor microenvironment that contribute to breast cancer metastases. Cancer Cell Int. (2023) 23:311. doi: 10.1186/s12935-023-03166-4
30. Wang J, Zuo Z, Yu Z, Chen Z, Meng X, Ma Z, et al. Single-cell transcriptome analysis revealing the intratumoral heterogeneity of ccRCC and validation of MT2A in pathogenesis. Funct Integr Genomics. (2023) 23:300. doi: 10.1007/s10142-023-01225-7
31. Ma B, Qin L, Sun Z, Wang J, Tran LJ, Zhang J, et al. The single-cell evolution trajectory presented different hypoxia heterogeneity to reveal the carcinogenesis of genes in clear cell renal cell carcinoma: Based on multiple omics and real experimental verification. Environ Toxicol. (2024) 39:869–81. doi: 10.1002/tox.24009
32. Ge Q, Zhao Z, Li X, Yang F, Zhang M, Hao Z, et al. Deciphering the suppressive immune microenvironment of prostate cancer based on CD4+ regulatory T cells: Implications for prognosis and therapy prediction. Clin Transl Med. (2024) 14:e1552. doi: 10.1002/ctm2.1552
33. Ding Y, Zhao Z, Cai H, Zhou Y, Chen H, Bai Y, et al. Single-cell sequencing analysis related to sphingolipid metabolism guides immunotherapy and prognosis of skin cutaneous melanoma. Front Immunol. (2023) 14:1304466. doi: 10.3389/fimmu.2023.1304466
34. Sun Z, Wang J, Fan Z, Yang Y, Meng X, Ma Z, et al. Investigating the prognostic role of lncRNAs associated with disulfidptosis-related genes in clear cell renal cell carcinoma. J Gene Med. (2024) 26:e3608. doi: 10.1002/jgm.3608
35. Sun Z, Wang J, Zhang Q, Meng X, Ma Z, Niu J, et al. Coordinating single-cell and bulk RNA-seq in deciphering the intratumoral immune landscape and prognostic stratification of prostate cancer patients. Environ Toxicol. (2024) 39:657–68. doi: 10.1002/tox.23928
36. Zhao J, Jiao W, Sui X, Zou J, Wang J, Lin Z. Construction of a prognostic model of luteolin for endometrial carcinoma. Am J Transl Res. (2023) 15:2122–39.
37. Lin Z, Sui X, Jiao W, Wang Y, Zhao J. Exploring the mechanism and experimental verification of puerarin in the treatment of endometrial carcinoma based on network pharmacology and bioinformatics analysis. BMC Complement Med Ther. (2022) 22:150. doi: 10.1186/s12906-022-03623-z
38. Lin Z, Fan W, Sui X, Wang J, Zhao J. Necroptosis-related lncRNA signatures for prognostic prediction in uterine corpora endometrial cancer. Reprod Sci. (2023) 30:576–89. doi: 10.1007/s43032-022-01023-9
39. Zhao J, Li D, Xie S, Deng X, Wen X, Li J, et al. Nomogram for predicting prognosis of patients with metastatic melanoma after immunotherapy: A Chinese population-based analysis. Front Immunol. (2022) 13:1083840. doi: 10.3389/fimmu.2022.1083840
40. Zou J, Lin Z, Jiao W, Chen J, Lin L, Zhang F, et al. A multi-omics-based investigation of the prognostic and immunological impact of necroptosis-related mRNA in patients with cervical squamous carcinoma and adenocarcinoma. Sci Rep. (2022) 12:16773. doi: 10.1038/s41598-022-20566-0
41. Wang Y, Zhao ZJ, Kang XR, Bian T, Shen ZM, Jiang Y, et al. lncRNA DLEU2 acts as a miR-181a sponge to regulate SEPP1 and inhibit skeletal muscle differentiation and regeneration. Aging (Albany NY). (2020) 12:24033–56. doi: 10.18632/aging.104095
42. Lin Z, Zou J, Sui X, Yao S, Lin L, Wang J, et al. Necroptosis-related lncRNA signature predicts prognosis and immune response for cervical squamous cell carcinoma and endocervical adenocarcinomas. Sci Rep. (2022) 12:16285. doi: 10.1038/s41598-022-20858-5
43. Zhao J, Zou J, Jiao W, Lin L, Wang J, Lin Z. Construction of N-7 methylguanine-related mRNA prognostic model in uterine corpus endometrial carcinoma based on multi-omics data and immune-related analysis. Sci Rep. (2022) 12:18813. doi: 10.1038/s41598-022-22879-6
44. Zhao ZJ, Zheng RZ, Wang XJ, Li TQ, Dong XH, Zhao CY, et al. Integrating lipidomics and transcriptomics reveals the crosstalk between oxidative stress and neuroinflammation in central nervous system demyelination. Front Aging Neurosci. (2022) 14:870957. doi: 10.3389/fnagi.2022.870957
45. Lin Z, Fan W, Yu X, Liu J, Liu P. Research into the mechanism of intervention of SanQi in endometriosis based on network pharmacology and molecular docking technology. Med (Baltimore). (2022) 101:e30021. doi: 10.1097/MD.0000000000030021
46. Smittenaar CR, Petersen KA, Stewart K, Moitt N. Cancer incidence and mortality projections in the UK until 2035. Br J Cancer. (2016) 115:1147–55. doi: 10.1038/bjc.2016.304
47. Turajlic S, Xu H, Litchfield K, Rowan A, Horswell S, Chambers T, et al. TRACERx renal consortium. Deterministic evolutionary trajectories influence primary tumor growth: TRACERx renal. Cell. (2018) 173:595–610.e11. doi: 10.1016/j.cell.2018.03.043
48. Zhang BY, Thompson RH, Lohse CM, Leibovich BC, Boorjian SA, Cheville JC, et al. A novel prognostic model for patients with sarcomatoid renal cell carcinoma. BJU Int. (2015) 115:405–11. doi: 10.1111/bju.12781
49. Gorrini C, Harris IS, Mak TW. Modulation of oxidative stress as an anticancer strategy. Nat Rev Drug Discov. (2013) 12:931–47. doi: 10.1038/nrd4002
50. Cairns RA, Harris IS, Mak TW. Regulation of cancer cell metabolism. Nat Rev Cancer. (2011) 11:85–95. doi: 10.1038/nrc2981
51. Rai P. Oxidation in the nucleotide pool, the DNA damage response and cellular senescence: Defective bricks build a defective house. Mutat Res. (2010) 703:71–81. doi: 10.1016/j.mrgentox.2010.07.010
52. Lee HC, Wei YH. Mitochondrial DNA instability and metabolic shift in human cancers. Int J Mol Sci. (2009) 10:674–701. doi: 10.3390/ijms10020674
53. Holmström KM, Finkel T. Cellular mechanisms and physiological consequences of redox-dependent signalling. Nat Rev Mol Cell Biol. (2014) 15:411–21. doi: 10.1038/nrm3801
54. D'Autréaux B, Toledano MB. ROS as signalling molecules: mechanisms that generate specificity in ROS homeostasis. Nat Rev Mol Cell Biol. (2007) 8:813–24. doi: 10.1038/nrm2256
55. Kuo CL, Chou HY, Chiu YC, Cheng AN, Fan CC, Chang YN, et al. Mitochondrial oxidative stress by Lon-PYCR1 maintains an immunosuppressive tumor microenvironment that promotes cancer progression and metastasis. Cancer Lett. (2020) 474:138–50. doi: 10.1016/j.canlet.2020.01.019
56. Kuo CL, Ponneri Babuharisankar A, Lin YC, Lien HW, Lo YK, Chou HY, et al. Mitochondrial oxidative stress in the tumor microenvironment and cancer immunoescape: foe or friend? J BioMed Sci. (2022) 29:74. doi: 10.1186/s12929-022-00859-2
57. Cheng AN, Cheng LC, Kuo CL, Lo YK, Chou HY, Chen CH, et al. Mitochondrial Lon-induced mtDNA leakage contributes to PD-L1-mediated immunoescape via STING-IFN signaling and extracellular vesicles. J Immunother Cancer. (2020) 8:e001372. doi: 10.1136/jitc-2020-001372
58. Chen G, Huang AC, Zhang W, Zhang G, Wu M, Xu W, et al. Exosomal PD-L1 contributes to immunosuppression and is associated with anti-PD-1 response. Nature. (2018) 560:382–6. doi: 10.1038/s41586-018-0392-8
59. Albiñana V, Recio-Poveda L, González-Peramato P, Martinez-Piñeiro L, Botella LM, Cuesta AM. Blockade of β2-adrenergic receptor reduces inflammation and oxidative stress in clear cell renal cell carcinoma. Int J Mol Sci. (2022) 23:1325. doi: 10.3390/ijms23031325
60. Costa A, Scholer-Dahirel A, Mechta-Grigoriou F. The role of reactive oxygen species and metabolism on cancer cells and their microenvironment. Semin Cancer Biol. (2014) 25:23–32. doi: 10.1016/j.semcancer.2013.12.007
61. Zhao ZJ, Zhao R. Purification and cloning of PZR, a binding protein and putative physiological substrate of tyrosine phosphatase SHP-2. J Biol Chem. (1998) 273:29367–72. doi: 10.1074/jbc.273.45.29367
62. Eminaga S, Bennett AM. Noonan syndrome-associated SHP-2/Ptpn11 mutants enhance SIRPalpha and PZR tyrosyl phosphorylation and promote adhesion-mediated ERK activation. J Biol Chem. (2008) 283:15328–38. doi: 10.1074/jbc.M801382200
63. Zhao R, Zhao ZJ. Identification of a variant form of PZR lacking immunoreceptor tyrosine-based inhibitory motifs. Biochem Biophys Res Commun. (2003) 303:1028–33. doi: 10.1016/s0006-291x(03)00484-4
64. Chen D, Cao L, Wang X. MPZL1 promotes tumor cell proliferation and migration via activation of Src kinase in ovarian cancer. Oncol Rep. (2019) 42:679–87. doi: 10.3892/or.2019.7199
65. Tan D, Zhang W, Tao Y, Galiya Y, Wang M. PZR promotes metastasis of colorectal cancer through increasing FAK and Src phosphorylation. Acta Biochim Biophys Sin (Shanghai). (2019) 51:356–64. doi: 10.1093/abbs/gmz019
66. Liu X, Huang J, Liu L, Liu R. MPZL1 is highly expressed in advanced gallbladder carcinoma and promotes the aggressive behavior of human gallbladder carcinoma GBC SD cells. Mol Med Rep. (2019) 20:2725–33. doi: 10.3892/mmr.2019.10506
67. Jia D, Jing Y, Zhang Z, Liu L, Ding J, Zhao F, et al. Amplification of MPZL1/PZR promotes tumor cell migration through Src-mediated phosphorylation of cortactin in hepatocellular carcinoma. Cell Res. (2014) 24:204–17. doi: 10.1038/cr.2013.158
68. Wang A, Fan F, Zang C, Li C, Hong H, Wang Z, et al. MPZL1 promotes lung adenocarcinoma progression by enhancing tumor proliferation, invasion, migration, and suppressing immune function via transforming growth factor-β1. Hum Gene Ther. (2023) 34:540–53. doi: 10.1089/hum.2022.234
69. Kavallaris M. Microtubules and resistance to tubulin-binding agents. Nat Rev Cancer. (2010) 10:194–204. doi: 10.1038/nrc2803
70. Mariani M, Karki R, Spennato M, Pandya D, He S, Andreoli M, et al. Class III β-tubulin in normal and cancer tissues. Gene. (2015) 563:109–14. doi: 10.1016/j.gene.2015.03.061
71. Borys F, Tobiasz P, Poterała M, Krawczyk H. Development of novel derivatives of stilbene and macrocyclic compounds as potent of anti-microtubule factors. BioMed Pharmacother. (2021) 133:110973. doi: 10.1016/j.biopha.2020.110973
72. Ferrandina G, Legge F, Salutari V, Paglia A, Testa A, Scambia G. Impact of pattern of recurrence on clinical outcome of ovarian cancer patients: clinical considerations. Eur J Cancer. (2006) 42:2296–302. doi: 10.1016/j.ejca.2006.03.025
73. Hetland TE, Hellesylt E, Flørenes VA, Tropé C, Davidson B, Kærn J. Class III β-tubulin expression in advanced-stage serous ovarian carcinoma effusions is associated with poor survival and primary chemoresistance. Hum Pathol. (2011) 42:1019–26. doi: 10.1016/j.humpath.2010.10.025
74. Koh Y, Jang B, Han SW, Kim TM, Oh DY, Lee SH, et al. Expression of class III beta-tubulin correlates with unfavorable survival outcome in patients with resected non-small cell lung cancer. J Thorac Oncol. (2010) 5:320–5. doi: 10.1097/JTO.0b013e3181ce684f
75. Duran GE, Wang YC, Francisco EB, Rose JC, Martinez FJ, Coller J, et al. Mechanisms of resistance to cabazitaxel. Mol Cancer Ther. (2015) 14:193–201. doi: 10.1158/1535-7163.MCT-14-0155
76. Smiyun G, Azarenko O, Miller H, Rifkind A, LaPointe NE, Wilson L, et al. βIII-tubulin enhances efficacy of cabazitaxel as compared with docetaxel. Cancer Chemother Pharmacol. (2017) 80:151–64. doi: 10.1007/s00280-017-3345-2
77. Reiman T, Lai R, Veillard AS, Paris E, Soria JC, Rosell R, et al. Cross-validation study of class III beta-tubulin as a predictive marker for benefit from adjuvant chemotherapy in resected non-small-cell lung cancer: analysis of four randomized trials. Ann Oncol. (2012) 23:86–93. doi: 10.1093/annonc/mdr033
78. Cicchillitti L, Penci R, Di Michele M, Filippetti F, Rotilio D, Donati MB, et al. Proteomic characterization of cytoskeletal and mitochondrial class III beta-tubulin. Mol Cancer Ther. (2008) 7:2070–9. doi: 10.1158/1535-7163.MCT-07-2370
79. Ploussard G, Terry S, Maillé P, Allory Y, Sirab N, Kheuang L, et al. Class III beta-tubulin expression predicts prostate tumor aggressiveness and patient response to docetaxel-based chemotherapy. Cancer Res. (2010) 70:9253–64. doi: 10.1158/0008-5472.CAN-10-1447
Keywords: HMGA1+ myofibroblasts, ccRCCs, ScRNA-seq, bulk RNA-seq, experiment validation
Citation: Zhou W, Lin Z and Tan W (2024) Deciphering the molecular landscape: integrating single-cell transcriptomics to unravel myofibroblast dynamics and therapeutic targets in clear cell renal cell carcinomas. Front. Immunol. 15:1374931. doi: 10.3389/fimmu.2024.1374931
Received: 23 January 2024; Accepted: 05 March 2024;
Published: 18 March 2024.
Edited by:
Chen Yang, German Cancer Research Center (DKFZ), GermanyReviewed by:
Zhou Sun, Jilin University, ChinaXinpei Deng, Sun Yat-sen University Cancer Center (SYSUCC), China
Copyright © 2024 Zhou, Lin and Tan. 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: Wang Tan, OTgxNjUzNTk4QHFxLmNvbQ==
†These authors have contributed equally to this work and share first authorship