The immune infiltration in clear cell Renal Cell Carcinoma and their clinical implications: A study based on TCGA and GEO databases

The tumor immune microenvironment in clear cell Renal Cell Carcinoma (ccRCC) still remains poorly understood. Previous methods to study the tumor immune microenvironment have a limitation when accounting for the functionally distinct cell types. In this study, we investigated the differently infiltrated immune cells and their clinical significance in ccRCC for the purpose of shedding some important light on the complex immune microenvironment in ccRCC. The devolution algorithm (CIBERSORT) was applied to infer the proportion of 22 immune infiltrating cells based on gene expression profiles of ccRCC bulk tissue, which were downloaded from TCGA and GEO databases. As a result, we observed considerable differences in immune cells percentage between ccRCC tumor tissue and paired normal tissue; meanwhile, we uncovered their internal correlations and associations with Fuhrman grade. Moreover, dendritic cells resting, dendritic cells activated, mast cells resting, mast cells activated and eosinophils were associated with favorable prognosis, whereas B cells memory, T cells follicular helper and T cells regulatory (Tregs) were correlated with poorer outcome.


Introduction
Clear cell Renal Cell Carcinoma (ccRCC) is the most common histology identified in renal carcinoma. Its manifestations are different both biologically and in clinic [1,2]. The cancer genome changes based on large scale sequencing researches contributed greatly to our understanding of underlying molecular mechanism of ccRCC [3,4]. Notably, a solid tumor is not only composed of cancer cells, but also by non-cancer cells, which can profoundly influence tumor progression in an elaborate and dynamic manner [5]. Among these non-cancer cells, the tumor infiltrating immune cells (TIICs) exert a central role in pro-and anti-tumorigenic processes; moreover, they have been found closely correlated with the clinical outcome and response to immunotherapy [6]. Compared to other carcinomas, the roles of immune cells recruited to microenvironment in ccRCC have not yet been elucidated and still remain an enigma, especially regarding the contradictory correlation of high CD8+ T cell infiltration with poor prognosis [7,8].
Previous traditional research techniques to evaluate TIICs include immunohistochemistry and flow cytometry, both of which inevitably are limited to a narrow viewpoint when analyzing the composition of immune cells comprehensively. Additionally, with the demanding sample processing, flow cytometry may result in cytolysis of certain cell Ivyspring International Publisher types [9,10]. However, the immune response in a tumor involves in plenty of specialized cell types [11]. To better understand the diversity and nature of infiltrating immune cells in ccRCC, it is a pre-requisite to enumerate the number of immune cells in an aggregative manner. CIBERSORT, a versatile gene expression-based devolution algorithm, can quantify cell fractions from gene expression profiles of bulk tissue [12]. Therefore, the different types of infiltrating immune cells can be quantified simultaneously, allowing this method to obviate the concern of various surface markers and possible cellular dissociation. Due to its superiority, we used CIBERSORT in this study to enumerate 22 distinct functional immune cell types in ccRCC to define the landscape of ccRCC tumor tissue and paired normal tissue; more importantly, we investigated its relationship with other immune cells, survival and pathological grade. We hoped this study will provide some important information regarding the complex immune microenvironment in ccRCC and help to reveal new therapeutic targets.

Data collection
This study made use of data from public datasets. Gene expression profiles and corresponding clinical information from primary ccRCC tumors, uploaded up to the 31st December 2018, were downloaded from Gene Expression Omnibus (GEO) and The Cancer Genome Atlas (TCGA) [13,14]. Duplicates and datasets with small sample sizes (N < 50) in GEO database were excluded. For TCGA datasets, preprocessing and aggregation of raw data were achieved by means of a robust multi-array average algorithm. Furthermore, voom (variance modelling at the observational level) was used to transform RNA sequencing data to values that are more similar to those from microarrays [15]. Furthermore, gene probe names must be transformed into gene names based on platform annotation flies for GEO datasets. Subsequently, we organized each sample and corresponding clinical data for further analysis; moreover, we manually identified and picked out tumor tissue and paired normal tissue to screen differentially infiltrated immune cells to investigate whether there is a difference in the infiltrated immune cells between different tissues. The overall study design and the different samples that were included at every stage of the analysis are illustrated as a flowchart in Figure 1.

Enumeration of tumor infiltrating immune cells
CIBERSORT is a robust analytic tool that uses gene expression signatures consisting of 547 genes. It characterizes each immune cell subtype and accurately quantifies distinct immune cell compositions using a deconvolution algorithm. The derived P-value reflects the statistical significance of deconvolution results and can help filter out samples with less significant accuracy.
Before running CIBERSORT, the original gene expression data downloaded from TCGA and GEO must be normalized as Binbin Chen et al. described previously [16]. Then, the data was uploaded to the CIBERSORT web portal (http://cibersort.stanford. edu) with a number of permutations being set to 100. Relative proportions of 22 infiltrating immune cells together with CIBERSORT metrics of CIBERSORT P-value, Pearson correlation coefficient and root mean squared error (RMSE) were evaluated for each sample simultaneously.

Statistical analysis
Only samples with a CIBERSORT P-value < 0.05 were regarded as statistically significant and included in further analysis. Correlations between different immune cell subtypes were established using the Pearson correlation coefficient. Associations between categorical and continuous variables were tested using the Kruskal-Wallis or Wilcoxon test. Survival analysis of specific immune cell subtype was conducted according to the median of the proportion of immune cell. Log-rank Mantel-Cox regression was applied to compare the survival curves between groups of patients using the Graphpad Prism 7.0 software. Additionally, multivariable analysis was adjusted according to age, gender, histological grade, T stage, lymph node metastasis, distant metastasis, and TNM stage using SPSS 24.0.
All analyses were conducted by R version 3.5.2 and all statistical tests performed were two-sided. A P-value < 0.05 was considered as statistically significant.

The landscape of immune infiltration in ccRCC
We first revealed the landscape of 22 immune cell subpopulations infiltration in ccRCC, and subsequently we investigated the difference between tumor tissue and paired normal tissue using the CIBERSORT algorithm. Detailed results are presented in Table 1. The fraction of immune cells varied distinctly between groups (Figure 2A, 2B). Compared with paired normal tumor tissue, ccRCC tissue contained a greater number of T cells CD8+, T cells follicular helper, T cells regulator (Tregs), Macrophages M0, Macrophages M1 and neutrophils. However, B cells naive, T cells CD4 naive, T cells CD4 memory resting, monocytes, dendritic cells resting and mast cells resting fractions were relative lower ( Figure 3A). The proportions of 22 TIICs were weakly-to-strongly correlated in tumor. T cells CD8+ and T cells follicular helper showed the strongest positive correlation (Pearson correlation = 0.54), while T cells CD8+ and T cells CD4+ memory resting showed the strongest negative correlation (Pearson correlation = 0.73); moreover, T cells CD8+ also indicated moderate negative correlation with Macrophages M2 (Pearson correlation = 0.56) ( Figure  3B). Altogether, these results revealed that the immune response of ccRCC acted as an intricate network and proceeded in a tightly regulated way.

Identification of clinical implications of TIICs subsets
Owning to the missing survival data in included GEO datasets, we investigated whether there was a statistical relationship between specific TIICs and ccRCC overall survival obtained from TCGA by univariate Cox regression through Graphpad Prism7.0. After a restriction of CIBERSORT filter to P < 0.05, there were 418 patients with available data on overall survival (142 events). The detailed 95% confidence intervals, unadjusted HRs and P-value for the median fractions of TIICs subtypes are presented in Table 2 Figure 4. Among those TIICs subtypes associated with overall survival, we evaluated their possibility of becoming independent prognostic factor using multivariable analysis adjusted for known prognostic factors. However, none represented an independent prognostic factor besides T stage (T1-T2/T3-T4), lymph node metastasis (N0/N1-N2), distant metastasis (M0/M1) and TNM stage.
Moreover, we first revealed the association between different immune cell subsets and ccRCC pathological grade by combing the clinical characteristics of TCGA and GEO databases.

Discussion
The tumor microenvironment, which comprises malignant tumor cells, various infiltrating immune cells, fibroblasts and numerous cytokines and chemokines, is well recognized as a complex biological process like an intricate and dynamic ecosystem. Among this ecosystem, immune response plays an important role in tumor growth, invasion and metastasis; therefore, it is treated as another therapeutic target beyond chemotherapy and radiation [11]. However, despite the astounding clinical successes in multiple tumor types, many more patients experienced minimal or no clinical response to the same immunotherapeutic intervention [17]. Until recently, the roles of TIICs have not been fully understood. Hence it is of vital importance to figure out the diversity and complexity of the tumor immune context to predict and guide immunotherapy, as well as to reveal novel biomarkers and targets for therapeutic modulation. With the advances in computational methods, a deconvolution algorithm called CIBERSORT was developed, which could infer the proportions of 22 TIICs subpopulations from tumor transcriptomes. This method was validated by FACS successfully and was conducted in breast cancer, and lung cancer patients [18,19]. In this study, we applied CIBERSORT to uncover distinct patterns of TIICs in ccRCC and associations of different immune cells subsets with clinical outcomes.
We observed significant differences in immune cell composition between ccRCC tumor tissue and paired normal tissue. Our data revealed the detailed profile of 22 TIICs subtypes infiltration in ccRCC that the proportions of total T cells accounted for more than 40%, in which the CD8+ T cells comprised of 14.8%. Secondly, the proportions of total macrophages accounted for more than 30%, in which 20.9% were M2 cells. Moreover, our work confirmed the findings that certain immune cells subsets can also predict clinical outcomes beyond the immunoscore. By univariate Cox regression analysis, we found that dendritic cells resting, dendritic cells activated, mast cells resting, mast cells activated and eosinophils are significantly associated with improved outcome, while B cells memory, T cells follicular helper and T cells regulatory (Tregs) indicate poorer outcome.
It is commonly believed that CD8+ T cells can recognize tumor specific antigens and play a role in tumor control [20]. High densities of tumor infiltrated CD8+ T cells are shown to be associated with favorable prognosis in vast majority of cancers [5]; therefore, many forms of immunotherapy aim at restoring T-cell mediated immune response [17]. Previous studies described ccRCC as a proinflammatory tumor where malignant cells and infiltrated neutrophils produce various kinds of cytokines that may help recruit and activate polyclonal CD8+ T cells [21][22][23]. Paradoxically, CD8+ T cells infiltration in ccRCC correlates with a poorer prognosis. Nicolas A. Giraldo et al. revealed that recruited CD8+ T cells correlated with favorable prognosis only with the present of fully functional mature dendritic cells (DC) [7]. Moreover, several studies indicated that dysfunctional DC maturation can be induced in the ccRCC microenvironment by down-regulating co-stimulatory molecules [24][25][26]. Recently, Kondou R et al. discovered that in PD-L1+ and CD8B+ patients, the gene expression profiling of fresh specimens exhibited an upregulation of dendritic cell maturation genes and T-cell activation genes. However, patients with PD-L1-and CD8Bexhibited a low expression of T-cell-activation genes [27]. Meanwhile, our results demonstrated that dendritic resting cells and dendritic activated cells are associated with favorable prognosis. These results revealed that malfunction of DC might involve in the process of T cells inhibition and become a potential combined therapeutic target.
Furthermore, research in human lung and colorectal cancer showed that CD8+ T cells are not only specific for tumor-derived antigens; but also recognize a variety of epitopes unrelated to cancer [28]. They demonstrated that abundant CD8+ T cells act as bystanders and are phenotypically heterogeneous within a tumor and across patients. From this perspective, the role of CD8+ T cells level in predicting prognosis differs between patients. Collectively, these explain the possible reasons for the correlation of CD8+ T cells infiltration with poor outcome. In our study, however, we did not find association between CD8+ T cells and overall survival; besides, the percentages of neutrophils in the entire TIICs was uncorrelated with clinical outcome, but Jensen HK et al. revealed that the presence of intramural neutrophils correlates with poor prognosis [29]. The discrepancy may be ascribed to the nature of algorithm, which lack information related to cellular heterogeneity and deeper spatial distribution.
Previous studies regarding regulatory T cells and Macrophages M2 showed that they all exert pro-tumorigenic function [30,31]. Our study revealed that T cells regulatory and T cells follicular helper are associated with poor prognosis. Interestingly, we also found the positive relationship between specific immune cells proportion and ccRCC Fuhrman grade including CD8+ T cells, T cells regulatory and T cells follicular helper. Moreover, corHeatmap indicated negative correlation between CD8+ T cells and Macrophages M2. Collectively, our work indicated that Macrophages M2, T cells regulatory and T cells follicular helper possibly play a role in T cell exhaustion/inhibition in an intricated "checks and balances" mannner [32].
Mast cells and eosinophils have been studied mainly in the allergic disease. It has been reported that mast cells and eosinophils can exhibit both tumor promoting and anti-tumor activity [33,34]. Underlying molecular mechanisms of the contradictory function have remained a conundrum. Increased mast cell infiltration predicts a poorer prognosis in many cancers including lung, colorectal, gastric, melanoma and cervical carcinoma [35]; nevertheless, mast cell infiltration is associated with improved prognosis for breast carcinoma and prostate cancer [36,37]. Our study revealed that infiltration of resting mast cells, activated mast cells and eosinophils are correlated with favorable prognosis in ccRCC. However, the fraction of resting mast cells decreased with the increased Fuhrman grade. With the merits of a low rate of cell division and long lifespan, mast cells are eligible alternative candidates for combined targeted immunotherapy. Moreover, Hollande et al. first demonstrated that using T cell-and eosinophils-targeted combination therapy yields increased tumor retardation, which advances our understanding of the anti-tumor role of eosinophils [38].
Nevertheless, our study has several limitations. We pooled together data from TCGA and GEO to enlarge our sample size, which may affect the repeatability of results on account of the heterogeneity. Secondly, our results are based on public database and computational algorithm. Although the accuracy of this technique has been testified using FACS, it is still required to further verify them by experiments in the future.
In conclusion, our analysis, based on a devolution algorithm, revealed significant differences in the cellular composition of infiltrated immune cells in ccRCC and associations between 22 immune cells subpopulations and clinical outcome. Particularly, DC, mast cells, eosinophils, T cells regulatory and T cells follicular helper emerge as potential immunotherapy targets. Moreover, the CIBERSORT algorithm makes the comprehensive analysis of ccRCC immune microenvironment possible by using gene expression data from bulk tissues.