Intrinsic Molecular Subtypes of Metastatic Castration-Resistant Prostate Cancer

Prostate cancer is a clinically and molecularly heterogenous disease. Outcomes for patients can vary substantially, as the presentation of prostate cancer can range from indolent tumors to highly aggressive and lethal disease. The molecular heterogeneity of the disease is reflected in the panoply of DNA alterations in key driver genes initially identified in localized prostate cancer (1–3). RNA-based gene expression profiling approaches further demonstrated a division between luminal-like or basal-like subtypes in prostate tumors (4–6) that have different prognoses and importantly may respond differently to anti-androgen therapy, the backbone of systemic therapy in prostate cancer. These clinically important subtypes are now being incorporated into national clinical trials to improve patient selection for various standard-of-care and experimental therapies.

Metastatic prostate cancer is a very different disease than localized prostate cancer; this is particularly true for metastatic castration-resistant prostate cancer (mCRPC), which represents end-stage disease with poor outcomes. The genomic landscape of mCRPC is characterized by a higher frequency of DNA alterations across many oncogenic drivers, especially in the androgen receptor (AR) due to the selective pressure of therapy (7–10). However, these DNA alterations have failed to capture the full heterogeneity of the disease and, with the exception of DNA-repair deficiencies and microsatellite instability in a minority of cases, cannot be used to guide clinical decision-making. In addition to DNA alterations, important transcriptional changes also take place during disease progression to mCRPC. Transcriptional changes have been identified on the basis of histologic differences between subsets of mCRPC; a number of studies have identified lineage plasticity, most commonly toward a neuroendocrine (NE) lineage, in response to prolonged AR-targeted therapy through which prostate cancer becomes independent of AR signaling for proliferation (10–14). Indeed, more recent studies focusing on differences in the NE and AR signaling axes between mCRPC samples have defined five different subgroups: AR-high prostate cancer, AR-low prostate cancer, amphicrine tumors composed of cells co-expressing AR and NE genes (AMPC), double-negative tumors (i.e., AR^–/NE^–), and tumors with small cell or NE gene expression without AR activity (SCNC, small-cell/NE prostate cancer; ref. 15). Lineage plasticity without NE differentiation may also be present toward a more stem cell phenotype through epithelial-to-mesenchymal transition or toward a gastrointestinal lineage (16). Additional RNA profiling studies have demonstrated that transcriptional subtypes that exist in localized prostate cancer may exist in mCRPC; a recent study demonstrated that luminal and basal subtypes in mCRPC predict response to AR signaling inhibitors (ASI), analogous to localized prostate cancer (17).

Although previous subtyping studies have revealed the biological heterogeneity of mCRPC, comprehensive unsupervised transcriptome-wide clustering approaches have seen limited use in metastatic samples to date, likely due to the lack of large cohorts of metastatic disease with molecular profiling, which requires metastatic tissue biopsies that are logistically difficult to obtain. Compared with localized prostate cancer, in which cohorts with thousands of patient samples have been published (18), the largest mCRPC cohorts with next-generation sequencing and clinical outcomes are comprised of <150 (7, 19, 20). To allow for unsupervised data-driven approaches that may reveal key information that have not yet been hypothesized and tested, we recently compiled publicly available RNA sequencing (RNA-seq) from mCRPC samples to form the largest dataset of its kind (17) that represents a unique opportunity for identifying the intrinsic molecular subtypes of mCRPC and assessing for molecular and clinical correlates.

We combined and batch-effect corrected data from four mCRPC cohorts to assess intrinsic subtypes, and analyzed RNA, DNA, and overall survival (OS) as previously published (17). The mCRPC cohorts used were from the Fred Hutchinson Cancer Research Center (FHCRC, N = 157; ref. 21), an SCNC-enriched cohort from Weill Cornell Medicine (WCM, N = 49; ref. 14), and cohorts from the Stand Up 2 Cancer/Prostate Cancer Foundation East Coast Dream Team (ECDT, N = 266; refs. 8, 19) and the West Coast Dream Team (WCDT, N = 162; refs. 7, 10, 22).

Normalized gene expression, mutation calls, and copy number for FHCRC, WCM, and the ECDT were obtained directly from cBioPortal (23). The same data from the WCDT were obtained from prior publications (7, 10, 17). RNA-seq batch correction was performed as previously described (17). In our primary analysis of the RNA-seq data, we included genes that had expression data in all samples. Hierarchical clustering was performed on the top 1,000 genes with the highest variance across all the samples. Spearman's correlation was used as the distance function, and default “ward.D” agglomeration method was used. Using this method, we identified four distinct subtypes based on the gene expression patterns. Clinical and pathologic variables were obtained from the original publications of these cohorts. Pathway scores using the Hallmark Pathways from MSigDb (24) were calculated using gene set variation analysis (25). Oncogene-activating alterations (amplification and/or mutation) and tumor-suppressor bi-allelic–inactivating alterations (copy-number loss and/or mutation) were defined as previously published (17). CIBERSORTx was used to calculate an absolute immune score (26).

Statistical testing between the subtypes and continuous variables was performed using a Wilcoxon rank-sum test. Statistical testing between the subtypes and categorical variables was performed using the Fisher's exact test. Survival analyses were performed using Cox regression and visualized using the Kaplan–Meier method. All statistical testing was performed in R version 4.0.4. All statistical testing was two-sided, and a P value of <0.05 was considered statistically significant.

Only previously published data were used for this study. FHCRC, WCM, and the ECDT were downloaded from www.cbioportal.org. WCDT genomics data are available from dbGaP (phs001648) and EGA (EGAD00001009065).

To identify intrinsic RNA subtypes within mCRPC, we performed hierarchical clustering on 634 total mCRPC samples with RNA-seq, using the top 1,000 genes with the highest variance (Fig. 1). This resulted in four distinct subtypes based on RNA expression patterns. These subtypes were not associated with the cohort of origin, supporting the effectiveness of the batch correction. However, there visually did appear to be associations with both biopsy site and SCNC (hereafter referring to the phenotypic definition as originally identified in each cohort). In the two subtypes associated with SCNC, the rate of SCNC was higher in one (78% SCNC, this RNA subtype hereafter referred to as NE-classic) versus the other (29% SCNC, with all of the samples originating from a liver biopsy, this RNA subtype hereafter referred to as NE-liver). This was reflected in gene expression patterns of the NE markers SYP and CHGA that were higher in the NE subtypes, more so in NE-classic (Fig. 2A and B; Supplementary Table S1). This was also reflected in AR and KLK3 (encoding PSA) gene expression, with the lowest expression in NE-classical, and intermediate expression in NE-liver subtype (Fig. 2C and D; Supplementary Table S1).

Next, we sought to understand the pathway level differences between the novel RNA subtypes using the MSigDb Hallmark pathways. We first examined the AR signaling pathway, where, as expected, the two adenocarcinoma subtypes had highest expression, followed by NE-liver, and then NE-classic, consistent with the AR and KLK3 expression above (Fig. 3A; Supplementary Table S2). SCNC is also known to be highly proliferative, and this characteristic was reflected in NE-classic having higher expression of proliferative pathways such as the G₂–M checkpoint and E2F targets compared with the adenocarcinoma subtypes, with NE-liver between the two (Fig. 3B and C; Supplementary Table S2).

There was a small but statistically significant difference between the two RNA-intrinsic adenocarcinoma subtypes for AR signaling (Fig. 3A; Supplementary Table S2), consistent with differences noted with previous hypothesis-driven approaches (15, 17). However, the most prominent difference was in immune pathways, with consistently higher expression of all eight immune-related hallmark pathways (Fig. 3D–K; Supplementary Table S2) in one RNA subtype (hereafter referred to adeno-immune) compared with the other (hereafter referred to adeno-classic). The results using the CIBERSORTx overall immune score were consistent with the pathway analyses, with overall immune scores higher in adeno-immune versus adeno-classic (Fig. 3L; Supplementary Table S2). These data suggest that there may be two immunogenically distinct subgroups of mCRPC adenocarcinoma. Comparison between the two RNA NE subtypes demonstrated that the NE-liver subtype exhibited higher expression of 7 out of the 8 immune-related hallmark pathways (all except for the allograft rejection signature; Fig. 3D–K; Supplementary Table S2). Although there was no significant difference in the CIBERSORTx overall immune score between the two NE subtypes (Fig. 3L; Supplementary Table S2), the pathway analysis suggests that there are qualitative differences in the tumor immunologic response between the two NE subtypes.

We next sought to investigate DNA alterations in key prostate cancer oncogenic-driver genes. Alterations in the tumor suppressors RB1, PTEN, and TP53 are all known to be associated with aggressive prostate cancer and SCNC (27). Of the two SCNC-enriched RNA subtypes, the NE-classic was more enriched for RB1 alterations (63% in NE-classic vs. 26% in NE-liver) and the NE-liver for PTEN alterations (21% in NE-classic vs. 39% in NE-liver), with similar TP53 alteration rates (43% in NE-classic vs. 37% in NE-liver; Fig. 4A). Although the rate of AR alterations was lower in both NE subtypes than in the adenocarcinoma subtypes, there were differences between the NE subtypes (33% in NE-classic vs. 63% in NE-liver) consistent with the AR and KLK3 expression data. There was also a difference in FOXA1 alteration rates (44% in NE-classic vs. 29% in NE-liver), but not MYC (70% in NE-classic vs. 71% in NE-liver; Fig. 4A; Supplementary Table S3). With regard to the adenocarcinoma subtypes, there were differences in PTEN alterations (43% in adeno-classic vs. 22% in adeno-immune). Mutation rates for APC and SPOP differed between both NE subtypes (APC 2.3% in NE-classic vs. 14% in NE-liver; SPOP 0% in NE-classic vs. 14% in NE-liver) and adenocarcinoma subtypes (APC 1% in adeno-classic vs. 5.2% in adeno-immune; SPOP 1% in adeno-classic vs. 11% in adeno-immune; Fig. 4B; Supplementary Table S3). We found no differences between the subtypes in tumor mutational burden (as a surrogate for MSI/dMMR) or BRCA2 deep deletion or mutations, which are the two genomic markers currently most used for treatment selection in mCRPC (for immunotherapy or PARP-inhibitor therapy, respectively; Supplementary Fig. S1).

The site of metastasis was different between the NE subtypes, with all of the NE-liver tumors coming from liver metastases, whereas the distribution was more evenly split in the NE-classic tumors between the different metastatic locations (Fig. 5; Supplementary Table S4). Interestingly, a difference between the two adenocarcinoma RNA subtypes was also observed, with adeno-immune being enriched for distant metastases whereas adeno-classic was enriched for lymph node (LN) metastases. The gene expression from bulk RNA-seq of tumor biopsies is a mix of the tumor cells and the surrounding non-tumor tissue that will influence comparison across tumor sites. However, immune-related measures should favor LNs preferentially, as the non-tumor background is mostly immune cells, which is in contrast with what we observe for adeno-immune versus adeno-classic, suggesting that there may indeed be a difference in immunogenicity between the two adenocarcinoma RNA-intrinsic subtypes.

The association between different RNA subtypes and clinical outcomes was also investigated. As expected, the two adenocarcinoma subtypes had significantly improved OS compared with the NE subtypes [Hazard Ratio (HR), 0.265; 95% confidence interval (CI), 0.159−0.44; P < 0.001; Fig. 6A], without a difference between them. However, when we account for ASI therapy, we observed that the adeno-classic subtype showed a significant benefit from ASI therapy (HR, 0.219; 95% CI, 0.1−0.476; P < 0.001) after the biopsy, whereas the adeno-immune subtype had a weaker trend in the same direction that was not significant (HR, 0.592; 95% CI, 0.337−1.04; P = 0.068; Fig. 6B). Although the NE-classic and NE-liver RNA-intrinsic subtypes also showed similar worse survival, there were different proportions of patients that were not originally identified as SCNC in these subtypes.

Data-driven and unsupervised approaches to derive intrinsic subtypes can reveal differences not immediately apparent based on hypothesis-driven approaches based on prior biological knowledge. Herein, we present the largest intrinsic RNA-based molecular subtyping of mCRPC in 634 samples with RNA-seq. We find two distinct adenocarcinoma subtypes, one of which (adeno-immune) is characterized by higher expression of immune pathways, higher CIBERSORTx immune scores, a lower HR for ASI benefit, and non-LN metastasis tropism compared with an adeno-classic subtype. We also identify two distinct NE subtypes, including an NE-liver subgroup characterized by liver metastasis tropism, PTEN loss, and APC and SPOP mutations compared with an NE-classic subgroup.

SCNC is an aggressive, androgen-independent subtype of mCRPC that can either be present de novo, or more commonly emerge during the course of treatment via lineage plasticity (12–14, 27, 28). SCNC has been shown to exhibit distinct transcriptomic and epigenomic patterns, and is genomically characterized by loss of the tumor suppressors PTEN, RB1, and TP53 (12–14, 27–31). When we compare our intrinsic RNA subtypes to the five subtypes previously described by AR and NE markers (15), the NE-classic subtype likely corresponds well to the SCNC (AR⁻, NE⁺) subgroup. The NE-liver subgroup shares similarities to some other subgroups described previously, such as AMPC (AR⁺, NE⁺) or intermediate atypical carcinoma (32), a histologic subgroup between adenocarcinoma and SCNC. It is also possible that it represents a transition state between adenocarcinoma and emerging SCNC, which on bulk sequencing would have gene expression between the two, or specific biology driving the aggressive tumors that metastasize to the liver. Indeed, both intrinsic RNA subtypes enriched for SCNC tumors identify tumors not originally diagnosed as SCNC, especially in the NE-liver group. However, the prognosis of both subtypes is equally poor, and the transcriptional similarities would suggest that similar therapeutic approaches could be tried. Interestingly, the difference in RB1/PTEN alterations between the RNA NE subtypes also suggest that these alterations are not equivalent. As PTEN-altered tumors begin to have new treatment options (33), the management of NE-liver tumors may diverge from NE-classic. We did not observe a distinct double-negative cluster in this unsupervised intrinsic clustering resulting in four main clusters (15), and such tumors may be grouped in with the other subtypes.

The differences between the intrinsic adenocarcinoma subtypes of mCRPC are also potentially clinically important. Interestingly, the two intrinsic adenocarcinoma subtypes do not differ via an AR high versus low dichotomy described in these previous subtyping efforts (15, 17). mCRPC has traditionally been thought of as not particularly immunogenic, with poor response rates to immune checkpoint blockade (34, 35). Our adeno-immune subtype suggests that not all tumors may be equally quiescent, and suggest potential biomarker-driven approaches for immunotherapies in mCRPC. Despite analyzing the largest harmonized cohort of mCRPC tumors with gene expression data, a limitation to this study is that limited clinical treatment annotation and follow-up data were only available for a subset of the patients, and did not allow for analysis of association between subtypes and treatments other than AR-targeted therapy.

Ultimately, our results emphasize the molecular heterogeneity of mCRPC. Our unbiased clustering approach confirms key findings identified through hypothesis-driven approaches in the field, such as the existence of distinct subgroups driven by AR signaling versus NE features; however, we also discover that AR-driven adenocarcinoma subgroups can be differentiated by signaling along immune pathways, and that NE subgroups may be grouped biologically, with distinctions in the site of metastases and DNA alterations. In total, our study suggests that future studies and clinical trials should consider transcriptome-wide profiling to better understand how these differences affect treatment and outcomes.

K.T. Helzer reports that their spouse is employed with Epic Systems. J.M. Lang reports personal fees from Janssen, Astellas, Arvinas, Gilead, Sanofi, Pfizer, 4D Pharma, and AstraZeneca outside the submitted work. R. Aggarwal reports personal fees from Dendreon, Pfizer, Axion Healthcare, Jubilant Therapeutics, Exelixis, Bayer, and Advanced Accelerator Applications; grants from Janssen, Amgen, and Clovis; and grants and personal fees from Merck and AstraZeneca outside the submitted work. E.J. Small reports personal fees and other support from Fortis, Teon, Janssen, Johnson and Johnson, and Ultragenyx, as well as other support from Harpoon during the conduct of the study. E.J. Small also reports personal fees and other support from Fortis, Teon, Janssen, Johnson and Johnson, and Ultragenyx, as well as other support from Harpoon outside the submitted work. M. Sjöström reports grants from the Swedish Research Council, the Swedish Society of Medicine, and the Prostate Cancer Foundation during the conduct of the study. S.G. Zhao reports grants from Department of Defense and National Institutes of Health during the conduct of the study; patent applications with Decipher Biosciences on molecular signatures in prostate cancer unrelated to this work; and their spouse is employed with Exact Sciences with RSUs. No disclosures were reported by the other authors.

E. Feng: Formal analysis, methodology, writing–review and editing. N.R. Rydzewski: Data curation, formal analysis, methodology, writing–review and editing. M. Zhang: Data curation, supervision, writing–review and editing. A. Lundberg: Data curation, writing–review and editing. M. Bootsma: Methodology, writing–review and editing. K.T. Helzer: Methodology, writing–review and editing. J.M. Lang: Supervision, writing–review and editing. R. Aggarwal: Data curation, project administration, writing–review and editing. E.J. Small: Data curation, project administration, writing–review and editing. D.A. Quigley: Data curation, supervision, project administration, writing–review and editing. M. Sjöström: Data curation, supervision, methodology, writing–original draft, project administration, writing–review and editing. S.G. Zhao: Conceptualization, data curation, supervision, funding acquisition, methodology, writing–original draft, project administration, writing–review and editing.

This research was supported by a Stand Up To Cancer-Prostate Cancer Foundation Prostate Cancer Dream Team Award (SU2C-AACR-DT0812; to E. J. Small) and by the Movember Foundation. Stand Up To Cancer is a division of the Entertainment Industry Foundation. This research grant was administered by the American Association for Cancer Research, the scientific partner of SU2C. M. Sjöström was supported by the Swedish Research Council (Vetenskapsrådet) with grant number 2018–00382, the Swedish Society of Medicine (Svenska Läkaresällskapet), and a Prostate Cancer Foundation Young Investigator Award. D.A. Quigley was funded by a Prostate Cancer Foundation Young Investigator Award and a BRCA Foundation Young Investigator Award. Additional funding was provided by a UCSF Benioff Initiative for Prostate Cancer Research award. S.G. Zhao was funded by the Department of Defense (PC190039) and the National Institutes of Health (1DP2CA271832). J.M. Lang and S.G. Zhao are supported by Department of Defense (PC200334) and the University of Wisconsin Carbone Cancer Center Support grant (P30 CA014520).

The publication costs of this article were defrayed in part by the payment of publication fees. Therefore, and solely to indicate this fact, this article is hereby marked “advertisement” in accordance with 18 USC Section 1734.