Prognostic Value of Genes and Immune Infiltration in Prostate Tumor Microenvironment

Background: Prostate cancer (PCa) is one of the most common cancers and the fifth leading cause of cancer-related death in men. Immune responses in the tumor microenvironment are hypothesized to be related to the prognosis of PCa patients; however, no studies are available to confirm the same. In this study, we aimed to explore the potential link between these two factors and identify new biomarkers to estimate the survival rate of PCa patients.

Methods: A total of 490 cases were obtained from The Cancer Genome Atlas (TCGA) database. The gene expression data were analyzed by the ESTIMATE algorithm to evaluate the immune and stromal scores. The survival rate was calculated according to the case-specific clinical data. Enrichment analysis was performed to discover the main biological processes and signaling pathways of immune responses. We further identified and analyzed hub genes in the protein-protein interaction (PPI) network and evaluated their prognostic values.

Results: Immune score significantly correlated with immune cell infiltration and overall survival of PCa patients. The genes CXCR4 and GPR183, identified as hub genes in the PPI network, correlated with immune cell infiltration and prognosis of PCa patients.

Conclusion: CXCR4 and GPR183 participate in immune cell infiltration and function in PCa patients. The immune score, as well as the expression of CXCR4 and GPR183 in prostate cancer tissues, could be potential indexes for the prognosis of prostate cancer.

Introduction

Prostate cancer (PCa) is one of the most common cancers of men and the fifth leading cause of cancer-related morbidity in men worldwide (1). PCa patients during the early stages are prescribed with radical prostatectomy and radiotherapy, which have a high rate of complete cure (2). For metastatic PCa patients, surgery or androgen deprivation is advised. However, androgen deprivation therapy or androgen receptor-targeted therapy can also induce tumor resistance (3). Therefore, new and promising therapeutic strategies for PCa are required. Anti-tumor immunotherapy with CTLA-4-targeted monoclonal antibodies is now under clinical investigation (4, 5). However, the precise relationship between the prognosis of PCa and the immune responses is still unknown.

Currently, prostate-specific antigen (PSA) values, histopathological scores such as Gleason score, are used as clinical parameters to diagnose PCa and assess the risk stratification (6). However, PSA is not a tumor-specific marker and the specificity of PSA is only 12.8% (7). This led to the introduction of prostate health index (PHI) to elevate the diagnostic capacity of PSA (8). The Oncotype DX Genomic Prostate Score (GPS) assay was used in clinical trials to make a risk assessment and predict the time of recurrence (9, 10). But this multi-gene assay may be unable to select the active surveillance candidates (11). Identification of specific markers with significant potential for the prognosis of PCa in terms of immune response is necessary.

The tumor microenvironment (TME) of the PCa patients is associated with inflammation (12). An increasing number of studies have reported that infiltrated immune cells play a pro-tumorigenic or anti-tumorigenic role in the TME of PCa patients (13–17). The Cancer Genome Atlas (TCGA) database has provided a series of global gene expression profiles and clinical data of the patients worldwide (18). Yoshihara et al. developed the ESTIMATE algorithm to assess the expression levels of specific molecular entities in stromal and immune cells of the TME (19), so that the non-tumor cell infiltrations in the TME can be predicted. Shah et al. used the ESTIMATE algorithm to evaluate the stromal score of prostate cancer; however, the investigators did not analyze the TCGA database or compare the survival rate in different groups based on estimated score (20).

In this study, we explored the prostate cancer data from TCGA databank, calculated the immune and stromal scores of every sample by ESTIMATE algorithm, and estimated their potential values of prognosis for PCa patients. To identify some specific genes to forecast the overall survival rate of the PCa patients, we performed enrichment analysis, and interaction analysis and identified two hub genes. These hub genes were found to be highly correlated with the prognosis of PCa and associated with tumor purity and immune infiltration. We hypothesized that these two genes can be the potential biomarkers for prognosis of PCa and guide the selection of immunotherapy strategy for the PCa patients. Further, these genes may also play a significant role in the underlying molecular mechanisms of PCa.

Materials and Methods

Raw Data From TCGA Database

We searched the TCGA database (https://portal.gdc.cancer.gov) and downloaded prostate cancer datasets after restricting the disease types to adenomas and adenocarcinomas. A total of 490 cases were included in this dataset that contained the gene expression files and clinical data of every patient. These gene expression files were processed by the ESTIMATE algorithm, (19) and the immune, stromal, and ESTIMATE scores of each sample were calculated. Since we only use the data from public online database in this study, so the ethic approval was not required.

Survival Rate Analysis

Included participants were divided into two groups according to their immune, stromal, and ESTIMATE scores, and the expression levels of the candidate genes (separated by median values). We used Kaplan-Meier correlation analysis (95% CI) to evaluate the association between the overall survival rate and the immune, stromal, and ESTIMATE scores, and the identified genes expression levels. A log-rank test was used to check the significance of the relationships. The analyses were carried out using R software (version 3.6.1), and P<0.05 was considered statistically significant.

Correlation Analysis

Included participants were divided into sub-groups to evaluate the relationships between different tumor stages including T1-T4, N0-N1, and M0-M1 stages and the estimated immune, stromal, and ESTIMATE scores. Wilcoxon signed-rank test was used to evaluate the relationship between the two groups; the Kruskal-Wallis test was used to evaluate the relationships among three or more groups. The analyses were conducted using R (version 3.6.1), and P<0.05 was considered statistically significant.

Identification and Cluster Analysis of Differentially Expressed Genes

Gene expression files were processed by limma package (21) of R and the differentially expressed genes (DEGs) in immune-score-high and stromal-score-high groups were identified. The cut-off values were set as |fold change|>2 and adjust. P<0.05. Cluster analysis and heatmaps were generated by pheatmap package (22) of R, and the upregulated genes in both the high-rank groups were calculated by the VennDiagram package (23). Further, we analysed the immune cell specific markers such as CD14 (monocytes) and CD3 (T cells) between the two groups.

Enrichment Analysis

Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Gene and Genomes (KEGG) pathway enrichment analysis of the upregulated genes were performed by R Packages clusterProfiler (24), enrichplot, and ggplot2 (25), and the significantly activated biological processes, molecular functions, cellular components, and signaling pathways were explored.

Protein-Protein Interaction (PPI) Network and Hub Genes Calculation

STRING tool (https://string-db.org/) was used to construct the Protein-Protein Interaction (PPI) network of the identified genes (CI=0.90). The network was reconstructed using CytoHubba, a plug-in of Cytoscape software (version 3.5). Using this software, we calculated the connection degree of every node in the network to identify the top 30 hub genes.

Specific Gene Expression Comparison and Immune Infiltration Analysis

TIMER tool (https://cistrome.shinyapps.io/timer/) was used to compare the specific gene expression levels between the tumor tissues and normal control tissues of different carcinomas. The correlation between specific gene expression and infiltration of different immune cells in PCa tissues was assessed by purity-corrected partial Spearman’s correlation analysis. Also, the relationship between the two particular genes was calculated by Spearman’s correlation analysis. P<0.05 was considered statistically significant.

Gene Set Enrichment Analysis (GSEA)

The gene expression files of the 490 cases downloaded from TCGA was separated into CXCR4 (or GPR183) high-express group and low-express group. Gene set enrichment analysis for KEGG pathway was conducted by using GSEA tool (https://www.gsea-msigdb.org/gsea/index.jsp). P<0.05 was considered as statistically significant.

GEO Datasets Analysis

We searched GEO database and downloaded the human PCa-related datasets. The transcriptome of CXCR4 and GPR183 were compared between tumor and paired normal tissues (GSE69223) (26), recurrent patients, and non-recurrent patients after surgical operation (GSE25163) (27), as well as castration-resistent patients and good prognosis patients (GSE37199) (28). P<0.05 was considered as statistically significant.

Results

Correlation Between Immune Score and Prognosis of Prostate Cancer Patients

Transcriptional expression files and clinical data of 490 prostate cancer (adenomas and adenocarcinomas) patients were downloaded from the TCGA databank. Out of these patients, one was American Indian (0.20%), 12 were Asian (2.45%), 58 were African American (11.84%), 407 patients were white (83.1%), and the race of rest of the 12 patients (2.41%) was unknown. All the patients whose data is included in this study were diagnosed with adenomas or adenocarcinomas.

Using the ESTIMATE algorithm, we evaluated the stromal, immune, and ESTIMATE scores of every sample. The stromal score of PCa patients varied from -1,925.26 to 1771.63, and the immune score varied from -1,339.23 to 1,646.97. The total ESTIMATE score, which is the combination of stromal score and immune score, ranged from -3,264.49 to 3,418.60.

To estimate the prospective relationship between the overall survival of PCa patients and their stromal, immune, and ESTIMATE scores, we separated the selected cases into a high-rank group and low-rank group based on these scores and compared the survival rate in the two groups. Intriguingly, the Kaplan-Meier survival curves revealed that PCa patients with higher immune scores live longer than those with lower scores (Figure 1A). Moreover, no significant differences were observed in the prognosis between stromal-scores-high and stromal-scores-low groups and between ESTIMATE-scores-high and ESTIMATE-scores-low groups (Figures 1B–C). So, we further compared the immune scores among the PCa patients under different stages, but no significant differences were observed in the immune scores among the different T stages, N stages, and M stages (Figures 1D–F). Further, we found that the immune score is significantly correlated with the overall survival of PCa patients, and it could be used as a potential index to forecast the prognosis of PCa.

FIGURE 1

Figure 1 Prostate cancer (PCa) patients with higher immune score have better prognosis. (A) PCa patients were divided into two groups according to the immune scores, and the survival rate of these two groups was shown as Kaplan-Meier survival curve, which indicated the PCa patients with higher immune score have better prognosis (P = 0.005). (B) PCa patients were divided into two groups according to the stromal scores, and their survival rates did not show any significant difference (P = 0.772). (C) PCa patients were divided into two groups according to the ESTIMATE scores, and there was no significant difference on survival rates between these two groups (P = 0.397). (D) Immune scores of PCs patients in each tumor stage were shown by box-plot, but no significant association between the immune score and tumor stage was found (P = 0.855). (E) Immune scores of PCs patients with or without lymph node metastasis were shown by box-plot, but there was no significant difference between the two groups (P = 0.128). (F) Immune scores of PCs patients with or without metastasis were shown by box-plot, and no significant difference was shown between the two groups (P = 0.128).

Immune-Related Genes Were Highly Expressed in High-Rank Groups of Immune Scores and Stromal Scores

DEGs analysis was performed using R software to compare the DEGs between the high rank and the low-rank groups of immune score and stromal score. A total of 1467 genes were found upregulated and 9 genes were found downregulated in the immune-scores-high group (Supplement Figure 1); 1,712 genes were found upregulated and 14 genes were found down-regulated in the stromal-scores-high group (Supplement Figure 2). The expressions of immune cell CD markers were also estimated. Expression levels of cell-specific markers such as CD14 (monocytes), CD3 (T cells), CD4, CD8, CD19 (B cells), and CD163 (macrophages); co-stimulatory factors such as CD28 and CD40; cell activation markers such as CD48 and CD79, were found to be higher in the immune-scores-high group as compared with an immune-scores-low group (Figure 2A). Moreover, 883 genes were found to be commonly upregulated in the immune-scores-high group and the stromal-scores-high group (Figure 2B).

FIGURE 2

Figure 2 PCa patients with higher immune score have more immune cell infiltration and more activated immune responses. (A) Heatmap of the CD markers gene expression levels between the immune score high group and low group (P < 0.05, Fold Change>2). (B) The number of commonly upregulated genes in immune score high group and stromal score high group shown by Venn diagram. (C) Top 10 GO terms, including BP, CC, and MF, respectively, enriched according to the commonly upregulated genes (P < 0.05). (D) Top 20 KEGG terms enriched according to the commonly upregulated genes (P < 0.05).

Next, we performed GO and KEGG enrichment analyses on these 883 commonly upregulated genes to identify their main functions. The top 10 of GO items including BP, CC, and MF items are listed in Figure 2C. The most significantly enriched GO items were adaptive immune response (GO: BP), plasma membrane (GO: CC) and antigen binding (GO: MF). Besides, these commonly upregulated genes were highly enriched in cytokine and cytokine receptor interaction, followed by items such as Th17 cell differentiation, Th1, and Th2 cell differentiation, and so forth (Figure 2D).

In conclusion, the PCa tumors with higher immune scores have more immune cell infiltration with monocytes, T-cells, B-cells, and macrophages. Also, the commonly upregulated genes pointed toward activated immune responses in the TME.

CXCR4, a Hub Gene in the PPI Network, Correlated With PCa Prognosis

PPI network was created through the STRING tool to analyze the connections of identified genes. GNG2, C3AR1, and C3 were located at the center of the network. Most chemokine and chemokine receptors such as CXCL3, CCR4, CXCR4, and CXCR3 were found to be tightly connected and located mainly at the lower-left region (Figure 3A); cytokine and cytokine receptors such as IL6, IL10, IL2RA, and IL2RB, co-stimulatory factors such as CD28, and CTLA4, and other important factors including IRF4, JAK3, and BTK were located at the upper-right region of the main network. To identify the critical genes in this network, we calculated the connection degree of each node in the network and identified the top 30 hub genes (Figure 3B). GNG2, with 55 connections with other nodes, was the most highly connected one, followed by C3, C3AR1, BDKRB2, ADCY7, FPR1, CCR5, PTAFR, FPR3, and CXCR4. We built a new network based on these 30 hub genes using CytoHubba plug-in of Cytoscape (version 3.5) (Figure 3C).

FIGURE 3

Figure 3 CXCR4 and GPR183 were hub genes from the Interaction Network of upregulated genes, and associated with PCa prognosis. (A) Protein-Protein Interaction Network (CI=0.90) based on the upregulated genes, and the size and color of node indicate the connection degree. (B) Top 30 hub genes calculated by the connection degree of each node. (C) Interaction Network of top 30 hub genes. (D) PCa patients were divided into two groups according to the CXCR4 expression level, and the survival rate of these two groups indicated the PCa patients with higher CXCR4 expression level have better prognosis (P = 0.035). (E) Similarly, the PCa patients with higher GPR183 expression level have better prognosis, although it was not significant (P = 0.062).

To identify the potential prognostic markers of PCa, we analyzed the relationship between the expression levels of these 30 hub genes and the prognosis of PCa (Figures 3D, E, and Supplement Figure 3). CXCR4 was found to be correlated with PCa prognosis (Figure 3D), and PCa patients with elevated expression of CXCR4 in tumor tissues have a higher survival rate than those with relatively low expression of CXCR4. Besides, GPR183 was another critical gene, which tended to be upregulated in PCa patients with a better prognosis (Figure 3E), but the difference in the GPR183 expression levels between PCa and normal control patients was not statistically significant (P=0.062). However, there is no significant differences of CXCR4 and GPR183 expression levels among different tumor stage patients (Supplement Figure 4).

Therefore, CXCR4 and GPR183 can be potential candidate biomarkers to forecast the prognosis of PCa patients.

CXCR4 and GPR183 Are Down-Regulated in PCa Tissue and Are Associated With Immune Cell Infiltration in PCa Tumors

We used the TIMER tool to elucidate the role of CXCR4 and GPR183 in PCa tumors. We compared the expression levels of the different genes from the tumor and the normal tissues (control) using the information from the TCGA databank. The expression level of CXCR4 was found to be similar between PCa tumor and normal tissues (Figure 4A); however, GPR183 was down-regulated in the PCa tumor tissues (Figure 4B) as compared with the normal tissues. To confirm the different expression pattern, we searched the GEO database and reanalyzed the GSE69223 dataset (26), which compare the transcriptome of PCa tissue and paired normal tissues. Consequently, CXCR4 and GPR183 are both found to be downregulated in PCa tissues (Supplement Figure 5).

FIGURE 4

Figure 4 CXCR4 and GPR183 correlated with immune cell infiltration in PCa tissue. (A) CXCR4 expression levels in different kinds of tumors and their control tissues according to the TCGA data, and there was no significant difference between prostate cancer tissue and its control normal tissue. (B) GPR183 expression levels in different kinds of tumors and their control tissues according to the TCGA data, and GPR183 was significantly decreased in prostate cancer tissue (P < 0.05). (C) Correlation analysis on CXCR4 and the tumor purity (Correlation Index=-0.412, P < 0.001) and immune cell infiltration of prostate cancer. The higher expression level of CXCR4 in tumor tissue was accompanied by increased infiltrated B cells, CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells (P < 0.001). (D) Correlation analysis on GPR183 and the tumor purity (Correlation Index=-0.427, P < 0.001) and immune cell infiltration of prostate cancer. The higher expression level of CXCR4 in tumor tissue was accompanied by increased infiltrated B cells, CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells (P < 0.001). (E) The expression level of CXCR4 and GPR183 in prostate cancer tissue was highly correlated (Correlation Index=0.842, P < 0.001).

According to the dataset of GES25136 (27), we analyzed the expression of CXCR4 and GPR183 in tumor of recurrent PCa patients and non-recurrent PCa patients. No significant differences of CXCR4 and GPR183 expression pattern was found between these two groups (Supplement Figure 6).

Further, neither CXCR4 nor GPR183 show a significantly different expression level in the blood of castration-resistant PCa patients and good prognosis PCa patients, when we reanalyzed the data from GSE37199 (28) (Supplement Figure 7).

Next, partial Spearman’s correlation analysis was used to assess the relationship among the expression of gene CXCR4, gene GPR183, and the immune cell infiltration levels in PCa tumors. The CXCR4 expression level was significantly associated with purity (correlation=-0.412, P<0.001) and positively correlated with the infiltration of B cells, CD4+ T-cells, CD8+ T-cells, macrophages, neutrophils, and dendritic cells (Figure 4C). Besides, GPR183 expression level had a negative relationship with purity (correlation=-0.427, P<0.001); infiltration of B cells, CD4+ T-cells, CD8+ T-cells, macrophages, neutrophils, and dendritic cells increased with the increase in expression of GPR183 (Figure 4D). Upon comparing the expressions of both CXCR4 and GPR183 in PCa tumors, we observed a high correlation (Figure 4E).

In addition, we separated the transcriptome data of the 490 cases into CXCR4-high (bigger than the median of all cases) and –low (smaller than the median of all cases) groups and GPR183-high and –low groups, respectively, to do gene set enrichment analysis (GSEA). Intriguingly, GSEA for KEGG pathway indicated that the CXCR4-high group are mostly enriched in Cytokine-Cytokine Receptor Interaction (M9809), Nature Killer Cell Mediated Signaling (M5669), Toll-Like Receptor Signaling (M3261), T Cell Receptor Signaling (M9904), and so on (Supplement Figure 8). Similarly, GPR183-high group are mainly enriched in MAPK Signaling (M10972), Cytokine-Cytokine Receptor Interaction (M9809), Toll-Like Receptor Signaling (M3261), and Nature Killer Cell Mediated Signaling (M5669), etc. (Supplement Figure 9). Therefore, we considered that CXCR4 and GPR183 have a deep relationship with immune reaction and immune cell function in PCa tissue.

In conclusion, CXCR4 and GPR183 were down-regulated in PCa tumors, and both of these genes were associated with tumor purity and immune cell infiltration.

Discussion

Prostate cancer is the second most common tumor among males worldwide (29). The incidence of PCa in the Chinese population is much lower than that of European and American countries, but the cases are substantially rising since the beginning of the 21^st century (30, 31). Since some immunotherapy-based therapies have started showing promise in PCa treatment (32, 33), the pro-tumorigenic or anti-tumorigenic role of infiltrated immune cells attracted more attention than ever before. However, there is a lack of valid prognostic biomarkers to evaluate the immune status of TME and predict the survival rate of PCa patients.

Using the TCGA databank, we analyzed the expression profiles of prostate cancer patients, explored the TME, and selected hub genes with significant prognostic value. Further, we assessed the immune, stromal, and ESTIMATE scores of the PCa patients using the ESTIMATE algorithm and divided them into two groups (high-rank group and low-rank group) according to the median values of the scores. We compared the survival rates of the high and low-rank groups and observed that the PCa patients with higher immune score live longer than those with a lower score, contrary to the results reported for the other tumors such as glioblastoma and breast cancer (34–36). These variations might be attributed to the type of immune cells in TME since the majority of the T-cells identified in glioblastoma and breast cancer tissues are regulatory T-cells or exhausted effect T-cells. Immune escape is one of the most critical reasons for tumorigenesis, and immune cell infiltration can support the chemotherapies (37). For instance, tumors with an increased number of CD8⁺ T-cells along with some Foxp3⁺ Treg cells infiltrated in the TME are more likely to respond to chemotherapies (38). Further, we compared the immune scores according to the tumor stages, but no significant differences were observed among the different T stages, N stages, and M stages. Since the development and transformation of tumor cells depend more on their characteristics rather than the immune cell infiltration and reactions, immune scores make less contribution in evaluating the tumor stages.

After comparing the expression profiles of the prostate tumor tissues, we found a total of 883 genes to be commonly upregulated in the stromal/immune score high-rank groups. A large number of immune-related CD markers were found to be highly expressed in the immune score high group. Therefore, we concluded that the monocytes, T-cells, B-cells, and macrophages may infiltrate the TME, consistent with a published study that used mass cytometry to identify the immune cells in the human prostate (39). As per the GO- and KEGG enrichment analysis, many of the upregulated genes were found to be involved in the immune responses such as adaptive immune response, immune-regulating cell surface receptor signaling pathway, and cytokine-cytokine receptor interaction. After selecting the 30 hub genes from the commonly upregulated genes by analyzing the PPI network, we performed survival rate analyses according to the expression levels of these 30 hub genes. CXCR4 was found to be correlated with PCa prognosis, and GPR183 was found to be upregulated in PCa patients with a better prognosis. Both TCGA dataset and GEO dataset confirmed the lower expression levels of CXCR4 and GPR183 in PCa tissue when compared with normal tissue. However, the mortality of PCa patients is very low, so that the number of samples in survival curve is a limitation in this research. We believe a follow-up study based on a large cohort of PCa patients is necessary to confirm the prognosic value of CXCR4 and GPR183.

CXCR4 is a chemokine receptor mainly expressed in most hematopoietic cells and is specific for stromal cell-derived factor-1 (SDF-1) (also known as CXCL12) (40, 41). SDF-1 activates cells through the receptors CXCR4 and CXCR7, and these two receptors are expressed in different tumors either individually or in combination (42). CXCR7 rather than CXCR4 is expressed on most of the human glioblastoma cells and small-cell lung cancer cells (43). CXCR4 is a critical receptor to modulate tumor-stromal interactions including cell invasion and migration, and therapeutic resistance (44), and blocking of CXCR4 by AMD3100 increases the anti-tumor effect of docetaxel in PCa patients with tumor metastasis (45). But, we should not ignore that the high expression of CXCR4 was only observed in bone metastasis lesions and that the difference between the total expression levels of CXCR4 in PCa tissues in situ and the normal control tissues was not statistically significant (Figure 4A). The CXCR4/CXCL12 axis plays a key role in immune surveillance of tissues (46), since CXCR4 is expressed in macrophages, monocytes, T lymphocytes, B lymphocytes, and neutrophils these cells can be recruited by the expression of CXCL12 by stromal cells of TME. Thus, while separating the PCa patients by the expression levels of CXCR4, we found that the patients with a higher level of CXCR4, which may recruit more immune cells, have a better prognosis. The CXCR4 expression level was found to be significantly associated with tumor purity and positively correlated with the infiltration of B-cells, CD4+ T-cells, CD8+ T-cells, macrophages, neutrophils, and dendritic cells (Figure 4C). However, the exact immune cell atlas in the TME of PCa patients with different CXCR4 expression profiles need to be further identified and compared to confirm the potential phenotypical and functional differences in the infiltrating immune cells.

GPR183 (or EBI2), upregulated in primary B lymphocytes after Epstein-Barr virus (EBV) infection, is predicted to encode a G-protein coupled receptor that is closely related to the thrombin receptor (47). GPR183 plays an important role in promoting B-cell localization to the outer follicle and mediating B-cell migration, together with CXCR5, to regulate the germinal center reactions (48, 49). Besides, GPR183 can promote follicular helper T (Tfh) cells differentiation through the positioning of the activated T cells at the follicle-T-zone interface and mediating the dendritic cells to induce T/B cell responses (50–52). However, no studies have reported any links between GPR183 and prostate cancer. This might be because the expression of GPR183 is mostly specific to the immune cells especially to the B cells, and most of the earlier research about PCa did not focus on the immune cells in TME. In the present study, GPR183 was found to be downregulated in PCa tissues as compared with normal control tissues and upregulated in both the high-rank groups (immune score- and stromal score high). Besides, GPR183 expression levels negatively correlated with purity and positively correlated with the infiltrated B-cells, CD4+ T-cells, CD8+ T-cells, macrophages, neutrophils, and dendritic cells (Figure 4D). These findings may indicate that the B cells in TME of the prostate tumors are more active in migration. The group with lower scores, as well as lower expression of GPR183, showed a significantly decreased survival rate that may be attributed to the higher level of migration. However, further studies are needed to elucidate the exact underlying mechanisms.

In conclusion, using the ESTIMATE algorithm, we estimated the immune and stromal scores of the TCGA PCa cohort and concluded that the patients with a higher immune score have a better survival rate. CXCR4 and GPR183 are the two hub genes with significant prognostic value for PCa patients, which may attribute to their contribution to the immune cell infiltration and immune reaction.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: https://portal.gdc.cancer.gov.

Author Contributions

Conceptualization: LJ. Data curation: WGS. Formal analysis: WGS and HS. Funding acquisition: LJ, WJS, and WGS. Investigation: HS and WGS. Methodology: HS, WGS, and ZY. Project administration: LJ and WJS. Software: ZY and XX. Supervision: LJ and WJS. Validation: LX, XQ, and WJS. Visualization: WGS, ZY, and XX. Writing—original draft: HS and WGS. Writing—review and editing: LJ and WJS. All authors contributed to the article and approved the submitted version.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2020.584055/full#supplementary-material

Supplementary Figure 1 | Heatmap of differently expressed genes of immune score high group and low group.

Supplementary Figure 2 | Heatmap of differently expressed genes of stromal score high group and low group.

Supplementary Figure 3 | Survival rate of the top 30 hub genes’ high-express group and low-express group.

Supplementary Figure 4 | Expression level of CXCR4 and GPR183 among different tumor status (T0, T1, T2, N0, N1 and M0, M1) of PCa patients.

Supplementary Figure 5 | Expression level of CXCR4 and GPR183 between tumor and paired normal tissue of prostate from PCa patients. This data is re-analysied by us according to the GEO dataset GSE69223.

Supplementary Figure 6 | Expression level of CXCR4 and GPR183 in PCa tissue between recurrent and non-recurrent patients. This data is re-analysied by us according to the GEO dataset GSE25136.

Supplementary Figure 7 | Expression level of CXCR4 and GPR183 in blood between castration resistant patients and non-resistent patients. This data is re-analysied by us according to the GEO dataset GSE37199.

Supplementary Figure 8 | Gene set enrichment analysis (GSEA) for the CXCR4 high-expressed group of PCa. Top10 enriched KEGG pathways are shown in this figure.

Supplementary Figure 9 | Gene set enrichment analysis (GSEA) for the GPR183 high-expressed group of PCa. Top10 enriched KEGG pathways are shown in this figure.

References

1. Wong MC, Goggins WB, Wang HH, Fung FD, Leung C, Wong SY, et al. Global Incidence and Mortality for Prostate Cancer: Analysis of Temporal Patterns and Trends in 36 Countries. Eur Urol (2016) 70:862–74. doi: 10.1016/j.eururo.2016.05.043

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Horwich A, Hugosson J, de Reijke T, Wiegel T, Fizazi K, Kataja V. Prostate cancer: ESMO Consensus Conference Guidelines 2012. Ann Oncol (2013) 24:1141–62. doi: 10.1093/annonc/mds624

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Watson PA, Arora VK, Sawyers CL. Emerging mechanisms of resistance to androgen receptor inhibitors in prostate cancer. Nat Rev Cancer (2015) 15:701–11. doi: 10.1038/nrc4016

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Kwon ED, Drake CG, Scher HI, Fizazi K, Bossi A, van den Eertwegh AJ, et al. Ipilimumab versus placebo after radiotherapy in patients with metastatic castration-resistant prostate cancer that had progressed after docetaxel chemotherapy (CA184-043): a multicentre, randomised, double-blind, phase 3 trial. Lancet Oncol (2014) 15:700–12. doi: 10.1016/S1470-2045(14)70189-5

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Beer TM, Kwon ED, Drake CG, Fizazi K, Logothetis C, Gravis G, et al. Randomized, Double-Blind, Phase III Trial of Ipilimumab Versus Placebo in Asymptomatic or Minimally Symptomatic Patients With Metastatic Chemotherapy-Naive Castration-Resistant Prostate Cancer. J Clin Oncol (2017) 35:40–7. doi: 10.1200/JCO.2016.69.1584.

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Kretschmer A, Tilki D. Biomarkers in prostate cancer - Current clinical utility and future perspectives. Crit Rev oncology/hematology (2017) 120:180–93. doi: 10.1016/j.critrevonc.2017.11.007

CrossRef Full Text | Google Scholar

7. Filella X, Foj L, Auge JM, Molina R, Alcover J. Clinical utility of %p2PSA and prostate health index in the detection of prostate cancer. Clin Chem Lab Med (2014) 52:1347–55. doi: 10.1515/cclm-2014-0027

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Le BV, Griffin CR, Loeb S, Carvalhal GF, Kan D, Baumann NA, et al. [-2]Proenzyme prostate specific antigen is more accurate than total and free prostate specific antigen in differentiating prostate cancer from benign disease in a prospective prostate cancer screening study. J Urol (2010) 183:1355–9. doi: 10.1016/j.juro.2009.12.056

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Klein EA, Cooperberg MR, Magi-Galluzzi C, Simko JP, Falzarano SM, Maddala T, et al. A 17-gene assay to predict prostate cancer aggressiveness in the context of Gleason grade heterogeneity, tumor multifocality, and biopsy undersampling. Eur Urol (2014) 66:550–60. doi: 10.1016/j.eururo.2014.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Cullen J, Rosner IL, Brand TC, Zhang N, Tsiatis AC, Moncur J, et al. A Biopsy-based 17-gene Genomic Prostate Score Predicts Recurrence After Radical Prostatectomy and Adverse Surgical Pathology in a Racially Diverse Population of Men with Clinically Low- and Intermediate-risk Prostate Cancer. Eur Urol (2015) 68:123–31. doi: 10.1016/j.eururo.2014.11.030

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Sanda MG, Cadeddu JA, Kirkby E, Chen RC, Crispino T, Fontanarosa J, et al. Clinically Localized Prostate Cancer: AUA/ASTRO/SUO Guideline. Part I: Risk Stratification, Shared Decision Making, and Care Options. J Urol (2018) 199:683–90. doi: 10.1016/j.juro.2017.11.095

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Sfanos KS, Yegnasubramanian S, Nelson WG, De Marzo AM. The inflammatory microenvironment and microbiome in prostate cancer development. Nat Rev Urol (2018) 15:11–24. doi: 10.1038/nrurol.2017.167

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Strasner A, Karin M. Immune Infiltration and Prostate Cancer. Front Oncol (2015) 5:128. doi: 10.3389/fonc.2015.00128

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Flammiger A, Bayer F, Cirugeda-Kuhnert A, Huland H, Tennstedt P, Simon R, et al. Intratumoral T but not B lymphocytes are related to clinical outcome in prostate cancer. APMIS Acta pathologica microbiologica immunologica Scandinavica (2012) 120:901–8. doi: 10.1111/j.1600-0463.2012.02924.x

CrossRef Full Text | Google Scholar

15. Ebelt K, Babaryka G, Frankenberger B, Stief CG, Eisenmenger W, Kirchner T, et al. Prostate cancer lesions are surrounded by FOXP3+, PD-1+ and B7-H1+ lymphocyte clusters. Eur J Cancer (Oxford Engl 1990) (2009) 45:1664–72. doi: 10.1016/j.ejca.2009.02.015

CrossRef Full Text | Google Scholar

16. Savage PA, Vosseller K, Kang C, Larimore K, Riedel E, Wojnoonski K, et al. Recognition of a ubiquitous self antigen by prostate cancer-infiltrating CD8+ T lymphocytes. Sci (New York N.Y.) (2008) 319:215–20. doi: 10.1126/science.1148886

CrossRef Full Text | Google Scholar

17. Lanciotti M, Masieri L, Raspollini MR, Minervini A, Mari A, Comito G, et al. The role of M1 and M2 macrophages in prostate cancer in relation to extracapsular tumor extension and biochemical recurrence after radical prostatectomy. BioMed Res Int (2014) 2014:486798. doi: 10.1155/2014/486798

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Cancer Genome Atlas Research Network. Comprehensive genomic characterization defines human glioblastoma genes and core pathways. Nature (2008) 455:1061–8. doi: 10.1038/nature07385

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Yoshihara K, Shahmoradgoli M, Martinez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun (2013) 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Shah N, Wang P, Wongvipat J, Karthaus WR, Abida W, Armenia J, et al. Regulation of the glucocorticoid receptor via a BET-dependent enhancer drives antiandrogen resistance in prostate cancer. eLife (2017) 6:e27861. doi: 10.7554/eLife.27861

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res (2015) 43:e47. doi: 10.1093/nar/gkv007

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Wang L, Wei B, Hu G, Wang L, Jin Y, Sun Z. Gene expression analyses to explore the biomarkers and therapeutic targets for gliomas. Neurol Sci (2015) 36:403–9. doi: 10.1007/s10072-014-1985-0

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Chen H, Boutros PC. VennDiagram: a package for the generation of highly-customizable Venn and Euler diagrams in R. BMC Bioinf (2011) 12:35. doi: 10.1186/1471-2105-12-35

CrossRef Full Text | Google Scholar

24. Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics J Integr Biol (2012) 16:284–7. doi: 10.1089/omi.2011.0118

CrossRef Full Text | Google Scholar

25. Callahan BJ, Sankaran K, Fukuyama JA, McMurdie PJ, Holmes SP. Bioconductor Workflow for Microbiome Data Analysis: from raw reads to community analyses. F1000Research (2016) 5:1492. doi: 10.12688/f1000research.8986.2

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Meller S, Meyer HA, Bethan B, Dietrich D, Maldonado SG, Lein M, et al. Integration of tissue metabolomics, transcriptomics and immunohistochemistry reveals ERG- and gleason score-specific metabolomic alterations in prostate cancer. Oncotarget (2016) 7:1421–38. doi: 10.18632/oncotarget.6370

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Sun Y, Goodison S. Optimizing molecular signatures for predicting prostate cancer recurrence. Prostate (2009) 69:1119–27. doi: 10.1002/pros.20961

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Olmos D, Brewer D, Clark J, Danila DC, Parker C, Attard G, et al. Prognostic value of blood mRNA expression signatures in castration-resistant prostate cancer: a prospective, two-stage study. Lancet Oncol (2012) 13:1114–24. doi: 10.1016/S1470-2045(12)70372-8

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Torre LA, Bray F, Siegel RL, Ferlay J, Lortet-Tieulent J, Jemal A. Global cancer statistics, 2012. CA: Cancer J clinicians (2015) 65:87–108. doi: 10.3322/caac.21262

CrossRef Full Text | Google Scholar

30. Jian L, Xie LP, Lee AH, Binns CW. Protective effect of green tea against prostate cancer: a case-control study in southeast China. Int J Cancer (2004) 108:130–5. doi: 10.1002/ijc.11550

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Kamangar F, Dores GM, Anderson WF. Patterns of cancer incidence, mortality, and prevalence across five continents: defining priorities to reduce cancer disparities in different geographic regions of the world. J Clin Oncol (2006) 24:2137–50. doi: 10.1200/JCO.2005.05.2308

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Patel VG, Oh WK. The evolving landscape of immunotherapy in advanced prostate cancer. Immunotherapy (2019) 11:903–12. doi: 10.2217/imt-2019-0019

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Redman JM, Gulley JL, Madan RA. Combining immunotherapies for the treatment of prostate cancer. Urologic Oncol (2017) 35:694–700. doi: 10.1016/j.urolonc.2017.09.024

CrossRef Full Text | Google Scholar

34. Jia D, Li S, Li D, Xue H, Yang D, Liu Y. Mining TCGA database for genes of prognostic value in glioblastoma microenvironment. Aging (2018) 10:592–605. doi: 10.18632/aging.101415

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Chen XY, Yeong J, Thike AA, Bay BH, Tan PH. Prognostic role of immune infiltrates in breast ductal carcinoma in situ. Breast Cancer Res Treat (2019) 177:17–27. doi: 10.1007/s10549-019-05272-2

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Ladoire S, Arnould L, Apetoh L, Coudert B, Martin F, Chauffert B, et al. Pathologic complete response to neoadjuvant chemotherapy of breast carcinoma is associated with the disappearance of tumor-infiltrating foxp3+ regulatory T cells. Clin Cancer Res (2008) 14:2413–20. doi: 10.1158/1078-0432.CCR-07-4491

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Ma Y, Pitt JM, Li Q, Yang H. The renaissance of anti-neoplastic immunity from tumor cell demise. Immunological Rev (2017) 280:194–206. doi: 10.1111/imr.12586

CrossRef Full Text | Google Scholar

38. Fritzsching B, Fellenberg J, Moskovszky L, Sapi Z, Krenacs T, Machado I, et al. CD8(+)/FOXP3(+)-ratio in osteosarcoma microenvironment separates survivors from non-survivors: a multicenter validated retrospective study. Oncoimmunology (2015) 4:e990800. doi: 10.4161/2162402X.2014.990800

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Fox JJ, Navarro HI, Hashimoto T, Garcia AJ, Goldstein AS. Mass cytometry reveals species-specific differences and a new level of complexity for immune cells in the prostate. Am J Clin Exp Urol (2019) 7:281–96.

PubMed Abstract | Google Scholar

40. Pozzobon T, Goldoni G, Viola A, Molon B. CXCR4 signaling in health and disease. Immunol Lett (2016) 177:6–15. doi: 10.1016/j.imlet.2016.06.006

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Bleul CC, Farzan M, Choe H, Parolin C, Clark-Lewis I, Sodroski J, et al. The lymphocyte chemoattractant SDF-1 is a ligand for LESTR/fusin and blocks HIV-1 entry. Nature (1996) 382:829–33. doi: 10.1038/382829a0

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Hattermann K, Mentlein R. An infernal trio: the chemokine CXCL12 and its receptors CXCR4 and CXCR7 in tumor biology. Ann Anat = Anatomischer Anzeiger Off Organ Anatomische Gesellschaft (2013) 195:103–10. doi: 10.1016/j.aanat.2012.10.013

CrossRef Full Text | Google Scholar

43. Hattermann K, Held-Feindt J, Lucius R, Muerkoster SS, Penfold ME, Schall TJ, et al. The chemokine receptor CXCR7 is highly expressed in human glioma cells and mediates antiapoptotic effects. Cancer Res (2010) 70:3299–308. doi: 10.1158/0008-5472.CAN-09-3642

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Domanska UM, Kruizinga RC, Nagengast WB, Timmer-Bosscha H, Huls G, de Vries EG, et al. A review on CXCR4/CXCL12 axis in oncology: no place to hide. Eur J Cancer (Oxford Engl 1990) (2013) 49:219–30. doi: 10.1016/j.ejca.2012.05.005

CrossRef Full Text | Google Scholar

45. Domanska UM, Timmer-Bosscha H, Nagengast WB, Oude Munnink TH, Kruizinga RC, Ananias HJ, et al. CXCR4 inhibition with AMD3100 sensitizes prostate cancer to docetaxel chemotherapy. Neoplasia (New York N.Y.) (2012) 14:709–18. doi: 10.1593/neo.12324

CrossRef Full Text | Google Scholar

46. Viola A, Luster AD. Chemokines and their receptors: drug targets in immunity and inflammation. Annu Rev Pharmacol Toxicol (2008) 48:171–97. doi: 10.1146/annurev.pharmtox.48.121806.154841

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Birkenbach M, Josefsen K, Yalamanchili R, Lenoir G, Kieff E. Epstein-Barr virus-induced genes: first lymphocyte-specific G protein-coupled peptide receptors. J Virol (1993) 67:2209–20. doi: 10.1128/JVI.67.4.2209-2220.1993

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Pereira JP, Kelly LM, Xu Y, Cyster JG. EBI2 mediates B cell segregation between the outer and centre follicle. Nature (2009) 460:1122–6. doi: 10.1038/nature08226

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Gatto D, Paus D, Basten A, Mackay CR, Brink R. Guidance of B cells by the orphan G protein-coupled receptor EBI2 shapes humoral immune responses. Immunity (2009) 31:259–69. doi: 10.1016/j.immuni.2009.06.016

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Li J, Lu E, Yi T, Cyster JG. EBI2 augments Tfh cell fate by promoting interaction with IL-2-quenching dendritic cells. Nature (2016) 533:110–4. doi: 10.1038/nature17947

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Suan D, Nguyen A, Moran I, Bourne K, Hermes JR, Arshi M, et al. T follicular helper cells have distinct modes of migration and molecular signatures in naive and memory immune responses. Immunity (2015) 42:704–18. doi: 10.1016/j.immuni.2015.03.002

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Gatto D, Wood K, Caminschi I, Murphy-Durland D, Schofield P, Christ D, et al. The chemotactic receptor EBI2 regulates the homeostasis, localization and immunological function of splenic dendritic cells. Nat Immunol (2013) 14:446–53. doi: 10.1038/ni.2555

PubMed Abstract | CrossRef Full Text | Google Scholar