Identification of an Independent Immune-genes Prognostic Index for Renal Cell Carcinoma

doi:10.21203/rs.3.rs-111386/v1

Download PDF

Research article

Identification of an Independent Immune-genes Prognostic Index for Renal Cell Carcinoma

https://doi.org/10.21203/rs.3.rs-111386/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Objective: Increasing evidence has indicated an association between immune micro-environment in clear cell renal cell carcinoma (ccRCC) and clinical outcomes. The aim of this research is to comprehensively investigate the effect of tumor immune genes on the prognosis of ccRCC patients.

Methods: 2498 immune genes were downloaded from ImmPort database. Additionally, we identified and downloaded the transcriptome data of patients with ccRCC from the TCGA database through the R package, as well as relevant clinical information. We apply certain survival R package to analyse the survival of hub-genes before analyzing the effect of immune genes on the prognosis of clear cell renal cell carcinoma (ccRCC) utilizing Cox regression analysis. Based on the statistical correlation between hub immune gene and survival ,immune risk score model was set up.We finally constructed a nomogram to predict the survival rate of ccRCC overally. In addition, whether the immune gene risk score model is an independent prognostic factor for ccRCC is comprehensively considered applying multivariate cox regression analysis. It is worth noting that throughout the data analysis, P< 0.05 was recognized to be of significance statistically.

Results: The results of the difference analysis showed that 556 immune genes exhibited differential expression between normal and ccRCC tissues (p<0. 05). Univariate cox regression analysis revealed 43 immune genes statistically correlated with ccRCC related survival risk (P<0.05). In addition, a 18-genes based immune genes risk scoring model was constructed through lasso COX regression analysis. KM curve indicated that patients in high-risk were associated with poor outcomes (p<0.001). ROC curve indicated that the immune risk score model was reliable in predicting survival risk (5-year OS, AUC=0.802). Our model showed satisfying AUC and survival correlation in the validation dataset ( 5-year OS AUC=0.705, P<0.001). Furthermore, multivariate cox regression analysis confirmed that the immune risk score model was an independent factor for predicting the prognosis of ccRCC. A nomogram was established to comprehensively predict the survival of ccRCC patients with the results of multivariate cox regression analysis. Finally, we found that 15 immune genes and risk scores were significantly associated with clinical factors and prognosis, and were involved in multiple oncogenic pathways.

Conclusion: Collectively, tumor immune genes played an essential role in the prognosis of ccRCC. Furthermore, immune risk score was an independent predictive factor of ccRCC, indicating a poor survival.

Cancer Biology

Oncology

clear cell renal cell carcinoma

immune genes

prognosis

immune risk score model

nomogram

Kidney cancer was one of the most prevalent malignant tumors in men and women all over the world. Its incidence has been increasing in the past decade, comprising up to 2% − 3% of all newly diagnosed tumor cases[1]. Histologically, clear cell renal cell carcinoma (ccRCC) is the prominent sub-types of kidney cancer, accounting for approximately 75% of kidney cancer cases[2]. Given that great development has undergone in screening, diagnosis, surgery and various tumor drug treatment[3–5], however, the clinical outcomes of ccRCC was still insufficient[2, 6]. Therefore, it was crucial to identify several prognostic factors and targets for improving the treatment and prognosis of ccRCC patients.

KIRC, an important cause of cancer-related death, is an immunogenic tumor, the growth of which is related to impaired anti-tumor immunity[7]. Despite the expression of activation markers by tumor infiltrating lymphocytes, which is puzzling through clinical observation, tumor immunity has been demonstrated to play an important role in controlling the disease[8]. For example, dysfunction of innate and adaptive immune cells has been shown to lead to failure of the immune system to control RCC[8]. Immunotherapy for ccRCC has been the subject of research for decades. Through in-depth studies of different immune cells, we have found that individual cell types do not act in isolation, but rather in complex networks of cell-cell interactions[9]. Thus improving our understanding of this immune network biology may help us to effectively utilize a variety of effector cells, enabling us to better develop new therapeutic strategies to successfully combat RCC. Considering that these interactions play a key role in the effective activation and function of effector cells, which is a prerequisite for the successful elimination of tumors[9]. Notably, they could be effectively targeted by drugs and affect the clinical outcomes of patients. In addition, immune genes have not been extensively studied in the context of ccRCC biology.

The purpose of this study is to describe the infiltration expression and lineage of immune genes in ccRCC and investigate the effects of immune genes on the prognosis of patients with ccRCC. Furthermore, an immune genes risk score model and a nomogram model is constructed to predict the survival of ccRCC.

Data acquisition

First of all, 2498 immune genes were downloaded from ImmPort database. Additionally, we identified and downloaded the transcriptome data of patients with ccRCC from the TCGA database through the R package, including 72 cases of paracancerous normal tissue and 507 cases of tumor tissue. Further, relevant clinical information of 507 ccRCC patients were obtained such as age, gender, stage, tumor&Lymph node&metastasis stage, survival status and survival duration (Table 1). Finally, "Limma" package in R software was utilized to correct the transcriptome data we have downloaded.

Table 1

Clinical characteristics of included patients in the study
Variables	Total (n = 507)	Training cohort (n = 252)	Validation cohort (n = 255)
Age (year)
< 40	17	9	8
40–59	225	110	115
60–79	253	130	123
80+	22	13	9
Gender
FEMALE	179	89	90
MALE	338	173	165
Grade
G1	13	6	7
G2	223	125	98
G3	203	91	112
G4	73	37	36
GX	5	3	2
Stage
I	257	136	121
II	55	25	30
III	123	63	60
IV	82	38	44
T stage
T1	263	137	126
T2	67	31	36
T3	176	87	89
T4	11	7	4
N stage
N0	236	124	112
N1	15	6	9
NX	266	132	134
M stage
M0	414	217	197
M1	77	35	42
MX	26	10	16

Gene function enrichment analysis

In order to explore the major biological process of selected hub-genes, methods were utilized to conduct the gene functional enrichment analyses including Kyoto Encyclopedia of Genes and Genomes (KEGG) and gene ontology (GO). We utilized Database for Annotation, Visualization, and Integrated Discovery (DAVID, https:// david.ncifcrf.gov/) to identify enriched KEGG and GO themes.

Survival analysis and expression comparison of hub-genes

Clinical information for TCGA-KIRC including survival time, survival state, and TNM staging were also downloaded from the TCGA database (samples with missing information were excluded). Survival R package was applied in survival analyses for hub-genes. For the overall survival rates, the logrank test was used to detect significant differences. The results were visualized using Kaplan-Meier survival curves, and P-value < 0.05 was considered as statistically significant.

Gene set enrichment analysis

Gene enrichment analysis (GSEA) (version 3.0, the broad institute of MIT and Harvard,http://software.broadinstitute.org/gsea/downloads.jsp) was conducted between ccRCC and paracancerous normal tissues to study the biological characteristics of renal carcinoma. In detail, the “collapse data set to gene symbols” was set to false, the number of marks was set to 1000, the “permutation type” was set to phenotype, the “enrichment statistic” was set to weighted, and the Signal2Noise metric was used for ranking genes. High expression group was used as experimental group and low expression group was used as reference group. “c2.cp.kegg.v7.0.symbols.gmt” gene sets database was used for enrichment analysis. Gene set size > 500 and < 15, FDR < 0.25, and nominal P-value < 0.05 were regarded as the cut-off criteria.

Statistical analysis

All analyses were performed using R 3.6.1. All statistical tests were two-sided, and P value < 0.05 was considered statistically significant. Continuous variables that conformed to the normal distribution were compared with the use of independent t test for comparison between groups, while continuous variables with skewed distribution were compared with the Mann-Whitney U test. The correlation matrix was constructed by R software based on Pearson Correlation Coefficient. The relationship between immune cell infiltration and overall survival was analyzed through the Kaplan-Meier curve which was evaluated by log-rank test. Time-dependent ROC curves were used to analyze the sensitivity and specificity of the recurrence prediction model. The univariate regression model was used to analyze the effects of individual variables on survival. The least absolute shrinkage and selection operator (LASSO) cox regression model was used to confirm independent impact factors associated with survival. The nomogram was constructed with the regression coefficients based on the cox analysis.

Differential expression screening of ccRCC

2498 immune genes were downloaded from ImmPort database. The transcriptome data of 507 cancer cases and 72 adjacent normal tissues cases was obtained from TCGA database for differential expression analysis. A total of 556 immune genes were identified as differentially expressed immune genes (DEIGs) between ccRCC and normal tissues, including 402 up-regulated and 154 down regulated (p < 0.05, Fig. 1A). The heatmap of the top 10 up-regulated and top 10 down regulated DEIGs was shown in Fig. 1B.

Functional annotation of these 556 DEIGs

In order to fully understand the biological attributes of these 556 DEIGs, we conducted Kyoto Encyclopedia of Genes and Genomes (KEGG) and gene ontology (GO) analysis. Based on the results of DAVID, the top three enriched GO terms for up-regulated genes and down-regulated genes were: cAMP − mediated signaling, humoral immune response and negative regulation of ERBB signaling pathway; regulation of lymphocyte activation, humoral immune response and regulation of leukocyte mediated immunity, respectively (Fig. 2A). The top biological pathway enriched for up-regulated genes and down-regulated genes were: JAK − STAT signaling pathway, PI3K − Akt signaling pathway and Rap1 signaling pathway; Cytokine − cytokine receptor interaction, Th1 and Th2 cell differentiation JAK − STAT signaling pathway, respectively (Fig. 2B).

Establishment of immune prognosis model

We constructed a PPI network based on differentially expressed genes and selected a total of 496 genes with more than 50 adjacent nodes (Figure S1). For the purpose of revealing the relationship between these 496 DEIGs and overall survival, 43 prognostic DEIGs were identified by utilizing univariate Cox regression analyses (Fig. 3A). TCGA ccRCC data were randomly divided into two sets (training set : validation set, 1 : 1). Then, lasso regression analysis was applied to increase the robustness and select the independent indicators for the overall survival based on training set and finally we got 18 DEIGs for the construction of prognostic index (Fig. 3B, Fig. 3C, Table 2). After the construction of prognostic index, patients were separated into high risk and low risk (Fig. 3D, Fig. 3F). Heap map was utilized to visualize the difference of gene expression profile in low- and high- risk patients in ccRCC training set (Fig. 2E). The results from K-M analysis indicated that high risk patients had lower overall survival than low risk patients in both training group and validation group (P < 0.001) (Fig. 4A, Fig. 4B, Fig. 4C). The ROC curve revealed that the risk model had a good sensitivity and specificity in predicting survival risk (AUC = 0.802, AUC = 0.705 for 5 years overall survival in training and validation group, respectively) (Fig. 4D, Fig. 4E, Fig. 4F). To explore whether the constructed immune risk scoring model was independent form age, gender, stage, and other clinical pathological parameters, we performed an univariate and multivariate cox regression analysis for age, gender, stage, grade, TNM and risk score. In univariate Cox model, age, pathological grade, pathological stage, pathological T, M stage and high risk score were associated with poor survival (Fig. 5A). In multivariate Cox model, only age, pathological stage and risk score worked as independent predicted factors (Fig. 5B). To better predict the prognosis of ccRCC patients at three and five years post-surgery, we constructed a new nomogram from the variables associated with OS (age, gender, histological grade, pathological stage TNM stage and risk score) (Fig. 5C, 5D, 5E).

Table 2

Multivariate cox regression analysis to establish immune genes risk score model
Gene	Coef
ICAM1	0.0058648074324021
IFNG	0.0336285898183753
CXCL5	0.00449394306613633
XCL1	0.178774321649646
TGFB1	0.0114011479579204
PDGFRA	0.0346969626392538
GNAI1	0.00617576278056856
TNFSF11	0.321235808441014
HMOX1	-0.00150622701372233
CCL22	-0.456892844189334
IL4	3.92829889144972
CRP	0.00127583274886622
EDN1	-0.00254186826990028
AVP	1.2503911973052
CSF2	0.873478979686199
GAL	0.0818075259757682
GNRH1	0.1244306056583
PPY	0.275246400052544

Clinical and prognostic correlation of 18 model genes and immune genes risk score

We further investigated the proportion of each model genes in different pathological stages. We demonstrated that EDN1, GNAL1 and ICAM1 were most significantly associated with development of ccRCC (Fig. 6). Meanwhile, we found that the expression of IFNG and XCL were associated with the infiltration of T cell CD4+, T cell CD8 + and Myeloid dendritic cell (Figure S2). Regard to the immune genes risk score, a strong correlation with grade, pathological stage and clinical TNM stage was identified (Fig. 7).

Gene set enrichment analysis of risk scores

To explore the biological relevance of risk scores involved in progression of ccRCC, we performed a gene set enrichment analysis of risk scores based on the TCGA breast cancer cohort. GSEA analysis indicated high risk scores was associated with IL6 JAK STAT3 SIGNALING, EPITHELIAL MESENCHYMAL TRANSITION and WNT BETA CATENIN SIGNALING pathway (Fig. 8).

ccRCC was a heterogeneous disease with different ethnic characteristics, which originated from renal epithelial cells[10]. It was estimated that ccRCC accounted for the majority of RCC related deaths [11]. Although radical nephrectomy has been proved to be an effective treatment for local renal cancer, many patients may experience development and metastasis after surgical resection. Given that targeted treatment for advanced and metastatic ccRCC has been fully developed, the response to treatment was diverse [12]. As was universally known, the determination of molecular mechanism and applicable prognostic factors may be the crucial lynchpin for the treatment of ccRCC [13]. Studies had shown that cancer prognosis was closely related to TME, especially cancer immune micro-environment [14, 15]. It is evident that different types of cancers had diverse immune genes sub-population. Therefore, it was crucial to investigate the immune genes subsets for the evaluation of risk and ccRCC prognosis.

In our study, we conducted a comprehensive and detailed assessment of immune genes in ccRCC, based on the data from a large set of samples. All gene expression data and patients clinical characteristics information were downloaded from TCGA dataset. We analyzed the 2498 immune genes between ccRCC and normal tissues from IMMPORT database, eventually, we verified 556 DEIGs. Moreover, we identified and constructed a 18 hub genes risk score model for breast cancer via univariate and lasso cox regression analysis, including TSLP, IL17B, NR3C2, RAC2, SERPINA3, HSPA2, CD79A, UNC93B1, NFKBIE, SDC1, IFNG, IRF7, GALP, TNFRSF18 and ULBP1. Furthermore, to investigate the prognostic value of the model, we performed the ROC curve and investigate the association between the model and clinical features. As expected, the high-risk group was correlated with worse overall survival and was inclined to have advanced stage and higher histological grade which might manifest poor outcome.

Previous studies have constructed ccRCC immune gene model by screening immune related lncRNAs. Moreover, a new prognostic gene marker based on immune lncrna in KIRC patients was found [16]. Zhao et al. integrated multiple levels of data to construct immune, inflammatory or KIRC oriented neighbor networks (IIKDN networks) and KIRC related gene directed networks (KIRCD networks). Their analysis showed that immune and inflammatory related genes have special topological characteristics and related KIRC expression patterns in the network. Further, they identified five core clusters to construct specific prognostic biomarkers for KIRC[17]. Another study evaluated the prognostic value of individual gene expression using TCGA data and ccRCC patient data[18]. A predictive nomogram was generated and independent prognostic factors were identified to assess overall survival (OS) and progression-free survival (PFS) of patients with ccRCC at 1, 5, and 8 years and to examine the functional involvement of individual genes in RCC in vitro and in vivo models[18]. The special feature of our study was that we used immune related genes to establish a prognosis model of renal cell carcinoma. Furthermore, we analyzed the expression profiles of these genes in the grading and staging of RCC. Finally, we analyzed the clinical and biological correlation of immune scores.

Moreover, among all of the genes metioned above, TSLP, IL17B, SERPINA3, HSPA2, UNC93B1, NFKBIE, IFNG, GALP have rarely been studied. But ULBP1 has been found to be a potential prognostic factor in patients with non-metastatic ccRCC after nephrectomy and was involved in RCC tumorigenesis[18, 19]. Rac family small GTPase 2 (RAC2) has also been studied in terms of RAC2 expression levels in human ccRCC tissues and cell lines, and it may serve as a promising prognostic and diagnostic biomarker for ccRCC[20]. CD79A [21], SDC1 [22] and NR3C2[23] were also proved to be a candidate gene associated with the prognosis and microenvironment of ccRCC in several researches.

However, there were some limitations in our research. Firstly, the sample size in our study was small and a larger cohort and more abundant sequencing results were needed. Secondly, we only focused on the gene expression, but ignored other events such as the gene methylation, and copy number amplification, which were also important in tumor progression.

In summary, our study sheds light on the utility of immune genes in the prognosis of ccRCC. The constructed immune genes risk scoring model is reliable in predicting the prognosis of ccRCC, and this risk scoring model is an independent influencing factor for the prognosis of ccRCC. With the rapid development of high-throughput technology, we have confidence to believe that our immune risk scoring model have great potential in clinical practice. And it may have critical value for exploring new anti-cancer immunodiagnosis and treatment strategies.

Authors contributions

Chao Qin and Ninghong Song designed this work. Guangyao Li and Xiyi Wei wrote the manuscript. Shifeng Su and Shangqian Wang performed the bioinformatics analysis. Shangqian Wang, Wei Wang, Yichun Wang, Xianghu Meng, Jiadong Xia performed the data review. All authors have read and approved the manuscript.

Acknowledgments

None.

Funding

This work was supported by the Jiangsu Province “Six Talent Peaks Project” (WSN-011), by the National Natural Science Foundation of China (grant number 81672531，81972386).

Ethical Approval and Consent to participate

Not applicable.

Competing interests

The authors have no conflicts of interest to declare.

Ljungberg B, Bensalah K, Canfield S, Dabestani S, Hofmann F, Hora M, Kuczyk MA, Lam T, Marconi L, Merseburger AS et al: EAU guidelines on renal cell carcinoma: 2014 update. Eur Urol 2015, 67(5):913-924.
Barata PC, Rini BI: Treatment of renal cell carcinoma: Current status and future directions. CA Cancer J Clin 2017, 67(6):507-524.
Dagher J, Kammerer-Jacquet SF, Dugay F, Beaumont M, Lespagnol A, Cornevin L, Verhoest G, Bensalah K, Rioux-Leclercq N, Belaud-Rotureau MA: Clear cell renal cell carcinoma: a comparative study of histological and chromosomal characteristics between primary tumors and their corresponding metastases. Virchows Arch 2017, 471(1):107-115.
Fernandez-Pello S, Hofmann F, Tahbaz R, Marconi L, Lam TB, Albiges L, Bensalah K, Canfield SE, Dabestani S, Giles RH et al: A Systematic Review and Meta-analysis Comparing the Effectiveness and Adverse Effects of Different Systemic Treatments for Non-clear Cell Renal Cell Carcinoma. Eur Urol 2017, 71(3):426-436.
Lin YW, Lee LM, Lee WJ, Chu CY, Tan P, Yang YC, Chen WY, Yang SF, Hsiao M, Chien MH: Melatonin inhibits MMP-9 transactivation and renal cell carcinoma metastasis by suppressing Akt-MAPKs pathway and NF-kappaB DNA-binding activity. J Pineal Res 2016, 60(3):277-290.
Garcia JA, Rini BI: Recent progress in the management of advanced renal cell carcinoma. CA Cancer J Clin 2007, 57(2):112-125.
Aparicio LMA, Fernandez IP, Cassinello J: Tyrosine kinase inhibitors reprogramming immunity in renal cell carcinoma: rethinking cancer immunotherapy. Clin Transl Oncol 2017, 19(10):1175-1182.
Yang J, Li H, Hu S, Zhou Y: ACE2 correlated with immune infiltration serves as a prognostic biomarker in endometrial carcinoma and renal papillary cell carcinoma: implication for COVID-19. Aging (Albany NY) 2020, 12(8):6518-6535.
Murphy KA, James BR, Guan Y, Torry DS, Wilber A, Griffith TS: Exploiting natural anti-tumor immunity for metastatic renal cell carcinoma. Hum Vaccin Immunother 2015, 11(7):1612-1620.
Wang XM, Lu Y, Song YM, Dong J, Li RY, Wang GL, Wang X, Zhang SD, Dong ZH, Lu M et al: Integrative genomic study of Chinese clear cell renal cell carcinoma reveals features associated with thrombus. Nat Commun 2020, 11(1):739.
Scelo G, Larose TL: Epidemiology and Risk Factors for Kidney Cancer. J Clin Oncol 2018:JCO2018791905.
Li JK, Chen C, Liu JY, Shi JZ, Liu SP, Liu B, Wu DS, Fang ZY, Bao Y, Jiang MM et al: Long noncoding RNA MRCCAT1 promotes metastasis of clear cell renal cell carcinoma via inhibiting NPR3 and activating p38-MAPK signaling. Mol Cancer 2017, 16(1):111.
Wang C, Wang Y, Hong T, Cheng B, Gan S, Chen L, Zhang J, Zuo L, Li J, Cui X: Blocking the autocrine regulatory loop of Gankyrin/STAT3/CCL24/CCR3 impairs the progression and pazopanib resistance of clear cell renal cell carcinoma. Cell Death Dis 2020, 11(2):117.
Bremnes RM, Al-Shibli K, Donnem T, Sirera R, Al-Saad S, Andersen S, Stenvold H, Camps C, Busund LT: The role of tumor-infiltrating immune cells and chronic inflammation at the tumor site on cancer development, progression, and prognosis: emphasis on non-small cell lung cancer. J Thorac Oncol 2011, 6(4):824-833.
Li T, Fan J, Wang B, Traugh N, Chen Q, Liu JS, Li B, Liu XS: TIMER: A Web Server for Comprehensive Analysis of Tumor-Infiltrating Immune Cells. Cancer Res 2017, 77(21):e108-e110.
Khadirnaikar S, Kumar P, Pandi SN, Malik R, Dhanasekaran SM, Shukla SK: Immune associated LncRNAs identify novel prognostic subtypes of renal clear cell carcinoma. Mol Carcinog 2019, 58(4):544-553.
Zhao E, Li L, Zhang W, Wang W, Chan Y, You B, Li X: Comprehensive characterization of immune- and inflammation-associated biomarkers based on multi-omics integration in kidney renal clear cell carcinoma. J Transl Med 2019, 17(1):177.
Cancer Genome Atlas Research N: Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature 2013, 499(7456):43-49.
Varela I, Tarpey P, Raine K, Huang D, Ong CK, Stephens P, Davies H, Jones D, Lin ML, Teague J et al: Exome sequencing identifies frequent mutation of the SWI/SNF complex gene PBRM1 in renal carcinoma. Nature 2011, 469(7331):539-542.
Liu Y, Cheng G, Song Z, Xu T, Ruan H, Cao Q, Wang K, Bao L, Liu J, Zhou L et al: RAC2 acts as a prognostic biomarker and promotes the progression of clear cell renal cell carcinoma. Int J Oncol 2019, 55(3):645-656.
Wan B, Liu B, Huang Y, Lv C: Identification of genes of prognostic value in the ccRCC microenvironment from TCGA database. Mol Genet Genomic Med 2020, 8(4):e1159.
Ozcan A, Celik E, Karslioglu Y, Basal S: CD138 expression in renal tumors and its usage in the differential diagnosis. Turk Patoloji Derg 2011, 27(2):110-115.
Zhao Z, Zhang M, Duan X, Deng T, Qiu H, Zeng G: Low NR3C2 levels correlate with aggressive features and poor prognosis in non-distant metastatic clear-cell renal cell carcinoma. J Cell Physiol 2018, 233(10):6825-6838.

FigureS1.pdf
Construction of PPI network with differentially expressed immune genes.
FigureS2.pdf
Correlation between 18 model immune genes and immune cell infiltration.

Download PDF

Editor assigned by journal
02 Dec, 2020
Submission checks completed at journal
18 Nov, 2020
Editor invited by journal
15 Nov, 2020
First submitted to journal
21 Oct, 2020

You are reading this latest preprint version

Identification of an Independent Immune-genes Prognostic Index for Renal Cell Carcinoma

Status:

Version 1

Abstract

Figures

Introduction

Materials And Methods

Data acquisition

Gene function enrichment analysis

Survival analysis and expression comparison of hub-genes

Gene set enrichment analysis

Statistical analysis

Results

Differential expression screening of ccRCC

Functional annotation of these 556 DEIGs

Establishment of immune prognosis model

Clinical and prognostic correlation of 18 model genes and immune genes risk score

Gene set enrichment analysis of risk scores

Discussion

Conclusion

Declarations

References

Supplementary Files

Status:

Version 1