Tumor immune contexture predicts recurrence after prostatectomy and efficacy of androgen deprivation and immunotherapy in prostate cancer
Journal of Translational Medicine volume 21, Article number: 194 (2023)
Prostate cancer is one of the most common cancers in men with notable interpatient heterogeneity. Implications of the immune microenvironment in predicting the biochemical recurrence-free survival (BCRFS) after radical prostatectomy and the efficacy of systemic therapies in prostate cancer remain ambiguous.
The tumor immune contexture score (TICS) involving eight immune contexture-related signatures was developed using seven cohorts of 1120 patients treated with radical prostatectomy (training: GSE46602, GSE54460, GSE70769, and GSE94767; validation: GSE70768, DKFZ2018, and TCGA). The association between the TICS and treatment efficacy was investigated in GSE111177 (androgen deprivation therapy [ADT]) and EGAS00001004050 (ipilimumab).
A high TICS was associated with prolonged BCRFS after radical prostatectomy in the training (HR = 0.32, 95% CI 0.24–0.45, P < 0.001) and the validation cohorts (HR = 0.45, 95% CI 0.32–0.62, P < 0.001). The TICS showed stable prognostic power independent of tumor stage, surgical margin, pre-treatment prostatic specific antigen (PSA), and Gleason score (multivariable HR = 0.50, 95% CI 0.39–0.63, P < 0.001). Adding the TICS into the prognostic model constructed using clinicopathological features significantly improved its 1/2/3/4/5-year area under curve (P < 0.05). A low TICS was associated with high homologous recombination deficiency scores, abnormally activated pathways concerning DNA replication, cell cycle, steroid hormone biosynthesis, and drug metabolism, and fewer tumor-infiltrating immune cells (P < 0.05). The patients with a high TICS had favorable BCRFS with ADT (HR = 0.25, 95% CI 0.06–0.99, P = 0.034) or ipilimumab monotherapy (HR = 0.23, 95% CI 0.06–0.81, P = 0.012).
Our study delineates the associations of tumor immune contexture with molecular features, recurrence after radical prostatectomy, and the efficacy of ADT and immunotherapy. The TICS may improve the existing risk stratification systems and serve as a patient-selection tool for ADT and immunotherapy in prostate cancer.
Prostate cancer is a malignant tumor with high incidence and mortality worldwide . It is characterized by a long and varied course of the disease and high intra- and inter-tumor heterogeneity . Prostate cancer is an androgen-dependent disease, and most prostate cancer patients initially respond to androgen deprivation therapy (ADT) . However, most androgen-sensitive prostate cancer patients eventually become resistant to ADT and develop castration-resistant prostate cancer (CRPC) .
The tumor immune microenvironment (TIME) is mainly composed of tumor cells, immune cells, various signaling molecules, and extracellular matrix . TIME was closely associated with survival outcome and treatment response . Several pieces of evidence have indicated that TIME plays a key role in the development of prostate cancer and the response to anti-tumor therapy [7, 8]. T-cell infiltration was related to tumor progression and cancer-specific survival in both localized and metastatic prostate cancer patients [9, 10]. For instance, T helper (Th) cells activated antigen-specific effector cells and recruited immune cells such as macrophages and mast cells, which could play crucial roles in the adaptive immune response to tumor cells [11, 12]. Th17 cells are a subgroup of CD4+ helper T cells, which can secrete interleukin (IL)-17 and IL-22, and play a critical role in the occurrence and development of prostate cancer . The preponderance of Th17-mediated inflammation has been negatively correlated with Gleason score in prostate cancer patients , and blocking the IL-17 pathway inhibited the formation of microinvasive prostate cancer in animal models . Given this, we speculated that TIME might be used to predict the biochemical recurrence-free survival (BCRFS) of prostate cancer.
As for immunotherapy, the anti-tumor effects of pembrolizumab and ipilimumab against prostate cancer have been demonstrated in clinical trials [16,17,18]. Immune checkpoint inhibitors (ICIs) such as the antibodies targeting cytotoxic T lymphocyte antigen-4 (CTLA-4) and programmed cell death 1 (PD-1)/programmed cell death-ligand 1 (PD-L1) could regulate the TIME and influence immune cells to promote anti-tumor immune responses [17, 19]. TIME might presumably also be a predictor for the ICI-induced antitumor activity in prostate cancer.
In this study, the multi-dimensional data of 1120 prostate cancer patients from The Cancer Genome Atlas-prostate adenocarcinoma (TCGA-PRAD) cohort and 6 Gene Expression Omnibus (GEO) cohorts were analyzed. In the training cohorts, we first evaluated the prognostic effects of 92 immune contexture-related signatures and constructed a tumor immune contexture score (TICS) based on the 8 signatures with prognostic effects. Then, the TICS was tested successfully in the validation cohorts. The TICS showed consistent prognostic effects across the prostate cancer patients with different clinicopathological characteristics (e.g., TNM stage, pre-treatment prostatic specific antigen [PSA], and Gleason score) and this score was an independent prognostic factor of BCRFS. Compared with the prognostic model constructed by the clinicopathological features without the TICS, the model with the TICS had a significantly higher power of prognostication. In addition, the correlations between the TICS with cancer-related signaling pathways and immunotherapy efficacy were also explored, which might provide promising insight into the individualized treatment of prostate cancer patients.
Prostate cancer datasets and pre-processing
Transcriptome profiles and clinicopathological characteristics of seven datasets of localized prostate cancer were downloaded from TCGA and GEO databases, including GSE46602, GSE54460, GSE70768, GSE70769, GSE94767, DKFZ2018, and TCGA-PRAD [20,21,22,23,24]. The Affymetrix microarray data from GSE46602 and GSE94767 were normalized using the robust multiarray analysis (RMA) algorithm. The Illumina microarray data from GSE70768 and GSE70769 were normalized using the quantile algorithm and log2 transformed. The RNA-seq data of GSE54460, DKFZ2018, TCGA-PRAD, GSE111177, and EGAS00001004050 produced by the Illumina HiSeq or TruSeq platform were normalized as Fragments Per Kilobase of exon model per Million mapped fragments (FPKM) or Reads Per Kilobase per Million mapped reads (RPKM) format and were transformed to log2 scale (for details, see Additional file 1: Table S1). The patients who had not received radical prostatectomy (RP) or had no BCRFS data were excluded. The seven PRAD datasets were split into the training set and the validation set based on techniques and platforms, ensuring that the training and validation sets each have datasets using microarray and RNA-seq (Additional file 1: Table S1).
To assess the utility of the TICS in predicting outcomes of patients receiving ADT, a cohort consisting of the pre-ADT samples from 20 prostate cancer patients was obtained from the GEO database (GSE111177) . For immunotherapy, the data of 18 metastatic castration-resistant prostate cancer (mCRPC) patients treated with ipilimumab were obtained from the work of Subudhi et al. (EGAS00001004050 at the European Genome-phenome Archive) .
In total, 1158 prostate cancer patients were included in this study. The detailed baseline characteristics are demonstrated in Additional file 1: Table S1. In microarray datasets, probes were replaced with matched gene symbols in later analysis (the first one was selected if multiple genes were matched), while the probes that did not match any gene were discarded. For all datasets, multiple rows with the same gene symbol were ordered by the median expression value, and the max one was selected for further analysis.
Prognosis-related markers selection and the development of the TICS algorithm
In total, 92 immune-related signatures were obtained from a previous work , and the single sample gene set enrichment analysis (ssGSEA) algorithm was applied to calculate the enrichment score (ES) of these immune-related signatures in the above cohorts. Based on the expression of all genes, the normalized Z-score value of mRNA expression profiles was used to calculate the score. Each signature was transformed into binary variables with the cutoff of median value in four training datasets (GSE46602, GSE54460, GSE70769, and GSE94767). The univariable Cox regression analysis was applied to identify the prognostic efficacy of ES stratification in each training set, then the meta-analysis was performed to estimate the hazard ratios (HRs) of four training sets for each signature. The signatures significantly associated with prognosis (P < 0.05) in the meta-analysis were selected to develop the TICS. The TICS was calculated by the following formula:
where HRi was the hazard ratio of ith signature from meta-analysis, and ssGSEA score represented the enrichment score of the corresponding signature.
To validate the prognostic performance of the TICS, three datasets (GSE70768, DKFZ2018, and TCGA-PRAD) were used as independent validation cohorts. The robustness and independence of the TICS as a predictor of BCRFS were evaluated by subgroup analysis and univariable/multivariable Cox regression. The added value of TICS to traditional prognostic factors such as pre-treatment PSA and Gleason score was evaluated by comparing the prognostic utility of two models constructed based on multivariable Cox regression. The first model consisted of pT, pN, surgical margin, PSA, Gleason score, and the TICS (continuous variable), and the second model involved all except the TICS.
Identification of biological characteristics related to the TICS
We compared genetic features between the high- and low-risk groups, including microsatellite instability (MSI) score, tumor mutation burden (TMB), tumor purity, and genomic scar signatures . Meanwhile, the Kyoto Encyclopedia of Genes and Genomes (KEGG) gene set (c2.cp.kegg.v7.5.1.symbols.gmt) was downloaded from MSigDB to investigate the potential difference in the biological function between the high- and low-TICS subgroups using the gene set enrichment analysis (GSEA) method. Additionally, the ssGSEA analysis based on published mitogenic-related gene sets was used to quantify the relative abundance in the TCGA-PRAD dataset .
All data analyses and graph plotting were conducted in R software (version 4.1.2), SPSS (version 20.0.0), and GraphPad Prism (version 8.0.2). The meta-analysis was performed by the “meta” R package. The Kaplan–Meier method was used to access the survival of different subgroups, and the log-rank (Mantel–Cox) test was used to measure the statistical significance. The univariable and multivariable Cox proportional hazards regression analyses were performed using the “ezcox” R package. Time-dependent receiver operating characteristic (tROC) curves, concordance index (C-index), and calibration curves were calculated by the “timeROC” package to evaluate the predictive power of the TICS. Mann–Whitney U test and Chi-square test were used to assess the statistical difference between the two groups. Spearman’s correlation analysis was used for correlational analyses. The “ComplexHeatmap” package was used to visualize the mutation landscape. Unless otherwise stated, P values of less than 0.05 were considered statistically significant.
Construction of the scoring algorithm based on recurrence-associated immune signatures
In total, 92 immune contexture-related signatures compiled in a previous study were investigated in our study and the genes in each signature are demonstrated in Additional file 1: Table S2. In the training cohorts (clinical information demonstrated in Additional file 1: Table S1), the association between BCRFS and the ssGSEA enrichment score of each signature as dichotomous variables according to the median value was analyzed. Then the mixed effects in the four training cohorts were calculated by meta-analysis. The prognostic effects of the 92 immune signatures in the training sets are demonstrated in Additional file 2: Fig. S1. Among these immune signatures, eight signatures showed a consistent association with BCRFS (Fig. 1A and Additional file 1: Table S3). In detail, the signatures related to type 17T helper cells, immature dendritic cells, cell death (Rotterdam_ERneg_PCA_15721472), and immunoglobulin (IGG_Cluster_21214954) heralded a long BCRFS (P < 0.05) [30, 31], while the signatures concerning CD56dim natural killer (NK) cell, wound healing (CHANG_CORE_SERUM_RESPONSE_UP), interferon signaling (Module3_IFN_score), and proliferation (Mudule11_Prolif_score) were associated with a poor BCRFS (P < 0.05) [32, 33]. The ssGSEA enrichment score for the eight signatures are demonstrated in Fig. 1B, showing comparable scores between the training cohorts. We also checked the overlap of the genes involved in these eight signatures and only 11 genes were repeated in the signatures of CHANG_CORE_SERUM_RESPONSE_UP (212 genes) and Mudule11_Prolif_score (120 genes). In addition, low multicollinearity was observed among the scores of these eight signatures (variance inflation factor [VIF] < 10 and kappa < 100, Additional file 1: Table S4), which infers the feasibility of constructing a prognostic model based on these eight signatures.
According to the HR values of the mixed effects, we calculated the coefficients of the eight signatures (shown on the left in Fig. 1B) and calculated the TICS by summing the products of each signature score and its coefficient. According to the TICS, samples were divided into two groups with the median value in each training cohort, and the patients in the high-TICS group had a more favorable BCRFS than the corresponding low-TICS group (GSE46602: P = 0.040, HR = 0.41, 95% CI 0.17–0.99; GSE54460: P < 0.001, HR = 0.31, 95% CI 0.17–0.55; GSE70769: P < 0.001, HR = 0.31, 95% CI 0.16–0.58; GSE94767: P < 0.001, HR = 0.33, 95% CI 0.18–0.59; Fig. 1C–F).
To further evaluate the prognostic power of the TICS, three additional independent cohorts (GSE70768, DKFZ2018, and TCGA-PRAD) were utilized for validation (clinical information demonstrated in Additional file 1: Table S1). The ssGSEA enrichment score for the eight signatures are also demonstrated in Fig. 1G. Among these cohorts, the two groups divided by the median value of the TICS also showed significant differences in BCRFS (GSE70768: P = 0.006, HR = 0.26, 95% CI 0.09–0.72; GSE54460: P = 0.019, HR = 0.28, 95% CI 0.09–0.87; TCGA: P < 0.001, HR = 0.51, 95% CI 0.35–0.73; Fig. 1H–J). The AUCs were 0.64–0.80 and the c-indices were 0.63–0.78 for all the included cohorts (Fig. 1K and Additional file 1: Table S5), revealing that the TICS was likely to be a robust prognostic marker for the BCRFS after RP. Pooled analysis revealed consistent results in the training cohorts (HR = 0.32, 95% CI 0.24–0.45, P < 0.001) and in the validation cohorts (HR = 0.45, 95% CI 0.32–0.62, P < 0.001, Fig. 1L). Altogether these results suggest the robust utility of the TICS in predicting biochemical recurrence in prostate cancers.
Furthermore, we compared the TICS with several previously published multi-gene signatures with exposed formulas [34,35,36,37]. In the three datasets in the validation set, the TICS showed a numerically higher C-index for predicting BCRFS than other signatures (Additional file 1: Table S6 and Additional file 2: Fig. S2), elucidating the superiority of the TICS over previously published signatures.
The added value of the TICS for prognostication in prostate cancer
We then assessed the associations between the TICS and clinicopathological features. The patients of all cohorts were grouped according to clinicopathological features, and the TICS were compared among different subgroups. A lower TICS was observed in elderly patients (≥ 60, P = 0.017), patients with advanced disease (pT3/4, P < 0.001; pN1, P < 0.001), patients with positive surgical margin (R1/2, P < 0.001), a higher pre-treatment prostate-specific antigen (PSA, P = 0.004), and a poorer Gleason score (P < 0.001, Fig. 2A).
Given the association between the TICS and prognostic indicators such as pathological stage, pre-treatment PSA, and Gleason score, we further performed a subgroup analysis of the association between the TICS and BCRFS. In the subgroups stratified by age, stage, surgical margin, PSA, and Gleason score, similar trends of the association between a high TICS and long BCRFS were observed (Fig. 2B), indicating the robustness of its prognostic utility. Furthermore, multivariable Cox regression analysis identified the TICS as an independent predictor of BCRFS, whether the TICS was treated as a categorical variable (grouped by the median value, multivariable HR = 0.50, 95% CI 0.39–0.63, P < 0.001) or as a continuous variable (multivariable HR = 0.61, 95% CI 0.53–0.70, P < 0.001, Fig. 2C), suggesting that the TICS was a prognostic factor independent of clinicopathological features.
We constructed two models by multivariable regression with or without the TICS to estimate the added value of the TICS. The first model consisting of pT, pN, surgical margin, PSA, Gleason score, and the TICS (continuous variable) achieved a bootstrap Dxy of 0.443. By contrast to this model, the second model involving all the above variables except the TICS yielded a relatively lower bootstrap Dxy of 0.387. The calibration curves of these two models are shown in Fig. 2D, in which the model with the TICS showed better agreement between the predicted and observed 1-, 3-, and 5-year BCRFS. Moreover, the addition of the TICS significantly improved the 1-year (0.790 vs. 0.757, P = 0.004), 2-year (0.759 vs. 0.716, P < 0.001), 3-year (0.745 vs. 0.709, P < 0.001), 4-year (0.753 vs. 0.729, P = 0.015), and 5-year AUC (0.752 vs. 0.730, P = 0.048, Fig. 2E). Adding the TICS into the model constructed by traditional prognostic markers markedly improved its prognostic utility, indicating the added value of the TICS for the prognostication based on the traditional indicators.
Genetic profiles of the prostate cancer in different TICS groups
Based on the TCGA dataset, we analyzed the somatic alterations of the cancerous samples in the high- and low-TICS groups (Fig. 3A). The frequency of the somatic mutations of homologous recombination (HR) genes and mismatch repair (MMR) genes trended higher in the low-TICS group (e.g., BRCA2 [P = 0.013], BRIP1 [P = 0.007], and MSH6 [P = 0.062], Fig. 3A). In addition, somatic alterations in other key genes were also investigated (Fig. 3A). TP53, FOXA1, and PTEN were reported to be involved in the regulation of cell survival and proliferation, and their somatic mutations were associated with cancer progression in prostate cancer [38,39,40]. In our study, mutation events of these three genes were significantly increased in the low-TICS group (TP53, P = 0.002; FOXA1, P = 0.016; PTEN, P < 0.001). Other genes associated with tumorigenesis of prostate cancer, such as ATM, MKI67, and SPOP [39, 41,42,43,44], showed non-significant increases in somatic mutation frequency in the low-TICS group.
As for the indices concerning genomic alterations, the low-TICS group had a high aneuploidy score (P < 0.001, Fig. 3B) and non-synonymous TMB (P < 0.001, Fig. 3C), indicating genomic instability. Based on the approval of poly-ADP ribose polymerase inhibitor (PARPi) in prostate cancer, we further analyzed the scores related to DNA damage repair. Consistent with the numerically more mutations of MMR and HR genes in the low-TICS group, the samples with a low TICS trended to have a high MSI score (P = 0.058, Fig. 3D) and high homologous recombination deficiency (HRD)-associated scores including telomeric allelic imbalance (NtAI, P < 0.001, Fig. 3E), large-scale state transitions (LST, P < 0.001, Fig. 3F), loss of heterozygosity (LOH, P < 0.001, Fig. 3G), and the total HRD score (P < 0.001, Fig. 3H).
The TICS represented alterations of cancer-related pathways
Based on the gene sets embodied in the KEGG database, we analyzed the changes in pathway activity between the high- and low-TICS groups by GSEA (Fig. 4A and Additional file 1: Table S7). Tumors in the high-TICS group had higher activity of the Hedgehog signaling, second messenger signaling, and cell–matrix adhesion-related pathways, while the high expression of the IL-17 family may indicate an intense inflammatory response (Fig. 4B). On the other side, the abnormal expression of genes involved in DNA replication and cell cycle in the low-TICS group suggested the endogenous aggressiveness of tumors, and the activation of steroid hormone biosynthesis and drug metabolism may lead to resistance to ADT and other treatments (Fig. 4B).
Further, we analyzed the alterations of the previously reported oncogenic signaling pathways by Sanchez-Vega et al. . The TICS was negatively correlated with the scores related to cell cycle (Rho = − 0.37, P < 0.001), p53 (Rho = − 0.37, P < 0.001), PI3K (Rho = − 0.25, P < 0.001), and Myc (Rho = − 0.12, P = 0.0058) pathways while positively correlated with the scores concerning Wnt (Rho = 0.21, P < 0.001), RTK/RAS (Rho = 0.15, P < 0.001), and Hippo (Rho = 0.15, P < 0.001) pathways (Fig. 4C, D). These results suggest that the TICS may represent the biological heterogeneity among different prostate cancers.
Associations between TICS and the efficacy of ADT and immunotherapy
Considering the association between the TICS and steroid hormone biosynthesis pathway, we hypothesized that the TICS might predict the response to ADT. To evaluate the relationship between the TICS and ADT efficacy, an independent cohort (GSE111177) of 20 patients with clinical stage Ic to IIIa prostate cancer was analyzed (Fig. 5A). Neither pre-treatment PSA (< 10 vs. ≥ 10, or < 20 vs. ≥ 20) nor Gleason score (≤ 6 vs. > 6) was able to distinguish BCRFS effectively, while the TICS-based dichotomy based on median value yielded two groups with significantly different BCRFS (P = 0.034, HR = 0.25, 95% CI 0.06–0.99; Fig. 5A). Taken together, the TICS rather than pre-treatment PSA and Gleason score might be an effective marker to identify patients susceptible to castration resistance to ADT.
As a score reflecting immune contexture, the TICS should be theoretically associated with immunotherapy efficacy, which was studied in a cohort consisting of 18 CRPC patients treated with ipilimumab monotherapy . First, based on the ROC curve illustrating the association between the TICS and disease control (Fig. 5B), an optimal cut-off was determined when the Youden’s index reached the maximum. According to this cut-off, the cohort was stratified into the high-TICS (n = 13) and the low-TICS groups (n = 5). Patients in the high-TICS group had significantly favorable BCRFS (P = 0.012, HR = 0.23, 95% CI 0.06–0.81), radiographic progression-free survival (rPFS, P = 0.003, HR = 0.17, 95% CI 0.04–0.65), and overall survival (OS, P = 0.018, HR = 0.23, 95% CI 0.06–0.87) compared to the low-TICS patients (Fig. 5C). We also analyzed the immune cell composition in the tumor tissues between both groups and found that the tumors in the high-TICS group were infiltrated with more T cells (P = 0.003, particularly cytotoxic T cells [P = 0.003] and memory T cells [P = 0.007]), B cells (P = 0.014), NK cells (P = 0.003) and dendritic cells (P = 0.046, Fig. 5D, E). These results indicate that the tumors with a low TICS may undergo the “immune desert” state, which is difficult to benefit from ICIs.
The recurrence of prostate cancer after RP is a serious threat, and the precise prediction of the risk of recurrence may facilitate appropriate personalized management. In the present study, we analyzed the expression profile of resected tumor samples and constructed a novel score to predict BCRFS based on the eight expression signatures related to immune contexture. The TICS was effective across the four training sets and the three validation sets and showed consistent prognostic effects in the subgroups stratified by key clinicopathological characteristics (e.g., TNM stage, pre-treatment PSA, and Gleason score). The TICS was further identified to be an independent predictor of BCRFS. Based on the prognostic model constructed by key clinicopathological features, adding the TICS significantly improved the power of prognostication. Moreover, the TICS may represent intrinsic discrepancies of prostate cancers at genomic and transcriptomic levels and possess the ability to predict the efficacy of ADT and ICI treatment.
In previous studies, the biological functions related to some of the signatures involved in the TICS calculation have been suggested to be associated with the prognosis. Inflammation mediated by Th17 was correlated with a low Gleason score , while the infiltration of dendritic cells conferred improved distant metastasis-free survival of patients [45, 46]. Prostate cancer at the onset of hormone resistance revealed an increased proliferative rate and a decreased apoptotic rate compared with hormone-responsive cancer , and the expression of genes involved in cell cycle progression was useful to predict biochemical recurrence . Nevertheless, there are also some signatures involved in the TICS that had inconsistent results with previous studies. Low NK cell activity was associated with positive surgical margins in 51 postoperative patients , while highly effective NK cells revealed a good prognosis in a study including 18 patients with metastatic cancer , but the high expression of CD56dim NK cell signature led to a decline of the TICS. Interferon-γ may suppress the invasive capacity of prostate cancer cells in vitro , and a low interferon-β level was an indicator to distinguish high-risk patients with localized cancer , but the Module3_IFN_score signature was negatively correlated with the TICS. These abnormalities may be related to the spatial and temporal distribution of cells with corresponding expression signatures in prostate cancer tissue, which will be elucidated in further studies.
Tumor stage, PSA level, and Gleason score have been used for prostate cancer risk stratification and are widely recommended in guidelines [53, 54]. However, each feature has its disadvantages, such as poor specificity , interobserver variation, and sampling problems , resulting in different manifestations among patients in the same risk group. The TICS in the present study predicted different outcomes in patients with similar TNM stage, PSA level, or Gleason score, and demonstrated its independence in multivariable analysis, indicating its ability as a novel prognostic biomarker. The prognostic model with the TICS performed better than the model constructed by clinicopathological features alone, suggesting the added value of the TICS to the existing risk stratification system. In addition, the TICS showed better performance compared with previously published models [34,35,36,37], indicating its added value.
Although existing studies have identified several factors associated with ADT response (such as the expression of PTEN or specific miRNAs) [57, 58], a clinical prognostic system to predict the outcomes of ADT is absent to date. ADT may cause remodeling of the immune microenvironment in prostate cancer and affect the patient's relapse-free survival , which inspired us to try to predict patients’ response to ADT with the TICS. The TICS demonstrated better prognostic efficacy than the PSA level or the Gleason score in a limited cohort, and this encouraging result needs to be validated on a larger scale. In previous clinical trials, few of the unselected patients benefited from treatment with ICIs, while several generally acknowledged immunotherapy biomarkers (i.e., high TMB, alterations in DNA damage repair [DDR] genes, or high-expression PD-L1) yielded controversial results in the different trials [26, 60,61,62]. Patients with low TICS revealed poor outcomes in the investigated cohort with CTLA-4 antibody treatment, while TMB or DDR gene mutations failed to make effective discrimination, suggesting that TICS may have a unique value in excluding patients who may not benefit from immunotherapy.
In our study, the transcriptomic data were from different platforms, which may lead to potential adverse effects on data analysis and the generalizability of results. Despite the lack of batch correction in our study, we did the following to minimize this disadvantage in developing a robust TICS. First, the seven PRAD datasets were split into the training set and the validation set based on platforms, ensuring that the training and validation sets each have datasets using microarray and RNA-seq. Second, the TICS was calculated by the rank-based ssGSEA method, suitable for cross-platform transcriptome measurements. In the xCell study using ssGSEA to estimate immune cell infiltration, the authors showed that their scores reliably predicted enrichment when using different RNA-seq techniques and different microarray platforms . The xCell scores were agnostic to normalization methods or concerns related to batch effects, making them robust to both technical and biological noise . Third, in the training set, the signature scores of different datasets were not directly combined. Instead, the prognostic value of each signature in each dataset independently was assessed by comparing the BCR of two subgroups classified by the median value of each signature in each dataset. Fourth, in model training, the coefficient of each signature score in the algorithm of the TICS was decided by the pooled HR by meta-analysis, reflecting the combined effect among the four datasets in the training set using different platforms. In validation, the TICS showed great performance in each of the three datasets in the validation set, including one using microarray (GSE70768) and two using RNA-seq (DKFZ2018 and TCGA-PRAD). These results indicate the potential generalizability of the TICS across platforms.
As for other limitations, first, the size of the cohorts used to analyze the associations of the TICS with the efficacy of ADT or immunotherapy efficacy was small, and castration-resistant patients in the immunotherapy cohort may have molecular characteristics different from the training cohorts. Our findings based on these two cohorts had hypothesis-generating rather than hypothesis-testing value. More data on prostate cancer patients with non-RNA-seq or microarray techniques (e.g., protein expression level) are needed to evaluate the precision and robustness of the TICS model in the future. Second, the TICS was calculated by the ssGSEA method which requires the whole microarray or RNA-seq data. Despite that RNA-seq may benefit PRAD patients by detecting fusions (e.g., NTRK), large genomic rearrangements (e.g., BRCA1/2), and splice variants (e.g., AR-V7) that may guide targeted therapies, its low cost-effectiveness may prevent it from being widely used in clinical settings. Third, the TICS consists of eight immune signatures/cells, while exploring the spatial relationship between immune cells and cancer cells must be warranted in future studies. Fourth, established with retrospective datasets, more prospective studies are necessary to examine the clinical utility and cost-effectiveness of the TICS.
In summary, our study identified the TICS as a robust biomarker based on the expression signatures in prostate cancer. In a retrospective study of multiple cohorts, TICS can optimize the accuracy of existing models for predicting BCR, while it is also associated with the outcomes of ADT or immunotherapy. In the PRAD patients after curative surgery with a normal PSA and/or a low Gleason score, those with a low TICS might benefit from radical treatment and close monitoring. In addition, the predictive value of the TICS for ADT and immunotherapy is worth further investigation. This work developed a robust and effective tool for improving the risk stratification scheme of prostate cancer patients, which may optimize post-operative management and personalized treatment for patients with prostate cancer in clinical practice.
Availability of data and materials
The datasets generated and analysed during the current study are available in Gene-Expression Omnibus (GEO) at GSE46602, GSE54460, GSE70768, GSE70769, GSE94767, and GSE111177, and in TCGA at TCGA-PRAD. The data from DKFZ2018 that was analyzed in this study was obtained from cBioPortal at http://www.cbioportal.org/study/summary?id=prostate_dkfz_2018. The data of mCRPC patients treated with ipilimumab were obtained from the additional data files of Subudhi et al.
Androgen deprivation therapy
Biochemical recurrence-free survival
Castration-resistant prostate cancer
Gene set enrichment analysis
Homologous recombination deficiency
Immune checkpoint inhibitor
Loss of heterozygosity
Large-scale state transitions
Number of telomeric allelic imbalance
Prostatic specific antigen
Single sample gene set enrichment analysis
Tumor immune contexture score
Tumor immune microenvironment
Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. 2021;71:209–49.
Dubost C, Blondeau P, d’Allaines C, Piwnica A, Guilmet D. Surgical treatment of tetralogy of fallot. (Complete correction under extracorporeal circulation). Arch Mal Coeur Vaiss. 1965;58:1552–63.
Li F, Mahato RI. MicroRNAs and drug resistance in prostate cancers. Mol Pharm. 2014;11:2539–52. https://doi.org/10.1021/mp500099g
Feng Q, He B. Androgen receptor signaling in the development of castration-resistant prostate cancer. Front Oncol. 2019;9:858.
Freeman JW. Structural biology of the tumor microenvironment. Adv Exp Med Biol. 2021;1350:91–100.
Hinshaw DC, Shevde LA. The tumor microenvironment innately modulates cancer progression. Cancer Res. 2019;79:4557–66.
Molavi O, Ma Z, Hamdy S, Lavasanifar A, Samuel J. Immunomodulatory and anticancer effects of intra-tumoral co-delivery of synthetic lipid A adjuvant and STAT3 inhibitor, JSI-124. Immunopharmacol Immunotoxicol. 2009;31:214–21.
Jafari S, Molavi O, Kahroba H, Hejazi MS, Maleki-Dizaji N, Barghi S, et al. Clinical application of immune checkpoints in targeted immunotherapy of prostate cancer. Cell Mol Life Sci. 2020;77:3693–710.
Flammiger A, Bayer F, Cirugeda-Kühnert A, Huland H, Tennstedt P, Simon R, et al. Intratumoral T but not B lymphocytes are related to clinical outcome in prostate cancer. APMIS. 2012;120:901–8.
McArdle PA, Canna K, McMillan DC, McNicol AH, Campbell R, Underwood MA. The relationship between T-lymphocyte subset infiltration and survival in patients with prostate cancer. Br J Cancer. 2004;91:541–3.
Knutson KL, Disis ML. Tumor antigen-specific T helper cells in cancer immunity and immunotherapy. Cancer Immunol Immunother. 2005;54:721–8.
Yang W, Chen X, Hu H. CD4+ T-cell differentiation in vitro. Methods Mol Biol. 2020;2111:91–9.
Zhu J, Yamane H, Paul WE. Differentiation of effector CD4+ T cell populations. Annu Rev Immunol. 2010;28:445–89.
Sfanos KS, Bruno TC, Maris CH, Xu L, Thoburn CJ, Demarzo AM, et al. Phenotypic analysis of prostate-infiltrating lymphocytes reveals T H17 and Treg skewing. Clin Cancer Res. 2008;14:3254–61.
Zhang Q, Liu S, Ge D, Cunningham DM, Huang F, Ma L, et al. Targeting Th17-IL-17 pathway in prevention of micro-invasive prostate cancer in a mouse model. Prostate. 2017;77:888–99.
Hellmann MD, Paz-Ares L, Bernabe Caro R, Zurawski B, Kim S-W, Carcereny Costa E, et al. Nivolumab plus ipilimumab in advanced non-small-cell lung cancer. N Engl J Med. 2019;381:2020–31.
Isaacsson Velho P, Antonarakis ES. PD-1/PD-L1 pathway inhibitors in advanced prostate cancer. Expert Rev Clin Pharmacol. 2018;11:475–86.
Antonarakis ES, Piulats JM, Gross-Goupil M, Goh J, Ojamaa K, Hoimes CJ, et al. Pembrolizumab for treatment-refractory metastatic castration-resistant prostate cancer: multicohort, open-label phase II KEYNOTE-199 study. J Clin Oncol. 2020;38:395.
Gao Z, Tao Y, Lai Y, Wang Q, Li Z, Peng S, et al. Immune cytolytic activity as an indicator of immune checkpoint inhibitors treatment for prostate cancer. Front Bioeng Biotechnol. 2020;8:930.
Mortensen MM, Høyer S, Lynnerup AS, Ørntoft TF, Sørensen KD, Borre M, et al. Expression profiling of prostate cancer tissue delineates genes associated with recurrence after prostatectomy. Sci Rep. 2015;5:16018.
Long Q, Xu J, Osunkoya AO, Sannigrahi S, Johnson BA, Zhou W, et al. Global transcriptome analysis of formalin-fixed prostate cancer specimens identifies biomarkers of disease recurrence. Cancer Res. 2014;74:3228–37.
Ross-Adams H, Lamb A, Dunning M, Halim S, Lindberg J, Massie C, et al. Integration of copy number and transcriptomics provides risk stratification in prostate cancer: a discovery and validation cohort study. EBioMedicine. 2015;2:1133–44.
Luca BA, Brewer DS, Edwards DR, Edwards S, Whitaker HC, Merson S, et al. DESNT: a poor prognosis category of human prostate cancer. Eur Urol Focus. 2018;4:842–50.
Gerhauser C, Favero F, Risch T, Simon R, Feuerbach L, Assenov Y, et al. Molecular evolution of early-onset prostate cancer identifies molecular risk markers and clinical trajectories. Cancer Cell. 2018;34:996–1011.
Sharma NV, Pellegrini KL, Ouellet V, Giuste FO, Ramalingam S, Watanabe K, et al. Identification of the transcription factor relationships associated with androgen deprivation therapy response and metastatic progression in prostate cancer. Cancers (Basel). 2018;10:379.
Subudhi SK, Vence L, Zhao H, Blando J, Yadav SS, Xiong Q, et al. Neoantigen responses, immune correlates, and favorable outcomes after ipilimumab treatment of patients with prostate cancer. Sci Transl Med. 2020;12:eaaz3577.
Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang TH, et al. The immune landscape of cancer. Immunity. 2018;48:812-830.e14.
Marquard AM, Eklund AC, Joshi T, Krzystanek M, Favero F, Wang ZC, et al. Pan-cancer analysis of genomic scar signatures associated with homologous recombination deficiency suggests novel indications for existing cancer drugs. Biomark Res. 2015;3:9.
Sanchez-Vega F, Mina M, Armenia J, Chatila WK, Luna A, La KC, et al. Oncogenic signaling pathways in The Cancer Genome Atlas. Cell. 2018;173:321-337.e10.
Wang Y, Klijn JG, Zhang Y, Sieuwerts AM, Look MP, Yang F, et al. Gene-expression profiles to predict distant metastasis of lymph-node-negative primary breast cancer. Lancet. 2005;365:671–9.
Fan C, Prat A, Parker JS, Liu Y, Carey LA, Troester MA, et al. Building prognostic models for breast cancer patients using clinical variables and hundreds of gene expression signatures. BMC Med Genomics. 2011;4:3.
Chang HY, Sneddon JB, Alizadeh AA, Sood R, West RB, Montgomery K, et al. Gene expression signature of fibroblast serum response predicts human cancer progression: similarities between tumors and wounds. PLoS Biol. 2004;2:E7.
Wolf DM, Lenburg ME, Yau C, Boudreau A, Van’t Veer LJ. Gene co-expression modules as clinically relevant hallmarks of breast cancer diversity. PLoS ONE. 2014;9:e88309.
Luan J, Zhang Q, Song L, Wang Y, Ji C, Cong R, et al. Identification and validation of a six immune-related gene signature for prediction of biochemical recurrence in localized prostate cancer following radical prostatectomy. Transl Androl Urol. 2021;10:1018–29.
Long G, Ouyang W, Zhang Y, Sun G, Gan J, Hu Z, et al. Identification of a DNA repair gene signature and establishment of a prognostic nomogram predicting biochemical-recurrence-free survival of prostate cancer. Front Mol Biosci. 2021;8:608369.
Wang X, Lv Z, Xia H, Guo X, Wang J, Wang J, et al. Biochemical recurrence related metabolic novel signature associates with immunity and ADT treatment responses in prostate cancer. Cancer Med. 2022. https://doi.org/10.1002/cam4.4856
Wu X, Lv D, Lei M, Cai C, Zhao Z, Eftekhar M, et al. A 10-gene signature as a predictor of biochemical recurrence after radical prostatectomy in patients with prostate cancer and a Gleason score ≥7. Oncol Lett. 2020;20:2906–18.
Wang S, Gao J, Lei Q, Rozengurt N, Pritchard C, Jiao J, et al. Prostate-specific deletion of the murine Pten tumor suppressor gene leads to metastatic prostate cancer. Cancer Cell. 2003;4:209–21.
Barbieri CE, Baca SC, Lawrence MS, Demichelis F, Blattner M, Theurillat JP, et al. Exome sequencing identifies recurrent SPOP, FOXA1 and MED12 mutations in prostate cancer. Nat Genet. 2012;44:685–9.
Mu P, Zhang Z, Benelli M, Karthaus WR, Hoover E, Chen CC, et al. SOX2 promotes lineage plasticity and antiandrogen resistance in TP53-and RB1-deficient prostate cancer. Science (80−). 2017;355:84–8.
Blattner M, Liu D, Robinson BD, Huang D, Poliakov A, Gao D, et al. SPOP mutation drives prostate tumorigenesis in vivo through coordinate regulation of PI3K/mTOR and AR SIGNALING. Cancer Cell. 2017;31:436–51.
Fraser M, Sabelnykova VY, Yamaguchi TN, Heisler LE, Livingstone J, Huang V, et al. Genomic hallmarks of localized, non-indolent prostate cancer. Nature. 2017;541:359–64.
Neeb A, Herranz N, Arce-Gallego S, Miranda S, Buroni L, Yuan W, et al. Advanced prostate cancer with ATM loss: PARP and ATR inhibitors. Eur Urol. 2021;79:200–11.
Hammarsten P, Josefsson A, Thysell E, Lundholm M, Hägglöf C, Iglesias-Gato D, et al. Immunoreactivity for prostate specific antigen and Ki67 differentiates subgroups of prostate cancer related to outcome. Mod Pathol. 2019;32:1310–9.
Zhao SG, Lehrer J, Chang SL, Das R, Erho N, Liu Y, et al. The immune landscape of prostate cancer and nomination of PD-L2 as a potential therapeutic target. J Natl Cancer Inst. 2019;111:301–10.
Zhang C, Chen T, Li Z, Liu A, Xu Y, Gao Y, et al. Depiction of tumor stemlike features and underlying relationships with hazard immune infiltrations based on large prostate cancer cohorts. Brief Bioinform. 2021;22:bbaa211.
Quinn DI, Henshall SM, Sutherland RL. Molecular markers of prostate cancer outcome. Eur J Cancer. 2005;41:858–87.
Cuzick J, Swanson GP, Fisher G, Brothman AR, Berney DM, Reid JE, et al. Prognostic value of an RNA expression signature derived from cell cycle proliferation genes in patients with prostate cancer: a retrospective study. Lancet Oncol. 2011;12:245–55.
Lu YC, Kuo MC, Hong JH, Jaw FS, Huang CY, Cheng JCH, et al. Lower postoperative natural killer cell activity is associated with positive surgical margins after radical prostatectomy. J Formos Med Assoc. 2020;119:1673–83.
Pasero C, Gravis G, Granjeaud S, Guerin M, Thomassin-Piana J, Rocchi P, et al. Highly effective NK cells are associated with good prognosis in patients with metastatic prostate cancer. Oncotarget. 2015;6:14360–73.
Creighton CJ. A gene transcription signature associated with hormone independence in a subset of both breast and prostate cancers. BMC Genomics. 2007;8:199.
Eiró N, Bermudez-Fernandez S, Fernandez-Garcia B, Atienza S, Beridze N, Escaf S, et al. Analysis of the expression of interleukins, interferon β, and nuclear factor-κ B in prostate cancer and their relationship with biochemical recurrence. J Immunother. 2014;37:366–73.
Mohler JL, Antonarakis ES, Armstrong AJ, D’Amico AV, Davis BJ, Dorff T, et al. Prostate cancer, Version 2. 2019, NCCN clinical practice guidelines in oncology. J Natl Compr Cancer Netw. 2019;17:479–505.
Schaeffer E, Srinivas S, Antonarakis ES, Armstrong AJ, Bekelman JE, Cheng H, et al. NCCN guidelines insights: prostate cancer, Version 1. 2021. J Natl Compr Cancer Netw. 2021;19:134–43.
Malik A, Srinivasan S, Batra J. A new era of prostate cancer precision medicine. Front Oncol. 2019;9:1263.
Hughes C, Murphy A, Martin C, Sheils O, O’Leary J. Molecular pathology of prostate cancer. J Clin Pathol. 2005;58:673–84.
Mithal P, Allott E, Gerber L, Reid J, Welbourn W, Tikishvili E, et al. PTEN loss in biopsy tissue predicts poor clinical outcomes in prostate cancer. Int J Urol. 2014;21:1209–14.
Konoshenko MY, Bryzgunova OE, Laktionov PP. miRNAs and androgen deprivation therapy for prostate cancer. Biochim. Biophys. Acta - Rev. Cancer. 2021.
Long X, Hou H, Wang X, Liu S, Diao T, Lai S, et al. Immune signature driven by ADT-induced immune microenvironment remodeling in prostate cancer is correlated with recurrence-free survival and immune infiltration. Cell Death Dis. 2020;11:779.
Cha HR, Lee JH, Ponnazhagan S. Revisiting immunotherapy: a focus on prostate cancer. Cancer Res. 2020;80:1615–23.
Graff JN, Beer TM, Alumkal JJ, Slottke RE, Redmond WL, Thomas GV, et al. A phase II single-arm study of pembrolizumab with enzalutamide in men with metastatic castration-resistant prostate cancer progressing on enzalutamide alone. J Immunother Cancer. 2020;8:e000642.
Fizazi K, González Mella P, Castellano D, Minatta JN, Rezazadeh Kalebasty A, Shaffer D, et al. Nivolumab plus docetaxel in patients with chemotherapy-naïve metastatic castration-resistant prostate cancer: results from the phase II CheckMate 9KD trial. Eur J Cancer. 2022;160:61–71.
Aran D, Hu Z, Butte AJ. xCell: digitally portraying the tissue cellular heterogeneity landscape. Genome Biol. 2017;18:1–14. https://doi.org/10.1186/s13059-017-1349-1
2018 Entrepreneurial Leading Talent of Guangzhou Huangpu District and Guangzhou Development District (2022-L023).
Ethics approval and consent to participate
Consent for publication
Y Liao, Y Xu, T Zhu, W Chen, G Wang, Y Han, and S Cai are employees of Burning Rock Biotech. The other authors declare no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Clinical characteristics of prostate cancer patients. Table S2. Genes involved in the 92 immune contexture-related signatures in the study. Table S3. Identified signatures correlated with BCRFS in meta-analysis of training cohorts. Table S4. Multicollinearity test used by variance inflation factors (VIF) and kappa test for selected signatures in all cohorts. Table S5. The C-index and AUC value of training and validation cohort. Table S6. Comparison of the TICS and other prognostic signatures. Table S7. GSEA results of the KEGG pathways in the high- and low-TICS groups in TCGA cohort.
The prognostic effect of immune signatures in the training sets. Fig. S2. Comparison of the TICS and other prognostic signatures.
About this article
Cite this article
Han, S., Shi, T., Liao, Y. et al. Tumor immune contexture predicts recurrence after prostatectomy and efficacy of androgen deprivation and immunotherapy in prostate cancer. J Transl Med 21, 194 (2023). https://doi.org/10.1186/s12967-022-03827-4