N-Myc–mediated epigenetic reprogramming drives lineage plasticity in advanced prostate cancer

Lineage plasticity, a process by which differentiated cells lose their identity and acquire an alternative lineage program, has recently been identified as an emerging mechanism of resistance to targeted therapies in several cancer types, including prostate cancer (1–3). Under pressure of androgen receptor–directed (AR-directed) therapies, a subset of patients with castration-resistant prostate cancer (CRPC) develop tumors that lose AR expression, AR-signaling dependence, and luminal prostate cancer markers and acquire alternative lineage programs, including the development of small-cell neuroendocrine carcinoma (neuroendocrine prostate cancer [NEPC]; see refs. 1, 4, 5). Few treatment options currently exist for patients developing NEPC and patient prognosis is poor. Although NEPC tumors arise clonally from prostate adenocarcinoma (PCa) and share genomic alterations, there is significant epigenetic deregulation during the plasticity process (1). In the last decade, we and others have identified several candidate drivers of plasticity including, but not limited to, loss of TP53 and RB1 (6–12), loss of REST (13), and upregulation of MYCN (1, 14–17), enhancer of zeste homolog 2 (EZH2) (7, 15), BRN2 (18), and SOX2 (7, 9, 18).

N-Myc (encoded by MYCN) is a transcription factor of the Myc family, essential for normal brain development during embryogenesis. Often amplified and overexpressed in neuroblastoma, medulloblastoma, retinoblastoma, astrocytomas, and glioblastoma multiforme (19–24), N-Myc has also been implicated in non-neuronal tumors, including prostate cancer, lymphomas, acute myeloid leukemia (AML), rhabdomyosarcomas, Wilms tumors, small-cell lung cancer, and pancreatic tumors (25–33) (reviewed in ref. 34). N-Myc is not normally expressed in the prostate epithelial lineage but is overexpressed in the majority of NEPC and a subset of CRPC adenocarcinomas (CRPC-Adenos) (14, 15). N-Myc forms a heterodimer with Max and specifically binds DNA at E-box motifs to activate or repress gene expression, depending on recruited cofactors (35). Genes involved in cell cycle, proliferation, and apoptosis are the most well-described targets of N-Myc across different cell types (36).

In this study, we describe what we believe is a novel role for N-Myc in prostate cancer, characterized by changes in the N-Myc cistrome, its interacting cofactors, as well as reprogramming of the epigenome in an androgen-dependent manner. This reprogramming is associated with an induction of lineage plasticity in prostate cancer cells and a switch towards a neural identity that favors the development of AR independence and NEPC.

Clinical NEPC is associated with neural lineage. Although there is a spectrum within the pathologic subtype of NEPC, we have found that NEPC tumor cells defined by morphologic features often lose AR expression and express neuroendocrine markers (37). This observation suggests a transition from epithelial to neural lineage, which may involve dedifferentiation. To fully appreciate the spectrum of lineage states, we performed whole transcriptome analyses on patient samples, including the largest cohort of NEPC patients to date. Gene expression was assessed from metastatic tumor biopsies of patients with pathology-confirmed NEPC (n = 36) and CRPC-Adeno (n = 73), as well as from localized PCa (n = 66). Relative to PCa, gene set enrichment analysis (GSEA) of our NEPC cohort data revealed a significant enrichment for stem cell genes associated not only with normal neuroendocrine cell precursors (neural crest stem cells [NCSCs]), but also with activated neural stem cells (aNSCs) (Figure 1A and Supplemental Table 1; supplemental material available online with this article; https://doi.org/10.1172/JCI127961DS1). Unsupervised clustering analysis of the NEPC versus PCa leading-edge genes from embryonic stem cell (ESC), NCSC, or NSC gene sets (38–40) segregated patients according to their tumor histological classification (PCa, CRPC-Adeno, or NEPC; Supplemental Figure 1A). Interestingly, the stratification of CRPC-Adeno and NEPC patients based on levels of MYCN expression correlated with expression of neural lineage genes and inversely with AR signaling (Figure 1A and Supplemental Figure 1B). To further validate the association with NSC genes, we purified murine NSCs from the subventricular zone and performed RNA sequencing (RNA-seq) and histone-3 lysine-4 trimethylation (H3K4me3) chromatin immunoprecipitation followed by DNA sequencing (ChIP-seq). From these studies, we combined the upregulated genes marked by H3K4me3 with publicly available single-cell RNA-seq data from activated NSCs (41) and revealed a similar enrichment of NSC pathways in NEPC patient samples versus PCa samples (Figure 1A).

Figure 1

Clinical NEPC is associated with neural lineage. (A) Top: Enrichment plots of the Neural Stem Cell Markers and Lee Neural Crest Stem Cell Up gene sets between indicated groups. Bottom: Targeted GSEA in the 5 NEPC samples with the highest (N-Myc^hi) or lowest (N-Myc^lo) level of MYCN expression versus PCa (n = 66) patient samples, NEPC N-Myc^hi versus NEPC N-Myc^lo, and on the 5 CRPC with the highest level of MYCN expression versus the 5 lowest. *FDR q value < 0.05, **FDR q value < 0.01, ***FDR q value < 0.001. (B) Kaplan-Meier plots of CRPC (n = 57) patients, NEPC (n = 24) patients, or CRPC plus NEPC (n = 81) patients, stratified into 2 categories according to the median value of MYCN mRNA expression. Survival analysis was performed using the Kaplan-Meier estimator (log-rank test). qNSC, quiescent NSC.

Moreover, we observed that high MYCN expression in patient CRPC-Adeno tumors was significantly associated with a worse overall survival compared with those with low MYCN expression (37.7 versus 80.3 months, hazard ratio [HR] = 1.95, 95% confidence interval [CI] 0.92–4.09, P value = 0.04). This was also true for patients with NEPC (25.9 versus 41.4 months, HR = 3.31, 95% CI 1.22–9.09, P value = 0.0064) as well as the combined cohort of all 81 patients (34.0 versus 76.5 months, HR = 2.27, 95% CI 1.24–4.11, P value = 0.002, Figure 1B). While these data reveal that MYCN expression correlates with neural lineage programs and poor clinical outcome, in both CRPC and NEPC patients, the underlying mechanism remains poorly described. We therefore sought to determine the precise mechanism that drives a prostate tumor epithelial cell to lose its luminal markers and gain a more neural-like lineage.

AR signaling alters the N-Myc transcriptome in vivo. While N-Myc cooperates with EZH2 to abrogate AR signaling (15), the role of N-Myc in driving lineage plasticity in an NEPC-like context of low-to-no AR signaling has not been characterized. Using our previously described Pb-Cre^+/– Pten^fl/fl LSL-MYCN^+/+ genetically engineered mouse (GEM) model, which overexpresses human MYCN (15), we found that N-Myc overexpression led to the formation of prostate tumors with 100% penetrance compared with Pb-Cre^+/– Pten^fl/fl LSL-MYCN^–/– mice (Supplemental Figure 2A). Interestingly, survival of Pb-Cre^+/– Pten^fl/fl LSL-MYCN^+/+ mice was increased by an average of 12 weeks in response to castration (Supplemental Figure 2A). This is likely due to the fact that the tumors are heterogeneous, with a component of AR-responsive adenocarcinoma, consistent with data from another GEM model of adenocarcinoma-to-NEPC transformation (7). Despite the prolonged survival, castrated mice developed invasive prostate tumors that metastasized to multiple locations, including the liver, at 6 months after castration. We also noted an increase in poorly differentiated foci that lost expression of luminal (AR and cytokeratin 8 [KRT8]) and basal (KRT5) markers and gained expression of the epithelial-mesenchymal transition marker vimentin (VIM) and the NEPC marker neural cell adhesion molecule 1 (NCAM1) (refs. 37, 42 and Figure 2A). Primary and metastatic lesions in castrated mice contained tumor foci with up to 90% of the total tumor area comprised of divergent differentiation (e.g., intestinal, squamous, and sarcomatoid as previously described in ref. 15, or chondroid differentiation) compared with up to 25% in the intact animals (Supplemental Figure 2, B and C). In addition, we observed large foci of neural differentiation characterized by ganglion-like cells. These cells had abundant eosinophilic cytoplasm and nuclei with prominent, centrally located nucleoli. Immunohistochemical staining of serial 4-μm tissue sections revealed populations of cells that were either positive for epithelial (AR and KRT8), mesenchymal (VIM), or neural markers (S100) but also contained cells that were positive for multiple lineages (AR and VIM) (Supplemental Figure 2D). This suggests that these double-positive cells may be transitioning from one lineage to another, similar to what has been observed in other mouse models (7, 12, 16). The ganglion differentiation is also consistent with previous observations made in N-Myc–driven neuroblastoma models (43). The loss of AR and gain of VIM and S100 suggest further differentiation towards the neural lineage, and some of the VIM- and S100-positive cells have invaded the local vasculature (Figure 2B), consistent with a prometastatic phenotype. These data suggest that the removal of androgen signaling enables a wider variety of N-Myc–induced differentiation programs.

Figure 2

AR signaling alters the N-Myc transcriptome in vivo. (A) Photomicrographs of H&E staining, vimentin (VIM), NCAM1, AR, cytokeratin 8 (KRT8), and 5 (KRT5) immunohistochemistry, and MYCN RNA in situ hybridization (RNAish) on primary prostate tumor region enriched with sarcomatoid differentiation (left) or liver metastatic lesion (right) from Pb-Cre^+/– Pten^fl/fl LSL-MYCN^+/+ mice 6 months after castration. Scale bars: 50 μm. (B) Photomicrographs of H&E staining or immunohistochemical staining for epithelial markers (AR and KRT8), a mesenchymal marker (VIM), or neural/ganglionic marker (S100) on 4-μm serial sections from mouse C1 (Supplemental Figure 2). Dotted lines indicate conventional adenocarcinoma adjacent to neural/ganglionic cells. Arrows indicate VIM⁺S100⁺ tumor cells that have invaded local vasculature. Scale bars: 50 μm. (C) Top: N-Myc signatures defined from 22Rv1-N-Myc xenografts versus 22Rv1 control (CTL) xenografts (–1 < log₂[fold change] < 1, adj. P value < 0.05, n = 3 biological replicates per condition). Bottom: GSEA analysis results comparing N-Myc castrated tumors versus the other 3 groups of tumors. FC, fold change.

To further define the transcriptional differences regulated by N-Myc in an androgen-dependent manner, we performed RNA-seq on N-Myc–expressing 22Rv1 xenografts grown in castrated or intact recipients. Transcriptome-wide analyses revealed that castration was associated with a significant increase in the number of N-Myc–deregulated target genes (Figure 2C). In castrated 22Rv1 xenograft mice, the N-Myc signature was enriched with genes associated with neural lineage pathways (Figure 2C), including neural progenitor cell (NPC) bivalent genes (H3K4me2/3 active and H3K27me3 repressive marks, see refs. 44, 45; normalized enrichment score = 2.00, FDR q value = 0.004), in addition to genes implicated in neural development (e.g., SOX11, SOX21, NTRK1, and NKX2-1), expressed in adult stem cells (e.g., HOXA2/A3/A9/A10 and WNT5A), ESCs (e.g., SOX2), or NEPC (e.g., CHGA), while epithelial lineage–associated genes were downregulated (Figure 2C). Similar observations were made in the GEM model (Supplemental Figure 2E). Importantly, both 22Rv1-CTL and 22Rv1-N-Myc cells express the AR variant ARv7 (Supplemental Figure 3, A and B). Despite ARv7 expression, androgen response gene sets were significantly downregulated in response to castration, as observed by GSEA (Supplemental Figure 3C), and N-Myc expression blocked the upregulation of ARv7 target genes (46) after castration (Supplemental Figure 3D). These data suggest that the removal of circulating androgen and AR signaling may impact N-Myc gene regulation. Since these data sets were generated from tumors following chronic androgen deprivation, we cannot rule out that these molecular changes did not evolve over time.

N-Myc expression leads to neural lineage gene expression and reduced androgen response. To directly assess the impact of acute removal of androgen on the N-Myc–induced molecular program, we performed, in replicate, RNA-seq in isogenic LNCaP cells (± MYCN, see ref. 15) following short-term androgen withdrawal. LNCaP cells, a hormone-naive prostate cancer cell line with a well-characterized luminal/epithelial phenotype, were starved in androgen-deprived media for 72 hours, and subsequently re-exposed to androgen or maintained in androgen-deprived conditions for an additional 24 hours (Figure 3A). RNA-seq data revealed a distinct N-Myc–driven transcriptional program with 8,585 genes differentially regulated in N-Myc cells compared with control (CTL) cells (Supplemental Figure 4A). Consistent with our observations in vivo, androgen withdrawal significantly altered the N-Myc signature in vitro (13% and 42% of the N-Myc target genes are differentially regulated specifically in the presence or absence of androgen, respectively; Supplemental Figure 4A). While, in the presence of androgen, AR signaling was enriched in LNCaP-CTL cells compared with N-Myc cells as expected (ref. 15 and Figure 3B), N-Myc–upregulated genes were enriched with stem cell signatures and markers of neural lineage differentiation in the absence of androgen (Figure 3B). These changes were shown to be N-Myc dependent, as shRNA-mediated knockdown of N-Myc reversed their expression (e.g., AUTS2 and NKX2-1; Supplemental Figure 4B). To determine if the transcriptional changes were stable, we maintained cells in an androgen-deprived condition for 41 days followed by a 24-hour androgen stimulation (Figure 3A). After long-term withdrawal, the response of AR target genes to androgen stimulation was dramatically reduced (over 80% for KLK4, Figure 3C). Moreover, expression of genes associated with the neural lineage was increased significantly in LNCaP-N-Myc cells (Figure 3, C and D, and Supplemental Figure 4, C and D) and N-Myc expression led to increased enrichment on day 42 compared with day 4 for neural lineage and stem cell gene sets, including an adult stem cell signature associated with small-cell neuroendocrine cancers from multiple epithelial tissues (47). N-Myc expression was also associated with reduced global responsiveness to androgen restimulation, with fewer genes being differentially expressed on day 42 in LNCaP-N-Myc cells compared with LNCaP-CTL cells (Figure 3E). Altogether, these data support the hypothesis that there is interplay between AR and N-Myc signaling that propagates lineage plasticity.

Figure 3

N-Myc expression leads to neural lineage gene expression and reduced androgen response. (A) Experimental schematic with LNCaP-N-Myc or CTL cells and corresponding time points for RNA-seq or ChIP-seq analyses (arrows) in the presence (+A, green) or absence (–A, red) of androgen. (B) Enrichment plots for the Androgen Response and the Neural Stem Cell Differentiation Pathways and Lineage-specific Markers gene sets from indicated conditions. (C) Gene expression (fragments per kilobase of transcript per million mapped reads, FPKM) of AR target genes (left) and neural lineage–associated genes (right) measured by RNA-seq in the indicated cells and conditions, on day 4 (D4, n = 3 biological replicates) and day 42 (D42, n = 2 biological replicates) of androgen withdrawal. *P < 0.05, **P < 0.01, Sidak-Bonferroni–adjusted 2-tailed t test. NS, not significant. (D) Targeted GSEA of RNA-seq data from LNCaP-N-Myc versus CTL cells, without androgen, at D4 or D42 as indicated. *FDR q value < 0.05, **FDR q value < 0.01, ***FDR q value < 0.001. (E) Number of genes differentially expressed (adj. P value < 0.05) in the indicated conditions. ASC, adult stem cell; qNSC, quiescent NSC.

The N-Myc cistrome is distinct from C-Myc and altered by AR signaling. As a transcription factor, N-Myc regulates gene expression by binding to DNA and modulating transcriptional activity. To fully describe the genomic loci where N-Myc is bound, we performed N-Myc ChIP-seq in LNCaP-N-Myc cells or isogenic LNCaP-CTL cells following acute androgen withdrawal (Figure 3A). We observed N-Myc binding as a narrow peak mainly within 2 kilobases (kb) of gene transcription start sites (TSSs) (Figure 4A) irrespective of androgen. Approximately 40%, 25%, and 25% of N-Myc peaks were within intronic, intergenic, and promoter/TSS regions, respectively. To define N-Myc–specific binding in prostate cancer cells, we compared the N-Myc cistrome (LNCaP-N-Myc cells with androgen) with the C-Myc cistrome (LNCaP cells with androgen) and with publicly available N-Myc ChIP-seq data from a MYCN-amplified neuroblastoma model (48). Approximately half of the N-Myc–bound sites were shared with C-Myc in LNCaP cells and/or N-Myc in neuroblastoma cells (Figure 4, B and C). The genes uniquely bound by C-Myc in LNCaP cells were enriched with cell cycle– and cell proliferation–related genes, corresponding to the most well-characterized functions of Myc family proteins (36), while the genes uniquely bound by N-Myc in neuroblastoma cells were neural lineage related. Despite these overlaps, we observed that 44% of the N-Myc–bound sites are specific to the prostate cancer cells and that these prostate-specific N-Myc–bound genes are enriched for neural lineage genes (Figure 4D). Moreover, while the C-Myc–specific and N-Myc/C-Myc–shared bound sites were enriched for E-box motifs, the N-Myc–specific sites in LNCaP cells were significantly enriched for forkhead motifs (Supplemental Figure 5A). The difference between N-Myc and C-Myc binding in prostate cancer cells is also in accordance with the different profiles of expression of MYC and MYCN mRNA in PCa (high MYC, low MYCN) versus NEPC (low MYC, high MYCN; Supplemental Figure 5B). Intriguingly, N-Myc binds near the TSS of the MYC locus and may help to enforce the selective pattern of N-Myc/C-Myc expression (Supplemental Figure 5C). Together these results suggest that, while Myc-family members share a number of common targets, they are not completely functionally redundant and have specific targets that are regulated in a cell-type-dependent manner.

Figure 4

The N-Myc cistrome is distinct from C-Myc and is altered by androgen signaling. (A) Left: Distributions and heatmaps of N-Myc ChIP-seq data generated from cells on day 4 (D4). Right: Proportion of N-Myc–bound sites at the indicated genomic annotation. (B) Overlap of Myc family member binding in prostate cancer cells (LNCaP, top) or N-Myc in prostate cancer cells versus published N-Myc in neuroblastoma (NB) cells (LNCaP/BE2C, bottom). (C) Examples of ChIP-seq tracks for indicated genes. (D) GSEA performed on the uniquely bound genes as identified in B. (E) Representation of N-Myc binding sites determined by ChIP-seq in LNCaP-N-Myc cells with and without androgen (left) or in 22Rv1 xenografts grown in castrated or intact recipients (right), and their distribution throughout the genome (n = 2 biological replicates per condition). com., common; enr., enriched; uni., unique.

To determine if the androgen-dependent changes in gene regulation could be explained by a change in N-Myc binding to chromatin, we compared the binding profiles obtained in the two conditions. We found that 42% of N-Myc binding was dynamic, changing in different contexts of AR signaling (Figure 4E). Upon androgen withdrawal, N-Myc binding was depleted or reduced in 36% of total bound sites, while it was enriched or newly bound to 6% (Figure 4E and Supplemental Figure 5D). These dynamic N-Myc peaks were less enriched at promoter/TSS regions and more enriched at intergenic and intronic regions (Figure 4E). The distribution of dynamic N-Myc binding was confirmed by performing ChIP-seq in N-Myc–expressing 22Rv1 xenografts from castrated or intact recipients (Figure 4E). We also found that the peaks shared with and without androgen were closer to the TSS than the peaks specific to one condition (Supplemental Figure 5E) and were enriched at genes involved in normal cellular homeostasis (Supplemental Figure 5F).

N-Myc interacts with known AR cofactors to alter DNA binding. Differential interactions between N-Myc and transcriptional cofactors or pioneering factors could explain the observed dynamic binding genome-wide. It has been shown that less than half of Myc binding sites on chromatin have consensus E-box motifs and many do not even show variant E-boxes (49, 50), suggesting that the E-box enhances chromatin association but is not a prerequisite for binding. Most likely, therefore, protein-protein interactions play a large role in targeting N-Myc to sites on chromatin. To reveal such cofactors, we performed a de novo motif discovery in 200-bp regions surrounding all N-Myc binding sites identified by ChIP-seq. Motif analysis of the unique and enriched N-Myc peaks revealed that, in addition to E-boxes, there was a significant enrichment of motifs consistent with forkhead box (FOX) and homeobox (HOX) family members (Figure 5A and Supplemental Figure 6A). Since FOXA1 and HOXB13 are critical AR coregulators (51–55), we performed HOXB13 and FOXA1 ChIP-seq in the LNCaP-N-Myc cells. Unlike for HOXB13, we found that androgen withdrawal dramatically altered the FOXA1 cistrome (31% unique peaks, Figure 5B). With respect to the unique N-Myc peaks in the absence of androgen, there was very little overlap (<20%) with the C-Myc, N-Myc in neuroblastoma, or AR cistromes, compared with the N-Myc common peaks (Figure 5C). In contrast, N-Myc unique peaks without androgen shared more binding sites with AR cofactor cistromes defined in the presence of androgen in N-Myc cells versus CTL cells (from 39% to 51% for FOXA1 and from 12% to 29% for HOXB13). This overlap was further increased for FOXA1 following androgen withdrawal (+10%, Figure 5C). Taken together, our data show that 57% of all enriched and unique N-Myc peaks without androgen were co-bound with AR, HOXB13, or FOXA1 either alone or in combination (Figure 5D). In agreement with the motif analysis, comparison of the genomic loci bound by N-Myc after androgen withdrawal revealed nearby binding of FOXA1 and/or HOXB13 (Figure 5E). Furthermore, while epithelial lineage genes were enriched at these co-bound sites as expected, there was also a significant enrichment of neural lineage gene sets (Supplemental Figure 6B). To directly assess the ability of a prostate lineage–defining factor, such as FOXA1, to modulate N-Myc binding at neural lineage genes, we performed siRNA-mediated knockdown of FOXA1 and performed N-Myc ChIP-qPCR at binding sites identified by ChIP-seq (Figure 5F). Despite a knockdown of only 50%, we observed a significant decrease in N-Myc binding on the target neural genes NKX2-1 and CHGA, suggesting a role for FOXA1 in regulating N-Myc binding (Figure 5F).

Figure 5

N-Myc interacts with known AR cofactors to alter DNA binding. (A) Motif analysis of unique N-Myc peaks with or without androgen obtained by ChIP-seq in LNCaP-N-Myc cells. Scores correspond to log₂(% target/% background). All motifs shown are enriched with a P value < 10^–5 and are listed with their best predicted match to a known protein family. *Percentage background of MEF2A motif was 0% and was subsequently adjusted to 0.001% to calculate a score. (B) Overlap between FOXA1 or HOXB13 ChIP-seq peaks in the LNCaP-N-Myc cells in the presence or absence of androgen on day 4. (C) Comparison of N-Myc binding with AR (GSE69045), FOXA1 (CTL cells: GSE69045) and HOXB13 binding with or without androgen, N-Myc in BE2C neuroblastoma cells (GSE80151), and C-Myc in LNCaP. Numbers represent the percentage of N-Myc peaks in each condition overlapping with the indicated cofactor. (D) Overlap of N-Myc peaks (enriched and unique in –A) with AR, HOXB13, or FOXA1 peaks in the indicated conditions. (E) Distribution of FOXA1 and HOXB13 binding at N-Myc–bound sites ± 4 kb in the indicated conditions. (F) Top: ChIP-seq tracks of genes cobound by N-Myc and FOXA1, independently of AR (in CTL cells), in the indicated conditions. Bottom: Effect of FOXA1 knockdown by siRNA (see Western blot inset) on N-Myc binding assessed by ChIP-qPCR. **P < 0.01 by Sidak-Bonferroni–adjusted 2-tailed t test. (G) Scatter plot of log₂(fold change of N-Myc–bound peptides versus IgG-bound peptides, identified by RIME) with (x axis) and without (y axis) androgen (n = 4 biological replicates per condition). Lines correspond to the regression line ± 1.7 Z. com., common; enr., enriched; uni., unique; FC, fold change; NB, neuroblastoma.

To identify other proteins that could regulate N-Myc binding, we performed rapid immunoprecipitation mass spectrometry of endogenous proteins (RIME) (56) in LNCaP-N-Myc cells. While a subset of common interactors demonstrated increased affinity between the conditions, no interactions specific to one condition versus another were observed (Figure 5G). The majority of interactions were shared between conditions and included well-known N-Myc–interacting proteins (e.g., MAX and TRRAP; see refs. 57–59), proteins associated with heterochromatin (e.g., chromobox homologs [CBX] 1, 3, and 5), as well as HOXB13 (Supplemental Table 2 and Figure 5G). These data suggest that, while the ability to physically interact with N-Myc is not dramatically affected by the presence or acute withdrawal of androgen, N-Myc cofactors may direct genomic binding to sites in chromatin that are accessible specifically in one condition either through competition at N-Myc binding sites or by altering the chromatin accessibility.

N-Myc promotes bivalency on neural lineage genes. On a global scale, gene transcription is also regulated by epigenetic modifications of histone proteins. Two of the major histone marks used to study the regulation of gene expression are H3K4me3 (associated with transcriptional activation) and H3K27me3 (associated with repression). Paradoxically, some genomic loci contain both marks (i.e., bivalent). Bivalently marked genes are typically transcriptionally poised, have been well characterized in ESCs and shown to be essential for development and lineage determination (60, 61). The bivalent mark changes during the differentiation of the cells, through the action of histone modifiers, and biases gene expression towards activation (H3K4me3 only) or repression (H3K27me3 only), depending on the identity of each cell. To interrogate a putative role of N-Myc in regulating bivalency, we assessed the overlap between H3K4me3 and H3K27me3 ChIP-seq profiles, following short-term or long-term androgen withdrawal (Figure 3A). Interestingly, the distribution and density of H3K27me3 histone marks near H3K4me3 marks varied dramatically, depending on N-Myc expression and androgen stimulation. While the H3K27me3 binding profile was narrow and centered on H3K4me3 peaks in LNCaP-CTL cells, the binding was redirected within a 3-kb range around the H3K4me3 peaks in LNCaP-N-Myc cells, yielding a multimodal distribution (Figure 6A). We also observed similar numbers of genomic loci marked by H3K4me3 alone, while only N-Myc–overexpressing cells showed a significant increase in H3K27me3 (8.6-fold) and bivalent peaks (18-fold, Figure 6B) after androgen withdrawal. Bivalent genes in N-Myc cells were enriched for polycomb repressor complex 2 (PRC2), neurogenesis, and neural lineage pathways, while bivalent genes in control cells were associated with gene regulation and stress response (Figure 6B). Integrating the RNA-seq data with H3K4me3/H3K27me3 ChIP-seq data revealed, as expected, low levels of gene expression associated with bivalent histone marks (Supplemental Figure 7A). Interestingly, we found an increase of nearly 2-fold in the number of N-Myc peaks that are bivalent in the absence versus presence of androgen, corresponding to 966 and 580 genes, respectively (Figure 6C and Supplemental Table 3). Moreover, androgen withdrawal dramatically increased the level of enrichment of PRC2 and neural-associated genes (Supplemental Figure 7B). Among the N-Myc–bound, bivalent genes, a subset showed a decrease of H3K4me3 and an increase of H3K27me3 levels in the absence of androgen. An example of this was observed for the desmocollin 3 gene (DSC3), which has been implicated in epithelial cell junctions (62). Inversely, genes such as NKX2-1, SOX2, SOX11, and SOX21 became bivalent by gaining H3K4me3 marks in the N-Myc cells, suggesting the activation of gene expression (Figure 6D and Supplemental Figure 7C). Binding at these regulated loci appeared to be specific for N-Myc in LNCaP cells compared with C-Myc in LNCaP and/or N-Myc in neuroblastoma (Supplemental Figure 7D). Similar results were obtained by performing N-Myc/H3K4me3/H3K27me3 ChIP-seq in an NEPC patient–derived organoid (PM154, Figure 6D and Supplemental Figure 7C). We thus hypothesized that the bivalent genes linked to neuronal identity would maintain H3K4me3 and be upregulated following chronic androgen withdrawal, while the genes associated with an epithelial lineage would maintain H3K27me3 and become downregulated over time. To assess this, we analyzed the expression levels of the 966 bivalently marked, N-Myc–bound genes by RNA-seq following acute or chronic androgen withdrawal. Genes that were identified as bivalently marked and N-Myc bound on day 4, were enriched in LNCaP-N-Myc versus LNCaP-CTL cells, on day 42 (Figure 6E). Unsupervised clustering of the differentially expressed genes in LNCaP N-Myc versus LNCaP-CTL cells on day 42 revealed distinct transcriptional modules (Figure 6F). Cluster 1 included repressed epithelium-related or AR target genes (ZBTB7A, ALDH1A3), consistent with an increase of the repressive mark H3K27me3 on AR target genes (Supplemental Figure 4D) and with previous findings (15). Conversely, cluster 4 contained many of the neural lineage genes, which are upregulated by N-Myc following long-term androgen withdrawal (Figure 6F). Furthermore, we confirmed that 57% of the 966 bivalently marked, N-Myc–bound genes on day 4 became marked only with H3K4me3 on day 42, with 30% remaining bivalently marked (Supplemental Figure 8, A and B), and the genes that became H3K4me3 marked only were significantly enriched for neural lineage gene sets (Supplemental Figure 8C).

Figure 6

N-Myc promotes bivalency on neural lineage genes. (A) H3K27me3 binding profiles within 8 kb centered at H3K4me3 peaks in LNCaP CTL and N-Myc cells, with and without androgen as specified. (B) Left: Number of H3K4me3, H3K27me3, or H3K4me3/H3K27me3 bivalent peaks in common (com.) or unique to the conditions on day 4 (D4) as indicated. Right: Top 5 gene sets from GSEA for uniquely bivalently marked genes in the absence of androgen from LNCaP-N-Myc cells (red) or CTL cells (black). (C) Left: Number of H3K4me3, H3K27me3, and bivalent peaks also bound by N-Myc with or without androgen. Right: H3K4me3, H3K27me3, and N-Myc binding profiles on bivalent peaks within 8 kb centered at H3K4me3 peaks, in LNCaP-N-Myc cells in the absence of androgen. (D) Examples of N-Myc and histone mark ChIP-seq tracks as indicated. (E) Enrichment plot of the bivalent genes identified in LNCaP-N-Myc –A cells on D4 measured in LNCaP-N-Myc –A cells versus LNCaP-CTL –A cells on D42. (F) Unsupervised clustering of the genes that were bivalent on D4 and differentially expressed (adj. P < 0.05) between LNCaP-N-Myc and LNCaP-CTL cells without androgen on D42.

To determine the clinical relevance of these findings, we queried the RNA-seq from our patient cohort and found that the 966 bivalent genes identified in androgen-deprived LNCaP-N-Myc cells successfully classified the PCa, CRPC-Adeno, and NEPC patient samples into groups based on NEPC score (Figure 7A and Supplemental Figure 9A) and were enriched in the NEPC samples (Figure 7B). The NEPC clinical samples were enriched for neural lineage–associated bivalent gene sets compared with PCa or benign samples (Supplemental Table 4). These findings also correlated with N-Myc expression level in both NEPC and CRPC-Adeno samples (Figure 7B). We observed a stepwise upregulation of the N-Myc–bound and bivalent neural-associated genes during the progression from PCa to NEPC, with upregulation of genes such as SOX11 and NKX2-1 and a corresponding downregulation of epithelial lineage–linked genes such as FGFRL1 and DSC3 (Figure 7C and Supplemental Figure 9B). This suggests that these genes associated with bivalent marks in N-Myc–expressing cells play a critical role in prostate cancer evolution and may help to identify patients who are most likely to develop NEPC. Lastly, we performed N-Myc ChIP-seq in the patient-derived NEPC organoid PM154 and found that, despite fewer numbers of total peaks, 79% of the peaks were in common with LNCaP-N-Myc and 22Rv1-N-Myc xenografts, including neural lineage–associated genes (NREP, ULK2, and RAB39A), bivalent genes (RAB39A), and the H3K27me3 writer EZH2 (Supplemental Figure 9, C and D). Moreover, we confirmed the N-Myc–dependent upregulation of these bivalent genes in 22Rv1 xenografts and NEPC patient–derived organoid models (Supplemental Figure 9E).

Figure 7

N-Myc–induced bivalent genes are clinically relevant. (A) Left: Unsupervised clustering of PCa (n = 66), CRPC (n = 73), and NEPC (n = 36) patient samples based on the expression level of the 966 bivalent and N-Myc–bound genes in LNCaP-N-Myc cells on day 4 (D4) without androgen. Right: NEPC score for each cluster group. Graph depicts the median value between the 25th and 75th percentiles, with whiskers indicating the range within 1.5 IQR, Student’s unpaired 2-tailed t test. *P < 1 × 10^–3, **P < 1 × 10^–5, ***P < 1 × 10^–10. (B) Targeted GSEA of bivalent-related gene sets in the 5 NEPC samples with the highest (N-Myc^hi) or lowest (N-Myc^lo) level of MYCN expression versus PCa (n = 66) patient samples, on the 5 CRPC with the highest level of MYCN expression versus the 5 lowest, and on PM154 cells treated with an EZH2 inhibitor versus vehicle. *FDR q value < 0.05, **FDR q value < 0.01, ***FDR q value < 0.001. (C) Heatmap of log₂(fold change) of genes in NEPC (n = 36) versus PCa (n = 66) or NEPC versus CRPC (n = 73) patient samples. Illustrated genes are bivalent and bound by N-Myc in LNCaP-N-Myc cells without androgen at D4. (D) Fold change expression of the indicated genes based on qRT-PCR data (n = 3 technical replicates) in PM154 following EZH2 knockdown (see Western blot inset). ***P < 0.001 by Sidak-Bonferroni–adjusted 2-tailed t test. FC, fold change.

We have previously shown that, in the context of active AR signaling, N-Myc cooperates with EZH2 and redirects its activity to downregulate AR target genes (15). Using an in situ proximity ligation assay (PLA) to monitor EZH2 and N-Myc complex formation, we confirmed the presence of the EZH2/N-Myc complexes in LNCaP and 22Rv1-N-Myc cells with and without active AR signaling and in an AR-negative, NEPC patient–derived organoid (PM154, Supplemental Figure 9F). This suggests that N-Myc and EZH2 maintain a protein-protein interaction in the absence of AR that may regulate the H3K27me3 status of the bivalent genes. To address this, we performed shRNA-mediated knockdown of EZH2 in PM154 cells. After validation of genes that have been previously shown to be regulated by EZH2, we observed that knockdown of EZH2 led to a dramatic upregulation of the bivalent genes that were previously downregulated by N-Myc and inversely (Figure 7D). In addition to shRNA strategies, PM154 cells were also treated with a pharmacologic EZH2 inhibitor (GSK503). Targeted GSEA demonstrated a significant de-enrichment of multiple neural lineage–associated bivalent gene sets (Figure 7B and Supplemental Figure 9G).

Despite recent advances in the development of highly effective AR-directed therapies for the treatment of prostate cancer, acquired resistance ultimately ensues. Since AR-indifferent prostate cancer is thought to arise and evolve in the face of hormonal therapy, there is concern that the development of more potent AR-targeted treatments will increase the incidence of this lethal form of prostate cancer. Lineage plasticity as a mechanism of resistance to targeted therapies has been implicated in both epithelial tumors and leukemias, suggesting that significant findings may successfully translate to additional tumor types (1, 3, 7, 9, 12, 63–65). A recently published study using a machine-learning approach has identified a link between stemness and what is termed oncogenic dedifferentiation (66).

Here, we present data that highlight how an oncogenic transcription factor, N-Myc, can promote the acquisition of alternative lineage identities, including stem-like precursor states of the neural lineage, and explain the progression of dedifferentiation steps required during disease progression. We show that N-Myc functions in the epigenomic and transcriptomic reprogramming of prostate cancer epithelial cells in an androgen-dependent manner. Data from in vivo models revealed that N-Myc–driven prostate tumors are more aggressive following castration, including the development of metastatic lesions and the aberrant acquisition of alternative lineage programs. Whether or not these changes involve a transient, undifferentiated stem-like state (7, 9, 47) or a direct transdifferentiation (12) and the role of the local microenvironment in dictating the choice of lineage remains to be elucidated.

The N-Myc–induced transcriptomic landscape also changed dramatically in vivo following castration. Integrative ChIP-seq/RNA-seq data revealed a rapid, androgen-dependent change in the N-Myc cistrome, transcriptome, and histone methylation. In the presence of androgen, N-Myc binds to regulatory sequences associated with AR binding, upstream of AR target genes, to suppress gene expression. Consistent with the overlap between AR and N-Myc binding, FOX and HOX consensus motifs were significantly enriched in chromatin regions directly bound by N-Myc. Moreover, FOXA1 binding was associated with androgen-dependent changes in the N-Myc binding at genes associated with the neural lineage and was essential for N-Myc binding at the loci tested, including NKX2-1, which has been implicated in neuroendocrine-like small-cell lung cancers (67). Proteomic data confirmed an interaction between chromatin-bound N-Myc and the AR cofactor HOXB13, as well as identified interactions with the heterochromatin-associated proteins CBX1/3/5. Recent studies have implicated CBX3 in neural lineage specification, which could cooperate with N-Myc to facilitate lineage plasticity (68). In the absence of androgen, the N-Myc cistrome is redirected towards promoters of genes expressed in neural precursors/NSCs and is associated with transcriptional activation and reprogramming of epigenetic H3K4me3/H3K27me3 bivalent marks at neural lineage–associated genes. Finally, using proteomic approaches, we identified a number of previously uncharacterized N-Myc–interacting proteins that will provide potential targets to inform future studies of Myc-driven cancers.

Another important finding from this study came from our direct cistromic comparison of N-Myc and its related family member C-Myc in prostate cancer or N-Myc in neuroblastoma. The C-Myc locus is commonly amplified in early prostate cancer development (69) and we show that its expression decreases during disease progression. This may be partially explained by the deposition of H3K27me3 marks in the gene body following long-term androgen withdrawal (Supplemental Figure 5C). Interestingly, a recent report suggested that, in a specific genetic context (Rb1/Trp53 loss; AKT; BCL2 induction), C-Myc expression in prostate basal cells leads to the development of NEPC (67). We found that in LNCaP cells (RB1 and TP53 wild type), 49% of N-Myc binding does not overlap with C-Myc binding. We also found that these N-Myc–specific sites are enriched for genes associated with the neural lineage. Altogether, these results suggest that N-Myc and C-Myc share common functions but N-Myc may also regulate molecular programs that are not driven by C-Myc in prostate cancer cells. This is consistent with other data from prostate cancer (69–72) as well as several other solid tumors, such as lung (73–79) and pancreas (26, 80, 81).

Our studies have identified a potential mechanism by which N-Myc overexpression and its subsequent DNA binding induce epigenomic and transcriptomic reprogramming, resulting in a castration-resistant, lineage-plastic phenotype that gives rise to NEPC. These data provide insights into the early events of lineage switching as a mechanism of escape from the effects of hormone therapy in prostate cancer. More intriguingly, the changes that occur in the face of AR-targeted therapy may also provide a molecular signature to classify prostate cancer patients and potentially help predict those approximately 20% of CRPC patients that may eventually develop NEPC (82, 83). Clinical trials have been developed for patients with NEPC, targeting Aurora kinase, Notch, and PD-1 pathways (NCT01799278, NCT02709889, and NCT03179410). Additionally, EZH2 inhibitors have entered early-stage clinical trials in patients with a wide range of tumor types (NCT03525795 and NCT02860286), including prostate (NCT03480646). Identifying the patients most at risk of developing NEPC could provide opportunities for earlier clinical intervention and allow for inclusion in future clinical trials of these novel targeted therapies in NEPC.

Cell culture

Human LNCaP and 22Rv1 male prostate cancer cells lines were obtained from ATCC and maintained as recommended. LNCaP cells overexpressing N-Myc or an empty vector control were cultured in RPMI 1640 medium (Gibco, 11875-093) supplemented with 10% FBS (Gemini, 900-108) and 1% penicillin/streptomycin (Gibco, 15140-122). For androgen-deprivation treatment, the cells were cultured in phenol red–free RPMI medium (Gibco, 11835-030) supplemented with 5% charcoal-stripped serum (Gibco, A33821-01) and 1% penicillin/streptomycin. The medium was then supplemented with 10 nM DHT (5α-dihydrotestosterone, Sigma-Aldrich, D-073) or the equivalent volume of ethanol for 24 hours. LNCaP cells maintained more than 7 days in androgen-deprived media detached from the flask and were then kept in culture as cells in suspension. PM154 and PM155 cells were cultured as previously described (84). PM154 were treated with GSK503 (5 μM) or vehicle for 6 days before RNA collection used for RNA-seq. CD133⁺EGFR⁺ NSCs were isolated from the subventricular zone of P1 neonatal mouse brain and cultured in N2 media with 20 ng/mL EGF (PeproTech, AF-100-15) and bFGF (PeproTech, 100-18B).

Immunoblot analysis

Protein lysates were collected in RIPA buffer (Thermo Fisher Scientific, 89901) supplemented with protease inhibitor cocktail and phosphatase inhibitors (Thermo Fisher Scientific, 78428/78430). Each protein sample was resolved by SDS–PAGE (Bio-Rad, 4561084), transferred onto a PVDF membrane (Thermo Fisher Scientific, IB24002), and incubated overnight at 4°C with primary antibodies against N-Myc (Santa Cruz Biotechnology, sc-53993), FOXA1 (Abcam, ab23738), EZH2 (Active Motif, 39901), AR (Millipore, 06-680), ARv7 (RevMab, 31-1109-00), and GAPDH (Millipore, AB2302). The membrane was then incubated for 1 hour at room temperature with horseradish peroxidase–conjugated secondary antibody (Cell Signaling Technology, 7074/7076; Abcam, ab97135) and immune complexes were visualized by enhanced chemiluminescence detection (Millipore, WBLUF0500).

Quantitative PCR

RNA extraction was performed using the NucleoSpin RNA Plus extraction kit (Macherey-Nagel, 740984) following the manufacturer’s recommendation. qRT-PCR was performed using a Power SYBR Green Kit (Applied Biosystems, 4389986) following the manufacturer’s recommendation and data were acquired with the QuantStudio 5 Real-Time PCR System (Applied Biosystems). Primer sequences can be found in Supplemental Table 5.

Clinical analysis and RNA extraction from prostate cancer samples

Patients were enrolled on an IRB-approved protocol with informed consent and data were analyzed as previously described (1, 85). For these studies, 81 patients (57 metastatic CRPC-Adeno, 24 metastatic NEPC) with outcomes data were analyzed for overall survival. Patients were divided into 2 categories according to the median value of MYCN expression. Univariate overall survival analysis was calculated from the time of initial diagnosis of metastatic disease to death from any cause. Patients still alive at time of last follow-up were censored. For transcriptomic analysis, 29 benign prostate tissues, 66 localized PCa, 73 metastatic CRPC-Adeno, and 36 metastatic NEPC were selected for further analysis. Tumor histology was confirmed in all cases by pathology review and morphologic classification (37). RNA was extracted from frozen material for RNA-seq using the Maxwell 16 LEV simplyRNA Tissue Kit (Promega, AS1280). Specimens were prepared for RNA-seq as previously described (86). RNA integrity was verified using the Agilent Bioanalyzer 2100 (Agilent Technologies). Paired-end sequencing was performed on Illumina GAII, HiSeq 2000, or HiSeq 2500 instruments.

RNA-seq analysis

Quality control of raw sequencing reads was performed using FastQC (Babraham Bioinformatics). Low-quality reads were removed using Trimmomatic (87) with a sliding window size of 4 bp and a quality threshold of 20. The resulting reads were aligned to GRCm38 or GRCh38 using STAR (88). Reads were sorted and indexed using SAMtools (89). Transcript abundance was calculated in FPKM using Cufflinks (90) and in gene counts using HTSeq (91). Differential gene expression was assessed using DESeq2 (92). Batch normalization of patient samples was done using ComBat (93) from the sva Bioconductor package(94). Rstudio (1.0.136) with R (v3.3.2) and ggplot2 (2.2.1) were used for the statistical analysis and the generation of figures depicting the expression levels in each of the classes.

ChIP

Cells were collected with trypsin (Gibco, 25300-120) and crosslinked in media containing 1% methanol-free formaldehyde (Thermo Fisher Scientific, 28908) and quenched for 8 minutes using 125 mM glycine. The cell pellets were centrifuged and washed twice in cold PBS. Each pellet was resuspended in 1 mL of lysis buffer (50 mM Tris HCl pH 8, 0.5% SDS, 10 mM EDTA with protease and phosphatase inhibitors; Thermo Fisher Scientific, 78428/78430) and lysed for 20 minutes at 4°C. The nuclei were collected by centrifugation and resuspended in a second lysis buffer (10 mM Tris HCl pH 7.5, 150 mM NaCl, 0.1% SDS, 1 mM EDTA, 1% NP-40, 1% sodium deoxycholate, with protease and phosphatase inhibitors; Thermo Fisher Scientific, 78428/78430). The protein-bound chromatin was sheared by sonication (Diagenode, Bioruptor Pico). Equal volumes of sheared chromatin were immunoprecipitated with antibodies against N-Myc (Santa Cruz Biotechnology, sc-53993), C-Myc (Santa Cruz Biotechnology, sc-764), H3K4me3 (Millipore, 17-678), H3K27me3 (Abcam, ab6002), FOXA1 (Abcam, ab23738), HOXB13 (Cell Signaling Technology, 90944), or IgG control (Santa Cruz Biotechnology, sc-2025). Following extensive washing, the immunoprecipitated sheared chromatin was eluted in a fresh buffer containing 100 mM NaHCO₃ and 1% SDS. The reverse crosslinking was performed with 45 mM Tris pH 7, 170 mM NaCl at 65°C overnight. After RNase and proteinase K treatment, the DNA fragments were purified using the NucleoSpin Gel and PCR Clean-up Kit (Macherey-Nagel, 740609) following the manufacturer’s recommendation. The libraries were generated using the Hyper Prep Kit (Kapa Biosystems, KK8502) and were assessed for quality, purity, and size using DNA High Sensitivity Bioanalyzer chips (Agilent, 5067-4626). Those passing quality control were quantified using the Illumina Library Quantification Kit (Kapa Biosystems, KK4854). The sequencing was performed by the Weill Cornell Medicine Genomics and Epigenomics Core using an Illumina HiSeq 2500 or 4000 instrument.

Public ChIP-seq data

Publicly available ChIP-seq data were obtained from the NCBI’s Gene Expression Omnibus (GEO) using the following accession numbers: AR and FOXA1 from parental LNCaP cells (GSE69045) and N-Myc from BE2C neuroblastoma cells (GSE80151).

ChIP-seq analysis

Quality control of raw sequencing reads was performed using FastQC (Babraham Bioinformatics). Low-quality reads were removed using Trimmomatic (87) with a sliding window size of 4 bp and a quality threshold of 20. The resulting reads were aligned to GRCh38 using Bowtie2 (95). PCR duplicates introduced during library creation were removed using SAMtools (89). ChIP-seq peaks were called using MACS2 (96) with default parameters and q-value threshold (0.05 for C-Myc ChIP-seq in LNCaP cells and N-Myc ChIP-seq in 22Rv1 xenografts and PM154 organoids; 0.01 for N-Myc ChIP-seq in BE2C cells and for H3K4me3, FOXA1, AR, or HOXB13 ChIP-seq in LNCaP cells; 0.0001 for N-Myc ChIP-seq in LNCaP cells), or with the broad-peak option enabled and a broad-cutoff q-value threshold (0.1 for H3K27me3 ChIP-seq in LNCaP cells). Sequencing reads from sonicated input chromatin derived from each individual cell line in each condition was used as a control for peak calling.

Downstream analysis of sequencing data

Enriched ChIP-seq peak regions were annotated to hg38 genomic features and assessed for the presence of transcription factor motifs ±100 bp from the ChIP-seq peak using HOMER (97). Regions of ChIP-seq overlap were defined using BEDTools (98). ChIP-seq enrichment profile plots and heatmaps were generated using deepTools (99). GSEA (39) was performed using gene sets included in the Molecular Signature Database.

RIME

RIME was performed as previously described (56) with modifications. Briefly, cells were collected with trypsin (Gibco, 25300-120) and crosslinked in media containing 1% methanol-free formaldehyde (Thermo Fisher Scientific, 28908) for 10 minutes at room temperature and quenched for 8 minutes using 125 mM glycine. The cell pellets were centrifuged and washed twice in cold PBS. Each pellet was resuspended in 1 mL of lysis buffer (50 mM Tris HCl pH 8, 0.5% SDS, 10 mM EDTA with protease and phosphatase inhibitors; Thermo Fisher Scientific, 78428/78430) and lysed for 20 minutes at 4°C. The nuclei were collected by centrifugation and resuspended in a second lysis buffer (10 mM Tris HCl pH 7.5, 150 mM NaCl, 0.1% SDS, 1 mM EDTA, 1% NP-40, 1% sodium deoxycholate, with protease and phosphatase inhibitors; Thermo Scientific, 78428/78430). The protein-bound chromatin was sheared by sonication for 10 minutes (Diagenode, Bioruptor Pico). Equal volumes of sheared chromatin were immunoprecipitated overnight with N-Myc antibody (Santa Cruz Biotechnology, sc-53993) or mouse IgG (Santa Cruz Biotechnology, sc-2025) bound to Protein G magnetic beads (Invitrogen, 10004D). Beads were washed 10 times in RIPA buffer and 2 times in 100 mM ammonium bicarbonate (AMBIC) solution (Alfa Aesar, A18566). On-bead enzymatic digestion was performed overnight at 37°C with 100 ng trypsin (Worthington Biochemicals, LS003740). An additional 100 ng was added for 4 hours prior to desalting.

Peptides were desalted on hand-packed C18 STAGE tip columns (100). Eluted peptides were dried down in a centrifugal evaporator, reconstituted in 5% formic acid, and analyzed by nanospray LC-MS/MS on an Orbitrap Fusion mass spectrometer (Thermo Fisher Scientific). Peptides were separated by reverse-phase HPLC on a hand-packed column (packed with 40 cm of 1.8-μm, 120-Å pores, Sepax GP-C18; Sepax Technologies) using an 85-minute gradient of 5%–27% buffer B (acetonitrile, 0.1% formic acid) at a 350 nL/minute. Peptides were detected using the Top20 method. For each cycle, one full MS scan of m/z = 375–1400 was acquired in the Orbitrap at a resolution of 120,000 at m/z with AGC target = 5 × 10⁵. Each full scan was followed by the selection of up to 20 of the most intense ions for collision-induced dissociation and MS/MS analysis in the linear ion trap. Selected ions were excluded from further analysis for 30 seconds. Ions with charge 1+ or unassigned were also rejected. Maximum ion accumulation times were 100 ms for each full MS scan and 35 ms for MS/MS scans.

MS2 spectra were searched using SEQUEST against a composite database containing the translated sequences of 20,193 reviewed human protein sequences in the Universal Protein Resource (UniProt) database (downloaded March 18, 2016) and their reverse complement, using the following parameters: a precursor mass tolerance of ±20 ppm, 1.0-Da product ion mass tolerance, tryptic digestion, up to 2 missed cleavages, static modifications of carbamidomethylation on cysteine (+57.0214), and a dynamic modification of methionine oxidation (+15.9949). Peptide spectral matches were filtered to 1% FDR using the target-decoy strategy (101) combined with linear discriminant analysis (LDA) using SEQUEST scoring parameters including Xcorr, ΔCn′, precursor mass error, and charge state (102). In-house software was used to extract the peptide peak signal-to-noise ratio and area under the curve from the precursor MS1 scans. Label-free protein abundance measurements were derived from the sum of these values for all unique peptides mapping to each protein.

GEM

Transgenic mice carrying an integrated CAG-LSL-MYCN gene at the Rosa26 (R26) locus (gift from Johannes Schulte, Charité – Universitätsmedizin, Berlin, Germany; ref. 103) were crossed with mice expressing Cre recombinase under the control of the probasin promoter and homozygously floxed Pten alleles (gift from Yu Chen, Memorial Sloan-Kettering Cancer Center, New York, New York, USA; ref. 104). Following removal of the LSL cassette by Cre, a chicken actin promoter drives human N-Myc expression in these models. All lines of mice were bred on the same mixed genetic background (C57BL6/129x1/SvJ) and have been previously described (15).

Xenografts

22Rv1-CTL or 22Rv1-N-Myc cells were engineered to express luciferase using a transposable element vector (gift from John Ohlfest, University of Minnesota Medical School, Minneapolis, Minnesota, USA) and 1 × 10⁶ cells were injected into NU/J mice (Jackson Laboratories). Mice were surgically castrated once the tumor volume reached 150 mm³. Mice were anesthetized with 4% isoflurane and the skin over the scrotum was disinfected with a 70% ethanol and Betadine solution. A 0.5-cm incision was made over the scrotum, the testes were exposed by pulling the adipose tissue, the blood vessels supplying the testis were cauterized, and the incision edges were closed with sterilized wound clips. Suture removal was performed between postoperative days 10 and 14.

Immunohistochemistry

Formalin-fixed, paraffin-embedded (FFPE) tissue sections were deparaffinized and endogenous peroxidase was inactivated. Antigen retrieval was accomplished using the Bond Epitope Retrieval Solution 1 (Leica Biosystems, AR9961) at 99°C to 100°C for 30 minutes. Following retrieval, the sections were incubated sequentially with the primary antibody for 25 minutes, post-primary for 15 minutes, and polymer for 25 minutes ending with colorimetric development with diaminobenzidine (DAB) for 10 minutes using the Bond Polymer Refine Detection Kit (Leica Biosystems, DS9800). All histological evaluations and quantifications (including hematoxylin and eosin–stained [H&E-stained] and immunohistochemistry images) were performed by a board-certified, genitourinary pathologist and follow criteria that have previously been described (105). Antibodies used were specific for keratin 5 (BioLegend, 905501), AR (Abcam, ab108341), vimentin (Santa Cruz Biotechnology, sc-7557), NCAM (Leica Biosystems, NCL-L-CD56-504), keratin 8 (Developmental Studies Hybridoma Bank, AB_531826), and S100 (Dako, Z0311).

In situ PLA

The in situ PLA was conducted using the Duolink In Situ Red Starter Kit (Sigma-Aldrich, DUO92101) as recommended by the manufacturer. Fixed and saturated cells were incubated with antibodies against N-Myc (Santa Cruz Biotechnology, sc53993) and EZH2 (Active Motif, 39901). For enzalutamide (provided by Astellas Pharma Inc. and Medivation Inc.) treatment, cells were treated with 1 μM drug for 72 hours prior to PLA.

Gene knockdown by interfering RNA

For lentiviral shRNA transduction, human organoids or LNCaP-N-Myc cells were transduced using lentiviruses containing shRNA constructs against EZH2 or MYCN, respectively. Organoids were dissociated with TrypLE (Gibco) and resuspended in organoid medium (84) containing Polybrene (Millipore) and Y27632 (Selleckchem, S1049). The dissociated organoid cells were combined with viral suspension and centrifugated at 600 g and 32°C for 60 minutes. The organoid/virus mix was then incubated at 37°C overnight. Organoid cells were subsequently collected, resuspended in 120 μL of Matrigel (Corning), and seeded in a 24-well plate. Antibiotic selection was performed using 1 μg/mL puromycin (Thermo Fisher Scientific) for 7 days. LNCaP-N-Myc cells were transduced in 6-well plates in medium containing Polybrene. Cells were selected with 2 μg/mL puromycin for 3 days.

For FOXA1 silencing, LNCaP cells were plated in T175 flasks and 24 hours later cultured in phenol red–free RPMI medium (Gibco, 11835-030) supplemented with 5% charcoal-stripped serum (Gibco, A33821-01) and 1% penicillin/streptomycin. Cells were transfected after 24 hours with siRNA FOXA1 or CTL (Supplemental Table 5). Proteins and cells were collected and fixed 48 hours later for immunoblot and ChIP-qPCR.

MYCN RNA in situ hybridization

RNA in situ hybridization was performed using the single-color chromogenic QuantiGene ViewRNA ISH Tissue Assay Kit (Thermo Fisher Scientific, QVT0012) with a human MYCN ViewRNA Type 1 probe (Thermo Fisher Scientific; assay ID, VA1-18174-VT; catalog VX-01) as previously described (15).

Data availability

The ChIP-seq and RNA-seq data generated in this study have been deposited in the NCBI’s GEO under accession numbers GSE117306, GSE117430, GSE117281, and GSE117282.

Statistics

Quantification. Quantitative data are presented as mean ± SEM unless otherwise indicated.

Statistical analysis. Statistical analysis was performed using GraphPad Prism as described in the figure legends. A P value less than 0.05 was considered significant (*P < 0.05, **P < 0.01, ***P < 0.001 unless otherwise indicated).

Study approval

Animals. Mice were maintained and all procedures were performed on male mice following protocols approved by the WCM-IACUC (protocol no. 2008-0019).

Human subjects. Male patients were enrolled on an IRB-approved protocol with informed consent (WCM IRB no. 1305013903 and 1210013164).