Identification of Survival Risk and Immune-Related Characteristics of Kidney Renal Clear Cell Carcinoma

Background Immunity exerts momentous functions in the progression and treatment of kidney renal clear cell carcinoma (KIRC). A better understanding of the relationship between KIRC and immunity may make a great contribution to evaluating the prognosis and immune-related therapeutic response of KIRC. Methods A series of information such as RNA sequence, clinical data, and tumor mutation burden (TMB) of KIRC patients were downloaded through The Cancer Genome Atlas (TCGA). Next, combining the survival information and gene expression data of TCGA and Gene Expression Omnibus (GEO), we established an immune gene-related prognosis model (IGRPM) and analyzed it. Then we constructed a nomogram which was convenient for clinicians to judge the prognosis of KIRC. Last but not the least, the expressions of some genes used to construct IGRPM in early KIRC, and adjacent normal tissues were verified through real-time fluorescence quantitative polymerase chain reaction (RT-qPCR). Perl (strawberry-perl-5.30.0.1-64bit), R software (4.0.3), and GraphPad Prism 7 were used to process the relevant data. Results The single-sample gene set enrichment analysis (ssGSEA) showed that there were significant differences in StromalScore, ImmuneScore, ESTIMATEScore, TumorPurity, 22 kinds of human immune cells infiltration, and HLA genes expression between high immunity group (Immunity_H) and low immunity group (Immunity_L). The Immunity_H expressed more immune-related genes and enriched more immune-related functions than the Immunity_L. In addition, compared with the low-risk group, the high-risk group had worse survival outcome and higher TMB. Combining IGRPM-based risk characteristic and TMB, we found that low-TMB + low-risk was the most beneficial to the survival outcome of KIRC patients. The risk characteristic based on IGRPM could be used as an independent prognostic factor for KIRC, and the nomogram constructed for evaluating the prognosis of KIRC showed excellent predictive potential. The RT-qPCR results suggested that not all the genes used to construct IGRPM showed differential expression in early KIRC compared with adjacent normal tissues, but all these genes had significant influence on the prognosis of KIRC. Conclusion These comprehensive immune assessments and survival predictions, integrating multiple aspects of data and clinical information, can provide additional value to the current Tumor Node Metastasis staging system for risk stratification of KIRC and may facilitate the development of KIRC immunotherapy.


Introduction
Renal cell carcinoma (RCC), which originates from nephron [1], is one of the 10 most common malignancies in humans [2]. However, kidney renal clear cell carcinoma (KIRC) is one of the most invasive pathological types of RCC and the main pathological type [3], accounting for about 75% of all RCC [4]. Moreover, KIRC is a crucial cause of cancer-related death worldwide [5]. It is important to analyze the clinicopathological features of KIRC patients according to the WHO classification and Tumor Node Metastasis (TNM) staging for selecting appropriate treatment methods [6]. However, there may be significant differences in prognosis even in patients with the same TNM staging and pathological grade after similar treatment [7,8], indicating that KIRC is a heterogeneous tumor, and the existing clinicopathological features and classification are not sufficient to evaluate the prognosis and risk stratification of KIRC patients. Therefore, finding a new method with high predictive value is of great significance to improve the prognosis of KIRC.
The development and occurrence of tumor and the effectiveness of immunotherapy are closely related to tumor immune microenvironment (TIM) [9,10]. More and more studies have emphasized that there is heterogeneity in the proportion of immune-related genes and immune cell populations in tumors [11,12], which may be one of the key reasons for different prognosis of patients. Oncologists are committed to exploring the immune responses and its regulatory mechanisms and screening immune-related biomarkers in order to coordinate the immune microenvironment of tumors and improve survival prognosis [13]. Immune cells are related to tumor metastasis and invasion and are potentially important forces to slow down or inhibit tumor growth [14]. However, although the research on the relationship between KIRC and immunity has increased in recent years, there is still a lack of research that can accurately predict the survival of KIRC patients in terms of immunity.
With the sharing of massive data, bioinformatics analysis has been widely used in the study of tumor prognosis [15]. In this study, we constructed an immune gene-related prognosis model (IGRPM) for patients with KIRC using data from The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO). In addition, we systematically understood and distinguished the immune-related characteristics of KIRC, and looked into the correlation of tumor immune-related genes with TIM in depth. Additionally, we constructed a nomogram which was convenient for clinical application and verified its similarity with the ideal model. Finally, we verified the expressions of some of these genes used to construct IGRPM in early KIRC and adjacent normal tissues through real-time fluorescence quantitative polymerase chain reaction (RT-qPCR). Moreover, we confirmed that all the genes used to construct IGRPM had a significant effect on the prognosis of KIRC through UALCAN database.

Materials and Methods
2.1. Data Source. We downloaded transcriptome and clinical data of KIRC from TCGA (https://portal.gdc.cancer.gov/; December, 2021) and GEO (https://www.ncbi.nlm.nih.gov/ geo/; December, 2021; GSE29609). The transcriptome data (HTseq-FPKM) of TCGA included 539 tumor samples and 72 paracancerous normal samples. 537 KIRC samples were downloaded from TCGA including clinical data. Masked somatic mutation data for 339 KIRC samples were downloaded from TCGA, of which tumor mutation burden (TMB) was successfully calculated for 336 samples. The GSE29609 cohort including 39 KIRC samples from GEO database and 530 KIRC samples from TCGA were used to construct prognostic models. During the data processing, we deleted some samples with incomplete information. The genetic data of related samples were converted into matrix files by Perl software (strawberry-perl-5. 30.0.1-64bit). We analyzed the data through R software (4.0.3).

Single-Sample Gene Set Enrichment Analysis (ssGSEA)
and Immune Classification Analysis. Using the ssGSEA algorithm of R software package ("GSEABase," "limma," and "GSVA"), the immune-related characteristics of each sample in TCGA-KIRC were evaluated comprehensively, based on 29 immune gene sets, including different immune cells, immune-related functions, and immune-related genes [16]. According to the immunological characteristics of the samples, 539 tumor samples from TCGA were analyzed by Ward's linkage and Euclidean distance, which were divided into high immunity group (Immunity_H) and low immunity group (Immunity_L) [17]. The accuracy of immunotyping of KIRC patients was verified by R package ("Rtsne").

Estimation of Tumor Purity and Immune
Characteristics. We estimated the infiltration of stromal cells and immune cells in 539 KIRC tumor samples from TCGA by R software package ("estimate"), thus determining the StromalScore, ImmuneScore, ESTIMATEScore, and Tumor-Purity of the samples [18]. We also calculated the abundance of 22 types of human immune cells in each sample by CIBERSORT method [19] and compared the cell composition of the two subtypes in 539 KIRC samples of TCGA.

Differentially Expressed
Immune Genes in Immunity_H and Immunity_L. We identified the differentially expressed genes (DEGs) between Immunity_H and Immunity_L by "limma" package [20]. Then, the DEGs were intersected with the immunologically relevant list of genes downloaded from ImmPort to obtain the differentially expressed immune genes (DEIGs) [21]. In addition, we also compared the differences in the expressions of HLA-related genes in Immu-nity_H and Immunity_L.
2.5. Functional Enrichment Analysis. We performed gene set enrichment analysis of 539 KIRC patients from TCGA by GSEA (Version4.1.0), in order to reveal the key signal pathways and functions involved in the two subtypes [22]. The screening thresholds of KEGG pathway and GO analysis in the Immunity_H and Immunity_L were p < 0:05. The enriched pathways were visualized with R packages ("reshape" and "ggplot2").

2
Journal of Immunology Research 2.6. Analyses of Prognosis Associated Immune Genes and Regulatory Networks. We obtained the DEIGs shared by TCGA and GEO through R packages ("limma" and "sva"). Next, we merged the shared DEIGs with survival data from TCGA and GEO into a single file and carried out univariate Cox proportional hazard regression analysis with R package ("survival"). It was considered that the survival risk gene of p < 0:001 was the prognosis associated immune genes (PIGs) of KIRC patients. Next, we used the CISTROME (http:// www.cistrome.org/) to find transcription factors (TFs) correlated to cancer occurrence and progression [23,24], retrieved the differentially expressed transcription factors (DETFs) from DEGs, and constructed the regulatory network of DETFs and PIGs by Pearson's correlation coefficient analysis. Correlation coefficient (|cor|) >0.4 and false discovery rate (FDR) <0.001 were the cut-off point of significant correlation. For coexpressed PIGs and DETFs, proteinprotein interaction (PPI) analysis was also performed using STRING (https://cn.string-db.org/).

Construction of Immune Gene-Related Prognostic Model
Based on TCGA and GEO Data. Using the LASSO Cox regression model based on R packages ("Matrix," "glmnet," and "survival") and taking PIGs as the research object, an optimal immune gene-related prognosis model (IGRPM) was constructed for KIRC. The formula for calculating the risk score is as follows: the risk score = Σ n k Expression ðkÞ × coef ðkÞ . In this formula, Expression (k) denotes the expression level of gene k in the TCGA and GEO samples, and coef (k) denotes the regression coefficient of characteristic gene k. All patients were separated into two groups based on IGRPM's median risk score: low-risk and high-risk; and R packages ("survival" and "survminer") were used to conduct the survival analysis. In order to reflect the predictive ability of risk characteristic based on IGRPM, we created a receiver operating characteristic curve (ROC) depending time through R packages ("survminer," "survival," "rms," and "timeROC") and calculated the area under the curve (AUC) of different survival times. Spearman correlation analysis of R software was used to explore the correlations of IGRPM-based risk characteristic with immune cells infiltration and TMB. Finally, we performed multivariate and univariate Cox regression analyses of IGRPM-based risk characteristic as well as other clinicopathological parameters to judge if IGRPM-based risk characteristic could be applied as an individual prognostic indicator for KIRC patients. Last but not the least, we constructed a nomogram through the risk characteristic based on IGRPM and the other clinicopathological features; what is more, the proposed nomogram was tested using ROC and calibration curve analyses to see if it was suitable for clinical application. PrimeScript RT reagent kit (RR037A, Takara) and Bio-Bad T100 Thermal Cycler according to the manufacturer's instructions were used to perform reverse transcription to synthesize cDNAs. Reverse transcription was performed with a final quality of 1000 ng of RNA. For quantification of gene expression levels, we performed RT-qPCR using FastStart Universal SYBR Green Master (04913914001, Roche) and Light Cycler96 Instrument (Roche). The comparative cycle threshold (Ct) method was used for relative quantification of each target [25][26][27]. △Ct = Ct ð gene of interestÞ − Ct ðGAPDHÞ. The GAPDH was used as an internal reference. The 2 -ΔΔCt method was applied to calculate the expression of target genes in tumor tissues compared with adjacent normal tissues. All groups set up three biological replicates, and experiments were repeated three times for each sample. Meanwhile, we explored the effect of all genes used to construct IGRPM on the overall survival (OS) of KIRC through the UALCAN database. The following primers were used for RT-qPCR.
GAPDH:   Journal of Immunology Research of TCGA were divided into two clusters: Immunity_L and Immunity_H (Figure 1(a)). Moreover, we further confirmed the cluster analysis of KIRC patients by using tSNE algorithm and obtained similar results (Figure 1(b)). In addition, we visually showed the scores of 29 immune cells or immunerelated functions as well as immune microenvironment scores in Immunity_H and Immunity_L through heat map and violin plot (Figures 1(c) and 1(d)). According to the two pictures, it can be concluded that the StromalScore, ImmuneScore, and ESTIMATEScore in Immunity_H were higher than those in Immunity_L, while the Immunity_H's TumorPurity was lower than the Immunity_L's TumorPurity (wilcox.test).      (Figure 2(a)). Next, we also explored the expressions of HLA-related genes between Immunity_H and Immunity_L. To our surprise, the expres-sions of all HLA-related genes in Immunity_H were significantly higher than that in Immunity_L (Figure 2(b)). All these results indicated that it was very reasonable and of great significance for us to classify KIRC patients into the Immunity_H and Immunity_L. Furthermore, we analyzed the DEGs between subtypes and obtained the DEIGs based on ImmPort database. Compared with Immunity_L, Immunity_H had more DEGs upregulation (Figure 2(c)). 1811 genes in the two subtypes were authenticated as DEGs, of which 402 were DEIGs (Figures 2(d)-2(f)).

The KEGG and GO Analyses of Two Subtypes.
We carried out KEGG and GO analyses according to DEGs in order Immunity_L Immunity_H  to further understand the functional differences between Immunity_H and Immunity_L (Immunity_H vs. Immu-nity_L). We found that the KEGG pathways enriched in Immunity_H were mainly antigen processing and presentation, cytokine-cytokine receptor interaction, T cell receptor signaling pathway, natural killer cell-mediated cytotoxicity, JAK-STAT signaling pathway, cell adhesion molecules cams, toll-like receptor signaling pathway, chemokine signaling pathway, B cell receptor signaling pathway, NOD-like receptor signaling pathway, complement and coagulation cascades, and leukocyte transendothelial migration; all these pathways are closely related to immunity (Figure 3(a)). However, there was no significant correlation of the enriched pathways with immunity in Immunity_L (Figure 3(a)). Biological process (BP), cellular component (CC), and molecular function (MF) of GO analysis also showed similar results as KEGG (Figures 3(b)-3(d)). These results demonstrated that the immunotyping of KIRC was accurate and that the immune activity was significantly improved in the Immunity_H.

Screening of PIGs and Construction of Regulatory
Networks. Next, combining the survival and gene expression data of TCGA and GEO, we identified 47 PIGs and calculated the hazard ratio (HR) by univariate Cox proportional hazard regression analysis (Figure 4(a)). There were 35 high-risk PIGs (HR >1) and 12 low-risk PIGs (HR <1). To better understand the multidimensional and complex regulatory network of PIGs involved in the onset and progression of KIRC, we first investigated the coexpression analysis of PIGs and DETFs. Surprisingly, all PIGs and DETFs that met our filter conditions showed a positive correlation (Supplementary table 2). The regulatory network between PIGs and DETFs was shown in Figure 4(b). In order to further confirm that there were significant correlations of DETFs with PIGs, we used STRING to analyze PPI and observed  OASL  TNFSF13B  OSM  NR3C2  THRB  IL17RD  SAA1  SAA2  MARCO  HNF4A  SLC11A1  LILRB3  IL6  KITLG  CD72  SEMA6D  IL2RA  CXCL5  PPARG  PLTP  TLR3  ISG20  RNASE2  IL1R2  TNFRSF18  LGR4  IL1RN  INHBE  CSF1  XCL1  ESRRG  CLDN4  CDH1  IGHG1  AQP9  PLA2G2A  SEMA7A  LBP  IL15RA  CRP  IKBKE  CXCL13  TNFSF14  PDGFRL  HCST  PAEP  TNFSF4   <0.  15 Journal of Immunology Research that there were complicated connections within them (Figure 4(c)).

Combining TCGA and GEO Samples to Construct an IGRPM of KIRC through Lasso Cox Regression Analysis.
To better utilize of the PIGs of KIRC, the dimension of PIGs was reduced by lasso Cox regression analysis. Combining with the KIRC cohorts of TCGA and GEO, we constructed an IGRPM containing 17 key genes to predict OS ( Figure 5(a)). The genes and coef used to determine the risk score for each subject were shown in Supplementary table 3. The corrected gene expression level, survival time, and survival status of TCGA and GEO samples were shown in Supplementary file1. We evaluated the risk score for KIRC patients of GSE29609 and TCGA and separated them into high-risk and low-risk groups based on the median risk score ( Figure 5(b)). It is not difficult to find that the high-risk group had considerably more deaths than the low-risk group ( Figure 5(c)). In addition, our further analysis showed that the risk score was negatively correlated with OS ( Figure 5(d)), suggesting that the OS of KIRC patients gradually decreased with the elevation of the risk score. We also drew the Kaplan-Meier curve of KIRC patients with R packages ("survival" and "survminer"). The results demonstrated that patients in the high-risk group had a significantly shorter survival time than those in the low-risk group (p < 0:001, Figure 5(e)). What is more, we conducted a time-dependent ROC analysis to better examine the prediction accuracy of IGRPM-based risk characteristic, and the predicted AUC at 1, 3, and 5 years was 0.751, 0.717, and 0.728, respectively ( Figure 5(f)). Additionally, the calibration diagram revealed the anticipated value of IGRPM-based risk characteristic was basically consistent with the real value ( Figure 5(g)).

The Correlations between IGRPM-Based Risk Characteristic and Immune Cells Infiltration and TMB.
Due to the absence of some data from GEO, we performed further analysis based on TCGA data. We analyzed the correlations between 17 genes used to construct IGRPM and 22 kinds of human immune cells by "spearman" method. Our results revealed that 17 key genes were closely associated with a variety of immune cells (Figure 6(a)), which once again demonstrated the importance of the immune microenvironment in the prognosis of KIRC patients. We further explored whether TMB differed in high-risk and low-risk groups. The results demonstrated that the TMB of highrisk group was substantially elevated compared with that of low-risk group (Figure 6(b)). Moreover, the survival time of patients with high TMB (H-TMB) was significantly shorter than that of patients with low TMB (L-TMB) (Figure 6(c)). Considering that both TMB and IGRPM had significant effects on the survival of patients with KIRC, we further investigated if the combination of IGRPM and TMB could better predict the survival outcome of KIRC patients. According to the risk score and TMB, we separated all samples from TCGA into 4 groups, including L-TMB + high-risk group, L-TMB + high-risk group, H-TMB + high-risk group, and H-TMB + low-risk group. The result revealed that the survival outcome of KIRC patients with L-TMB + low-risk group was the best in four groups (Figure 6(d)).   Journal of Immunology Research

Independent Prognostic Analysis and Establishment and
Verification of Nomogram. Since the risk score was substantially associated with the progression of KIRC, univariate and multivariate Cox regression analyses were performed to verify whether IGRPM-based risk characteristic could be used as an independent prognostic factor for KIRC patients. The results showed that the p value of the risk score was less than 0.001 in both univariate and multivariate regression analyses (Figures 7(a) and 7(b)), which confirmed that the risk score based on IGRPM was an independent variable affecting the survival of KIRC patients.
Combining with the patient's age, sex, TNM staging, grade, and risk score, we designed a nomogram to expand clinical applicability of the risk characteristic based on IGRPM (Figure 7(c)). By scoring the comprehensive characteristics of each patient, in a certain range, the higher the total scored, the worse the prognosis of the patient would be. In addition, the ROC proved that nomogram could accurately judge the survival status of KIRC patients (Figure 7(d)). Moreover, the calibration diagram showed a good similarity between the performance of the nomogram and the ideal model (Figure 7(e)).

Validation of the Expression of Genes Used to Construct
IGRPM in Early KIRC and Adjacent Normal Tissues. In order to further deepen the understanding of the genes used to construct IGRPM, we verified the expressions of some of these genes in early KIRC and adjacent normal tissues. The results of RT-qPCR indicated that the expressions of CLDN4, NR3C2, OASL, SAA1, and THRB in early KIRC were considerably suppressed than those in adjacent normal tissues (Figures 8(a)-8(e)); however, the expression of LGR4 in early KIRC was no statistical significance compared with that in adjacent normal tissues (Figure 8(f)). These results demonstrated that not all genes used to construct IGRPM were significantly differentially expressed in early KIRC compared with adjacent normal tissues. Last but not least, we explored the impact of all genes used to construct IGRPM on the prognosis of KIRC through the UALCAN database. To our surprising, the expressions of all these genes significantly affected the OS of KIRC (Figures 9(a)-9(q)), which showed that all these genes used to construct IGRPM had a significant impact on the prognosis of KIRC, although not all these genes showed differential expression in early KIRC compared with adjacent normal tissues. Therefore, it was more accurate to judge the prognosis of KIRC patients by integrating the expressions of these genes, which indicated that the IGRPM constructed by us was also more reliable and meaningful than the expression of a single or few genes to judge the prognosis of KIRC patients.

Discussion
As the prognosis of KIRC patients varies greatly, it is particularly important to establish a stable classifier for risk stratification of different KIRC patients, which can maximize the individualized treatment and timely follow-up of KIRC patients. Therefore, scientists have invested a lot of time and energy to explore the complex mechanism of KIRC, but the current understanding of KIRC, especially the    Journal of Immunology Research understanding of TIM and immune-related prognostic factors of KIRC, is far from satisfactory. In this study, the multifaceted mining of transcriptome data, tumor immunerelated characteristics, and TMB was designed to comprehensively understand the immune characteristics of KIRC and to construct a simple and practical method to help clinicians determine the prognosis of patients with KIRC, so as to develop personalized treatment options for KIRC patients.
First of all, we divided KIRC into two immune-related subtypes, Immunity_L and Immunity_H, and analyzed the

23
Journal of Immunology Research molecular characteristics of subtype specificity, such as immune-related gene expression, function, and pathway enrichment. Immunity_H had lower tumor purity and higher immune score than Immunity_L, indicating that immune cells exerted critical role in inhibiting tumor progression. These findings supported previous researches that the activation of immune-related responses can effectively control the progression of KIRC [28,29]. Additionally, Immunity_H had higher levels of HLA genes expression and immune cells infiltration than Immunity_L, indicating that Immunity_H had stronger immunogenicity. It has been reported that the expressions of HLA genes play an important role in antitumor activity and in the recognition and elimination of cancer cells [30,31]. Then we visualized the DEGs and DEIGs in Immunity_H and Immu-nity_L, and it could be seen that the expressions of immune-related genes were found to be significantly different in the two groups. GO and KEGG analyses showed that Immunity_H had more abundant immune-related signals than Immunity_L. These results demonstrated the accuracy and necessity of dividing KIRC into two subtypes associated with immunity. Then we combined the gene expression data and survival data of TCGA and GEO to identify the genes that may have a significant impact on the prognosis of KIRC patients. We screened out 47 genes that had a significant impact on the prognosis of KIRC patients. Most of these genes have been reported in previous studies whose results were similar to ours. For example, it has been reported that NR3C2 expression could be employed as an independent prognostic factor of KIRC [32], and overexpression of NR3C2 have potential to restrict the proliferation and invasion of KIRC [33]. For HNF4A (Hepatocyte Nuclear Factor 4, Alpha), a tumor suppressor, its expression could affect the prognosis of KIRC [34]. Overexpression of IKBKE (Inhibitor of Kappa Light Polypeptide Gene Enhancer in B-cells, Kinase Epsilon) is associated with an increased risk of death in KIRC [35]. These results fully showed that the PIGs screened out by us had great reliability and that our method of screening prognostic genes had certain superiority. Then coexpression analysis of KIRC-related transcription factors and PIGs was performed, and their interactions were again demonstrated by PPI. It has been reported that PPARG (peroxisome proliferation-activated receptor gamma) promotes apoptosis and inhibits the migration and proliferation of KIRC by inhibiting SIX2 (SIX homeobox 2) [36]. BATF (Basic Leucine Zipper ATF-Like Transcription Factor) affects the prognosis of many kinds of cancer by combining with HHLA2 (HERV-H LTR-associating 2) DNA, especially the prognosis of KIRC [37]. CIITA (Class II Major Histocompatibility Complex Transactivator) methylation is associated with immune heterogeneity in hematological tumors [38]. These transcription factors are associated with tumors through different mechanisms. All these results provided a basis for us to further understand the occurrence and development mechanism of KIRC. Based on the above information, we could have a more comprehensive understanding of the relationship between KIRC and immunity, and dividing KIRC patients into different subtypes through immune-related genes would lay a foundation for us to estimate the prognosis of KIRC patients from a new perspective.
Based on the expressions of 17 immune genes, we constructed an IGRPM to assess the prognosis of patients with   [39], NR3C2 [32,33], SAA1 [40], CXCL5 (C-X-C Motif Chemokine Ligand 5) [41], TLR3(toll-like receptor 3) [42], CLDN4 [43,44], THRB [45], and LGR4 [46]. However, the specific molecular mechanisms of these genes affecting the prognosis of tumor patients are not clear, so the molecular mechanisms of the effect of these genes on KIRC patients need to be further studied in the future. Nevertheless, we speculated that the 17 genes affecting the       For the sake of further confirming the scientificity of classifying KIRC patients into the high-and low-risk groups, we compared the TMB of the two groups and found that the TMB of the high-risk group was substantially elevated compared with that of the low-risk group, and the survival time of KIRC patients with L-TMB was significantly longer than that of KIRC patients with H-TMB. This was consistent with the conclusion that TMB could be used as a biomarker to judge the prognosis of tumor patients reported in other researches [47,48]. Similar results have been obtained from the research of lung squamous cell carcinoma [49]. Moreover, combining the risk score and TMB, we found that L-TMB + low-risk was most beneficial to the survival outcome of patients with KIRC, which provided a more accurate method for us to evaluate the prognosis of KIRC patients.
Finally, through univariate and multivariate Cox regression analyses, we proved that the risk characteristic based on IGRPM could be employed as an independent prognostic factor for patients with KIRC. Meanwhile, in order to make it convenient for clinicians to judge the prognosis of KIRC patients simply and effectively, we constructed a nomogram combining the clinicopathological information and the risk score of patients and verified its reliability by ROC and calibration diagram. Last but not the least, we verified the expressions of some of these genes used to construct IGRPM in early KIRC and adjacent normal tissue by RT-qPCR. The RT-qPCR results demonstrated that not all these genes used to construct IGRPM were significantly differentially expressed in early KIRC compared with adjacent normal tissue, but all these genes had significant influence on the prognosis of KIRC. Because of this, it further indicated that it was more accurate and effective to judge the prognosis of patients with KIRC by risk stratification based on the expression of these genes than by the expression of a single or few genes.
Although our study could make an important contribution to differentiating the immune status and predicting the prognosis of patients with KIRC, admittedly, there are some limitations to our study. On the one hand, the immune-related features of KIRC are constructed through public datasets, but these features have not yet been verified by clinical cohorts. On the other hand, due to the different information contained in TCGA and GEO databases, some data can only be obtained from TCGA, which makes some of the data less representative. Nevertheless, our results still have important guiding significance for survival prediction and immune assessment of patients with KIRC and may promote the development of immunotherapy for KIRC.

Conclusion
In a word, we elucidated the information of the interaction between tumor and immunity in KIRC from different aspects and used a variety of algorithms to establish IGRPM-based risk characteristic to evaluate the prognosis of KIRC patients. In addition, we constructed a nomogram which was

34
Journal of Immunology Research convenient for clinicians to combine the clinicopathological information of patients to judge the prognosis of patients with KIRC. These comprehensive immune assessments and survival predictions, integrating multiple aspects of data and clinical information, can provide additional value to the current Tumor Node Metastasis staging system for risk stratification of KIRC and may facilitate the development of KIRC immunotherapy.

Data Availability
The data used to support the findings of this study are available in TCGA (https://gdc.cancer.gov/) and GEO databases.
The raw data of RT-qPCR will be made available by the authors, without undue reservation, to any qualified researcher.

Ethical Approval
This research was approved by the Ethics Committee of the Second Affiliated Hospital, Hengyang Medical School, University of South China. The patients provided their written informed consent to participate in this study.

Conflicts of Interest
The authors declare that they have no conflicts of interest.