Yanlin Tang1,2†
Chujin Ye2†
Jiayi Zeng2,3Ping Zhu4Shouyu Cheng2,5Weinan Zeng1,2Bowen Yang2,6
Yanjun Liu4*
Yuming Yu2,3*- 1Shantou University Medical College, Shantou, China
- 2Department of Urology, Guangdong Provincial People's Hospital (Guangdong Academy of Medical Sciences), Southern Medical University, Guangzhou, China
- 3The Second School of Clinical Medicine, Southern Medical University, Guangzhou, China
- 4Department of Immunology, School of Basic Medical Science, Southern Medical University, Guangzhou, China
- 5School of Medicine, South China University of Technology, Guangzhou, China
- 6Guangdong Cardiovascular Institute, Guangdong Provincial People's Hospital, Guangdong Academy of Medical Sciences, Guangzhou, China
Clear cell renal cell carcinoma (ccRCC) belongs to one of the 10 most frequently diagnosed cancers worldwide and has a poor prognosis at the advanced stage. Although multiple therapeutic agents have been proven to be curative in ccRCC, their clinical application was limited due to the lack of reliable biomarkers. Considering the important role of basement membrane (BM) in tumor metastasis and TME regulation, we investigated the expression of BM-related genes in ccRCC and identified prognostic BM genes through differentially expression analysis and univariate cox regression analysis. Then, BM-related ccRCC subtypes were recognized through consensus non-negative matrix factorization based on the prognostic BM genes and evaluated with regard to clinical and TME features. Next, utilizing the differentially expressed genes between the BM-related subtypes, a risk scoring system BMRS was established after serial analysis of univariate cox regression analysis, lasso regression analysis, and multivariate cox regression analysis. Time-dependent ROC curve revealed the satisfactory prognosis predictive capacity of BMRS with internal, and external validation. Multivariate analysis proved the independent predictive ability of BMRS and a BMRS-based nomogram was constructed for clinical application. Some featured mutants were discovered through genomic analysis of the BMRS risk groups. Meanwhile, the BMRS groups were found to have distinct immune scores, immune cell infiltration levels, and immune-related functions. Moreover, with the help of data from The Cancer Immunome Atlas (TCIA) and Genomics of Drug Sensitivity in Cancer (GDSC), the potential of BMRS in predicting therapeutic response was evaluated and some possible therapeutic compounds were proposed through ConnectivityMap (CMap). For the practicability of BMRS, we validated the expression of BMRS-related genes in clinical samples. After all, we identified BM-related ccRCC subtypes with distinct clinical and TME features and constructed a risk scoring system for the prediction of prognosis, therapeutic responses, and potential therapeutic agents of ccRCC. As ccRCC systemic therapy continues to evolve, the risk scoring system BMRS we reported may assist in individualized medication administration.
1 Introduction
Globally, more than 430,000 individuals suffering from kidney cancers were newly diagnosed in 2020 and approximately 180,000 people died from this type of cancer (Sung et al., 2021). As one of the major subtypes, clear cell renal cell carcinoma (ccRCC) could earn a good prognosis when treated at an early stage. However, around one-third of ccRCC patients were found to be in the metastatic stage, requiring systemic therapy other than radical surgery (Jonasch et al., 2014). Although novel treatments including immunotherapies and targeted therapies were demonstrated to be curative in this chemoresistant cancer type, their clinical effects were uncontrollable due to the lack of predictive biomarkers for the therapeutic response and adverse events (Jhaveri and Perazella, 2018; Motzer et al., 2020a). Besides, it was demonstrated that combined therapy of immune checkpoint inhibitors (ICIs) and vascular endothelial growth factor (VEGF) tyrosine kinase inhibitors (TKIs) could exert better curative effects than monotherapy, leading to a requirement for more individualized markers for treatment selection (Amin and Hammers, 2018). Therefore, reliable predictive biomarkers ought to be developed for the prognosis and therapeutic response of ccRCC.
The tumor microenvironment (TME) is a complex ecosystem including immune cells, stromal cells, and extracellular matrix, surrounding and interacting with tumor cells (Blankenstein et al., 2012). As a highly immune-infiltrated cancer type, ccRCC cells were able to modulate the TME including immune cells for evasion of anti-cancer immunity through multiple mechanisms (Díaz-Montero et al., 2020). The understanding of these mechanisms could facilitate the application of cancer-specific therapies, such as ICIs, which restored anti-cancer immunity through interrupting the suppressive signals from the ccRCC cells (Motzer et al., 2015). Besides, accumulating evidence indicated that molecular classification of ccRCC into groups with distinct TME features could distinguish their prognosis and therapeutic response (de Velasco et al., 2017). Thus, it would be valuable to investigate the TME in ccRCC for the discovery of novel predictive biomarkers.
Basement membrane (BM) is a thin sheet of extracellular matrix (ECM) lining beneath endothelial and epithelial tissues, mainly composed of collagen IV and laminin (Yurchenco, 2011). It serves as one of the barriers preventing cancer cells from invasion, but its remodeling and stiffness would contribute to the metastasis of tumor (Chang and Chaudhuri, 2019). Some BM-related genes were revealed to be associated with the prognosis of RCC. Wragg et al. (2016) demonstrated that the high expression of LAMA4, a laminin component, was correlated with the poor prognosis of RCC (Wragg et al., 2016; Ho et al., 2017). Moreover, BM could mediate the signal transduction between the microenvironment and cells. As a major component of BM, laminin was revealed to have the ability to modulate the migration, activation and functionality of T lymphocytes within tumors (Liu et al., 2022). With these concerns, investigating the BM-related genes in ccRCC may assist in understanding the relationship between ccRCC and TME and developing predictive biomarkers.
In the current study (Figure 1), we investigated the expression of BM-related genes in ccRCC and used them to classify ccRCC patients into distinct subtypes, based on which a risk scoring system, BMRS, was established. Comprehensive analyses were conducted to evaluate the capacity of BMRS in distinguishing the prognosis, TME features, and therapeutic response of ccRCC. In this way, we constructed a gene-based BMRS for prognosis and treatment prediction of ccRCC and provided molecular candidates as novel therapeutic targets.
FIGURE 1. The flow chart of the current study.
2 Materials and methods
2.1 Data acquisition
A ccRCC cohort including 539 tumor samples and 72 normal samples from the KIRC project of The Cancer Genome Atlas (TCGA) was selected and its RNA expression data, somatic mutation data, and the corresponding clinical data were extracted from the Genomic Data Commons Data Portal (https://portal.gdc.cancer.gov). ArrayExpress (https://www.ebi.ac.uk/arrayexpress) is another public database containing high-throughput genomic data from more than 75,000 experiments. Gene expression data of a ccRCC cohort, E-MTAB-1980, with 101 tumor samples were acquired from ArrayExpress and the survival information was obtained from previous research (Li et al., 2018). The expression data of both cohorts were transformed into a data format of transcripts per million (TPM) for better analysis. All the data were publicly available and no ethical consent was required.
2.2 Investigation of the basement membrane genes in ccRCC
Jayadev et al. discovered 222 protein-coding genes that were related to BM and human health (Jayadev et al., 2022). The expression of these genes was extracted from the TCGA cohort and differentially analyzed between tumor and normal samples with the help of the R package “limma” to identify the differentially expressed BM genes (Ritchie et al., 2015). Thereafter, univariate cox regression analysis facilitated the selection of the prognostic BM genes from those differentially expressed genes. For a better understanding of the expression of the prognostic BM genes in ccRCC, correlation analyses were performed to reveal their interrelationship. Besides, utilizing the R package “maftool” (Mayakonda et al., 2018), the variation in these genes was depicted including both somatic mutation status and copy number variations (CNV).
2.3 Discovery and investigation of basement membrane-related clusters in ccRCC
Consensus non-negative matrix factorization (CNMF) is a powerful method for the dimension reduction of genomic data to discover distinguished molecular patterns. Through the R package ‘CancerSubtypes’ (Xu et al., 2017), the expression data of the prognostic BM genes were used to classify the ccRCC samples in the TCGA cohort into 2, 3, and 4 clusters. Silhouette width is a measurement for the evaluation of the classifications and a larger average silhouette width means more valuable subtypes. Clusters with the highest average silhouette width were selected for the following analyses. To further investigate the value of these BM-related clusters, Kaplan-Meier survival analysis and chi-square test were conducted to visualize the clinical difference between them. Meanwhile, though GSEA 4.1.0, gene sets enrichment analysis (GSEA) together with gene ontology (GO) gene sets assisted the understanding of the molecular functions enriched in each cluster. Furthermore, Estimation of STromal and Immune cells in MAlignant Tumor tissues using Expression data (ESTIMATE) (Yoshihara et al., 2013) was utilized to quantify the immune and stromal status of ccRCC samples and the resulting scores were compared between each cluster to discover whether there were differences in their TME. Besides, we calculated the infiltrating levels of 23 tumor infiltrating immune cells (TIIC) through single sample GSEA (ssGSEA) and 14 stoma cells by xCell (https://xcell.ucsf.edu/). By comparing the infiltrating level of these cells between different clusters, their TME differences could be better understood.
2.4 Construction of a gene expression-based risk score system
After discovering the distinct features between the BM-related subtypes, we tended to construct a risk score system based on these clusters. First, the differentially expressed genes between the clusters were identified. Then, the TCGA cohort was randomly divided into a training cohort and a testing cohort in a ratio of 7 to 3. In the training cohort, expression of the differentially expressed genes was extracted and used for univariate cox regression analysis to select the prognostic genes. Through R package ‘glmnet’ (Friedman et al., 2010), Least Absolute Selection and Shrinkage Operator (LASSO) regression analysis facilitated identifying those genes with higher association with ccRCC prognosis, and multivariate cox regression analysis was utilized to choose the best genes for construction of the BM-related risk score system (BMRS): risk score = coefficient 1*gene expression 1 + … + coefficient n * gene expression n (1 to n represent each prognostic genes).
2.5 Analysis of the prognostic predictive capacity of BMRS
To investigate the predictive ability of BMRS, the TCGA training cohort, TCGA testing cohort, TCGA cohort, and E-MTAB-1980 cohort were respectively divided into high and low risk groups according to the median risk scores. Kaplan-Meier analysis was applied in each pair of high and low risk groups to reveal the relationship between risk score and the overall survival (OS) of ccRCC patients. Subsequently, with the help of time-dependent Receiver Operating Characteristic (tROC) curves, the predictive power of the risk scoring system for the 1-, 3-, and 5-year OS of each cohort were illustrated and the corresponding values of area under curve (AUC) were calculated. Additionally, the clinical features (age, gender, Fuhrmann grade, and AJCC stage) in high and low risk groups of the TCGA cohort were compared through chi-square test. Meanwhile, the ability of BMRS in differentiating ccRCC OS was tested under different clinical statuses by Kaplan-Meier survival analysis.
2.6 Establishment of a clinical predictive nomogram
For better application in the clinic, BMRS (high risk vs. low risk) combined with age (<65-year-old vs. >=65-year-old), gender (female vs. male), Fuhrmann grade (Grade 1/2 vs. Grade 3/4), AJCC stage (Stage I/II vs. Stage III/IV), T stage (T1/2 vs. T3/4), N stage (N0 vs. N1), and M stage (M0 vs. M1) were incorporated into univariate and multivariate cox regression analyses to demonstrate its predictive value. Furthermore, the resulting values of the multivariate cox regression analysis were used to construct a nomogram for the prediction of ccRCC survival. The predictive capability of the nomogram for 3-and 5-year survival was assessed by calibration plot. Besides, decision curve analysis (DCA) was utilized to compare the net benefit of BMRS, AJCC stage and two previous BM-related gene signatures (Zhou et al., 2022a; Xiong et al., 2022) in predicting 3- and 5-year survival of ccRCC. TCGA cohort was separated into high and low risk groups according to the median of nomogram scores and was analyzed to identify their survival difference.
2.7 Investigation of the genomic variation in different BMRS groups
The genomic features in high and low BMRS groups were illustrated and compared through the R package ‘maftool’. The variation status of the top 20 mutation genes in both groups was depicted and all the mutated genes were compared between groups to identify the group-specific mutation. At the same time, by applying pairwise Fisher’s exact test between every two genes, we wanted to discover whether there were some exclusive or co-occurrence gene pairs in the high BMRS group. Furthermore, some tumor heterogeneity-related features, including single nucleotide variation (SNV), homologous recombination defects (HRD), cancer testis antigen (CTA), and intratumor heterogeneity (ITH) were introduced from previous research (Thorsson et al., 2018) and compared between different BMRS groups.
2.8 Analysis of the immune landscape related to BMRS
In order to investigate the BMRS-related biological processes, cellular components, and molecular functions, the genes that were differentially upregulated in the high BMRS group were chosen and incorporated into GO functional analysis through the R package ‘clusterProfiler’ (Wu et al., 2021). Concerning the immune microenvironment in ccRCC, the immune scores calculated from ESTIMATE and the infiltrating scores of 23 TIICs were compared between BMRS groups. Each TIIC was analyzed through Kaplan-Meier survival analysis to discover its relationship with ccRCC prognosis. Meanwhile, ssGSEA algorithm was utilized to induce scores representing immune suppression and some immune-related functions with the help of previously published gene signatures (Yi et al., 2020). In addition, a publicly accessible website called Tumor Immune Dysfunction and Exclusion (TIDE, http://tide.dfci.harvard.edu) which provides a platform for estimation of T cells dysfunction scores based on gene expression data, was adopted to compare the T cells status between high and low BMRS clusters.
2.9 Exploration of the therapeutic predictive potential of BMRS
The gene expression of several popular inhibitory immune checkpoints was extracted and differentially analyzed between high and low BMRS groups. Additionally, the immunophenoscores, which are estimated scores of immunotherapeutic response to ICIs for TCGA samples, were obtained from The Cancer Immunome Atlas (TCIA, https://tcia.at/home) and compared between BMRS groups. Except for immunotherapeutic response, sensitivity toward some therapeutic drugs (Axitinib, Pazopanib, Sorafenib, Sunitinib, and Mitomycin C) clinically used for ccRCC patients were predicted through R package ‘pRRophetic’ (Barbour et al., 2014), which connects a large amount of gene expression and drug sensitivity data in Genomics of Drug Sensitivity in Cancer (GDSC, https://www.cancerrxgene.org). Furthermore, the gene expression of immune chemokines, immune receptors, and major histocompatibility complex (MHC) molecules was correlated with BMRS with the concern of discovering potential therapeutic targets. Besides, we uploaded the top 150 upregulated genes in high and low BMRS groups respectively to ConnectivityMap (CMap, https://clue.io) for exploration of potential therapeutic compounds for patients in each group.
2.10 Transcriptome sequencing analysis
To further prove the applicability of BMRS, we collected clinical ccRCC and adjacent normal samples (18 ccRCC and 6 adjacent normal samples) for transcriptome sequencing analysis. This project was supported by the hospital ethics committee and consent was acquired from all the patients. In line with the protocol of the manufacturer, each ccRCC and adjacent normal samples underwent paired-end sequencing on the NovaSeq 6000 high-throughput sequencing platform (Illumina, United States) to remove sequencing reads containing aptamer sequences and low-quality reads as well as bases. Then, high-quality pairwise reads were aligned to the human genome GRCh38 through HISAT2 (v2.1.1), generating BAM files. BAM files were arranged by samtools (v1.15.1) and then counted with the help of Subread (v2.0.1). Raw counts of transcripts per gene were converted to the format of TPM, allowing better analysis of gene expression between samples. Thereafter, the genes enrolled in BMRS were differentially analyzed between the normal and ccRCC samples.
2.11 Statistical analysis
All the analyses process in the current study were achieved through the usage of R 4.1.0 and R studio Desktop 2022.07.1 + 554. The graphs displayed were drawn by R studio Desktop 2022.07.1 + 554 and Adobe Illustrator CS6 (64 Bit). During differentially expression analysis, genes with the absolute value of log fold change (logFC) more than were selected. p-value less than 0.05 was regarded as significant for all analyses.
3 Results
3.1 Identification of the differentially expressed and prognostic BM genes
Differential expression analysis of the 222 BM genes between 539 ccRCC samples and 72 normal samples revealed that there were 106 differentially expressed BM genes, in which 39 genes were downregulated and 67 genes were upregulated (Supplementary Table S1). Then, after univariate cox regression analysis, 49 BM genes were demonstrated to be prognostic including 26 protective genes and 23 risk genes (Supplementary Table S2). A substantial positive correlation existed between these prognostic BM genes indicating that they were highly interconnected (Figure 2A). As for the genomic variation that happened in these genes, 25.6% of the ccRCC samples possessed prognostic BM gene alterations and more of these alterations were missense mutations (Figure 2B). The top 3 mutated genes were HMCN1 (5%, risk gene), COL6A3 (3%, protective gene), and COL4A5 (2%, protective gene). In addition, the CNV of these genes was analyzed and a relatively low frequency of CNV was discovered, except for SPARC (21.0%, gain of function), and TGFBI (19.8%, gain of function) (Figures 2C, D).
FIGURE 2. Genetic analyses of the BM-related genes. (A) The inter-relationship between each BM-related gene and the role of these genes in the prognosis of ccRCC. (B) The waterfall plot depicting the mutation of the prognostic BM genes. (C) The CNV frequency of the prognostic BM genes. (D) The distribution of the prognostic BM genes in each chromosome with CNV gain colored red, CNV loss colored blue and no CNV colored black.
3.2 Recognition of BM-based clusters of ccRCC
The 49 differentially expressed and prognostic BM genes were further used for the classification of ccRCC patients through CNMF algorithm. As the results show, the samples could be classified into two separate clusters A and B, in which cluster A had significantly lower OS than cluster B and the average silhouette width was 0.94 (Figures 3A–C). For the reason that the average silhouette widths were relatively low, classifications of 3, 4, and 5 clusters were not under consideration (Supplementary Figure S1). By comparing the clinical features in both clusters, we discovered that the ccRCC samples in cluster A had more advanced features including Fuhrmann grade, AJCC stage, T stage, N stage, and M stage. There were more males in cluster A than in cluster B while the age distributions in both clusters were similar (Figure 3D). The differentially expressed analysis also indicated that some metastasis-related genes, such as MMP13 and ROS1 were upregulated in cluster A. The following GSEA analysis demonstrated that basement membrane-related functions like collagen fibril organization and collagen catabolic process were more prominent in cluster A. Meanwhile, some immune-related biological processes were notably enriched in cluster A, such as positive regulation of T helper 1 type immune response and positive regulation of interleukin 17 production (Figures 3E, F). Considering these enriched immune functions, we further investigated the microenvironment components in both clusters through ESTIMATE and ssGSEA. It was demonstrated that the stromal score, immune score and ESTIMATE score in cluster A were higher than those in cluster B, indicating the rich stromal, immune components and low tumor purity in cluster A (Figure 3G). At the same time, cluster A had a higher level of nearly every type of infiltrating immune cell than cluster B except for eosinophils, neutrophils, and plasmacytoid dendritic cell (Figure 3H). As for stromal cells, more adipocytes, chondrocytes, fibroblasts, and mesangial cells were found in cluster A while cluster B possessed more endothelial cells, lymphatic endothelial cells and microvascular endothelial cells (Figure 3I). Therefore, the prognostic BM genes could divide ccRCC samples into two clusters with distinct clinical and TME features.
FIGURE 3. BM gene-based classification of ccRCC patients. (A) Survival analysis of the two ccRCC clusters classified based on BM genes. (B) The sample similarity matrix plot of the two identified clusters. (C) The Silhouette width plot of the two BM-related subtypes. (D) The heatmap depicting the distinct clinical features and gene expression between cluster (A) and cluster (B). (E, F). The biological processes that enriched in cluster (A). (G) Comparison of the stromal score, immune score and ESTIMATE score between cluster (A) and cluster (B). (H, I) The different infiltration levels of immune cells and stromal cells between the two BM-related subtypes.
3.3 Construction of BMRS risk scoring system for prognosis prediction
Since the BM gene-based ccRCC subtypes could discriminate prognostic and clinical features, we would like to derive a more applicable risk scoring system from these subtypes. In the TCGA training cohort, the differentially expressed genes between clusters A and B were identified (Supplementary Figure S2A). Based on these genes, univariate cox regression analysis was conducted and induced 414 prognostic genes (Supplementary Table S3). Subsequently, LASSO algorithm was applied to obtain 12 genes that were significantly associated with the prognosis of ccRCC (Supplementary Figure S2B, C). Then, these 12 genes were incorporated into multivariate cox regression analysis (Figure 4A), inducing a risk scoring system, BMRS, based on 7 genes, CCDC85A, AJAP1, ANK3, P4HA3, C8G, ADAM8, HJURP.
FIGURE 4. Construction and evaluation of the risk scoring system BMRS. (A) The forest plot illustrating the 7 genes resulted from multivariate cox regression analysis. (B, C) Distribution of the risk scores and survival status of ccRCC patients in the TCGA training cohort. (D, F) The Kaplan Meier survival curve comparing the overall survival between high and low risk groups in the TCGA training cohort and E-MTAB-1980 cohort respectively. (E, G) The tROC curve of the risk score for 1-, 3-, and 5-year survival in the TCGA training cohort and E-MTAB-1980 cohort. (H) Comparison of the clinical features in high and low risk groups.
3.4 Evaluation of the predictive ability of BMRS
For the evaluation of BMRS system, we calculated risk scores for each sample in TCGA training cohort and divided them into high and low risk groups according to the median risk scores (Figure 4B). It could be recognized that samples with high risk were more frequently dead than those with low risk (Figure 4C). At the same time, survival analysis demonstrated that the high risk group had longer OS than the low risk group (Figure 4D). The results of tROC curve proved that BMRS system possessed a great prognostic capacity for 1-year (AUC = 0.813), 3-year (AUC = 0.753), and 5-year (AUC = 0.777) survival (Figure 4E). The same analyses were applied in the TCGA testing cohort and TCGA cohort as an internal validation and delivered similar results (Supplementary Figure S2D–K). As an external validation, the survival analysis result and tROC curve derived from E-MTAB-1980 also supported that BMRS could discriminate the prognosis of ccRCC and had high accuracy (Figures 4F, G, 1-year AUC = 0.867, 3-year AUC = 0.832, and 5-year AUC = 0.843). Furthermore, survival analyses of high and low risk groups in the TCGA cohort under various clinical situations were conducted and revealed that BMRS could discriminate ccRCC survival in most situations (Supplementary Figure S3). Meanwhile, the Fuhrmann grade, AJCC stage, T stage, N stage, and M stage in the high risk group were notably higher than those in the low risk group, which was in line with the poor prognosis in the high risk group (Figure 4H). Moreover, we evaluated the gene expression of the genes used for BMRS and the results demonstrated that most of the genes were differentially expressed in ccRCC (Supplementary Table S4), consistent with the above analysis. In this way, a reliable risk scoring system BMRS was constructed with considerable predictive capacity.
3.5 Establishment of a clinical nomogram based on BMRS
For better application of the BMRS in clinic, we evaluated its predictive ability with the consideration of the clinical variables. Analyzing risk score and clinical variables cooperatively, univariate cox regression analysis revealed that only age and gender could not predict ccRCC prognosis (Figure 5A). The following multivariate cox regression analysis uncovered the independent prognostic predictive ability of BMRS, T stage and M stage (Figure 5B). Thereafter, the Fuhrmann grade, T stage, M stage and BMRS risk were incorporated to establish a nomogram (Figure 5C). The calibration curve depicted the 3-year and 5-year OS predicted by the nomogram had a satisfactory consistency with those observed in clinic (Figure 5D). While compared with the clinically popular TNM stage, Zhou’s and Xi’s gene signatures, the nomogram could deliver higher net benefit in 3-year (Figure 5E) and 5-year OS prediction (Figure 5F). The survival analysis also demonstrated that samples with high points had significantly lower OS than those with low points (Figure 5G). Thus, BMRS could not only independently predict ccRCC survival but also assist in establishing more competitive predictive methods.
FIGURE 5. Establishment and evaluation of the clinical predictive nomogram. (A, B) The forest plots depicting the results of univariate cox regression analysis and multivariate cox regression analysis respectively. (C) The nomogram based on grade, T stage, M stage, and risk score. The red dot lines represented an example of the total point calculation. (D) The calibration curve testing predictability of the nomogram for 3- and 5-year survival. (E, F) The DCA curves comparing the net benefit of nomogram, TNM stage, Zhou’s and Xi’s gene signatures for prediction of 3- and 5-year survival respectively. (G) The Kaplan Meier survival curve comparing the overall survival of ccRCC patients with high and low points.
3.6 The genomic variation in different risk groups
After discovering the clinical significance of BMRS, further analyses of the underlying genomic alterations were under concern. In both risk groups, the mutation of VHL (44%), PBRM1 (38%), TTN (13%), and SETD2 (11%) comprised the major part of all gene alterations (Figure 6A). The following comparison revealed that 10 mutants were occurring more frequently in the high risk group (Figure 6B). Among these mutants, SETD2 was the most frequent and significant one and most of its mutation happened in the non-coding area (Figure 6C). Besides, the discovery of the relationship between gene variations uncovered a notable co-occurance of SETD2 and PBRM1 mutation. At the same time, the mutation of MUC16 and BAP1 also exhibited a significant positive inter-relationship while MUC17 and VHL mutation were mutually exclusive (Figure 6D). Apart from these mutational differences, the high risk group was demonstrated to have higher scores of SNV antigens and HRD, indicating the high genomic heterogeneity in patients with high BMRS risk (Figures 6E, F). Meanwhile, the CTA and ITH scores in the high risk group were significantly higher than those in the low risk group (Figures 6G, H).
FIGURE 6. Genomic analysis in accordance with the BMRS risk groups. (A) The waterfall plot representing the top 20 mutated genes in both risk groups. (B) The mutants that significantly differentiated between high and low risk groups. (C) Exhibition of the mutation sites and types of SETD2 in both risk groups. (D) The correlation of the mutants in the high risk group. (E–H) Comparison of SNV neoantigens, HRD, CTA score, and ITH between high and low risk groups.
3.7 The association between BMRS and the immune landscape of ccRCC
Considering that BMRS was derived from the BM-related clusters with distinct functional and immunal features, we also investigated the functions associated with BMRS. Similar to the above results, the high risk group was enriched with BM-related functions, such as extracellular structure organization and extracellular matrix organization. Meanwhile, there were some immune involvements recognized in the high risk group including humoral immune response and acute inflammatory response (Figure 7A). The comparison of immune scores between the two risk groups also indicated that the high risk group had a higher immune level than the low risk group (Figure 7B). For a better understanding of this immune involvement, we further evaluated the infiltrating immune cells in each group. As the results showed, most of the immune cells were more frequently infiltrated in the high risk group than in the low risk group except for eosinophil and neutrophil whose infiltrating levels were lower in the high risk group (Figure 7C). Survival analyses were conducted for each immune cell and it was revealed that some innate immune cells like neutrophil, and mast cell were protective cells. Most of the immune cells enriched in the high risk group were associated with poor survival, including suppressive immune cells (e.g., MDSC) and effector immune cells (e.g., activated CD8 T cell). For an exploration of this connection between abundant immune infiltration and poor survival, the functionality of the immune environment was assessed, delivering a result that ccRCC samples with high risk scored higher on immune suppression and T cell dysfunction than those with low risk (Figures 7D, E). Moreover, although high-risk samples had a higher level of most immune functions like APC co-stimulation and T cell co-stimulation, their scores of some negative immune functions such as APC co-inhibition and T cell co-inhibition were also greater than low-risk samples (Figure 7F). Therefore, BMRS may have a relationship with the suppressive immune microenvironment of ccRCC.
FIGURE 7. The immune landscape in different BMRS risk groups. (A) The enriched GO functions in the high risk group. (B) The different cumulative distributions of the immune score in different risk groups. (C) The distinct infiltration level of the immune cells in high and low risk groups. (D, E) Comparison of the scores of immune suppression and T cell dysfunction between the two BMRS groups. (F) The difference in immune functions between high and low risk groups.
3.8 The therapeutic predictive potential of BMRS
After the investigation of BMRS and the immune microenvironment, the potential relationship between BMRS and immunotherapies was under consideration. Differential analysis of the gene expression of inhibitory immune checkpoints showed that PD1, CTLA4, LAG3, and TIGIT were significantly upregulated in the high risk group (Figure 8A). Data on immunotherapeutic response predicted by TCIA demonstrated that high-risk ccRCC patients possessed a higher level of IPS to CTLA4 inhibitor and combination of PD1 and CTLA4 inhibitors than low-risk patients, indicating high-risk patients may respond better toward these two strategies of immunotherapies (Figure 8B). Apart from immunotherapies, we also correlated BMRS with some therapeutic drugs used in clinic. The half maximal inhibitory concentration (IC50) of sorafenib, sunitinib, and mitomycin C. were notably lower in the high risk group than in low risk group, suggesting that ccRCC patients with high risk may have better outcomes receiving these drugs (Figure 8C).
FIGURE 8. The potential of BMRS in guiding ccRCC treatment. (A) Comparison of the inhibitory immune checkpoint expression in different risk groups. (B) The IPS difference of4 types of ICI therapies in high and low BMRS groups. (C) Differential analysis of the IC50 of 5 therapeutic agents between BMRS risk groups. (D) Correlation between BMRS risk scores and the expression level of immune checkpoints, immune receptors and MHC molecules. (E) Potential therapeutic compounds predicted for high and low risk groups respectively.
To discover some additional therapeutic targets, we analyzed the correlation between BMRS and some immune-related molecules (Figure 8D). There were multiple immune checkpoints (e.g., CXCL13, and CCL25), immune receptors (e.g., CXCR5, and CCR6) and MHC molecules (e.g., HLA-E) significantly related to BMRS, providing potential targets for therapeutic investigation. Moreover, as predicted by CMap (Figure 8E), molecular compounds like BRD-K23875128, which is Rho associated kinase inhibitor, could serve as therapeutic drugs for high-risk ccRCC patients. Some other compounds were also estimated to be effective in low-risk ccRCC patients, for example, PF-04217903, which is a kind of arginase inhibitor.
4 Discussion
With the progression of ccRCC research, more treatment options were available to improve the prognosis of ccRCC patients (Choueiri and Motzer, 2017). However, the existing predictive markers for the prognosis and therapeutic response of ccRCC had limited functions due to the heterogeneity of this type of patient (Greef and Eisen, 2016). As a component of ECM, BM serves as not only a physical barrier against tumor invasion and metastasis but also a mediator of signals between microenvironment and cells (Bissell and Hines, 2011), indicating its potential as a target for the investigation of novel biomarkers. Through analyses of the expression of BM-related genes, we discovered two subtypes with distinct clinical and microenvironmental features and established a 7-gene risk scoring system for the prediction of prognosis and treatment response of ccRCC patients.
According to statistics, the prognosis of ccRCC was largely influenced by the occurrence of distant metastasis, which could decrease the 5-year survival rate to less than 20% (Siegel et al., 2022). In line with this finding, our results demonstrated that BM-related genes could divide ccRCC patients into two clusters with distinct prognosis features. This could be explained by the defensive role of BM in preventing the tumor cells from invading the stroma (Valastyan and Weinberg, 2011). Matrix metalloproteinases (MMPs), a series of zinc-dependent proteinases, played an important role in the degradation of BM to weaken its barrier capacity. Without MMPs, the tumor cells could hardly squeeze through the nanosized pores in the BM (Wolf et al., 2013; Eatemadi et al., 2017). It was demonstrated that the increased expression of MMPs was correlated with the invasion and poor prognosis of carcinomas (Winer et al., 2018), and in our study, the upregulation of MMP13 was identified in cluster A. Moreover, the high infiltration of fibroblasts in cluster A could also be associated with the expression of MMPs. It was reported that tumor cells induced the secretion of MMPs from cancer-associated fibroblasts (CAFs) to modulate BM (Kessenbrock et al., 2010). Conversely, MMPs such as MMP1 secreted from tumor cells promote the transdifferentiation of fibroblasts into CAFs (Heneberg, 2016). In addition to the chemical changes, BM could also be modulated into an invasion-favored status. A recent study revealed that in those BM with high plasticity, cells could mechanically enlarge the nanosized pores and migrate through BM without the help of MMPs (Wisdom et al., 2018). Meanwhile, the stiffness of BM could influence tumor invasion through a protein called netrin-4 (Net4). Net4 in BM mechanistically bound to laminin and diluted laminin ternary node complex, softening the BM and making it more resistant to tumor cell invasion (Reuten et al., 2021). Therefore, the different statuses of BM could influence the prognosis of ccRCC and may facilitate the survival prediction.
Based on the BM-related clusters, we constructed a risk scoring system called BMRS. Genes included in BMRS were found to be correlated with the progression and metastasis of cancers including both risk and protective genes. Holliday junction recognition protein (HJURP) is a kind of centromeric protein, being essential for the stimulation of chromosome division and cell mitosis (Zhang et al., 2021). Its upregulation was correlated with an increased invasion and migration capacity of cancer (Chen et al., 2019). Prolyl 4-hydroxylase alpha subunit 3 (P4HA3) is also a risk gene that could strengthen the motility and invasiveness of tumors (Song et al., 2018). As an enzymic subunit of prolyl 4-hydroxylase, P4HA3 was critical for the stability of collagen and its dysregulation would directly activate the invasiveness potential of tumor cells (Nakasuka et al., 2021). Meanwhile, a disintegrin and metalloproteinase 8 (ADAM8) was another gene related to enzyme production and could promote tumor metastasis probably through the degradation of ECM components (Conrad et al., 2019). Apart from risk genes, BMRS also contained protective genes including adherens junctional associated protein-1 (AJAP1) and ankyrin G (ANK3). AJAP1 belonged to multi-protein complexes named adherens junction, which is critical for cell adhesion and growth inhibition (Xu et al., 2019; Zhou et al., 2022b). ANK3, a family member of ankyrins, assisted in maintaining cell stability through anchoring cytoskeleton to the cell membrane (Wang et al., 2016). The downregulation of both genes could lead to the proliferation and invasion of cancer cells. However, little was known about coiled-coil domain-containing protein 85A (CCDC85A) and complement component 8 gamma (C8G) in tumorigenesis and further investigations were required.
Further analysis of the TME in ccRCC revealed that BM-based clusters and risk groups had distinct immune landscape. Higher immune scores and infiltration levels of immune cells were discovered in the cluster or group that possessed more enriched BM-related functions. It was reported that BM was one of the barriers through which the extravasation of lymphocytes into tumor sites should be overcome (Marchand et al., 2019). During this process, laminin, the major component of BM, played an important role in mediating the functions of lymphocyte trafficking (Pozzi et al., 2017). Evidence suggested that laminins containing LAMA4 favored the transmigration of T cells through providing some permissive signals while those with LAMA5 tended to oppose this process (Sixt et al., 2001). Besides, T cells might adhere less strongly to LAMA4-bearing laminins but migrate faster across these lamins than those with LAMA5 (Zhang et al., 2020). BM distributed with a higher level of LAMA4-containing laminins would have a higher potential for lymphocyte extravasation (Wu et al., 2009). However, we discovered that cluster with high immune cell infiltration was related to poor prognosis and most of the infiltrated lymphocytes in ccRCC including CD8 T lymphocyte and CD4 T lymphocyte were associated with low survival time. This could be due to the dysfunction of lymphocytes by the suppressive lymphocytes and tumor cells. Myeloid-derived suppressor cells (MDSCs) were a cluster of immune suppressive myeloid cells frequently found in cancers (Gabrilovich et al., 2007). As their name suggested, MDSCs exerted suppressive functions on various cells, especially T lymphocytes (Gabrilovich et al., 2012). They could induce tolerance of antigen-specific T lymphocytes mainly through the production of reactive oxygen species (ROS) to nitrate receptors on T lymphocytes and reduce their responsiveness, inhibiting the anti-tumor function (Nagaraj et al., 2007). Apart from suppressive immune cells, tumor cells themselves could transduce inhibitory signals to lymphocytes for immune evasion (Muenst et al., 2016). Inhibitory immune checkpoints such as PD1 were targeted by tumor cells to transform T lymphocytes into suppressive status (Daassi et al., 2020). Thus, BM remodeling may assist in the development of the rich but suppressive immune microenvironment in ccRCC.
With an increasing understanding of the mechanisms underlying tumorigenesis, therapeutic strategies developed against them were shown to be effective in cancer therapy. ICIs, functioning through interrupting the suppressive signals transduced by CTLA4 or PD1 to reactivate the anti-tumor immunity and prevent tumor immune evasion, was approved to be successful in the treatment of multiple cancers (Sharma and Allison, 2015). It was believed that normalization of the suppressive environment and restoration of the anti-tumor immunity would be more effective than directly enhancing the immune function (Sanmamed and Chen, 2018). In line with this statement, our finding suggested that clusters with high immunogenicity potential but suppressed responded better to ICIs. The higher level of SNV, HRD, and CTA in the high risk group indicated it possessed a higher potential to generate tumor-associated neoantigens, which were critical for the immune system to exert anti-tumor function (Meng et al., 2021; van Wilpe et al., 2021; Wang et al., 2021). Meanwhile, once the suppressive environment was removed, the relatively high infiltration level of immune cells served as an immune reservoir that supported a powerful immune reaction (Waldman et al., 2020). Additionally, the increased expression of immune checkpoints in the high risk group was thought to be predictive for immunotherapy (Thommen et al., 2018). Apart from ICI monotherapy, recent studies had domonstrated the benefits of treatment combined ICI and TKI, which could have the better therapeutic capacity (Quhal et al., 2021). Our results also demonstrated that the high risk group was more sensitive to TKIs including sorafenib and sunitinib, as well as chemotherapy like mitomycin C. The mutation analysis helped discover a co-occurance mutant pair, SETD2 and PBRM1, in which SETD2 mutation was correlated with a favorable outcome of ICI-treated patients (Lu et al., 2021) and PBRM1 mutation was associated with high angiogenesis and TKIs therapeutic outcomes (Motzer et al., 2020b), indicating the potential of the high risk group for combination therapy. It was demonstrated that TKIs therapy was associated with improved vessel extravasation and enhanced drug delivery to tumors (Zhou et al., 2008; Zhou and Gallo, 2009). Moreover, after TKI treatment, the amount of suppressive immune cells including Treg and MDSCs was found to be decreased in TME (Finke et al., 2008; Ko et al., 2009). Therefore, BMRS constructed in the current study could not only help individualized administration of immunotherapy but also assist in the combination therapy of ICIs and TKIs.
In addition to the existing therapeutic agents, the investigation was also focus on novel therapeutic targets. In patients with high BMRS, we identified the upregulation of chemokine (C-X-C motif) ligand 13 (CXCL13) and its receptor CXCR5. Their interaction could impede the tumor-specific cytotoxic function of CD8 T lymphocytes and be related to the recruitment of suppressive immune cells including MDSCs and Treg (Ammirante et al., 2014). Inhibitor targeting the CXCL13/CXCR5 axis was demonstrated to have an encouraging effect on cancer treatment (Hussain et al., 2019). A Rho associated kinase inhibitor, BRD-K23875128, may be a potential therapeutic agent for high risk patients. Rho kinase pathway took part in multiple cell functions and was implicated in tumor metastasis as well as ECM remodeling. It may be a candidate for combination treatment due to its ability to increase sensitivity to other therapeutic drugs (Kim et al., 2021). As for patients with low BMRS, there may also be some curative agents. Histocompatibility leucocyte antigen E (HLA-E) and its cognate inhibitory receptor NKG2A could serve as a novel immune checkpoint to be targeted for inducing anti-tumor immunity (Borst et al., 2020). As predicted by CMap, arginase inhibitor, could also potentially treat patients with high BMRS by promoting the T cells activation and proliferation to exert an anti-tumor immune response (Borek et al., 2020).
Overall, the current study identified BM-related subtypes of ccRCC and constructed a risk scoring system BMRS for prognosis and therapeutic prediction based on both public and clinical data. The mutation, TME and treatment analysis also provided potential novel therapeutic agents for further investigations. Comparing with the existing similar gene signatures, BMRS possessed higher net benefit in predicting the prognosis of ccRCC patients. Besides, we had validated the scoring system and its gene expression using external dataset and clinical samples. However, the lack of accessibility to clinical treatment-related data restricted our analysis of BMRS in predicting therapeutic response. More in-depth study could be conducted to validate the proposed potential therapeutic targets for ccRCC.
5 Conclusion
In summary, we recognized BM-related subtypes of ccRCC with distinct survival and TME features. A risk scoring system BMRS was established for prognosis prediction and individualized treatment instruction. Mechanistic investigations based on BM-related clusters and risk groups helped identified some therapeutic candidates. With more and more studies focusing on combination therapy, our results may provide certain practical instructions for clinical application and future research.
Data availability statement
The data and information demonstrated and analyzed throughout the present study were obtained from the Genomics Data Commons Data Portal (https://portal.gdc.cancer.gov/), ArrayExpress (https://www.ebi.ac.uk/arrayexpress/), Gene Set Enrichment Analysis (http://www.gsea-msigdb.org/gsea/index.jsp), xCell (https://xcell.ucsf.edu/), The Cancer Immunome Database (https://tcia.at/home), Genomics of Drug Sensitivity in Cancer (https://www.cancerrxgene.org), ConnectivityMap (CMap, https://clue.io). The patient information of the clinical transcriptome data could be accessed from the supplementary material (Supplementary Table S1).
Ethics statement
The studies involving human participants were reviewed and approved by Ethics Committee of Guangdong Provincial People’s Hospital. The patients/participants provided their written informed consent to participate in this study.
Author contributions
YT: Conception and design, data analysis and interpretation, collection and assembly of data; CY: Conception and design, data analysis and interpretation; JZ and PZ: Data analysis and interpretation; SC, WZ, and BY: Collection and assembly of data; YY and YL: Administrative support. All authors wrote and approved the final version of the manuscript.
Funding
This work was supported by Guangdong Basic and Applied Basic Research Foundation (Grant/Award Number: 2020A1515011473) and Science and Technology Program of Guangzhou (Grant/Award Number: 202002030060).
Acknowledgments
The information utilized in this study was retrieved from The Cancer Genome Atlas (TCGA), ArrayExpress, xCell, cBioPortal for Cancer Genomics, The Cancer Immunome Database, Genomics of Drug Sensitivity in Cancer and Gene Set Enrichment Analysis. We sincerely appreciate their sharing of the data and information.
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/fgene.2023.1038924/full#supplementary-material
References
Amin, A., and Hammers, H. (2018). The evolving landscape of immunotherapy-based combinations for frontline treatment of advanced renal cell carcinoma. Front. Immunol. 9, 3120. doi:10.3389/fimmu.2018.03120
Ammirante, M., Shalapour, S., Kang, Y., Jamieson, C. A. M., and Karin, M. (2014). Tissue injury and hypoxia promote malignant progression of prostate cancer by inducing CXCL13 expression in tumor myofibroblasts. Proc. Natl. Acad. Sci. U. S. A. 111 (41), 14776–14781. doi:10.1073/pnas.1416498111
Barbour, J. D., Cox, N., and Huang, R. S. (2014). pRRophetic: An R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS ONE 9 (9), e107468. doi:10.1371/journal.pone.0107468
Bissell, M. J., and Hines, W. C. (2011). Why don't we get more cancer? A proposed role of the microenvironment in restraining cancer progression. Nat. Med. 17 (3), 320–329. doi:10.1038/nm.2328
Blankenstein, T., Coulie, P. G., Gilboa, E., and Jaffee, E. M. (2012). The determinants of tumour immunogenicity. Nat. Rev. Cancer 12 (4), 307–313. doi:10.1038/nrc3246
Borek, B., Gajda, T., Golebiowski, A., and Blaszczyk, R. (2020). Boronic acid-based arginase inhibitors in cancer immunotherapy. Bioorg Med. Chem. 28 (18), 115658. doi:10.1016/j.bmc.2020.115658
Borst, L., van der Burg, S. H., and van Hall, T. (2020). The nkg2a-HLA-E Axis as a novel checkpoint in the tumor microenvironment. Clin. Cancer Res. 26 (21), 5549–5556. doi:10.1158/1078-0432.CCR-19-2095
Chang, J., and Chaudhuri, O. (2019). Beyond proteases: Basement membrane mechanics and cancer invasion. J. Cell Biol. 218 (8), 2456–2469. doi:10.1083/jcb.201903066
Chen, T., Zhou, L., Zhou, Y., Zhou, W., Huang, H., Yin, S., et al. (2019). HJURP promotes epithelial-to-mesenchymal transition via upregulating SPHK1 in hepatocellular carcinoma. Int. J. Biol. Sci. 15 (6), 1139–1147. doi:10.7150/ijbs.30904
Choueiri, T. K., and Motzer, R. J. (2017). Systemic therapy for metastatic renal-cell carcinoma. N. Engl. J. Med. 376 (4), 354–366. doi:10.1056/NEJMra1601333
Conrad, C., Benzel, J., Dorzweiler, K., Cook, L., Schlomann, U., Zarbock, A., et al. (2019). ADAM8 in invasive cancers: Links to tumor progression, metastasis, and chemoresistance. Clin. Sci. (Lond) 133 (1), 83–99. doi:10.1042/CS20180906
Daassi, D., Mahoney, K. M., and Freeman, G. J. (2020). The importance of exosomal PDL1 in tumour immune evasion. Nat. Rev. Immunol. 20 (4), 209–215. doi:10.1038/s41577-019-0264-y
de Velasco, G., Culhane, A. C., Fay, A. P., Hakimi, A. A., Voss, M. H., Tannir, N. M., et al. (2017). Molecular subtypes improve prognostic value of international metastatic renal cell carcinoma database consortium prognostic model. Oncologist 22 (3), 286–292. doi:10.1634/theoncologist.2016-0078
Díaz-Montero, C. M., Rini, B. I., and Finke, J. H. (2020). The immunology of renal cell carcinoma. Nat. Rev. Nephrol. 16, 721–735. doi:10.1038/s41581-020-0316-3
Eatemadi, A., Aiyelabegan, H. T., Negahdari, B., Mazlomi, M. A., Daraee, H., Daraee, N., et al. (2017). Role of protease and protease inhibitors in cancer pathogenesis and treatment. Biomed. Pharmacother. 86, 221–231. doi:10.1016/j.biopha.2016.12.021
Finke, J. H., Rini, B., Ireland, J., Rayman, P., Richmond, A., Golshayan, A., et al. (2008). Sunitinib reverses type-1 immune suppression and decreases T-regulatory cells in renal cell carcinoma patients. Clin. Cancer Res. 14 (20), 6674–6682. doi:10.1158/1078-0432.CCR-07-5212
Friedman, J., Hastie, T., and Tibshirani, R. (2010). Regularization paths for generalized linear models via coordinate descent. J. Stat. Softw. 33 (1), 1–22. doi:10.18637/jss.v033.i01
Gabrilovich, D. I., Bronte, V., Chen, S. H., Colombo, M. P., Ochoa, A., Ostrand-Rosenberg, S., et al. (2007). The terminology issue for myeloid-derived suppressor cells. Cancer Res. 67 (1), 425. doi:10.1158/0008-5472.CAN-06-3037
Gabrilovich, D. I., Ostrand-Rosenberg, S., and Bronte, V. (2012). Coordinated regulation of myeloid cells by tumours. Nat. Rev. Immunol. 12 (4), 253–268. doi:10.1038/nri3175
Greef, B., and Eisen, T. (2016). Medical treatment of renal cancer: New horizons. Br. J. Cancer 115 (5), 505–516. doi:10.1038/bjc.2016.230
Heneberg, P. (2016). Paracrine tumor signaling induces transdifferentiation of surrounding fibroblasts. Crit. Rev. Oncol. Hematol. 97, 303–311. doi:10.1016/j.critrevonc.2015.09.008
Ho, T. H., Serie, D. J., ParasraMka, M., Cheville, J. C., Bot, B. M., Tan, W., et al. (2017). Differential gene expression profiling of matched primary renal cell carcinoma and metastases reveals upregulation of extracellular matrix genes. Ann. Oncol. 28 (3), 604–610. doi:10.1093/annonc/mdw652
Hussain, M., Adah, D., Tariq, M., Lu, Y., Zhang, J., and Liu, J. (2019). CXCL13/CXCR5 signaling axis in cancer. Life Sci. 227, 175–186. doi:10.1016/j.lfs.2019.04.053
Jayadev, R., Morais, M. R. P. T., Ellingford, J. M., Srinivasan, S., Naylor, R. W., Lawless, C., et al. (2022). A basement membrane discovery pipeline uncovers network complexity, regulators, and human disease associations. Sci. Adv. 8 (20), eabn2265. doi:10.1126/sciadv.abn2265
Jhaveri, K. D., and Perazella, M. A. (2018). Adverse events associated with immune checkpoint blockade. N. Engl. J. Med. 378 (12), 1163. doi:10.1056/NEJMc1801663
Jonasch, E., Gao, J., and Rathmell, W. K. (2014). Renal cell carcinoma. BMJ 349, g4797. doi:10.1136/bmj.g4797
Kessenbrock, K., Plaks, V., and Werb, Z. (2010). Matrix metalloproteinases: Regulators of the tumor microenvironment. Cell 141 (1), 52–67. doi:10.1016/j.cell.2010.03.015
Kim, S., Han, J., and Kim, I. S. (2021). Rho-kinase as a target for cancer therapy and its immunotherapeutic potential. Int. J. Mol. Sci. 22 (23), 12916. doi:10.3390/ijms222312916
Ko, J. S., Zea, A. H., Rini, B. I., Ireland, J. L., Elson, P., Cohen, P., et al. (2009). Sunitinib mediates reversal of myeloid-derived suppressor cell accumulation in renal cell carcinoma patients. Clin. Cancer Res. 15 (6), 2148–2157. doi:10.1158/1078-0432.CCR-08-1332
Li, P., Ren, H., Zhang, Y., and Zhou, Z. (2018). Fifteen-gene expression based model predicts the survival of clear cell renal cell carcinoma. Med. Baltim. 97 (33), e11839. doi:10.1097/MD.0000000000011839
Liu, X., Qiao, Y., Chen, J., and Ge, G. (2022). Basement membrane promotes tumor development by attenuating T cell activation. J. Mol. Cell Biol. 14 (2), mjac006. doi:10.1093/jmcb/mjac006
Lu, M., Zhao, B., Liu, M., Wu, L., Li, Y., Zhai, Y., et al. (2021). Pan-cancer analysis of SETD2 mutation and its association with the efficacy of immunotherapy. NPJ Precis. Oncol. 5 (1), 51. doi:10.1038/s41698-021-00193-0
Marchand, M., Monnot, C., Muller, L., and Germain, S. (2019). Extracellular matrix scaffolding in angiogenesis and capillary homeostasis. Semin. Cell Dev. Biol. 89, 147–156. doi:10.1016/j.semcdb.2018.08.007
Mayakonda, A., Lin, D. C., Assenov, Y., Plass, C., and Koeffler, H. P. (2018). Maftools: Efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 28 (11), 1747–1756. doi:10.1101/gr.239244.118
Meng, X., Sun, X., Liu, Z., and He, Y. (2021). A novel era of cancer/testis antigen in cancer immunotherapy. Int. Immunopharmacol. 98, 107889. doi:10.1016/j.intimp.2021.107889
Motzer, R. J., Banchereau, R., Hamidi, H., Powles, T., McDermott, D., Atkins, M. B., et al. (2020). Molecular subsets in renal cancer determine outcome to checkpoint and angiogenesis blockade. Cancer Cell 38 (6), 803–817.e4. doi:10.1016/j.ccell.2020.10.011
Motzer, R. J., Escudier, B., McDermott, D. F., George, S., Hammers, H. J., Srinivas, S., et al. (2015). Nivolumab versus everolimus in advanced renal-cell carcinoma. N. Engl. J. Med. 373 (19), 1803–1813. doi:10.1056/NEJMoa1510665
Motzer, R. J., Jonasch, E., Boyle, S., Carlo, M. I., Manley, B., Agarwal, N., et al. (2020). NCCN guidelines insights: Kidney cancer, version 1.2021. J. Natl. Compr. Canc Netw. 18 (9), 1160–1170. doi:10.6004/jnccn.2020.0043
Muenst, S., Laubli, H., Soysal, S. D., Zippelius, A., TzAnkov, A., and Hoeller, S. (2016). The immune system and cancer evasion strategies: Therapeutic concepts. J. Intern Med. 279 (6), 541–562. doi:10.1111/joim.12470
Nagaraj, S., Gupta, K., Pisarev, V., Kinarsky, L., Sherman, S., Kang, L., et al. (2007). Altered recognition of antigen is a mechanism of CD8+ T cell tolerance in cancer. Nat. Med. 13 (7), 828–835. doi:10.1038/nm1609
Nakasuka, F., Tabata, S., Sakamoto, T., Hirayama, A., Ebi, H., Yamada, T., et al. (2021). TGF-beta-dependent reprogramming of amino acid metabolism induces epithelial-mesenchymal transition in non-small cell lung cancers. Commun. Biol. 4 (1), 782. doi:10.1038/s42003-021-02323-7
Pozzi, A., Yurchenco, P. D., and Iozzo, R. V. (2017). The nature and biology of basement membranes. Matrix Biol. 57-58, 1–11. doi:10.1016/j.matbio.2016.12.009
Quhal, F., Mori, K., Bruchbacher, A., Resch, I., Mostafaei, H., Pradere, B., et al. (2021). First-line immunotherapy-based combinations for metastatic renal cell carcinoma: A systematic review and network meta-analysis. Eur. Urol. Oncol. 4 (5), 755–765. doi:10.1016/j.euo.2021.03.001
Reuten, R., Zendehroud, S., Nicolau, M., Fleischhauer, L., Laitala, A., Kiderlen, S., et al. (2021). Basement membrane stiffness determines metastases formation. Nat. Mater 20 (6), 892–903. doi:10.1038/s41563-020-00894-0
Ritchie, M. E., Phipson, B., Wu, D., Hu, Y., Law, C. W., Shi, W., et al. (2015). Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43 (7), e47. doi:10.1093/nar/gkv007
Sanmamed, M. F., and Chen, L. (2018). A paradigm shift in cancer immunotherapy: From enhancement to normalization. Cell 175 (2), 313–326. doi:10.1016/j.cell.2018.09.035
Sharma, P., and Allison, J. P. (2015). The future of immune checkpoint therapy. Science 348 (6230), 56–61. doi:10.1126/science.aaa8172
Siegel, R. L., Miller, K. D., Fuchs, H. E., and Jemal, A. (2022). Cancer statistics. CA Cancer J. Clin. 72 (1), 7–33. doi:10.3322/caac.21708
Sixt, M., Engelhardt, B., Pausch, F., Hallmann, R., Wendler, O., and Sorokin, L. M. (2001). Endothelial cell laminin isoforms, laminins 8 and 10, play decisive roles in T cell recruitment across the blood–brain barrier in experimental autoimmune encephalomyelitis. J. Cell Biol. 153 (5), 933–946. doi:10.1083/jcb.153.5.933
Song, H., Liu, L., Song, Z., Ren, Y., Li, C., and Huo, J. (2018). P4HA3 is epigenetically activated by slug in gastric cancer and its deregulation is associated with enhanced metastasis and poor survival. Technol. Cancer Res. Treat. 17, 1533033818796485. doi:10.1177/1533033818796485
Sung, H., Ferlay, J., Siegel, R. L., Laversanne, M., Soerjomataram, I., Jemal, A., et al. (2021). Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin. 71 (3), 209–249. doi:10.3322/caac.21660
Thommen, D. S., Koelzer, V. H., Herzig, P., Roller, A., Trefny, M., Dimeloe, S., et al. (2018). A transcriptionally and functionally distinct PD-1(+) CD8(+) T cell pool with predictive potential in non-small-cell lung cancer treated with PD-1 blockade. Nat. Med. 24 (7), 994–1004. doi:10.1038/s41591-018-0057-z
Thorsson, V., Gibbs, D. L., Brown, S. D., Wolf, D., Bortone, D. S., Ou Yang, T. H., et al. (2018). The immune landscape of cancer. Immunity 48 (4), 812–830 e14. doi:10.1016/j.immuni.2018.03.023
Valastyan, S., and Weinberg, R. A. (2011). Tumor metastasis: Molecular insights and evolving paradigms. Cell 147 (2), 275–292. doi:10.1016/j.cell.2011.09.024
van Wilpe, S., Tolmeijer, S. H., Koornstra, R. H. T., de Vries, I. J. M., Gerritsen, W. R., Ligtenberg, M., et al. (2021). Homologous recombination repair deficiency and implications for tumor immunogenicity. Cancers (Basel) 13 (9), 2249. doi:10.3390/cancers13092249
Waldman, A. D., Fritz, J. M., and Lenardo, M. J. (2020). A guide to cancer immunotherapy: from T cell basic science to clinical practice. Nat. Rev. Immunol. 20 (11), 651–668. doi:10.1038/s41577-020-0306-5
Wang, T., Abou-Ouf, H., Hegazy, S. A., Alshalalfa, M., Stoletov, K., Lewis, J., et al. (2016). Ankyrin G expression is associated with androgen receptor stability, invasiveness, and lethal outcome in prostate cancer patients. J. Mol. Med. Berl. 94 (12), 1411–1422. doi:10.1007/s00109-016-1458-4
Wang, Y., Shi, T., Song, X., Liu, B., and Wei, J. (2021). Gene fusion neoantigens: Emerging targets for cancer immunotherapy. Cancer Lett. 506, 45–54. doi:10.1016/j.canlet.2021.02.023
Winer, A., Adams, S., and Mignatti, P. (2018). Matrix metalloproteinase inhibitors in cancer therapy: Turning past failures into future successes. Mol. Cancer Ther. 17 (6), 1147–1155. doi:10.1158/1535-7163.MCT-17-0646
Wisdom, K. M., Adebowale, K., Chang, J., Lee, J. Y., Nam, S., Desai, R., et al. (2018). Matrix mechanical plasticity regulates cancer cell migration through confining microenvironments. Nat. Commun. 9 (1), 4144. doi:10.1038/s41467-018-06641-z
Wolf, K., Te Lindert, M., Krause, M., Alexander, S., Te Riet, J., Willis, A. L., et al. (2013). Physical limits of cell migration: Control by ECM space and nuclear deformation and tuning by proteolysis and traction force. J. Cell Biol. 201 (7), 1069–1084. doi:10.1083/jcb.201210152
Wragg, J. W., Finnity, J. P., Anderson, J. A., Ferguson, H. J. M., Porfiri, E., Bhatt, R. I., et al. (2016). MCAM and LAMA4 are highly enriched in tumor blood vessels of renal cell carcinoma and predict patient outcome. Cancer Res. 76 (8), 2314–2326. doi:10.1158/0008-5472.CAN-15-1364
Wu, C., Ivars, F., Anderson, P., Hallmann, R., Vestweber, D., Nilsson, P., et al. (2009). Endothelial basement membrane laminin alpha5 selectively inhibits T lymphocyte extravasation into the brain. Nat. Med. 15 (5), 519–527. doi:10.1038/nm.1957
Wu, T., Hu, E., Xu, S., Chen, M., Guo, P., Dai, Z., et al. (2021). clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. Innov. (N Y) 2 (3), 100141. doi:10.1016/j.xinn.2021.100141
Xiong, X., Chen, C., Yang, J., Ma, L., Wang, X., Zhang, W., et al. (2022). Characterization of the basement membrane in kidney renal clear cell carcinoma to guide clinical therapy. Front. Oncol. 12, 1024956. doi:10.3389/fonc.2022.1024956
Xu, C., Liu, F., Xiang, G., Cao, L., Wang, S., Liu, J., et al. (2019). β-Catenin nuclear localization positively feeds back on EGF/EGFR-attenuated AJAP1 expression in breast cancer. J. Exp. Clin. Cancer Res. 38 (1), 238. doi:10.1186/s13046-019-1252-6
Xu, T., Le, T. D., Liu, L., Su, N., Wang, R., Sun, B., et al. (2017). CancerSubtypes: An R/bioconductor package for molecular cancer subtype identification, validation and visualization. Bioinformatics 33 (19), 3131–3133. doi:10.1093/bioinformatics/btx378
Yi, M., Nissley, D. V., McCormick, F., and Stephens, R. M. (2020). ssGSEA score-based Ras dependency indexes derived from gene expression data reveal potential Ras addiction mechanisms with possible clinical implications. Sci. Rep. 10 (1), 10258. doi:10.1038/s41598-020-66986-8
Yoshihara, K., Shahmoradgoli, M., Martinez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring tumour purity and stromal and immune cell admixture from expression data. Nat. Commun. 4, 2612. doi:10.1038/ncomms3612
Yurchenco, P. D. (2011). Basement membranes: Cell scaffoldings and signaling platforms. Cold Spring Harb. Perspect. Biol. 3 (2), a004911. doi:10.1101/cshperspect.a004911
Zhang, F., Yuan, D., Song, J., Chen, W., Wang, W., Zhu, G., et al. (2021). HJURP is a prognostic biomarker for clear cell renal cell carcinoma and is linked to immune infiltration. Int. Immunopharmacol. 99, 107899. doi:10.1016/j.intimp.2021.107899
Zhang, X., Wang, Y., Song, J., Gerwien, H., Chuquisana, O., Chashchina, A., et al. (2020). The endothelial basement membrane acts as a checkpoint for entry of pathogenic T cells into the brain. J. Exp. Med. 217 (7), e20191339. doi:10.1084/jem.20191339
Zhou, H., Huang, S., Shao, C., Zou, J., Zhou, A., Yu, J., et al. (2022). miR-1266-3p suppresses epithelial-mesenchymal transition in colon cancer by targeting P4HA3. Anal. Cell Pathol. (Amst) 2022, 1542117. doi:10.1155/2022/1542117
Zhou, Q., and Gallo, J. M. (2009). Differential effect of sunitinib on the distribution of temozolomide in an orthotopic glioma model. Neuro Oncol. 11 (3), 301–310. doi:10.1215/15228517-2008-088
Zhou, Q., Guo, P., and Gallo, J. M. (2008). Impact of angiogenesis inhibition by sunitinib on tumor distribution of temozolomide. Clin. Cancer Res. 14 (5), 1540–1549. doi:10.1158/1078-0432.CCR-07-4544
Zhou, T., Chen, W., Wu, Z., Cai, J., and Zhou, C. (2022). A newly defined basement membrane-related gene signature for the prognosis of clear-cell renal cell carcinoma. Front. Genet. 13, 994208. doi:10.3389/fgene.2022.994208
Abbreviation
ccRCC, Clear cell renal cell carcinoma; ICIs, Immune checkpoint inhibitors; VEGF, Vascular endothelial growth factor; TKIs, Tyrosine kinase inhibitors; TME, Tumor microenvironment; BM, Basement membrane; ECM, Extracellular matrix; TCGA, The cancer genome atlas; TPM, Transcripts per million; CNV, Copy number variations; LASSO, Least absolute selection and shrinkage operator; BMRS, BM-related risk score system; CTA, Cancer testis antigens; HRD, Homologous recombination repair deficiency; ITH, Intratumoral heterogeneity; SNV, Single-nucleotide variation; TIDE, Tumor immune dysfunction and evasion; TMB, Tumor mutation burden; MDSC, Myeloid-derived suppressor cell; Treg, Regulatory T cell; WGCNA, Weight gene correlation network analysis; CNMF, Non-negative matrix factorization; ESTIMATE, Estimation of stromal and immune cells in malignant tumor tissues using expression data; ssGSEA, Single-sample gene set enrichment analysis.
Keywords: clear cell renal cell carcinoma, basement membrane, gene expression, tumor microenvironment, individualized therapy
Citation: Tang Y, Ye C, Zeng J, Zhu P, Cheng S, Zeng W, Yang B, Liu Y and Yu Y (2023) Identification of a basement membrane-based risk scoring system for prognosis prediction and individualized therapy in clear cell renal cell carcinoma. Front. Genet. 14:1038924. doi: 10.3389/fgene.2023.1038924
Received: 07 September 2022; Accepted: 20 January 2023;
Published: 03 February 2023.
Edited by:
Shun Lu, University of Electronic Science and Technology of China, ChinaReviewed by:
Tianhang Li, Nanjing University, ChinaYoutao Lu, University of Pennsylvania, United States
Copyright © 2023 Tang, Ye, Zeng, Zhu, Cheng, Zeng, Yang, Liu and Yu. 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: Yuming Yu, yuym72@163.com; Yanjun Liu, yanjun@smu.edu.cn
†These authors have contributed equally to this work