Microenvironment-associated gene HSD11B1 may serve as a prognostic biomarker in clear cell renal cell carcinoma: a study based on TCGA, RT‑qPCR, Western blotting, and immunohistochemistry

ABSTRACT Clear cell renal cell carcinoma (ccRCC) is one of the most common malignant tumors worldwide. The clinical treatment of ccRCC is strongly associated with the tumor microenvironment (TME). Identifying potential markers of ccRCC is important to improve prognosis. Therefore, in the present study, the levels of immune/stromal components and the proportion of tumor-infiltrating immune cells (TIICs) were determined in 611 ccRCC samples using the ESTIMATE and CIBERSORT analytical tools. Subsequently, hydroxysteroid 11-beta dehydrogenase-1 (HSD11B1) was identified by univariate Cox regression analysis, protein-protein interaction (PPI) networks and clinical survival analysis to be associated with ccRCC prognosis. At the same time, the abundance of HSD11B1 increased significantly in ccRCC was verified by western blotting, RT‑qPCR and immunostaining analysis. Furthermore, Gene Set Enrichment Analysis (GSEA) and TME suggested that HSD11B1 was involved in TME immune-related status. Taken together, the results of the present study demonstrated that HSD11B1 is a potential prognostic biomarker associated with immune cell infiltration in ccRCC.


Introduction
Clear cell renal cell carcinoma (ccRCC) is the third most common type of cancer of the urinary system with a high mortality rate [1,2], accounting for 2-3% of all adult cancer cases and 70% of all renal cancer cases [3][4][5]. Currently, surgical resection, radiotherapy and chemotherapy are the main treatment approaches to ccRCC. However, these treatment strategies do not significantly prolong the survival time of patients [6,7]. Several studies have suggested that changes of specific genes expression are highly related to the tumorigenesis and development of ccRCC [8]. Hence, investigating key gene correlated with ccRCC prognosis and identifying optimal immune-related biomarkers is crucial [9].
The tumor microenvironment (TME) is extremely complex, and consists of immune cells, various types of stromal cells, alongside with tumor cells [10][11][12]. Previous studies have also suggested that immune and stromal cells are the main components of TME, which are closely related to tumor progression and clinical prognosis [13,14]. Several clinical and genomic studies have also reported that ccRCC is a highly immune infiltration tumor [15]. The components of TME have an significant influence on the occurrence and progression of ccRCC [16,17]. Emerging evidence has indicated that the tumor-infiltrating immune cells (TIICs) within TME may affect therapeutic efficacy [18].
The Cancer Genome Atlas (TCGA) data of patients with ccRCC were downloaded, and the proportions of immune/stromal cells and TIICs in samples were calculated using the ESTIMATE [19] (Estimation of Stromal and Immune cells in Malignant tumor tissues using Expression data) and CIBERSORT (A general computational method for estimating the composition of TIICs populations from gene expression data) algorithms [20]. Eventually, a predictive biomarker, hydroxysteroid 11-beta dehydrogenase-1 (HSD11B1) was selected. The HSD11B1 gene encodes the type 1 isoform of 11-β-hydroxysteroid dehydrogenase, which converts glucocorticoids into the active form, which plays an important role in the regulation of cell proliferation, differentiation and immune response [21,22]. Proliferative phenotypes induced by high expression of HSD11B1 are associated with poor cancer prognosis [23]. Previous studies suggested that increased HSD11B1 expression was interconnected with poor survival of ccRCC patients, supporting that HSD11B1 have the potential to be a prognostic biomarker in patients with ccRCC receiving immunotherapy.
Our work aims to further evaluate significance of HSD11B1 expression in ccRCC prognosis and the correlation with clinicopathological characteristics using TCGA data and in vitro experiments, thereby providing additional evidence of HSD11B1 as a prognostic biomarker associated with immune cell infiltration in ccRCC.

Data collection
Transcriptome single-cell RNA sequencing profiling data (611 cases: 539 tumor cases and 72 normal cases, workflow type: HTseq-FPKM) and related clinicopathological characteristics (Table 1), including age, gender, pathological stage, grade, T/M/N classification of ccRCC samples were downloaded from TCGA database [24] (https://portal.gdc.cancer.gov/). The ESTIMATE [19] algorithm, an open-source web tool, was used to calculate the proportions of immune/stromal cells and scores in ccRCC samples using 'estimate' package (HTTP: // R-forge.R-project. org;Repos = Rforge, dependencies = TRUE). The scores were positively correlated with the proportion of the immune/ stromal/ ESTIMATE components in TME.

Identification of differentially expressed genes (DEGs)
Subsequently, the ccRCC samples were divided into high/low-score groups according to the median score. R (version 4.0.3) and R language packages 'limma Bioconductor' [25] were applied for DEGs between the high/low-score group. FDR<0.05 and |logFC|>1 were considered significant for screening DEGs.

PPI network construction and univariate Cox regression analysis
The PPI network among the significantly enriched DEGs was constructed using the Search Tool for the Retrieval of Interacting Genes (STRING) [28] online database (https://string-db.org/) and the Cytoscape [29] (version 3.7.2) platform was then utilized to visualize the interactive network among DEGs with a confidence level > 0.95 were used for building network. The top 30 DEGs with the largest number of nodes were filtered by interactions between genes. Univariate Cox proportional hazards regression was used to further screen for meaningful DEGs with prognostic value, and forest plots were generated using R software package 'survival' (FDR<0.05).
Among the first 30 core DEGs of the PPI network with predictive value based on the univariate Cox regression analysis, five hub DEGs with strong predictive value were selected.

Association of HSD11B1 expression with survival and clinicopathological characteristics
The R packages 'survival' (https://CRAN. R-project.org/package=survival, version = 3.1-8) was applied for Kaplan-Meier survival analysis. In addition, Clinical parameters, including age, gender, pathological stage, grade, T/M/N classification, were analyzed by R software and Wilcoxon rank-sum test. P < 0.05 was statistically significant.

Gene set enrichment analysis (GSEA)
GSEA [30] [31] to identify the enriched pathways. Each enrichment analysis carried out one thousand times gene set permutations. Pathways with the false discovery rate (FDR) <0.05 and NOM P < 0.05 were considered to be significantly enriched.

Cell culture
ccRCC cell lines (HK-2, 786-O) was obtained from the Department of Nephrology, the Affiliated Hospital of Qingdao University. HK-2 and 786-O cell lines were cultured in DMEM-F12 and RPMI 1640 with 10% FBS (Gibco) and 1% streptomycinpenicillin at 37°C.

Western blotting
RIPA (Meilunbio, MA0151) buffer supplemented with protease inhibitors and phosphatase inhibitors was used to extract total protein from cells. The protein samples were electrophoresed on SDS-PAGE gels and transferred to PVDF membranes before being blocked with 5% skimmed milk for 1 h. The membrane was incubated with primary antibodys, including HSD11B1 (1:1500 dilution; Abbkine) and β-actin (1:50,000 dilution; Abcam) at 4°C overnight. Secondary goat antimouse IgG-HRP (1:5000 dilution; Abcam) antibody was incubated for 1 hour at room temperature. The proteins were visualized using ECL detection reagents (Meilunbio, MA0186) and quantitatively analyzed by Image J. Western blotting was carried out according to the experimental process of our laboratory [32].

Immunostaining analysis
Immunohistochemistry staining was performed on tissues. Renal clear cell carcinoma and adjacent para-carcinoma tissues (N = 10) were obtained from the Department of Pathology, Affiliated Hospital of Medical College Qingdao University, Qingdao, China. Immunochemical staining was performed according to our previous experiment procedure [32]. Sections were incubated with antibodies against HSD11B1 (1:200 dilution; Abbkine) and secondary antibody (1:200 dilution; Absin).

Analysis of TIICs
CIBERSORT [34] analysis tool is a deconvolution algorithm, using R software to run CIBERSORT algorithm. The CIBERSORT analytical tool can accurately estimate the proportion of TIICs from the expression profiles of complex samples. The algorithm calculates the proportion of 22 kinds of TIICs for each ccRCC sample based on LM22 signatures and 1000 matrix permutations. The proportion of immune cells is presented as a bar graph. P < 0.05 can be further correlation analysis.
The Wilcoxon rank-sum test and Spearman's correlation coefficient were applied to analyze differential and immune cell correlation analyses by using R packages 'Pheatmap', 'Corrplot' and 'Vioplot' to evaluate the association between TIICs and HSD11B1 expression. P < 0.05 was considered statistically significant.

Statistical analysis
Statistical analysis was performed using GraphPad Prism 5. The unpaired t-test (two-tailed) with confidence intervals of 95% was used to detect the differences in HSD11B1 expression between two independent groups and P < 0 .05 were considered statistically significant.

Results
In the present research, we aim to establish a predictive model to identify the key genes affecting the prognosis of ccRCC, so as to provide a theoretical basis for predicting the prognosis of ccRCC patients and further seeking new treatment options. In this study, we first downloaded 611 ccRCC gene-expressed data and 537 related clinical information from the Cancer Genome Atlas (TCGA) database. Then we screened HSD11B1 as a possible prognostic biomarker of ccRCC. In addition, we performed association analysis between HSD11B1 and immune cells to evaluate its correlation with tumor microenvironment. Finally, the bioinformatics results were verified by in vitro experiments. We determined that HSD11B1 may be a novel biomarker for the diagnosis and prognosis of ccRCC. This finding is expected to benefit ccRCC patients.

Correlation analysis between immune/stromal scores and patient survival and clinicopathological characteristics
To assess the association between survival rate and immune/stromal scores, 611 ccRCC samples (539 tumor and 72 normal cases) were divided according to scores, and Overall survival analysis was performed. Kaplan-Meier survival curve showed that the proportion of immune cells was notably associated with patients survival (P = 0.033; Figure 1(a)). Although patients with low stromal score had a higher median Overall survival rate than those with higher stromal score, no significant correlation was observed between stromal score and survival (P = 0.316; Figure 1(b)). The aforementioned findings suggested that the proportion of immune cells was closely related to the Overall survival in ccRCC patients. Furthermore, the paired clinicopathological characteristics of ccRCC patients were collected from TCGA database and their association with the proportion of immune and stromal cells was determined. The analysis revealed that the immune score was closely associated with several clinicopathological characteristics, including pathological stage (stage III vs. I, P = 0.00043; stage IV vs. I, P = 0.00022; Figure 1(c)), grade (G4 vs. G1, P = 0.01; G3 vs. G2, P = 0.00028; G4 vs. G2, P = 3.8E−06; G4 vs. G3, P = 0.027; Figure 1 (e)), T classification (T2 vs. T1, P = 0.03; T3 vs. T1, P = 9.2E−05; Figure 1(g)), and M classification (P = 0.0046; Figure 1(i)). However, stromal score was only associated with pathological stage (stage III vs. II, P = 0.026; Figure 1(d)), grade (G4 vs. G2, P = 0.034; G4 vs. G3, P = 0.036; figure 1(f)), and T classification of the tumor-node-metastasis staging system (T2 vs. T1, P = 0.04; T3 vs. T2, P = 0.0075; Figure 1(h)). Overall, the results demonstrated that the proportion of immune/stromal cells was closely related to factors promoting ccRCC progression, such as immune cell infiltration, metastasis and prognosis.

Correlation between ccRCC gene expression profiles and immune/stromal scores
ccRCC samples were divided into two groups based on their respective median immune/stromal scores. A total of 656 DEGs were selected according to scores (high vs. low scores). Among them, 510 upregulated genes and 146 downregulated genes were obtained from ImmuneScore ( Supplementary Figure 1(a) and Figure 2(a)). Similarly, among the 411 DEGs from StromalScore, 259 and 152 DEGs were up-and downregulated, respectively (Supplementary Figure S1b and Figure 2(b)). Venn plot were used to determine the intersection between the two types of cells. Therefore, 93 common DEGs were overlapping genes in immune and stromal groups, including 44 upregulated and 49 downregulated genes. The aforementioned findings indicated that DEGs may determine the TME status. Subsequently, the top 10 DEGs, significantly enriched in biological processes (BP), cellular component (CC) and molecular function (MF) were selected from GO enrichment analysis (Figure 2(c)). The functional clusters of these DEGs corresponded to immune-related GO terms, such as 'humoral immune response' and 'B-cell proliferation'. KEGG analysis revealed that DEGs were also significantly enriched in immune response-related signaling pathways, including 'cytokine-cytokine receptor interaction', 'primary immunodeficiency' and 'cytokine activity' (Figure 2(d)). Overall, the aforementioned results indicated that DEGs obtained by the intersection of immune/stromal scores were correlated with immune response.

Identification of DEGs by PPI networks and univariate Cox regression analysis
Subsequently, two methods were used to screen for DEGs. Therefore, to study interactions among DEGs, we constructed a PPI network using STRING and Cytoscape software. The top 30 DEGs with the largest number of adjacent nodes were selected for further analysis (Figure 3(a)). Furthermore, univariate Cox regression analysis was used to identify the DEGs significantly associated with prognosis in ccRCC patients (Figure 3  (b)). The results of both methods were combined into a Venn plots, and five overlapping genes, namely HSD11B1, TNFSF13B, MZB1, IGLL5 and PPARGC1A, were obtained (Figure 3(c)). The details of the five overlapping genes are presented in Table 2.

Association between HSD11B1 expression levels with survival and clinicopathological characteristics
The HSD11B1 gene encodes the type 1 isoform of 11-β-hydroxysteroid dehydrogenase, which converts glucocorticoids into the active form, which plays an important role in the regulation of metabolic syndrome and immune response. This finding suggested that HSD11B1 may play a key role in TME and be a key factor in ccRCC. Therefore, the present study indicated that HSD11B1 potential be a prognostic marker in ccRCC. Herein, ccRCC samples were divided into two groups using the median HSD11B1 expression level as cutoff values. Wilcoxon rank-sum test showed that HSD11B1 expression in tumor group was significantly higher than that in normal group (Figure 4(a)). Furthermore, increased HSD11B1 expression was positively correlation with unfavorable prognosis (Figure 4(b)). In addition, HSD11B1 expression was closely related to several clinicopathological characteristics with ccRCC patients. Therefore, pathological stage (III vs. I, P = 0.041; IV vs. I, P = 0.036), grade (G4 vs. G1, P = 0.00021; G4 vs. G2, P = 6.6E−08; G4 vs. G3, P = 1.6E−06), T classification (T3 vs. T1, P = 0.022, T4 vs. T1, P = 0.0061, T4 vs. T2, P = 0.0094) and N classification (P = 0.0043) were strongly correlated with HSD11B1 expression (Figure 4(c-g)).
Overall, the aforementioned results suggested that HSD11B1 expression levels were negatively related to the prognosis in ccRCC patients.

GSEA of HSD11B1
As shown in Figure 4(h-i), GSEA was used to evaluate the correlation between HSD11B1 expression with immune response. The bioinformatics analysis demonstrated that several immunerelated signaling pathways, including 'B-cell receptor signaling pathway', 'cytokine-cytokine receptor interaction' and 'chemokine signaling pathway' were primarily enriched in C2 KEGG gene sets  in HSD11B1 high expression group (Figure 4(h)). In addition, the C7 immunological gene sets (Figure 4(i)) were enriched in the high expression group of HSD11B1, including 'CD4 T cell vs. CD8 T cell up', 'CD4 T cell vs. alphabeta CD8 T cell up', 'naive vs. secondary memory CD8 T cell up'. The result showed that HSD11B1 is not only associated with immune-related activities, but may also be a potential biomarker of the TME status.

Positive expression of HSD11B1 in 786-O cell lines and ccRCC clinical specimens
The high expression of HSD11B1 in cell lines was verified in protein level compared to controls, which is consistent with our bioinformatic analyses ( Figure 5(c)). Furthermore, HSD11B1 expression in clinical tissues ( Figure 5(a-b)) were observed by immunochemical staining, which showed that significantly higher HSD11B1 expression in ccRCC clinical tissues than that in paracarcinoma specimens. These results confirmed that the abundance of HSD11B1 increased significantly in ccRCC. RT-qPCR (P&lt;0.0001) showed that HSD11B1 was up-regulated in ccRCC ( Figure 5(d)). This result is consistent with the result in Figure 4.

Association between the expression of HSD11B1 and the proportion of TIICs
To further elucidate the association between HSD11B1 expression and the TME status in ccRCC, the proportion of TIICs in ccRCC was analyzed using the CIBERSORT algorithm ( Figure 6(a)). The violin plot (Figure 6 (Figure 7(a)). Subsequently, the results of the differential analysis and those from the immune cell association analysis were intersected (Figure 7(b), Table 3). The violin plot showed the correlation analysis between the infltration level of immune cells and low (green) and high (red) HSD11B1 expression groups relative to the median of HSD11B1 expression. Wilcoxon rank sum test was commonly believed to be significantly.

Discussion
In the present study, ccRCC samples from TCGA database were used to elucidate the genes associated with prognosis. Bioinformatics analysis results indicated that HSD11B1 could be a novel prognostic biomarker of ccRCC. TME serves an important part in tumor development and progression. Nowadays, the research into potential prognostic targets for advanced renal cancer is growing rapidly. Compared with traditional surgical resection, patients undergoing targete.d therapy have a significantly longer survival time [7,35,36]. Identifying potential prognostic biomarkers is of great importance for the diagnosis, treatment and prognosis of cancer targeted immunotherapy. Transcriptome analysis of ccRCC in TCGA database indicated a close association between ImmuneScore with survival of ccRCC patients. Furthermore, it was significantly correlated with multiple clinicopathologic characteristics of TCGA-ccRCC samples, including grade, pathological stage, T classification, and M classification. It further demonstrated that the components of the immune system within TME could serve a vital role in the development and prognosis of ccRCC. GO analysis confirmed a significant enrichment of DEGs in immune-  related terms. Emerging evidence has suggested that immune responses promote the development of several types of cancer [37][38][39]. Herein, KEGG enrichment analysis showed that DEGs were closely linked to immune-related activities. These findings indicated that DEGs may have a significant impact on regulating TME and influencing the prognosis of ccRCC. Since The HSD11B1 gene encodes the type 1 isoform of 11β-hydroxysteroid dehydrogenase, which converts glucocorticoids into the active form, which plays an important role in the regulation of metabolic syndrome and immune response [21]. Further analyses were carried out on HSD11B1 expression profile. These studies revealed that HSD11B1 was upregulated in ccRCC samples and can be considered as prognostic indicator of ccRCC, as predicted by univariate Cox regression analysis. HSD11B1 elevated expression was significantly related to poor survival rate and several clinicopathological characteristics. The results obtained by western blotting, RT-qPCR and immunohistochemistry were consistent with the results of this study. Overall, these findings demonstrated convincingly that HSD11B1 is a potential valuable prognostic marker. HSD11B1 is a member of the short-chain dehydrogenase/ reductase (SDR) superfamily, which is the only enzyme in the body that converts inactive cortisone into active glucocorticoid hormones [40]. The enzyme is abundant in liver and can be induced to express in immune cells. Previous studies have shown that glucocorticoids encoded by HSD11B1 contribute to the regulation of cell proliferation and differentiation, and recognized the important role HSD11B1 play in promoting the development and progression of colon tumors and adenocarcinoma cells [22,23,41,42]. Proliferative phenotypes induced by high expression of HSD11B1 are associated with poor cancer prognosis [23]. Herein, increased HSD11B1 expression was interconnected with poor survival of ccRCC patients, supporting a potential association between HSD11B1 expression and tumor progression and clinical outcome. These findings indicated that HSD11B1 could act as oncogene in ccRCC, which was consistent with the tumorpromoting effect of HSD11B1 on other types of cancer. GSEA showed that HSD11B1 expression was associated with several immune-related signaling pathways, suggesting that the immune responses may be activated during the progression of ccRCC. In addition, the analysis revealed that the up-regulation of chemokine signaling pathway also promoted the development of CCRCC. Previous studies have shown that chemokines, such as CXCL13 and CXCR5, are crucial for the progression and poor prognosis of ccRCC [43,44].
The results also demonstrated that seven kinds of TIICs were significantly associated with the TME status in ccRCC. In the present study, the proportion of Macrophages M0 was positively associated with HSD11B1 expression. Compared with the low HSD11B1 expression group, the number of Macrophages M0 was obvious abundant in HSD11B1 high expression group. Clinical trials indicated that Macrophages M0, as inhibitors of the antitumor immune responses, were related to poor prognosis in sizable proportion of all cancers, such as lung [45,46], prostate [47], colorectal [48], breast [49,50], hepatocellular [51] and head and neck [52] cancer. Numerous studies have verified that macrophages, as an important member of myeloid origin cells in TME, can participate in the occurrence, development and immunosuppression of tumors [53]. Another immunosuppressive cell type in cancer, Plasma cells, may significantly affect the survival of patients with ccRCC. Plasma cells are known to be products of B-cell differentiation that play a key role in humoral immune system by producing and secreting antibodies [54]. More and more studies have shown that plasma cells have become a significantly survival biomarkers for a variety of solid tumors [55,56]. Recent studies have also confirmed that the results of the present study, demonstrating that increased numbers of Plasma cells were closely associated with poor survival. In the majority of solid tumors, CD8 T-cell infiltration is often associated with good prognosis [57]. Consistently, the results of this study showed a negative correlation between CD8 T cell infiltration and HSD11B1 expression. Our results support B cells provided the robust anti-tumor immunity by served as antigenpresenting cells [58,59]. Other studies have reported that immune cells, for instance, NK cells activated and Dendritic cells activated are involved in tumor immune surveillance and exert antitumor effects. A study showed that NK cells could mediate cytolysis to induce apoptosis of target cells [60]. In addition, NK cells can enhance the antitumor activity via secreting proinflammatory chemokines and cytokines [61]. A recent study demonstrated that NK cells, as immune effector cells, were potentially involved in the treatment of kidney cancer via enhancing chimeric antigen receptor (CAR) modification [62]. As mentioned above, active Dendritic cells can activate tumor immunity and aggregate immune effector cells to the tumor site [63]. Therefore, the negative association of the proportion of NK cells activated, Dendritic cells activated with the expression levels of HSD11B1 in ccRCC patients indicated that HSD11B1 potential play a pivotal role in tumorigenesis and development.

Conclusion
In conclusion, the results of the present study supported that the interaction between ccRCC and TME affected cancer progression, thereby affecting the overall prognosis of patients with ccRCC. The analysis of TCGA-ccRCC patient samples and in-vitro experimental verification revealed that HSD11B1 as a potential prognostic indicator of ccRCC. Finally, further studies on HSD11B1 could provide a foundation for understanding the complex association between TME and the prognosis of patients with ccRCC.

Highlights
• The prognosis of ccRCC was correlated with immune regulation of TME. • The expression level of HSD11B1 was correlated with the prognosis of ccRCC patients. • HSD11B1 can be used as a prognostic marker in TME of ccRCC.