A Novel Prognostic Signature Associated With the Tumor Microenvironment in Kidney Renal Clear Cell Carcinoma

Background: The tumor microenvironment (TME) is a complex and evolving environment, and the tumor immune microenvironment in kidney renal clear cell carcinoma (KIRC) has a strong suppressive profile. This study investigates the potential prognostic role and value of genes of the tumor microenvironment in KIRC.

Methods: The transcriptome sequencing data of 530 cases and 39 cases of KIRC and the corresponding clinical prognosis information were downloaded from TCGA data and GEO data, respectively, and TME-related gene expression profiles were extracted. A prognostic signature was constructed and evaluated using univariate Cox regression analysis and LASSO regression analysis. Gene set enrichment analysis (GSEA) was used to obtain the biological process of gene enrichment in patients with high and low-risk groups.

Results: A prognostic signature consisting of eight TME-related genes (LRFN1, CSF1, UCN, TUBB2B, SERPINF1, ADAM8, ABCB4, CCL22) was constructed. Kaplan-Meier survival analysis yielded significantly lower survival times for patients in the high-risk group than in the low-risk group, and the AUC values for the ROC curves of this prognostic signature were essentially greater than 0.7, and univariate and multifactorial Cox regression analyses indicated that the risk score was independent risk factors for KIRC prognosis. GSEA analysis showed that immune-related biological processes were enriched in the high-risk group and that risk values were strongly associated with multiple immune cell scores and immune checkpoint-related genes (PDCD1, CTLA4).

Conclusions: The prognostic signature can accurately predict the prognosis of KIRC patients, which may provide new ideas for future precision immunotherapy of KIRC.

Introduction

Renal cell carcinoma (RCC) is a common genitourinary malignancy that causes nearly 170,000 deaths each year (1). Kidney renal clear cell carcinoma (KIRC) is the predominant histological type, accounting for approximately 75% of all RCC cases (2). The main treatment modality for KIRC is currently surgical resection, and early-stage patients can achieve good results with surgery. However, in patients with advanced KIRC with recurrence and metastases, surgery is difficult, recurrence rates are high, conventional radiotherapy is not sensitive, and the prognosis is usually very poor (3). However, no reliable biomarkers have been identified to predict the prognosis of KIRC patients (4).

The tumor microenvironment (TME) is a complex and evolving environment that varies in composition depending on the tumor type and consists mainly of stromal cells, immune cells, and the extracellular matrix (ECM), with immune cells being a key component of the TME (5). Renal cancer TME is also a dynamic system that plays a key role in driving immune escape (6). Unlike other tumor types, the tumor immune microenvironment of KIRC is characterized by a high degree of immune cell infiltration, with the highest degree of T-cell infiltration (7). Chevrier et al. also found that T cells and tumour-associated macrophages (TAMs) were the major immune cell populations in KIRC, accounting for an average of 51% and 31% respectively, with CD8+ T cells being mostly depleted and functionally deficient (8). A growing number of studies have shown that the therapeutic efficacy of immune checkpoint inhibition (ICIs) is closely related to the components of the TME. Recent clinical trial results have shown that ICIs combined with anti-angiogenic agents, or a combination of different ICIs are superior to monotherapy, making it the most effective treatment strategy for advanced KIRC today (9). However, immunosuppressive cell subsets and molecules in the TME can lead to the insensitivity of KIRC to immunotherapy. Therefore, it is essential to elucidate the occurrence and development of KIRC and TME-related genes, and transcriptome sequencing can be used to screen and identify potential targets for disease treatment (10).

It also has been shown that cancer cells and tumor-infiltrating immune cells in TME play an important role in regulating cancer progression. They play an important role in determining the type of malignancy. Tumor-infiltrating lymphocytes, including T cells and B cells, are an important class of cells in TME. CD4+ helper T cells and cytotoxic CD8+ T cells play an important role in tumor prevention by targeting antigenic tumor cells, and CD8+ T cells are associated with better clinical outcomes and immunotherapeutic responses in many cancers. In addition, it was recently observed that tumor-associated B cells play an important role in the immune system by producing antibodies and presenting antigens that can predict survival and response to immune checkpoint blockade therapy. In addition, an association between TME genetic signature and lower survival was observed in KIRC patients, and tumor-associated macrophage and T cell phenotypes were found to correlate with clinical outcomes. These observations highlight the importance of analyzing TME, including immune cell variability, to identify target tumors for each specific treatment and to design new effective cancer therapies.

Given the high heterogeneity of TME in KIRC (11), understanding TME-related gene expression changes has clinical implications for clarifying the prominent molecular characterization and prognosis of KIRC. Therefore, we identified two subtypes with different TME characteristics based on TME-related genes and constructed a prognosis-related predictive scoring signature to provide a molecular basis for the pathogenesis and treatment of KIRC.

Materials and Methods

Data Sources and Processing

Gene expression data and clinical information on KIRC were downloaded from The Cancer Genome Atlas (TCGA) database (72 normal samples and 539 tumor samples), and after excluding patients with no recorded survival time, 530 KIRC tumor samples were finally included, defined as an entire TCGA cohort. The dataset GSE29609 (39 tumor samples) containing the KIRC for survival time was also downloaded from the Gene Expression Omnibus (GEO) as a model test set to validate the predictive power of the signature. We eliminate the non-conforming data through quality control. Then all data expression values are transformed to a comparable level for subsequent analysis by standardization of the data. The main standardization method is quantile standardization. Based on previous studies (12–14), 4061 TME-related genes were obtained. TME-related genes were extracted, and differential genes were analyzed using the “limma” package in R. A threshold value of P < 0.05 and |log2FC| > 1 was defined to identify differentially expressed genes (DEGs).

Identification of Molecular Subtypes of TME-Related DEGs

The differentially expressed TME-associated genes in KRIC tissues and adjacent tissues were screened, and the screened differential genes and their expression were organized into a gene expression matrix with a corrected p<0.05 and the absolute value of differential expression multiplicity >1 (FDR<0.05 and | log2Fold Change|>1) was set as the threshold value, and the “NMF” package (15) was used to extract the biological correlation coefficients of the data in the above-mentioned differential gene expression matrix. The samples were grouped by organizing the genes and samples to capture the internal structural features of the data. When clusters k = 2, the clusters showed appropriate performance and stability, resulting in two subtypes (C1 and C2). And survival prognosis analysis, including overall survival (OS) and progression-free survival (PFS), was performed using ‘survival’.

Differences in Immunological Characteristics of Different Subtypes

A comparative analysis of the immune profile of patients with different subtypes was performed to clarify the differences in the immune profile of each subtype. Immune cells are an important component of the tumor immune microenvironment, and we used the “MCPcounter” R package (13) to calculate the immune scores of 10 immune cells, including CD8+ T cells, cytotoxic lymphocytes, fibroblasts, monocytic lineage, myeloid dendritic cells, NK cells, T cells, neutrophils, endothelial cells, B lineage. The immune scores of immune cells were compared between C1 and C2.

Prognostic Signature Construction and Validation

Using the ‘caret’ package (16), the 530 TCGA samples were randomly divided into the most appropriate TCGA training cohort and TCGA testing cohort using a 3:7 subgroup sampling ratio, with both groups being similar in terms of clinical characteristics. The TCGA training cohort was used for signature construction, and the TCGA testing cohort and the entire TCGA cohort were used to internally validate the predictive power of the signature, while external validation was performed in the GSE29609. In the TCGA training cohort, the differential genes were subjected to univariate Cox regression, set at P ≤ 0.05, to screen for the corresponding TME-related prognostic genes. LASSO regression and cross-validation were performed using the “glmnet” package to obtain the optimal gene set, and 10-fold cross-validation was used to construct a prognostic risk score signature with the formula: risk score = risk gene expression₁ × coef₁ + risk gene expression₂ × coef₂ +… + risk geneexpression × coef_n (coef was the risk coefficient). The prognostic signature was used to calculate the risk values for each sample of TCGA training cohort, TCGA testing cohort, entire TCGA cohort, and GSE29609, respectively. The risk scores of the TCGA training cohort groups were ranked from lowest to highest and patients were divided into low- and high-risk groups based on the median. Kaplan-Meier analysis and log-rank tests were used to assess the survival of the two groups of patients, and receiver operating characteristic (ROC) curve analysis and decision curve analysis (DCA) were used to evaluate the predictive ability of the signature. We also established a nomogram to better predict the prognosis of KIRC patients, and the concordance index (C-index) and DCA were used to evaluate the performance of the nomogram. To demonstrate the superior performance of the signature developed in this study, we compared the signature with three recently published prediction models (17–19).

Independence Analysis of Prognostic Signature and Clinical Characteristics

In addition, the predictive power of the predictive signature and other clinical characteristics (Age, gender, grade,stage, TNM staging) were compared using univariate Cox regression analysis and multivariate Cox regression analysis. In the entire TCGA cohort, Kaplan-Meier survival analysis stratified by TNM staging and stage staging was performed on KIRC patients to further validate the performance of the prognostic signature.

Gene Set Enrichment Analysis for the Prognostic Signature and its Significance in Clinical Treatment

GSEA was used to analyze the enrichment of biological processes (BP) in the high-risk and low-risk groups, using GSEA 4.0.3 software for GSEA analysis. There is a strong relationship between the effectiveness of immunotherapy and the tumor immune microenvironment. Therefore, we calculated the correlation between the risk scores and the MCP algorithm for the immune score of immune cells using the Spearman method. We further explored the correlation of risk scores with the expression of representative genes of immune checkpoints, DNA replication, mismatch repair, and epithelial-mesenchymal transition (EMT).

Results

Correlation Between Immune Scores and KIRC Subtypes

A total of 530 KIRC samples from TCGA and 39 samples from GSE29609 were included after pre-treatment (Table 1). The TME-related transcriptomic data downloaded from the TCGA database for KIRC were analyzed for differential expression in cancer versus para cancer groupings, and the results are shown in (Figures 1A, B and Table S1). In these TME related DEGs, the NMF algorithm is used to find that the clustering result is the best when k=2, so we build two clusters: Cluster1 and Cluster2 (Figures 1C and S1). Cophenet index is a measure of correlation between the distance of points in the feature space and the distance on the tree graph. Usually, it obtains all possible pairs of points in the data and calculates the Euclidean distance between these points. With this analysis, we determined that the C1 and C2 groupings are two independent clusters (Figures 1D, E). A comparison of OS and PFS between the two clusters showed significant differences (Figures 1F, G), with patients in the Cluster2 subtype having a longer survival time and a better prognosis. Our molecular typing results were compared with six immune subtypes of international translational immune typology of solid tumors. The results are shown in Figure 1H. To explore the relationship between the immune infiltration of KIRC and each subtype, the immune scores of individual immune cells were calculated for each sample by the MCP-counter algorithm and were compared.

TABLE 1

Table 1 Clinicopathological characteristics of KIRC patients from the TCGA and GEO databases.

FIGURE 1

Figure 1 Two clusters were identified based on TME-related differentially expressed genes. (A, B) Volcano and heatmap of the distribution of differentially expressed genes, with P < 0.05 and |log2 FC| > 1.0 as cut-off values. (C) Consensus map of 530 cases of KIRC via the NMF algorithm. (D, E) Factorization rank for k = 2–10. (F, G) The Kaplan‐Meier survival curve showed the OS and PFS of the two subtypes of patients. (H) Percentage of the four immune subtypes accounting for each of the two clusters. TME, tumor microenvironment; KIRC, kidney renal clear cell carcinoma; NMF, non-negative matrix factorization; OS, overall survival; PFS, progression-free survival.

The results showed that patients in the Cluster1 group had higher immune scores of CD8+ T cells, cytotoxic lymphocytes, fibroblasts, monocytic lineage, myeloid dendritic cells, T cells (P<0.05, Figures 2A–G), and lower immune scores for neutrophils and endothelial cells (P<0.05, Figures 2H–I).

FIGURE 2

Figure 2 Differences in the distribution of immune cell content via the MCP-counter algorithm between the two clusters of KIRC. (A–I) Patients in the C1 group had higher immune scores of CD8+ T cells, cytotoxic lymphocytes, fibroblasts, monocytic lineage, myeloid dendritic cells, T cells, and lower immune scores for neutrophils and endothelial cells (P < 0.05). MCP-counter, microenvironment cell populations-counter; KIRC, kidney renal clear cell carcinoma.

Establishment and Validation of the TME-Related Signature

The entire TCGA cohort of 530 cases was randomly divided into the most appropriate training cohort (374 samples) and the test cohort (156 samples) in a 3:7 subgroup sampling ratio and the two groups were similar in terms of clinicopathological characteristics (Age, gender, grade, stage, TNM staging), as shown in Table 2. In the TCGA training cohort, the above DEGs were subjected to univariate Cox regression, and P ≤ 0.05 was set to obtain genes significantly associated with the prognosis of KIRC patients. A total of 8 genes were identified for the construction of the prognostic signature (Figures 3A, B) by removing over-fitted genes through LASSO regression analysis. These genes include LRFN1, CSF1, UCN, TUBB2B, SERPINF1, ADAM8, ABCB4, CCL22. The main biological functions involved include chemotaxis of some non-characterized immune cells as well as lymphocytes. Most tissue macrophages and osteoclasts are regulated by colony-stimulating factor-1 (CSF-1, also known as macrophage CSF). LRFN (leucine-rich repeat and fibronectin type-III domain-containing protein) recognizes bacteria and promotes hemocytic phagocytosis. CCL22 expression by dendritic cells (DCs) promotes the formation of cell–cell contacts and interaction with regulatory T cells (T reg) through their CCR4 receptor. These genes also include some of the genes that make up the basic cellular units, including Tubulin and the extracellular matrix. The TUBB2B gene provides instructions for making one version of a protein called beta-tubulin (β-tubulin). This protein is part of the tubulin family of proteins that form and organize cell structures called microtubules.

TABLE 2

Table 2 Comparison of TCGA training and testing cohort.

FIGURE 3

Figure 3 Construction of the TME-related signature for KIRC. (A) LASSO coefficients of the prognostic TME-related gene. (B) The LASSO model parameters (λ) were selected using a minimum criterion of ten cross-validations. (C–F) ROC curves for prognostic models of the entire TCGA cohort, TCGA training cohort, TCGA testing cohort, and GSE29609. TME, tumor microenvironment; KIRC, kidney renal clear cell carcinoma; LASSO, least absolute1 shrinkage and selection operator; ROC, receiver operating characteristic.

Based on the expression of each gene and the corresponding coefficients, the following formula was constructed to calculate the risk score: Risk score = (LRFN1 * 0.433) + (CSF1 * 0.448) + (UCN * 0.376) + (TUBB2B * 0.190) + (SERPINF1 * 0.201) + (ADAM8 * 0.271) + (ABCB4 * -0.579) + (CCL22 * -1.055). To assess the accuracy of our signature predictions, we introduced a time-dependent ROC curve analysis, which showed that the prognostic signature showed low accuracy in GSE29609 (Figure 3F). However, the entire TCGA cohort, TCGA training cohort, and TCGA testing cohort showed higher accuracy, with AUC values almost all greater than 0.7 in predicting 1-, 3-, and 5-year survival (Figures 3C–E). Prognostically, as shown in Figures 4C–E, patients at higher risk had a worse prognosis compared with those in the lower risk group (P<0.001), and we obtained the same result in the GSE29609 dataset (Figure 4F). To provide clinicians with a quantitative method to predict the prognosis of patients with KIRC, we further created a nomogram (Figure 4A), and we assessed the predictive efficacy of the nomogram by calibration plot (Figure 4B). Based on the calibration plot, we can find that the 1-, 3-, and 5-year survival rates predicted by the nomogram overlap well with the actual survival rates of KIRC patients, which indicates that the nomogram can better predict the survival rates of KIRC patients.

FIGURE 4

Figure 4 Establishment and evaluation of nomogram. (A) Nomogram was constructed by the prognostic model risk score, patient age, gender, stage, grade, and TNM Staging based on multivariate regression analysis. (B) 1-year, 3-year and 5-year calibration curves. (C-F) Kaplan-Meier survival curves for prognostic models of entire TCGA cohort, TCGA training cohort, TCGA testing cohort, and GSE29609.

The TME-Related Signature Was an Independent Predictor of Prognosis

The prognostic significance of different clinical characteristics was assessed using univariate and multivariate Cox regression analysis, and the results showed that risk score was an independent predictor of prognosis (Table 3). The ROC and DCA analysis found that the signature and nomogram were better predictors of patient prognosis than traditional clinical characteristics (Figure 5). Further stratification showed that risk scores increased with disease progression (Figures 6A–D) and were not affected by Stage (I-II and III-IV) (Figures 6E, F), and these results are sufficient to demonstrate the accuracy and stability of our signature.

TABLE 3

Table 3 Univariable and multivariable Cox regression to analyze the relationship between the RS and clinical prognosis.

FIGURE 5

Figure 5 Independent predictive analysis of TME related characteristics. In the whole TCGA cohort, ROC curves of our characteristics and other clinical indicators were analyzed for 1 year, 3 years and 5 years respectively (Figures A–C). DCA curve analysis of our characteristics and other clinical indicators, the time is also one year, three years and five years (Figures D, E). Tumor microenvironment; ROC, receiver operating characteristics; Decision curve analysis.

FIGURE 6

Figure 6 Analysis of the TME-related signature stratified by TNM staging and Stage staging. (A–D) Risk scores also tend to be higher in KIRC patients with combined T3-4, N1, M1, Stage III-IV (P < 0.05). (E, F) Patients with low-risk scores all had a better prognosis in the StageI-II and Stage III-IV groups (P < 0.05). TME, tumor microenvironment; KIRC, Kidney renal clear cell carcinoma.

Comparison of Predictive Performance Between the TME-Related Signature and Previously Published Signatures

To highlight the predictive performance of our 8-gene prognostic signature, three other established risk signatures were selected for comparison. The results showed that the AUC value of our signature was greater than those of the other three established prediction signatures we selected (Figures 7A–D) in predicting 1-, 3- and 5-year OS and that high-risk patients had a worse prognosis than low-risk patients (P<0.05, Figures 7E–H). The C-index of each model was calculated to evaluate the predictive power of the signature, and our signature had the highest C-index of 0.744 (Figure 7I). The restricted mean survival (RMS) results show that the model we have built performs best over a time of more than 8 years compared to the signatures built by Chen/Zhu (Figure 7J). The results demonstrate the robust predictive performance of our new signature.

FIGURE 7

Figure 7 Performance of the TME-related signature compared with other prognostic signatures in predicting OS in KIRC patients. (A–D) ROC curves were used to predict 1-, 3-, and 5-year OS in KIRC patients with our signature and three other published gene signatures. (E–H) Kaplan-Meier survival curves of KIRC patients with our signature and three other published gene signatures. (I) The C-index of the four prognostic signatures including our signature, has the highest C-index. (J) The RMS time curves for our prognostic signature and the Chen/Zhu prognostic signatures showed an overlap of 8 years. TME, tumor microenvironment; OS, overall survival; KIRC, kidney renal clear cell carcinoma; C-index, concordance index; RMS, restricted mean survival.

GSEA and Immune Correlation Analyses of the TME-Related Signature

GSEA analysis showed that high-risk group genes were mainly enriched in immune-related biological processes, such as activation of the immune response, adaptive immune response based on somatic ecombination of immune receptors built from immunoglobulin superfamily domains, antigen receptor mediated signalling pathway, B cell activation, B cell-mediated immunity (Figure 8A), the low-risk group genes were mainly enriched in fatty acid metabolic process, spliceosomal snrnp assembly, spliceosomal tri snrnp complex assembly, apical part of cell, anion transmembrane transporter activity (Figure 8B). Immunotherapy is now a treatment strategy with great potential after targeted therapy for patients with KIRC, and the tumor immune microenvironment is an important factor influencing the response to immunotherapy (20, 21). We used the MCP-counter algorithm to estimate immune cell scores for each patient with KIRC and then calculated their correlation with risk scores by Pearson’s method. The results showed that the risk score was significantly positively correlated with B lineage, Monocytic lineage, and Fibroblasts, while significantly negatively correlated with Neutrophils and Endothelial cells (Figures 8C–E). Interestingly, the risk score was also positively correlated with immune checkpoint-related genes (PDCD1, CTLA4), DNA replication-related genes (POLE2, FEN1, MCM6), and EMT-related genes (FAP, LOXL2) (Figure 8D).

FIGURE 8

Figure 8 GSEA and immune correlation analyses of the TME-related signature. (A, B) GO-BP analyses for high- and low-risk groups. (C) Correlation between the risk score for our signature, immune cell score via the MCP-counter algorithm, and TMB (Pearson method). (D) The relationship between the risk score and expression of representative genes for immune checkpoints, DNA replication, mismatch repair, and epithelial-mesenchymal transition (Pearson method). (E) The relationship between the risk score and tumor immune infiltrating cells (Pearson method). GSEA, gene set enrichment analysis; TME, tumor microenvironment; GO, gene ontology; BP, biological process; MCP-counter, microenvironment cell populations-counter; TMB, tumor mutation burden.

Discussion

Currently, the insidious onset of RCC and the lack of specific clinical signs in the early stages result in 20-30% of patients having metastases by the time of presentation (22). As KIRC has been studied more intensively, immunotherapy has been found to significantly improve the prognosis of patients (23–25). Results of completed clinical trials have shown various other immunomodulatory pathways in the tumor immune microenvironment, which can affect tumor cell survival and significantly influence immunotherapeutic response (26–28). Due to the heterogeneity of KIRC (29, 30) and immune relevance (31), the genetic information of TME is crucial to optimize the therapeutic approach of KIRC. In this study, we screened the differential genes based on the TME gene profile of KIRC And a prognostic index was constructed by association of these differential genes. This indicator is superior to those commonly used in the clinic and provides a prognostic method with lower sampling cost and better reproducibility for personalized treatment of patients.

In this study, we systematically analyzed 4061 TME-related genomic data from 537 KIRC cases in TCGA data to construct two clusters. The two clusters differed significantly in immune score comparisons for most tumor immune cells, as well as differences in survival prognostic analysis, which are consistent with the heterogeneity of the TME. The validated prognostic signature was constructed from eight TME-related genes (LRFN1, CSF1, UCN, TUBB2B, SERPINF1, ADAM8, ABCB4, CCL22) identified by univariate regression analysis and LASSO analysis. Patients were divided into a low-risk group and a high-risk group according to their risk scores. The difference in survival rates between the two groups was statistically significant, with the higher the risk score, the worse the prognosis of the patients. The stratified analysis revealed that the risk values increased with disease progression. The AUC values of our prognostic signature were almost >0.7, confirming its effectiveness in predicting the prognosis of KIRC patients. CSF1 acts as a tumour-promoting cytokine by recruiting macrophages to the tumour area, which in turn leads to the release of various tumour-promoting growth factors in the microenvironment (32). In vitro and in vivo studies have confirmed that CSF1 promotes proliferation and reduces apoptosis in renal cell carcinoma (33). TUBB2B alterations can lead to increased incidence of kidney disease and KIRC (34). Studies have shown that ADAM8 is highly expressed in a variety of malignant tumours and is strongly associated with tumour metastasis and poor patient prognosis (35, 36).

The nomogram is a practical and intuitive evaluation tool. The establishment of the nomogram can more accurately predict the prognosis of KIRC patients, and the ROC, DCA, and calibration curves showed the validity of the nomogram. Univariate and multifactorial Cox regression analyses combined with clinicopathological parameters showed this prognostic signature to be an independent prognostic factor for KIRC patients. When we compared our risk prognostic signature with the prognostic models developed by Yang et al (17), Chen et al (18), and Zhu et al (19), the overall performance of our signature was better than these three signatures.

Further biological process GSEA enrichment analysis was performed on patients in the high-risk and low-risk groups. The high-risk group was mainly enriched for several immune response-related biological processes. Interestingly, immune correlation analysis showed that risk values were positively correlated with B lineage, Monocytic lineage, Fibroblasts, and negatively correlated with Neutrophils, Endothelial cells. Monocytes are also involved in the immune response and will lose control of the tumour as it progresses, and the body’s defence system becomes less functional. Moreover, as a precursor to macrophages, high-risk patients have significantly higher levels of M0 macrophages relative to low-risk patients with KIRC (37). In recent years, immune checkpoints have received increasing attention in the study of tumour immunotherapy (38, 39). Kidney cancer is a highly immunogenic tumor in which the tumor cells produce an immunosuppressive environment through a variety of mechanisms, such as increased expression of immunosuppressive molecules(PDCD1, CD274, CTLA4) (40–42) in the TME and the occurrence of immune escape (43). Our analysis also confirmed that PDCD1 and CTLA4 expression levels were upregulated with increasing risk values, which may cause poor prognosis in patients at high risk of KIRC. In addition, our risk values were strongly and positively correlated with DNA replication-related genes (POLE2, FEN1, MCM6) and EMT-related genes (FAP, LOXL2). DNA replication-related genes play a crucial role in DNA replication and cell cycle regulation (44), and malfunctioning expression of these proteins has the potential to promote cell proliferation and tumorigenesis (45, 46). Previous studies have also shown that POLE2 and MCM6 are overexpressed in KIEC tissue and that high levels of POLE2 and MCM6 expression are associated with poor prognosis (47, 48). These immunosuppressive molecules and cell cycle regulatory genes may be involved in the development and progression of KIRC, informing its diagnosis and regression.

The Stage parameters need to be obtained from multiple imaging results and pathological staining and scored by a professional, our gene signature only requires the scoring of a specific gene after sequencing. This scoring process does not require the involvement of a professional and requires only simple calculations, thus providing better reproducibility and accessibility. And the gene information can be obtained from a small amount of tissue or even exosomes, thus greatly reducing the cost of obtaining these prognostic indicators. From the view in Figure 5D, the area under the AUC curve of our prognostic model is superior to all other clinically used parameters. The net benefit is also better than the stage indicator in the case of low risk. Therefore, our model is more advantageous and worthy of physicians’ choice. Patients in group C1 had higher immune scores for CD8+ T cells, cytotoxic lymphocytes, fibroblasts, monocyte lineage, myeloid dendritic cells, and T cells compared to group C2, while neutrophils and endothelial cells had lower immune scores instead, suggesting that the two subgroups C1 and C2 correspond to different TME immunophenotypes, respectively, and are potential indicators of immunophenotyping.

Admittedly, the current study still has some limitations that must be considered. This study is a retrospective study based on a public database and we need to validate these results in prospective cohorts, for animal models and at the cellular level, which will be the focus of our follow-up studies.

Conclusion

In conclusion, this study successfully constructed a TME-related prognostic signature that accurately predicted the survival prognosis of KIRC patients. The nomogram established in combination with the risk score and other clinicopathological parameters can also individually predict patient survival. These results could provide the basis for future studies on potential individualized treatments for KIRC patients.

Data Availability Statement

The origina contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

Author Contributions

DP and CX contributed equally to this work. DP, CX and HZ designed this work. DW, XS and YL integrated and analyzed the data. DP, YZ, JG and NL wrote this manuscript. HZ and CX edited and revised the manuscript. All authors approved this manuscript.

Funding

This study was funded by Henan Provincial Health Commission (LHGJ20190422), Henan Provincial Health Commission (LHGJ20210486) and the Key Scientific Item of Henan Province Education Department (21A320037).

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.

The reviewer XJ declared a shared parent affiliation with the authors to the handling editor at the time of review.

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/fonc.2022.912155/full#supplementary-material

Abbreviations

TME, Tumor microenvironment; KIRC, Renal clear cell carcinoma; NMF, non-negative matrix factorization; GSEA, Gene set enrichment analysis; RCC, Renal cell carcinoma; ECM, extracellular matrix; TAMs, tumour-associated macrophages; ICIs, immune checkpoint inhibition; TCGA, Cancer Genome Atlas; GEO, Gene-Expression Omnibus; OS, Overall survival; MCP-counter, microenvironment cell populations-counter; PFS, Progression-free survival; ROC, receiver operating characteristic; DCA, decision curve analysis; C-index, concordance index; BP, Biological process; EMT, epithelial-mesenchymal transition; GO, Gene Ontology; RMS, restricted mean survival; TMB, Tumor mutation burden

References

1. 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(6):394–424. doi: 10.3322/caac.21492

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Hsieh JJ, Purdue MP, Signoretti S, Swanton C, Albiges L, Schmidinger M, et al. Renal Cell Carcinoma. Nat Rev Dis Primers (2017) 3:17009. doi: 10.1038/nrdp.2017.9

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Siegel RL, Miller KD, Jemal A. Cancer Statistics, 2019. CA Cancer J Clin (2019) 69(1):7–34. doi: 10.3322/caac.21551

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Zhu J, Armstrong AJ, Friedlander TW, Kim W, Pal SK, George DJ, et al. Biomarkers of Immunotherapy in Urothelial and Renal Cell Carcinoma: PD-L1, Tumor Mutational Burden, and Beyond. J Immunother Cancer (2018) 6(1):4. doi: 10.1186/s40425-018-0314-1

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Anderson NM, Simon MC. The Tumor Microenvironment. Curr Biol (2020) 30(16):R921–r925. doi: 10.1016/j.cub.2020.06.081

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Díaz-Montero CM, Rini BI, Finke JH. The Immunology of Renal Cell Carcinoma. Nat Rev Nephrol (2020) 16(12):721–35. doi: 10.1038/s41581-020-0316-3

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring Tumour Purity and Stromal and Immune Cell Admixture From Expression Data. Nat Commun (2013) 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Chevrier S, Levine JH, Zanotelli VRT, Silina K, Schulz D, Bacac M, et al. An Immune Atlas of Clear Cell Renal Cell Carcinoma. Cell (2017) 169(4):736–749.e718. doi: 10.1016/j.cell.2017.04.016

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Michot JM, Bigenwald C, Champiat S, Collins M, Carbonnel F, Postel-Vinay S, et al. Immune-Related Adverse Events With Immune Checkpoint Blockade: A Comprehensive Review. Eur J Cancer (2016) 54:139–48. doi: 10.1016/j.ejca.2015.11.016

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Li S, Xu W. Mining TCGA Database for Screening and Identification of Hub Genes in Kidney Renal Clear Cell Carcinoma Microenvironment. J Cell Biochem (2019). doi: 10.1002/jcb.29511

CrossRef Full Text | Google Scholar

11. Chen WJ, Cao H, Cao JW, Zuo L, Qu FJ, Xu D, et al. Heterogeneity of Tumor Microenvironment is Associated With Clinical Prognosis of non-Clear Cell Renal Cell Carcinoma: A Single-Cell Genomics Study. Cell Death Dis (2022) 13(1):50. doi: 10.1038/s41419-022-04501-9

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Zheng M, Long J, Chelariu-Raicu A, Mullikin H, Vilsmaier T, Vattai A, et al. Identification of a Novel Tumor Microenvironment Prognostic Signature for Advanced-Stage Serous Ovarian Cancer. Cancers (Basel) (2021) 13(13):3343. doi: 10.3390/cancers13133343

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Becht E, Giraldo NA, Lacroix L, Buttard B, Elarouci N, Petitprez F, et al. Estimating the Population Abundance of Tissue-Infiltrating Immune and Stromal Cell Populations Using Gene Expression. Genome Biol (2016) 17(1):218. doi: 10.1186/s13059-016-1070-5

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Rooney MS, Shukla SA, Wu CJ, Getz G, Hacohen N. Molecular and Genetic Properties of Tumors Associated With Local Immune Cytolytic Activity. Cell (2015) 160(1-2):48–61. doi: 10.1016/j.cell.2014.12.033

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Gaujoux R, Seoighe C. A Flexible R Package for Nonnegative Matrix Factorization. BMC Bioinf (2010) 11:367. doi: 10.1186/1471-2105-11-367

CrossRef Full Text | Google Scholar

16. Zheng H, Li Z, Li J, Zheng S, Zhao E. Construction, Validation, and Visualization of Two Web-Based Nomograms to Predict Overall and Cancer-Specific Survival in Patients With Gastric Cancer and Lung Metastases. J Oncol (2021) 2021:5495267. doi: 10.1155/2021/5495267

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Yang X, Miao Y, Wang J, Mi D. A Pan-Cancer Analysis of the HER Family Gene and Their Association With Prognosis, Tumor Microenvironment, and Therapeutic Targets. Life Sci (2021) 273:119307. doi: 10.1016/j.lfs.2021.119307

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Chen S, Wang C, Su X, Dai X, Li S, Mo Z. KCNN4 is a Potential Prognostic Marker and Critical Factor Affecting the Immune Status of the Tumor Microenvironment in Kidney Renal Clear Cell Carcinoma. Transl Androl Urol (2021) 10(6):2454–70. doi: 10.21037/tau-21-332

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Zhu J, Xiao J, Wang M, Hu D. Pan-Cancer Molecular Characterization of M(6)A Regulators and Immunogenomic Perspective on the Tumor Microenvironment. Front Oncol (2020) 10:618374. doi: 10.3389/fonc.2020.618374

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Peng YL, Xiong LB, Zhou ZH, Ning K, Li Z, Wu ZS, et al. Single-Cell Transcriptomics Reveals a Low CD8(+) T Cell Infiltrating State Mediated by Fibroblasts in Recurrent Renal Cell Carcinoma. J Immunother Cancer (2022) 10(2):e004206. doi: 10.1136/jitc-2021-004206

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Ged Y, Voss MH. Novel Emerging Biomarkers to Immunotherapy in Kidney Cancer. Ther Adv Med Oncol (2021) 13:17588359211059367. doi: 10.1177/17588359211059367

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Hahn AW, Klaassen Z, Agarwal N, Haaland B, Esther J, Ye XY, et al. First-Line Treatment of Metastatic Renal Cell Carcinoma: A Systematic Review and Network Meta-Analysis. Eur Urol Oncol (2019) 2(6):708–15. doi: 10.1016/j.euo.2019.09.002

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Tannir NM, Papadopoulos KP, Wong DJ, Aljumaily R, Hung A, Afable M, et al. Pegilodecakin as Monotherapy or in Combination With Anti-PD-1 or Tyrosine Kinase Inhibitor in Heavily Pretreated Patients With Advanced Renal Cell Carcinoma: Final Results of Cohorts A, G, H and I of IVY Phase I Study. Int J Cancer (2021) 149(2):403–8. doi: 10.1002/ijc.33556

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Labadie BW, Balar AV, Luke JJ. Immune Checkpoint Inhibitors for Genitourinary Cancers: Treatment Indications, Investigational Approaches and Biomarkers. Cancers (Basel) (2021) 13(21):5415. doi: 10.3390/cancers13215415

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Jang A, Adler DM, Rauterkus GP, Bilen MA, Barata PC. Immunotherapies in Genitourinary Oncology: Where Are We Now? Where Are We Going? Cancers (Basel) (2021) 13(20):5065. doi: 10.3390/cancers13205065

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Lei X, Lei Y, Li JK, Du WX, Li RG, Yang J, et al. Immune Cells Within the Tumor Microenvironment: Biological Functions and Roles in Cancer Immunotherapy. Cancer Lett (2020) 470:126–33. doi: 10.1016/j.canlet.2019.11.009

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Jian Y, Yang K, Sun X, Zhao J, Huang K, Aldanakh A, et al. Current Advance of Immune Evasion Mechanisms and Emerging Immunotherapies in Renal Cell Carcinoma. Front Immunol (2021) 12:639636. doi: 10.3389/fimmu.2021.639636

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Ho WJ, Jaffee EM, Zheng L. The Tumour Microenvironment in Pancreatic Cancer - Clinical Challenges and Opportunities. Nat Rev Clin Oncol (2020) 17(9):527–40. doi: 10.1038/s41571-020-0363-5

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Gerlinger M, Rowan AJ, Horswell S, Math M, Larkin J, Endesfelder D, et al. Intratumor Heterogeneity and Branched Evolution Revealed by Multiregion Sequencing. N Engl J Med (2012) 366(10):883–92. doi: 10.1056/NEJMoa1113205

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Hanahan D, Weinberg RA. Hallmarks of Cancer: The Next Generation. Cell (2011) 144(5):646–74. doi: 10.1016/j.cell.2011.02.013

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Ricketts CJ, De Cubas AA, Fan H, Smith CC, Lang M, Reznik E, et al. The Cancer Genome Atlas Comprehensive Molecular Characterization of Renal Cell Carcinoma. Cell Rep (2018) 23(1):313–326.e315. doi: 10.1016/j.celrep.2018.03.075

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Lin W, Xu D, Austin CD, Caplazi P, Senger K, Sun Y, et al. Function of CSF1 and IL34 in Macrophage Homeostasis, Inflammation, and Cancer. Front Immunol (2019) 10:2019. doi: 10.3389/fimmu.2019.02019

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Wang S, Yang X, Xie W, Fu S, Chen Q, Li Z, et al. LncRNA GAPLINC Promotes Renal Cell Cancer Tumorigenesis by Targeting the miR-135b-5p/CSF1 Axis. Front Oncol (2021) 11:718532. doi: 10.3389/fonc.2021.718532

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Jeruschke S, Jeruschke K, DiStasio A, Karaterzi S, Büscher AK, Nalbant P, et al. Everolimus Stabilizes Podocyte Microtubules via Enhancing TUBB2B and DCDC2 Expression. PloS One (2015) 10(9):e0137043. doi: 10.1371/journal.pone.0137043

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Schlomann U, Koller G, Conrad C, Ferdous T, Golfi P, Garcia AM, et al. ADAM8 as a Drug Target in Pancreatic Cancer. Nat Commun (2015) 6:6175. doi: 10.1038/ncomms7175

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Conrad C, Benzel J, Dorzweiler K, Cook L, Schlomann U, Zarbock A, et al. ADAM8 in Invasive Cancers: Links to Tumor Progression, Metastasis, and Chemoresistance. Clin Sci (Lond) (2019) 133(1):83–99. doi: 10.1042/CS20180906

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Hua X, Chen J, Su Y, Liang C. Identification of an Immune-Related Risk Signature for Predicting Prognosis in Clear Cell Renal Cell Carcinoma. Aging (Albany NY) (2020) 12(3):2302–32. doi: 10.18632/aging.102746

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Saleh R, Toor SM, Sasidharan Nair V, Elkord E. Role of Epigenetic Modifications in Inhibitory Immune Checkpoints in Cancer Development and Progression. Front Immunol (2020) 11:1469. doi: 10.3389/fimmu.2020.01469

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Rotte A, Jin JY, Lemaire V. Mechanistic Overview of Immune Checkpoints to Support the Rational Design of Their Combinations in Cancer Immunotherapy. Ann Oncol (2018) 29(1):71–83. doi: 10.1093/annonc/mdx686

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Qin S, Xu L, Yi M, Yu S, Wu K, Luo S. Novel Immune Checkpoint Targets: Moving Beyond PD-1 and CTLA-4. Mol Cancer (2019) 18(1):155. doi: 10.1186/s12943-019-1091-2

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Woo SR, Turnis ME, Goldberg MV, Bankoti J, Selby M, Nirschl CJ, et al. Immune Inhibitory Molecules LAG-3 and PD-1 Synergistically Regulate T-Cell Function to Promote Tumoral Immune Escape. Cancer Res (2012) 72(4):917–27. doi: 10.1158/0008-5472.CAN-11-1620

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Casey SC, Tong L, Li Y, Do R, Walz S, Fitzgerald KN, et al. MYC Regulates the Antitumor Immune Response Through CD47 and PD-L1. Science (2016) 352(6282):227–31. doi: 10.1126/science.aac9935

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Dunn GP, Bruce AT, Ikeda H, Old LJ, Schreiber RD. Cancer Immunoediting: From Immunosurveillance to Tumor Escape. Nat Immunol (2002) 3(11):991–8. doi: 10.1038/ni1102-991

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Lei Y, Wang S, Liu J, Yan W, Han P, Tian D. Identification of MCM Family as Potential Therapeutic and Prognostic Targets for Hepatocellular Carcinoma Based on Bioinformatics and Experiments. Life Sci (2021) 272:119227. doi: 10.1016/j.lfs.2021.119227

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Issac MSM, Yousef E, Tahir MR, Gaboury LA. MCM2, MCM4, and MCM6 in Breast Cancer: Clinical Utility in Diagnosis and Prognosis. Neoplasia (2019) 21(10):1015–35. doi: 10.1016/j.neo.2019.07.011

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Zhang P, Chen X, Zhang L, Cao D, Chen Y, Guo Z, et al. POLE2 Facilitates the Malignant Phenotypes of Glioblastoma Through Promoting AURKA-Mediated Stabilization of FOXM1. Cell Death Dis (2022) 13(1):61. doi: 10.1038/s41419-021-04498-7

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Zhang C, Shen Y, Gao L, Wang X, Huang D, Xie X, et al. Targeting POLE2 Creates a Novel Vulnerability in Renal Cell Carcinoma via Modulating Stanniocalcin 1. Front Cell Dev Biol (2021) 9:622344. doi: 10.3389/fcell.2021.622344

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Jang NR, Baek J, Ko Y, Song PH, Gu MJ. High MCM6 Expression as a Potential Prognostic Marker in Clear-Cell Renal Cell Carcinoma. In Vivo (2021) 35(1):299–306. doi: 10.21873/invivo.12259

PubMed Abstract | CrossRef Full Text | Google Scholar