The long noncoding RNA landscape of neuroendocrine prostate cancer and its clinical implications

The long noncoding RNA landscape of neuroendocrine prostate cancer and its clinical implications Background: Treatment-induced neuroendocrine prostate cancer (tNEPC) is an aggressive variant of late-stage metastatic castrate-resistant prostate cancer that commonly arises through neuroendocrine transdifferentiation (NEtD). Treatment options are limited, ineffective, and, for most patients, result in death in less than a year. We previously developed a first-in-field patient-derived xenograft (PDX) model of NEtD. Longitudinal deep transcriptome profiling of this model enabled monitoring of dynamic transcriptional changes during NEtD and in the context of androgen deprivation. Long non-coding RNA (lncRNA) are implicated in cancer where they can control gene regulation. Until now, the expression of lncRNAs during NEtD and their clinical associations were unexplored. Results: We implemented a next-generation sequence analysis pipeline that can detect transcripts at low expression levels and built a genome-wide catalogue (n = 37,749) of lncRNAs. We applied this pipeline to 927 clinical samples and our high-fidelity NEtD model LTL331 and identified 821 lncRNAs in NEPC. Among these are 122 lncRNAs that robustly distinguish NEPC from prostate adenocarcinoma (AD) patient tumours. The highest expressed lncRNAs within this signature are H19, LINC00617, and SSTR5-AS1. Another 742 are Received: 28 April 2017; Revised: 15 December 2017; Accepted: 1 May 2018 The Author(s) 2018. Published by Oxford University Press. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 2 NEPC lncRNAs associated with the NEtD process and fall into four distinct patterns of expression (NEtD lncRNA Class I, II, III, and IV) in our PDX model and clinical samples. Each class has significant (z-scores >2) and unique enrichment for transcription factor binding site (TFBS) motifs in their sequences. Enriched TFBS include (1) TP53 and BRN1 in Class I, (2) ELF5, SPIC, and HOXD1 in Class II, (3) SPDEF in Class III, (4) HSF1 and FOXA1 in Class IV, and (5) TWIST1 when merging Class III with IV. Common TFBS in all NEtD lncRNA were also identified and include E2F, REST, PAX5, PAX9, and STAF. Interrogation of the top deregulated candidates (n = 100) in radical prostatectomy adenocarcinoma samples with long-term follow-up (median 18 years) revealed significant clinicopathological associations. Specifically, we identified 25 that are associated with rapid metastasis following androgen deprivation therapy (ADT). Two of these lncRNAs (SSTR5-AS1 and LINC00514) stratified patients undergoing ADT based on patient outcome. Discussion: To date, a comprehensive characterization of the dynamic landscape of lncRNAs during the NEtD process has not been performed. A temporal analysis of the PDX-based NEtD model has for the first time provided this dynamic landscape. TFBS analysis identified NEPC-related TF motifs present within the NEtD lncRNA sequences, suggesting functional roles for these lncRNAs in NEPC pathogenesis. Furthermore, select NEtD lncRNAs appear to be associated with metastasis and patients receiving ADT. Treatment-related metastasis is a clinical consequence of NEPC tumours. Top candidate lncRNAs FENDRR, H19, LINC00514, LINC00617, and SSTR5-AS1 identified in this study are implicated in the development of NEPC. We present here for the first time a genome-wide catalogue of NEtD lncRNAs that characterize the transdifferentiation process and a robust NEPC lncRNA patient expression signature. To accomplish this, we carried out the largest integrative study that applied a PDX NEtD model to clinical samples. These NEtD and NEPC lncRNAs are strong candidates for clinical biomarkers and therapeutic targets and warrant further investigation. Keywords: neuroendocrine prostate cancer; transdifferentiation; small cell carcinoma; long non-coding RNA ing data suggest drivers include splice factor SRRM4 [19–21], Introduction master neural transcription factor BRN2 [22], and forkhead box Prostate cancer (PCa) is the most common cancer affecting men A1 (FOXA1) [23]. NEPC tumours have been characterized with (1) and is the third highest cause of cancer death in developed gains in MYCN and AURKA [5]; (2) overexpression of PEG10 [24], countries globally [1]. Advances in detection and treatment for HP1α [25], N-Myc [26, 27], SOX2 [28], and SOX11 [18]; (3) down- PCa have translated to many men being successfully treated regulation of PHF8, KDM3A [29, 30], REST [31], and SPEDF [32]; by surgery and/or radiation. Concomitantly, androgen depriva- and lastly (4) disease dependency on GPX4 [33]. Discoveries such tion therapy (ADT) has resulted in significant survival gains for as these continue to define the protein-coding transcriptome men with metastatic PCa. Commonly administered therapeu- of NEPC. The process of transdifferentiation, however, is highly tics include Enzalutamide, Bicalutamide, and Abiraterone [2]. complex and likely involves multiple layers of genetic and epi- These drugs inhibit the androgen signaling axis, a growth and genetic regulation. differentiation-inducing pathway mediated by the androgen re- Dysregulation of long non-coding RNAs (lncRNAs) could pro- ceptor (AR). Despite these successes, with the steady accumula- vide an additional mechanism for the gene expression alter- tion of facilitating genomic and epigenomic aberrations, a more ations that occur during NEtD. LncRNAs are broadly defined as aggressive tumour capable of growing in castrate levels of testos- large (>200 nucleotides/nt) RNA transcripts, with the most abun- terone can develop [3], termed castration-resistant prostate can- dant subtypes classified as antisense RNAs, pseudogenes, and cer (CRPC). Three main classes of treatment resistance to AR- long intergenic noncoding RNAs (lincRNA) [34]. They are impli- targeted therapies exist, falling into two broad categories asso- cated in a variety of diseases, and their association with cancer ciated with AR signaling [4]. The majority of CRPC reactivate the progression is reported through mechanisms such as remod- AR signaling axis (AR CRPC). However, some tumour cells lever- eling of chromatin, transcriptional co-activation or repression, age their inherent plasticity and progress to an AR-negative state modulation of protein activity, post-transcriptional regulation, − − (AR CRPC), circumventing AR dependence. AR CRPC is highly or as decoy elements [35–37]. LncRNAs form an important regu- heterogeneous, but a major established aggressive subtype is latory layer in global gene expression and, as such, alterations of neuroendocrine prostate cancer (NEPC) [5]. NEPC is pathologi- lncRNAs in cancer are identified as one of the driving forces for cally and clinically similar to small cell carcinoma of the prostate tumorigenesis [38, 39], cancer progression, and metastasis [40, (SCPC), which has been defined as a distinct morphological sub- 41]. More specifically in PCa, lncRNAs have been reported to play type of PCa with neuroendocrine differentiation [6]. Xenograft critical roles at every stage, including the transformation of nor- NEPC models have shown expression of a dominant and irre- mal prostate cells to prostate intraepithelial neoplastic cells, the versible neuronal-like phenotype [7] where conventional CRPC development of localized tumours, and finally progression to ad- therapies are ineffective. Platinum-based chemotherapy is only vanced metastatic disease [42]. These roles in initiation and pro- transiently effective, resulting in poor overall survival [8]with gression are due to aberrant lncRNA expression, which changes most patients surviving ∼7months[9]. Molecular pathology the balance of protein-coding genes involved in processes such markers include expression of chromogranin A (CHGA), synap- as proliferation and apoptosis, thereby facilitating cellular trans- tophysin (SYP), neuron-specific enolase (NSE) [ 10], cell-surface formation. marker CEACAM5 [11], and negative (or low) levels of AR and AR- We recently developed a first-in-field transplantable patient- regulated genes such as PSA [7]. NEPC can arise de novo but much derived xenograft (PDX) model of NEtD: a treatment-na¨ıve ade- more commonly occurs as a consequence of ADT via an adap- nocarcinoma (LTL331) that upon host castration initially re- tive process termed neuroendocrine transdifferentiation (NEtD) gresses (LTL331–8 and 12 week) but then rapidly relapses as ter- [7, 12] and frequently metastasizes to visceral organs [13]. Pre- minally differentiated NEPC (LTL331R) [7]. In our previous study disposing aberrations for NEtD include loss of RB1 [14], TP53 [15], using this model, we demonstrated a lack of evidence for NEPC mutation of Trp53 [16], and/or PTEN inactivation [17, 18]. Emerg- cells before host castration and the conservation of genome Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Table 1: Xenograft model systems used in the study Resistance Molecular Characteristics Model System TMPRSS2- PTEN Name Name Source Phenotype AN TE EZ BI AR PSA SYP SPINK1 ERG ERG PTEN GENE p53 RB LTL313B 313 Primary PCa AD ---- + + - - + + - -/- M WT LTL313BR 313 LTL313B CRPC - ++++ + - - + + - -/- M WT LTL418 418 Primary PCa AD ---- + + - + - - + +/+ WT WT LTL418BR 418 LTL418B CRPC - + - - + + - ————————–Same as parental xenograft above————– LTL331-3 331 Primary PCa AD ---- + + - - + + - -/- M M LTL331-7 331 Primary PCa AD ---- + + - - + + - -/- M M LTL331-5-8 week 331 LTL331-5 AD ---- + - - ————————Same as parental xenograft above————— LTL331-5-12 week 331 LTL331-5 AD ---- + - - LTL331-3-R 331 LTL331-3 NEPC - + - - - - + - - + - -/- M M LTL331-7-R 331 LTL313-3-R NEPC - + - - - - + - - + - -/- M M + − AR and AR CRPC xenograft model samples and their associated molecular characteristics. AD, adenocarcinoma; AN, grows in the absence of androgen; BI, resistant to Bicalutamide; CRPC, castration-resistant prostate cancer; EZ, resistant to Enzalutamide; M, mutation present; NEPC, neuroendocrine prostate cancer; TE, grows in the absence of supplemented testosterone; WT, wild type. Ramnarine et al. 3 characteristics pre- and post-castration, strongly suggesting a phase transition or state change from adenocarcinoma to NEPC [24]. With this model we have identified protein-coding tran- scripts such as PEG10 [24], SRRM4 [19], and HP1α [25]thatare active in the phase transition and validated the discovery of BRN2 [22]. In addition to these, our model has led to the iden- tification of candidate biomarkers and therapeutic targets for NEPC, including the DEK proto-oncogene [43] and epigenetic reg- ulators CBX2 and EZH2 [44] (members of the polycomb group family of transcriptional repressors). We now report the com- prehensive characterization of lncRNAs in our NEtD model. In the current study, we used the longitudinal genomic profiling of our PDX-based NEtD model focusing on lncRNA transcripts. We hypothesized that lncRNA expression across the “time se- ries” would associate with the development of NEPC. Our objec- tive was to comprehensively characterize the dynamic lncRNA landscape of NEtD and NEPC, identify putative functional motifs within lncRNA sequences, determine the clinical relevance of lncRNA expression, and identify associated clinicopathological features. To accomplish this, we implemented a sequence anal- ysis pipeline optimized for the detection of lncRNAs, identified a clinical signature that can robustly distinguish NEPC from ade- nocarcinoma (AD) tumours, and identified four NEtD-associated lncRNA expression profiles. We also identified significant en- richment of well-known transcription factor motifs within the lncRNA sequences. Lastly, we observed that a subset of these lncRNAs is associated with rapid metastasis in treated patients and can stratify tumours based on patient outcome. We present here for the first time a comprehensive landscape of NEPC lncR- NAs and their clinical associations. Results Comprehensive catalogue of lncRNAs in NEPC To identify lncRNAs involved in NEPC, we performed next- generation polyadenylated RNA sequencing on the PDX CRPC models and patient samples. We implemented a sequence anal- ysis pipeline composed primarily of algorithms from the Tuxedo suite of analysis tools [45]. Typically, lncRNAs are expressed at low levels, so the pipeline was augmented to include windowed- adaptive quality control corrections (see Methods and Supple- mentary Figs S1 and S2) that increase the ability to detect low abundance transcripts. We applied this pipeline to all PDX (n = 10) and clinical specimens (n = 117) acquired from the Vancou- ver Prostate Centre (VPC) and Weill Cornell Medicine (WCM) (Ta- bles 1 and 2). Using a quasi de novo mapping strategy combined with amalgamating all sample transcriptome assemblies, we identified 210,999 annotated transcripts spanning 38 Ensembl transcript classes. Defined by Ensembl’s core biotypes, tran- scripts are classified as either protein-coding RNAs, lncRNAs, short ncRNAs, or pseudogenes, which totaled 102,334 (48%), 82,846 (39%), 9,803 (5%),and 16,016 (8%), respectively (Fig. 1A, pie chart 1–2). Within lncRNA, seven subclasses exist: processed transcripts (n = 31,142), retained intron (n = 28,455), lincRNA (n = 12,047), antisense (10,012), sense intronic (n = 821), sense overlapping (n = 340), and three prime overlapping ncRNA (n = 29) (Fig. 1A, pie chart 3; due to their small totals, sense in- tronic, sense overlapping, and three prime overlapping ncRNAs are labeled as “Other”). Despite pseudogenes not being included within Ensembl’s lncRNA classes (listed above), they are by def- inition considered under the umbrella of lncRNA [34]. For each of the eight lncRNA subclasses and their corre- sponding transcripts, we performed unsupervised hierarchi- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 4 NEPC lncRNAs Figure 1: Transcriptome composition and study design. A) Proportions and totals of transcripts detected using our sequence analysis pipeline. Transcripts were separated into protein coding (mRNA) or non-coding RNA (ncRNA) and as defined by Ensembl’s core biotypes as either mRNA, lncRNA, short ncRNA, or pseudo gene. Within lncRNA, there exist seven classes, including processed transcripts, retained intron, lincRNA, antisense, sense intronic, sense overlapping, and three prime overlapping ncRNA (the last three labelled as “other”). Transcript totals are denoted around each pie slice. B) The three transcript classes used in this study due to their ability to separate AD and NEPC tumours, which collectively totalled 37,749 lncRNAs. ∗The pseudogene total was the combination of eight pseudogene subclasses and collectively referred to as pseudogene here. These subclasses include processed pseudogene, unprocessed pseudogene, transcribed unprocessed pseudogene, transcribed processed pseudogene, translated processed pseudogene, polymorphic pseudogene, unitary pseudogene, and pseudogene. These lncRNAs formed the basis for all down-stream analysis and C) the studies project workflow and study design. AUC, area under the curve; GRID, GenomeDx Decipher GRID databa se; JHSM, Johns Hopkins School of Medicine; PDX, patient derived xenograft; ROC, receiver-operating characteristic. See Table 2 for cohort clinical features and compositions. cal clustering (UHC) and principle component analysis (PCA) (n = 12,047),are collectively referred to as lncRNAs here on in (n on the VPC and WCM cohorts (see Methods: Statistical analy- = 37,749 transcripts in total; Fig. 1B). It should be noted that im- sis). Five were incapable of distinctly separating NEPC and AD munoglobulin and T cell receptor genes (n = 326) were removed clinical samples due to insufficient transcript counts, incorrect from the pseudogene transcript total. We explored these lncR- transcript classification, or in general poor transcript annota- NAs in our samples through two analytical workflows (model- tion. The remaining three subclasses were capable of separating based discovery and patient-based discovery), which we later NEPC and AD (Supplementary Figs S3 and S4) and became the fo- merged for clinicopathological analysis. The outline presented cus of all downstream analysis. These three lncRNA subclasses, in this figure represents the study’s overall workflow (Fig. 1C). antisense (n = 10,012), pseudogenes (n = 15,690), and lincRNAs Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Table 2: Clinical samples used in the study Treatment Status Gleason Grade Clinical End Points Cohort +RMET -RMET Institute Name Clinical Group Total Naive NHT ADT RT -6 7 8+ BCR MET PCSM +ADT +ADT VPC VPC AD-NAIVE 56 56 0 0 0 23 0 33 VPCVPC AD-NHT 14 0140 0 0 014 VPC VPC NEPC 50150005 VPC VPC CRPC 53250113 WCM RUBIN NEPC 7 NA NA NA WCM RUBIN AD 30 2235 JHSM LOTAN AD2 17 00 12 JHSM LOTAN NEPC 16 NA NA NA GRID MCI AD 545 0 0 124 54 63 271 211 388 212 132 11 113 GRID MCII AD 232 0 0 77 24 18 117 97 124 75 34 19 24 Total 927 59 17 211 78 107 412 380 512 287 166 30 137 Patient samples and their associated clinical variables including treatment status, Gleason grading, and clinical endpoints. AD, adenocarcinoma; ADT, androgen deprivation therapy; BCR, biochemical recurrence; CRPC, castration- resistant prostate cancer; MET, adenocarcinoma metastasis; NAIVE, adenocarcinoma na¨ıvely treated; NEPC, neuroendocrine prostate cancer; NHT, adenocarcinoma with neoadjuvant treatment; PCSM, prostate cancer specific mortality; +RMET+ADT, ADT-treated rapid metastasis (<36 months) with at least 10 years of clinical follow-up; -RMET+ADT, ADT treated no metastasis with at least 10 years of clinical follow-up. a b Patient overlap exists across different sites of metastasis. Contains a subset of mixed histology tumours (see Methods for breakdown). Cells with a dash are clinical features that are unknown for cohort. Ramnarine et al. 5 A lncRNA expression signature for NEPC − + Recently it has been shown that AR and AR CRPC share sub- stantial genomic overlap yet display significant epigenetic dif- ferences [32]. Here, we hypothesized that the lncRNA transcrip- tome would similarly show unique and common expression al- + − terations between AR and AR CRPC (unexplored to date). To +/− investigate this, we used the AR CRPC xenograft models (Ta- ble 1) to identify changes occurring temporally within the same tumour pre- and post-castration. Once castrated the three AD models (LTL313, LTL418, and LTL331) progress to either AR CRPC (LTL313BR and LTL418BR) or AR CRPC/NEPC (LTL331R). This al- lowed for the identification and quantification of differentially expressed transcripts between pre- and post-CRPC. We inte- grated this data with patient tumour data having matched clini- cal information to ensure the results were clinically relevant and to remove any model-based bias. As we suspected, of all lncR- NAs altered between pre- and post-CRPC (>2 fold, P < 0.05), only 8% (n = 300) were commonly deregulated in both CPRC subtypes. The remaining transcripts (n = 2,669) displayed unique changes + − in the AR or AR CRPC subtype (Supplementary Table S2 and Supplementary Fig. S5). These data support the notion that AR and AR CRPC contain largely distinct lncRNA landscapes. LncRNA expression may be useful as additional biomark- ers beyond those currently used in the diagnosis of NEPC (i.e., CGHA, SYP, and NSE). Moreover, an lncRNA expression signa- ture would strongly support the involvement of lncRNAs in NEPC at a molecular and cellular level. These lncRNAs would be candidates for mechanisms in the activation of a develop- mental pathway and/or plasticity involving previously identified protein-coding genes (PEG10, HP1α, NMYC, SOX2, SRRM4, REST, BRN2, etc.) in NEPC/NEtD. Conversely, since some of these genes (NMYC, SOX2, BRN2, and SRRM4) are well-studied transcription or splicing factors, NEPC lncRNA could be under their regulation. To build an lncRNA expression signature for NEPC, we selected the top fifth percentile of transcripts based on standard devia- tions of expression for the VPC and WCM cohorts independently and performed UHC. All uncharacterized transcripts (i.e., RP##- ######.#, AC######.#, etc.) were removed from the analysis at this point. This produced 265 and 490 NEPC lncRNAs in the VPC and WCM cohorts, respectively. Taking the intersection of these lists and then repeating UHC generated an expression signature of 122 lncRNAs (Supplementary Table S5) that distinctly segre- gated NEPC from AD tumours (Fig. 2A and B). To assess the ro- bustness of this signature, we validated it on an external clin- ical cohort of tumours (n = 33 ; Table 2) from Johns Hopkins School of Medicine (JHSM). These tumours contained 17 AD and 16 NEPC samples and were profiled on the Human Exon array 1.0 ST platform (see Methods) compared with the sequenced dis- covery cohorts. Using the same approach (UHC), a clear sepa- ration of NEPC and AD was observed (Fig. 2C). Observing consis- tent results across different technologies, institutes, and clinical samples further strengthens the robust nature of the NEPC sig- nature. To our knowledge, this is the first report of lncRNAs ex- hibiting a unique, unbiased expression classifier capable of seg- regating NEPC and AD patient samples. Some lncRNA from the patient-derived signature have been previously reported as altered in other cancer types. These lncR- NAs include MALAT1 (alias NEAT2), PCA4 (aliases GDEP, PCAN1, or PCAT4), DSCAM-AS1, and SNHG12. MALAT1 is one of the most well-characterized and studied lncRNAs in cancer and has been identified as a regulator of metastasis and cell migration, a prognostic marker, and a transcriptional regulator of alter- native splicing in lung cancer [46]. PCA4 has been identified Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 6 NEPC lncRNAs Colour Key (A) (B) (C) VPC WCM JHSM −2 −1 0 1 2 Row Z−Score Figure 2: NEPC lncRNA expression signature and clinical classifier. Unsupervised hierarchical clustering of the 122 identified lncRNAs from (A) VPC and (B) WCM cohorts. Validation of this signature was shown in the (C) JHSM cohort. Samples (columns) are labelled as adenocarcinomas (blue) or neuroendocrine (yellow) tumours. See Supplementary Fig. S3 for row/lncRNA labels for each plot. as a prostate and retinal specific transcript [ 47] and frequently uate SSTR-targeted therapy for neuroendocrine tumours in cir- mutated in PCa [48]. DSCAM-AS1 mediates tumour progression culating tumour cells [55], and its use in patient management and tamoxifen resistance in breast cancer through interacting is being tested in a Phase IV clinical trial (NCT02075606). Over- protein hnRNPL [49]. SNHG12 is induced by c-MYC and regu- all, the identification of the NEPC lncRNA expression signature lates cell proliferation, apoptosis, and migration in triple neg- has provided a previously unexplored component of the NEPC ative breast cancer [50]. We were interested in identifying the transcriptome, revealed candidate NEPC biomarkers, and asso- most highly expressed lncRNAs in the signature. Therefore we ciations to NEPC biology. ranked each according to their fold changes when compared with AD samples, required concordance in fold changes across both of the cohorts, and >10-fold change in magnitude. H19, Distinct expression profiles of lncRNAs are associated LINC00617 (alias TUNA/TUNAR), NKX2–1-AS1, and SSTR5-AS1 with NEtD were the only four that fit these thresholds and each with pre- A major goal of this study was to characterize the lncRNA land- vious reports in cancer. Of note, H19 is the most studied among scape during the dynamic phase transition from AD to NEPC the four lncRNAs and is implicated in numerous cancer types using our unique PDX model LTL331 [7](Fig. 3A). To accom- [51]. It is involved in proliferation and both differentiation pro- plish this, we sequenced six samples of our PDX NEtD model cesses related to metastasis, epithelial-to-mesenchymal tran- representing three primary time points along disease progres- sition (EMT), and mesenchymal-to-epithelial transition (MET) sion: two samples from each terminal point (AD and NEPC) and [52]. LINC00617 in breast cancer regulates EMT, cancer progres- two samples post-castration (postTX). Time points 8- and 12- sion, and metastasis through activation of the transcription of week post-castration were selected to represent postTX due to SOX2 [53]. SSTR5-AS1 has not been functionally characterized, tumour volume and serum PSA levels reaching nadir (Table 1 but its sense form SSTR5 has and is a biomarker for neuroen- and Fig. 3A). We identified and quantified all lncRNA transcripts docrine tumours [54]. In fact, recently it has been used to eval- that were altered across the time series and defined four pat- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 AD NEPC Ramnarine et al. 7 Figure 3: Xenograft model of neuroendocrine transdifferentiation, phenotype-driven data integration, and NEtD-associated lncRNAs. A) Schematic depicting the time points at which xenograft tumours were collected along the transdifferentiation of AD to NEPC (adapted from Akamatsu et al., 2015 [24]). B) Phenotypes for clinical samples (coloured circles) that align to various time points from above xenograft model and group-wise comparisons (black connector bars) analyzed for clinical samples. C) Four isolated expression profiles (grey triangles) from select time points in A (light grey circles) with appropriate clinical group-wis e comparisons overlaid and integrated. Unsupervised hierarchical clustering with NEtD lncRNAs (Class I, Deactivated: black bars; Class II, Activated: orange bars; and Class III, Persistent: red bars) identified from integration outlined in (C). Distinct clusters of AD and NEPC clinical samples are observed in the (D) VPC and (E) WCM cohorts. Cla ss IV, Transient lncRNAs were excluded from the clustering due to the lack of clinical samples that would represent this intermediate state. terns of transcript expression: (a) continuous decline in expres- tegration, the following patient group-wise comparisons were sion (Class I, Deactivated, n = 1,613); (b) increasing expression performed: (a) NEPC vs AD, (b) NEPC vs NHT, (c) CRPC vs AD, (d) ¨ ¨ from either AD to postTX or postTX to NEPC (Class II, Activated, NHT vs NAIVE (untreated AD), and (e) NHT vs NAIVE in combi- n = 4,281); (c) continuous increased expression (Class III, Persis- nation with NEPC vs NHT. This produced 1,927;713; 975; 1,045; tent, n = 1,054); and (d) maximum expression at postTX (Class and 117 transcripts, respectively (>2 fold with P < 0.05; total n = IV, Transient, n = 2,668)(total n = 7,627;Fig. 3A). The NEtD model 3,154;Fig. 3B; Supplementary Table S4). Integrating these results and postTX state represents a biological process that currently with the PDX NEtD model transcripts above led to 475, 222, 84, is not characterized as a clinical entity but offers invaluable in- and 45 lncRNAs identified within Class I (Deactivated), Class II sight into the transcriptome of transdifferentiating AD cells. (Activated), Class III (Persistent), and Class IV (Transient), respec- To determine the clinical relevance of Class (I–IV) lncRNAs, tively (total n = 742; Fig. 3C). Unsupervised hierarchical cluster- we integrated patient samples (VPC and WCM, Table 2 - col- ing of Class I-III within WCM (Fig. 3D) and VPC (Fig. 3E) cohorts umn “Clinical group”) with time points in our model (see Fig. exhibited (as expected) a distinct separation of AD and NEPC tu- 3B for alignment of time points to patient groups). Terminal mours and a distinct separation between lncRNAs in Class I-III time points were appropriately aligned to AD and NEPC sam- (rows of heat map). Class IV transcripts were excluded from this ples. However, due to the lack of clinical specimens undergo- illustration due to their lack of altered expression between AD ing NEtD, we hypothesized that neoadjuvant hormone ther- and NEPC clinical samples. Collectively these 742 NEtD lncRNAs apy (NHT) given to AD patients may exhibit characteristics of are associated with the pathogenesis of tNEPC. the postTX state. The transcriptomes from these patients have Prominent examples identified by this biological integration been shown to display the effects of therapy response and more of our NEtD model (Fig. 4A), WCM cohort (Fig. 4B), and VPC co- specifically androgen depletion [ 56]. In fact, neuroendocrine dif- hort (Fig. 4C) illustrate each of these NEtD-defining transcript ferentiation has been shown to increase after only three months classes. PCA3, PCAT1, and PCGEM1 were selected as controls of NHT in a retrospective analysis of 103 radical prostatectomy for this study due to their elevated expression in PCa and high specimens [57]. These early events are the specific alterations we level of characterization. As expected, their expression patterns sought to isolate from the postTX time points of our PDX model. followed the trend in the PCa and NEPC samples (Fig. 4A–C, We also postulated that a subset of Class I (down-regulated in NEtD controls; P < 0.001). SOCS2-AS1 and HOXA11-AS are se- our PDX model) would be up-regulated in the (AR ) CRPC clini- lect examples that characterize the NEtD lncRNA Class I De- cal samples due to reactivation of the AR signalling axis in clas- activated (Fig. 4A–C, Deactivated; P < 0.01). HOXA11-AS, asso- sical CRPC [56, 58–60]. Based on this model-to-patient data in- ciated with the cell cycle through E2F1 [61], has been seen to Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 8 NEPC lncRNAs (A) (B) (C) Figure 4: Select NEtD lncRNAs that exemplify each expression pattern are shown from the (A) NEtD PDX LTL331 model, (B) WCM cohort, and (C) VPC cohort. The expression for NEtD lncRNAs within Class IV, Transient was only identified through the VPC cohort due to the presence of NHT samples, which were not pres ent within the WCM cohort. All boxplots showed significant separation ( P < 0.05) between groups based on a standard Student’s t-test with the exception of ∗ lncRNAs. promote gastric cancer proliferation and invasion (with EZH2) and a focus of our future research and functionalization. Taken and can act as a “molecular sponge” for EZH2 by absorbing (via together, these NEtD lncRNAs (n = 742) characterize the transd- direct interaction) miR-1297 [62]. SOCS2-AS1 is another lncRNA ifferentiation that occurs post-castration and is associated with in this class that has been identified as an AR-regulated tran- tNEPC. script [63] and further supports our hypothesis of AR-regulated lncRNAs within NEtD Class I Deactivated. NKX2–1-AS1 exempli- fies the NEtD Class II Activated (Fig. 4A–C, Activated; P < 0.05) NEtD lncRNAs are enriched with distinct transcription and has been previously seen to characterize lung cancer sub- factor binding motifs types AD and squamous [64]. CDKN2B-AS1 (alias ANRIL) and H19 LncRNAs are not translated and carry out their functions post- are prime illustrations for persistently expressed NEtD Class III transcription in their secondary or tertiary RNA form. This is (Fig. 4A–C, Persistent; P < 0.05). Both of these lncRNAs have been unlike protein-coding transcripts that function in their post- identified across a number of cancer studies (H19 [ 51], ANRIL [65, translational form. Thus, identifying sequence motifs within 66]);however, depending on the cancer type, each has functioned lncRNAs could identify interacting transcripts or proteins that as a tumour suppressor (i.e., ANRIL deactivating tumour sup- provide clues to function. Enrichment of transcription factor pressors CDKN2A/B in cis by three different epigenetic mech- (TF) binding sites (TFBS) was determined by calculating z-scores anisms [67–69]) and as an oncogene (i.e., H19 acts as a sponge for overrepresentation of motifs present in the NEtD lncRNA for FOXM1 by absorbing miR-342–3p [70]). Two demonstrations Classes (I–IV) against their genomic background (Supplementary for transiently expressed NEtD lncRNA Class IV are FENDRR Tables S6–S17 and Methods: Genomatix overrepresented TFBS). and CASC15 (Fig. 4A–C, Transient; P < 0.01). These lncRNAs are We also integrated each of these class-specific enrichment re- well studied in cancer: FENDRR for its prognostic value and its sults to identify unique TFBS for each NEtD Class (Supplemen- involvement in gastric cancer metastasis [71] and CASC15 for tary Tables S18 and S19). All TF and TFBS descriptions, family its regulation of SOX4 in RUNX1-rearranged leukemia [72]and classifications, and annotation is described in Supplementary harboring a risk SNP for susceptibility of neuroblastoma [73]. Table 25. CASC15 has also been identified as a mediator of neural growth In NEtD Class I we identified 33 significant and uniquely en- and differentiation [74], which we believe could be occurring in riched TFBS (Supplementary Tables S6, S18, and S20). Interesting our NEtD model based on the data presented here. Each of these results included binding motifs for TP53, scratch family tran- lncRNAs are among the top candidates identified in this study scriptional repressor 2 (SCRT2), and POU Class III homeobox 3 Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 9 (POU3F3) (z-scores = 4.02, 4.27, and 2.41, respectively). TP53 of- With FOXA1 as one of the characterizing TFBS in Class IV, we ten absent in NEPC, could be an activating TF for many of these sought to explore the persistently expressed (Class III) in con- Deactivated lncRNAs, and suggests an apoptosis or cell cycle ar- junction with the transiently expressed transcripts (Class IV). rest role is present here. Scratch family transcriptional repres- We hypothesized that subsets of these lncRNAs have mechanis- sor 2 has been linked as a neural-specific Snail family transcrip- tic involvement in the transdifferentiation process. To investi- tional repressor and critical for neuronal differentiation [75]. gate this we repeated the TFBS enrichment analysis on Class III Similar to REST, this TF is likely causing the down-regulation of a and IV together and identified six significantly enriched TFBS subset of these lncRNAs. Lastly, POU3F3/BRN1 (a member of the (Supplementary Tables S10, S18, and S20). Confirming our hy- POU family of TFs, as is BRN2) is involved in the development pothesis was the presence of TWIST1 (z-score = 4.01), an es- of the nervous system, expressed in small cell lung cancer cells sential member of the EMT transcriptional reprogramming fac- (which has pathological overlaps with NEPC), and involved in tors [87]. Interestingly, TWIST1 and AURKA have very recently proneural/neuroendocrine differentiation [76]. Considering this been seen to form a feedback loop promoting metastasis and and the significant enrichment of these TF motifs, this suggests highly aggressive phenotypes in pancreatic carcinoma [88], and a role in proliferation and differentiation in NEtD Class I. TWIST1 is a marker for EMT in neuroendocrine tumours [89, Performing TFBS enrichment analysis in NEtD Class II and III 90]. Concerning PCa, it has been identified as AR regulated (and identified 12 and 15 significant and distinct TFBS motifs, respec- repressed via NKX3–1), whereas in the absence of AR it is up- tively (Supplementary Tables S7, S8, S18, and S20). Interestingly, regulated and present in metastatic disease [91]. both classes had significant enrichment for at least one ETS and We further investigated global functional characteristics HOX family member, suggesting overlapping functional roles for across all NEtD lncRNAs. Specifically, we wanted to identify TFBS their respective lncRNAs. For Class II this included ELF5, SPIC, that were significantly enriched and common across all classes. and HOXD1 (z-scores = 2.17, 2.43, and 2.22, respectively) and Due to the high number of lncRNAs (n = 2,147), we decided to for Class III, PDEF (alias SPDEF) and HOX/PBX (z-scores = 3.03 perform this analysis at the TF family level; therefore for each and 2.71, respectively). Members of the ETS family fused to TM- class and the full lncRNA set, we repeated the motif enrichment PRSS2 is the most frequent genomic alteration in PCa; therefore analysis and integrated all of their results (Supplementary Ta- the prevalence of their motifs in these classes is not surprising. bles S12–S17, S19, and S20). We identified 62 significantly com- While the ETS fusion transcript is relatively more specific to PCa mon TFBS families (z-score = >2;Supplementary Table S20). Not vs NEPC, ETS TFs on their own are involved in a wide variety surprising were families involving cell cycle regulation, cyclin of functions, including cellular differentiation and angiogene- B2/CCNB2 and the E2F family (z-scores = 40.66 and 67.36, re- sis. In fact, recently SPDEF was found to be down-regulated in spectively). We also observed both the ETS (z-score = 9.9) and metastatic NEPC due to DNA methylation [32] and was also sig- REST (z-score = 20.52) families of TFs, which reaffirmed our hy- nificantly down-regulated in treated vs untreated high-risk PCa pothesis that these lncRNAs are involved in tumour progres- patients [77]. Conversely, the HOX family has never been linked sion and neuronal pathways. Surprising was the presence of two to PCa or NEPC for that matter, and so this result was unex- PAX families, PAX5 (z-score = 10.15) and PAX9 (z-score = 18.18). pected. In neuroblastoma, however, the HOX genes have been The PAX family is known to regulate lineage specification and linked to differentiating cells [78] and specifically HOXD1 iden- progenitor cell maintenance. In developmental biology, PAX5 is tified here (as well as HOXC6 and HOXD8) are associated with involved in B-cell differentiation and PAX9 in neural crest de- differentiation towards a neuronal phenotype [79]. velopment. PAX5 has been observed as overexpressed in other Performing TFBS enrichment analysis in NEtD Class IV iden- neuroendocrine tumours [92, 93], overexpressed in neuroblas- tified enrichment of 17 distinct TFBS motifs (Supplementary Ta- toma [94], and shown to positively regulate c-Met transcription bles S9, S18, and S20). Class IV transcripts are only expressed in small cell lung cancer [95]. In lung NETs, PAX6 is prognostic during treatment (castration) response. Interestingly, heat shock for aggressiveness [96]. Their role in NEPC is yet to be charac- TFs HSF1 (z-score = 2.06) and HSF2 (z-score = 3.87) were within terized; however, evidence here supports their global involve- these results. Heat shock proteins (HSPs) are expressed at low ment in NEtD and lncRNA function. Lastly, the selenocysteine levels under normal conditions, up-regulated by cellular stress, tRNA activating factor (STAF, z-score = 15.18) was very intrigu- and function as molecular chaperones to control client protein ing to us. A recent Nature study by Schreiber et al. suggested stability and function. Their candidacy as therapeutic targets that treatment resistance in NEtD of PCa depends on a drug- has been well studied in PCa [80]and AR CRPC [81]. In breast gable lipid-peroxidase pathway that protects against ferropto- cancer, HSF1 specifically induces a cancer stem cell phenotype sis (a non-apoptotic form of cell death) [33]. The increased lipid in vitro [82]. In PCa, HSPs bind dihydrotestosterone to the AR metabolism creates a dependency on GPX4, which prevents fer- and enhance AR-mediated transcription. One of the functions roptosis. GPX4 is a selenocysteine-containing enzyme and 1 of of lncRNAs is to facilitate this type of mechanism. For exam- only 25 proteins with this rare amino acid in the entire human ple, LINC00152/CYTOR (identified within this class) binds and re- genome. The data presented here suggest that some of these cruits EZH2 to its target promoters p15 and p21 in gastric cancer lncRNAs may be involved in the selenocysteine pathway via [83] and IL24 in lung cancer [84] and thereby causes repression STAF and in selenoprotein biosynthesis of molecules (i.e., GPX4). of their expression. Considering the transient expression of the Identifying and targeting these lncRNAs could be a path for up- lncRNAs in this class, this data suggest a subset may be stress stream inhibition of GPX4 up-regulation and therefore allow cell response mediators via HSPs. Lastly, FOXA1 showed a significant death in these resistant cells to occur naturally by ferroptosis. enrichment (z-score = 2.8) in this class. Recently, FOXA1 loss was Comprehensive in vitro experimentation would need to be car- identified as a driver of NEtD [ 23], which leads to AR reprogram- ried out to confirm this therapeutic avenue. ming [85] and EMT through direct regulation of SLUG expression [86]. This suggests that some of the lncRNAs in this class could NEtD lncRNAs contain NEPC-related TFBS have a functional role in maintaining cellular identity when un- der the control of FOXA1. It is now well established that complex cellular reprogramming occurs during NEtD, and master regulators such as REST [31], Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 10 NEPC lncRNAs BRN2 [22], SOX2 [28], and SOX11 [18] have been identified as key AS1 regulates EMT and Notch signaling to promote colorectal TFs involved in this process. Identification of well-known TFBS cancer [110]. H19 has been identified as a mediator of breast such as these would test our current hypotheses on the func- cancer plasticity during EMT and its reverse process MET [111], tional involvement of individual lncRNAs in NEtD pathogenesis as well as having a role in stemness in prostate cells [112]. (Supplementary Tables S21–S24). TFBS identification was carried LINC00152 is involved in EMT (combined with cell migration and out using MatInspector [97–99] (see Methods: Genomatix Mat- invasion) in gastric cancer [113]. LINC00478 (alias MONC) inter- base and MatInspector) on each of the NEtD Class for select TFs. feres with hematopoietic lineage decisions and enhances pro- As before, all TF and TFBS descriptions, family classifications, liferation of immature progenitor cells in acute megakaryoblas- and annotation is described in Supplementary Table 25. tic leukemia [114]. Lastly, again in gastric cancer, lncRNA SNHG6 With the dominance of AR-regulated genes in AD, the lack of has been seen to promote cell proliferation and EMT [115]. Based expression that defines Class I NEtD lncRNAs is likely caused by on these data, Class III and IV lncRNAs could have a role in de- the absence of androgen (post-castration), and therefore these veloping a cellular “plastic” state during NEtD. lncRNAs are putatively AR-regulated. To test this, we searched To test the involvement of known NEPC-involved TFs in all for androgen and glucocorticoid response elements (ARE and NEtD lncRNA, we searched for BRN2, ARE/GRE, REST, SOX11, GRE, respectively). The results showed that 107 lncRNAs con- SOX2, NMYC, ETVI, ETS, and NKX3 motifs (Supplementary Ta- tained ARE and/or GRE motifs, of which 16 contained only an ble S24). Since each of the classes had different sizes, this would ARE motif, 49 contained only a GRE motif, and 21 contained both influence the distribution and presence/absence of these mo- ARE and GRE motifs (Supplementary Table S21). To further test tifs, so we extracted the top 25 lncRNA within each class (n = and support our AR-regulated lncRNA hypothesis, we explored 100 NEtD lncRNA), ranked by their magnitude of fold change. all previously reported AR-regulated lncRNAs. Currently, the fol- Observing the distribution of these TFs separated by NEtD class lowing 17 lncRNAs have been identified with experimental evi- revealed an interesting pattern (Fig. 5B). TFs SOX2, SOX11, and dence: PCGEM1 [100], PlncRNA-1/CBR3-AS1 [101], PCAT-18 [102], REST had a relatively more balanced distribution across each PCAT29 [103], SOCS2-AS1 [63], RP1–45I4.2 [104], SUZ12P1 [104], class compared to NKX3, ETSF, ETVI, and NMYC, which showed SNHG5 [104], LINC01138 [104], SNHG1 [104], KLKP1 [104, 105], a preference to binding persistent and transiently expressed LINC00969 [104], LINC-PINT [104], TUG1 [104], MIR17HG [104], lncRNA. Interestingly, more than 50% of ARE/GRE motifs were POTEF-AS1 [106], and CTBP1-AS1 [107]). Of these, 4 (PCAT29, present in transiently expressed lncRNA vs relatively few in SUZ12P1, SNHG1, and CTBP1-AS1) were not within our pipeline Class I Deactivated and Class II Activated. Conversely, BRN2 mo- lncRNA class annotation, and 8 of 13 (61%) were represented in tifs were relatively more present in Class I and II. These patterns NEtD Class I Deactivated lncRNAs (PCGEM1, PlincRNA-1, PCAT- suggest a time-dependent or cellular phase-dependent usage of 18, SOCS-AS1, KLKP1, LINC00969, LINC-PINT, and POTEF-AS1). TFs post-castration and during the NEtD process. Due to our integrative study design (Fig. 3A–C), the remaining 5 did not move forward in the analysis. However, removing the NEPC and NEtD lncRNAs identify putative NEPC integrative steps, down-regulation of these lncRNAs did occur subtypes in either our model or patient samples independently. Overall, of the 13 lncRNA annotated by our pipeline and reported as AR To corroborate the lncRNA expression in an external NEPC regulated, all overlapped in this study. (extNEPC) cohort [32], we visualized NEtD lncRNA Classes II-IV Due to the elevated pattern of expression that defines Class and the up-regulated NEPC lncRNA expression, including ge- II and III NEtD lncRNAs, we hypothesized that a subset of these nomic profiles (copy number and mutation), through an Onco- lncRNAs are constituents of the neuronal phenotype present in Print schematic. The cohort consisted of 44 NEPC specimens NEPC. To test this hypothesis, we analyzed these lncRNAs for (largest published to date) from 30 patients that were classi- the presence/absence of the following select TFs known to ac- fied based on their histomorphology [ 6]. Duetothe exomese- tivate this cell type: APOU Class III homeobox 2 (POU3F2), also quencing performed on this cohort, not all lncRNAs were rep- known as BRN2, and RE1 silencing transcription factor (REST). resented/detectable in this sequencing profiling. We also plot- Activation of BRN2 and deactivation of REST are involved in ted previously reported NEPC predisposing genes, oncogenes, neuronal differentiation and regulation of neurogenesis, respec- drivers, and the TFs we identify above to provide “transcriptome tively. Again, using the MatInspector algorithm, we identified 11, context” for the altered lncRNAs (Supplementary Figs S6 and 22, and 21 lncRNAs in Class II or III with TFBS for BRN2, REST, S7). From the 58 NEPC and 243 NEtD lncRNAs represented in or both, respectively (Supplementary Table S22). Taken together, the extNEPC exome sequencing profiling, 43% (25/58) and 27% this evidence supports involvement for a subset of these lncR- (66/243) showed altered expression in 2–34% of NEPC patients, NAs to neuronal function/pathways in NEtD. respectively (Supplementary Figs S6, S8–S11). To further support the hypothesis of mechanistic involve- Surprisingly, these testable lncRNAs (58 from the NEPC lncR- ment for the NEtD process in Classes III and IV, we expected TFBS NAs and 66 from the NEtD lncRNA) in combination with known related to plasticity and stemness to be present. Therefore, we oncogenes/tumour suppressors/transcription factors (Supple- used MatInspector to identify binding motifs for members of the mentary Fig. S7) resulted in identifying three distinct subsets of following well-studied cellular differentiation TF families: HOX NEPC patients within the extNEPC cohort. Group 1 had relatively [108], SOX, STAT3 [109], and “STEM” (STEM members are defined higher mutation frequencies, higher ploidy, mixed tumour sites, by Matbase and include POU5F1/OCT4, SALL4B, SOX2, NANOG, and mixed pathological classifications. Group 2 had a relatively and TCF7L1). We observed 42, 49, 30, and 33 lncRNAs with TFBS low mutation frequency and low ploidy, derived mostly from for HOX, SOX, STAT3, and STEM genes, respectively (Supplemen- pelvic masses and with pathological classification D (large-cell tary Table S23). In fact, some lncRNA had TFBS within more neuroendocrine carcinoma). Group 3 tumours, however, were than one of these TFs (Fig. 5A). Previous studies have linked 6/7 mostly derived from the prostate with a pathological classifica- of these (highlighted in Fig. 5A) to various components of EMT tion B, and likely primary (de novo) NEPC samples where NEtD and/or cellular plasticity. FENDRR (antisense lncRNA to FOXF1) has not occurred. Of note, copy number loss or mutations in regulates gastric cancer metastasis via fibronectin1 [ 71]. FOXD2- TP53 and RB1 were present in 60% of patients (26/44), spread Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 11 (A) SOX STAT3 HOX STEM FOXD2-AS1 5 1 1 0 7 1 SNHG6 TDRG1 LINC00152 H19 0 2 FENDRR LINC00478 0 10 (B) NKX3 ETSF ETVI Class I NMYC Class II SOX2 SOX11 Class III NRSF Class IV ARE+GRE BRN2 0% 20% 40% 60% 80% 100% Figure 5: TFBS Venn diagram and distribution plots. A) Common and unique TFBS for HOX, SOX, STAT3, and STEM families of transcription factors within Class III and IV of the NEtD lncRNAs. B) Distribution of TFBS for known NEPC-involved TFs within NEtD Class lncRNAs. across the cohort, and did not appear to be associated with a marily with adverse pathology (i.e., high grade/stage) and long- particular group (Supplementary Fig. S7). The three groups could term follow-up for treatment and outcomes (median 18 years). be revealing an lncRNA expression signature that is specific for From these cohorts, a subset (n = 211) received adjuvant ADT tumour site, degree of genomic mutations (SNPs or CNVs), and post-radical prostatectomy (RP). To determine the most clini- pathological classification. However, it is important to note this cally relevant lncRNA transcripts, we first ranked the NEtD/NEPC is an observational result requiring statistical validation in a lncRNAs within their respective classes and selected the top larger cohort. The specificity of these genomic and lncRNA tran- deregulated from each. The ranking was performed based on scriptome profiles would need to be explored across a variety fold changes observed within the clinical groups (see Methods). of metastatic sites and NEPC pathologies to validate these three This produced 100 top-ranking NEtD/NEPC lncRNAs that we in- novel NEPC molecular subtypes. vestigated within the GRID cohorts (Fig. 1C and Supplementary Table S26). We validated 11 of these (2 from each NEtD Class and 3 from the NEPC lncRNA signature) by quantitative real-time PCR to confirm expression changes identified in the model and clin- NEPC and NEtD lncRNAs are associated with ical samples (Supplementary Fig. S12). Due to the difference in treatment-related metastasis profiling platforms between GRID (Affymetrix microarray) and Prognostic and predictive biomarkers for NEtD and NEPC are in VPC/WCM cohorts (Illumina Sequencing), it was necessary to dire need since ADT is not effective for a cancer that has un- remap the GRID microarray probes (see Methods) that aligned dergone NEtD and thus circumvents the AR signalling axis. We within NEtD/NEPC lncRNA sequenced regions. This resulted in examined if the NEtD (n = 742) and NEPC (n = 122) lncRNAs are 81/100 being present and quantifiable on the microarray plat- associated with NEPC related clinical outcomes in patients with form. primary prostatic adenocarcinoma. To accomplish this, we ex- A characteristic of NEPC patients in the clinic is the occur- plored the candidates in two cohorts from the Mayo Clinic (MCI rence of rapid metastasis following treatment [118], and so we [116]and MCII [117]) from the Decipher GRID database (GRID) first tested the lncRNAs’ ability to predict rapid metastasis post- (n = 777, Table 2). We could not perform this analysis within ADT. We performed receiver-operating characteristic (ROC) anal- VPC/WCM cohorts due to their small sample sizes and short- ysis to compare the sensitivity and specificity of predicting rapid term clinical follow-up. The GRID cohorts represent tumours pri- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 12 NEPC lncRNAs (A) 0 (events = 15) 0 (events = 13) 1 (events = 22) 1 (events = 20) p−value: 0.005 p−value: 0.010 0 20 40 60 80 100 120 140 020 40 60 80 100 120 140 Time(months) RP−Met Time(months) RP−Met Group=0 21 20 15 15 11 9988 2 Group=0 21 21 16 16 12 10 9 7 6 Group=1 22 18 11 86 5 3 1 Group=1 22 17 10 7 5 4 3 222 SSTR5-AS1 LINC00514 (B) 0 (events = 20) 0 (events = 19) 1 (events = 20) 1 (events = 21) p−value: 0.900 p−value: 0.832 0 20 40 60 80 100 120 140 020 40 60 80 100 120 140 Time(months) RP−Met Time(months) RP−Met Group=0 28 23 20 19 14 11 10 986 Group=0 28 23 20 19 14 11 10 9 9 4 Group=1 28 24 23 17 15 12 10 983 Group=1 28 24 23 17 15 12 10 9 7 5 VPC WCMC (C) SSTR5-AS1 LINC00514 SSTR5-AS1 LINC00514 Figure 6: Kaplan-Meier estimates and expression for SSTR5-AS1 and LINC00514. Kaplan-Meier estimates for metastasis-free survival in the MCII cohort comparing low (blue lines) and high (yellow lines) expression (split by median) in (A) treated patients that received post-prostatectomy adjuvant ADT for SSTR5-AS1 (left) and LINC00514 (right) and (B) patients not receiving ADT treatment. C) Box plot expression for the top two NEPC lncRNA candidates (SSTR5-AS1 and LINC00514) within the VPC and WCM cohorts. metastasis (within 36 months) for each lncRNA. We then cal- The expression of two NEtD/NEPC lncRNA transcripts (SSTR5- culated the area under the curve (AUC) for each lncRNA ROC AS1 and LINC00514) was able to separate patients more likely in both cohorts using probe set region expression summarized to develop metastatic disease from those that did not (P = 0.005 across the full lncRNA transcript (Supplementary Table S26). and P = 0.010, respectively; Fig. 6A). To increase our confidence This identified eight lncRNAs with the highest scores: NR2F1- that the results were associated with treatment status, we gen- AS1, LINC00654, FENDRR, PCAT2, and NKX2–1-AS1 in MCI (AUC erated Kaplan-Meier estimates for these transcripts in untreated >0.70) and LINC00478, LINC00173, and LINC00514 in MCII (AUC patients from the same cohort; neither showed significant sep- > 0.70). These lncRNAs serve as candidates for predicting rapid aration in their performance (P = 0.905 and P = 0.832, respec- metastasis in patients receiving ADT. Selecting all NEtD/NEPC tively; Fig. 6B). Expression for SSTR5-AS1 and LINC00514 in the lncRNAs with AUC >0.65 (n = 25), we performed survival analy- VPC and WCM cohorts illustrate their distinct and elevated ex- sis to ascertain their ability to separate patients for metastasis as pression in NEPC vs. AD patient samples (Fig. 6C). These results an outcome and end-point. Specifically, we calculated Kaplan- suggest a strong association between treatment status and in- Meier estimates for metastatic disease progression stratified by creased probability of metastatic disease in patients with dif- median expression in ADT-treated samples of the MCII cohort. ferential expression of these lncRNAs. This, together with re- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 AD NEPC EXPRESSION (LOG2) Probability of Mets Free Survival Probability of Mets Free Survival 02468 10 12 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 Probability of Mets Free Survival Probability of Mets Free Survival 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 Ramnarine et al. 13 sults from the NEtD model and NEPC clinical samples, implicate In this study, we characterized the unexplored global lncRNA SSTR5-AS1 and LINC00514 in NEPC and these lncRNAs serve as landscape during NEtD to provide insights into the NEPC non- strong candidates as predictive biomarkers for metastatic dis- coding milieu of this lethal and treatment-induced process. This ease post-RP following ADT. required the implementation of a sequence analysis pipeline One of the mechanisms observed with lncRNAs is direct with increased sensitivity towards lower expressed transcripts, RNA-RNA interaction with mRNA, resulting in regulation of their characteristic of lncRNAs. The pipeline was able to detect 37,749 expression (activation or repression). This type of investigation lncRNA transcripts (subclassified as either lincRNA, antisense, is computationally intensive, and there are limited algorithms or pseudogene) and quantify them in the two clinical cohorts available to identify putative mRNA targets genome-wide. How- (VPC and WCM). The novelty of this study lies in the use of pa- ever, a method was recently published to predict lncRNA-mRNA tient samples integrated with the NEtD PDX model to detect clin- interactions genome-wide [119], and so we sought to identify ically relevant lncRNAs involved in the NEtD/phase transition candidate mRNA transcripts interacting with SSTR5-AS1 and process. In this study, we identified 742 lncRNAs associated with LINC00514. The pipeline’s three core algorithms include Rac- NEtD and identified a robust 122 NEPC lncRNA patient signature cess [120] for the identification of accessible regions within the capable of classifying NEPC from AD patient samples. The motif lncRNA, IntaRNA [121] to calculate nucleotide interaction en- analysis identified significantly enriched TFBS that were unique ergies, and RactIP [122] to predict joint secondary structures. to NEtD Classes I Deactivated (TP53 and BRN1), II Activated (ELF5, Applying this methodology to SSTR5-AS1 and LINC00514 pro- SPIC, and HOXD1), III Persistent (SPDEF and HOX), IV Transient duced a list of predicted interacting partners for these lncRNAs (TP53, HSF1, HSF2, and FOXA1), and III Persistent combined with (Supplementary Tables S27–S28). The top-ranked mRNAs were IV Transient (TWIST1). Through similar analysis, we also iden- KDM4B and TADA3 that are predicted to hybridize and form joint tified common TFBS (CCNB2, E2F, ETS, REST, PAX5, PAX9, and structures independently with SSTR5-AS1 and LINC00514, re- STAF) enriched across all of the NEtD lncRNAs. From among the spectively (Supplementary Figs S12 and S13). In the clinical co- 100 top-ranking lncRNA, we observe that a subset have strong horts, TADA3 is down-regulated in NEPC vs AD (>2-fold), while clinical associations with metastatic PCa patients after receiv- KDM4B is up-regulated (>5 fold); however, only the deregulation ing ADT. In previous lncRNA studies in cancer, several have been of TADA3 is statistically significant (VPC P = 0.003 and WCM P = linked to malignant transformation with key roles affecting vari- 0.017). Both genes have NEPC associations (see Discussion), and ous aspects of cellular homeostasis, including proliferation, sur- our data suggest they are being regulated by these lncNRAs. vival, migration, and genomic instability [141]. Similarly, lncR- NAs identified in this study, including SSTR5-AS1 and LINC00514 with their association with poor outcome, FENDRR for its asso- ciation with rapid metastasis, and H19 and LINC00617 for their Discussion concordantly high expression across both of the discovery co- Primary NEPC arises de novo in 0.5% to 2% of all prostate can- horts, could be the missing links in the mechanisms causing cer patients [123]. However, tNEPC can develop in 20–30% of NEtD. These five represent the top candidates discovered in this metastatic castrate-resistant prostrate cancer tumours [124]and study due to this evidence but also for their characterization in increases with disease progression [125]. The real incidence of other cancer types. tNEPC may be higher because of under-recognition due to tu- FENDRR is a top deregulated lncRNA in NEtD Class IV Tran- mour heterogeneity, the limited number of metastatic tumour sient and may have a role in the NEtD process. It is implicated biopsies performed, lack of uniform consensus definition based in a lethal lung development disorder [142], lung cancer [143], within a mutational hotspot that is copy number lost in PCa on histology or biomarker expression, and frequent misclassifi- cation as high-grade PCa (most notable in tumours with mixed [144], and can bind to PRC2 [145, 146]. PRC2 plays a significant histologies) [126]. NEPC can be induced in vitro in AR LNCaP role in tumour progression through binding of HOTAIR (a very cells in androgen-depleted culture conditions [127, 128], sim- well-studied lncRNA). Together, HOTAIR and PRC2 are involved ilarly in vivo [7, 129], and in patient tumours long-term ADT in the control of chromatin structure and associated gene activ- has increased neuroendocrine differentiation [118, 124, 130]. It ity [147]. FENDRR may be involved in tumorigenesis like HOTAIR is now common to observe treatment-resistant tumours with due to its known interaction with PRC2. A recent study showed neuroendocrine features upon metastatic biopsy, and the pre- down-regulation of FENDRR is associated with poor prognosis vailing consensus is that epithelial plasticity enables tumour in gastric cancer and regulates cancer cell metastasis through adaptation in response to AR-targeted therapies [7, 9, 118, 126, fibronectin [ 71]. Functionally, this could be occurring in NEPC as 131–134]. This evidence supports the notion that tNEPC inci- well due to FENDRR’s transient expression in the NEtD model dence through NEtD will increase as new powerful ADTs en- and its association to rapid metastasis in ADT-treated PCa pa- ter the clinic. There is an urgency for therapeutic strategies and tients from the GRID (MCI) cohort. Another putative function of this transcript is through upregulating FOXF1, which is a clinical biomarkers defining NEtD/NEPC. Currently, the only op- tion for patients is the short-lived effects of platinum-based protein-coding gene and the sense form for the antisense tran- chemotherapy. Optimism is on the rise, as there is an AURKA in- script FENDRR. Antisense transcripts are known to regulate their hibitor (MLN8237) in a Phase 2 clinical trial (NCT01799278), com- sense forms (positively or negatively). Using TANRIC, an interac- binatorial approaches using AURKA with PARP inhibitors un- tive resource for the exploration of lncRNAs in large patient co- der investigation [135], indirect methods that resensitize the tu- horts within 20 TCGA cancer types [148], we see that FENDRR ex- mour to Enzalutamide [136] or platinum-based chemotherapy pression is positively correlated to FOXF1 in 16 of 20 cancer types −9 [137] (Phase 2 clinical trial NCT02489903 with a Phase 3 clinical (P < 3.71 × 10 ; Supplementary Table S29). In fact, FOXF1 dele- trial being planned), a SSTR4/5 analogue (Pasireotide/SOM230) in tion has been seen to significantly reduce FENDRR in endothe- four independent clinical trials at various Phases (NCT01646684, lial cells [149]. FOXF1 is also a target gene of p53 and is seen to NCT01313559, NCT01468532, and NCT01794793) with one al- regulate cancer cell migration and invasiveness [150]. Together these transcripts may play a transient coordinated role in NEtD ready reporting promising clinical efficacy [ 138], and increased study of NEPC/NEtD in general [134, 139, 140]. through PRC2 or fibronectin. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 14 NEPC lncRNAs LINC00514 is amongst the highest expressed lncRNAs in functions (SSTR5-AS1: SSTR5 or SSTR5-AS1: KDM4B: N-Myc) re- NEtD Class III Persistent. It has not been characterized. It is pre- quire thorough in vitro and in vivo exploration to ascertain their dicted to bind to TADA3 (Supplementary Fig. S14), potentially validity. causing a reduction of its activity. This is intriguing because Although multiple layers of genetic and epigenetic deregula- TADA3 is involved in the stabilization and activation of p53 [151, tion likely cooperate to facilitate NEtD, understanding the non- 152], and this putative interaction (LINC00514: TADA3) could be coding contribution to this multifarious process is necessary to an alternative mechanism for loss of p53 activity, already known design effective novel therapeutics. Using the five independent to be frequently lost in NEPC [15]. H19 and LINC00617 were two patient cohorts and our proven NEtD PDX LTL331 model, lncR- of the four highest (>10-fold) NEPC-expressed lncRNAs in this NAs such as FENDRR, LINC00514, LINC00617, H19, SSTR5-AS1, study and fortunately (unlike most lncRNAs) have both been and others identified in this study may provide more in-depth thoroughly characterized functionally. LINC00617 is highly con- insights into NEtD and NEPC. Research identifying the relation- served across vertebrate genomes and is required for mainte- ship of these lncRNAs to other known drivers, oncogenes, and nance of pluripotency and neural differentiation in embryonic Activated pathways in NEtD is now required. This study is the stem cells [153]. It controls this lineage commitment through first to report the lncRNA landscape of NEtD, a robust NEPC RNA-binding proteins PTBP1, hnRNP-K, and Nucleolin. These lncRNA expression clinical classifier, and provides numerous RNA-binding protein complexes have been detected at promot- candidate biomarkers and therapeutic targets. ers of NANOG, SOX2 (promoter of lineage plasticity in NEPC [28]), and FGF4 [153]. H19 has also been identified in neural differen- Methods tiation of pluripotent stem cells [154] but with unknown mech- anisms. With such an elevated level of expression in the clini- PDXs cal cohorts (∼30- to 40-fold and ∼20- to 30-fold in VPC/WCM for Animal ethics, care, experiments, xenograft generation, and all LINC00617 and H19, respectively), these lncRNA could be respon- protocols were carried out in accordance with the guidelines of sible for maintaining the neuronal component of NEPC through the Canadian Council of Animal care as previously described [7]. epigenetic regulation. Specific xenograft models used in this study have been previ- SSTR5-AS1 is the highest expressed lncRNA in the NEPC clin- ously published (protein-coding transcriptomes) by Akamatsu ical samples when requiring expression concordance in VPC et al. [24] and Mo et al. [162]. In brief, six LTL331, two LTL313, and WCM cohorts. It is an antisense transcript of SSTR5, which and two LTL418 PDXs were raised in NOD-SCID mice (NOD.CB17- is a member of the superfamily of somatostatin receptors. So- Prkdcscid/J) at the Living Tumor Laboratory [163]. Xenograft tis- matostatins are peptide hormones that regulate diverse cellu- sue was harvested after fixed lengths of time post host castra- lar functions such as neurotransmission, cell proliferation, and tion, tissue was measured, fixed for histopathological analysis, endocrine signalling, as well as inhibiting the release of many and processed for RNA analysis. hormones and other secretory proteins. The SSTR family (1–5) are markers for neuroendocrine tumours of the lung [155], with SSTR1 and SSTR5 the most dominant forms of SSTR in neu- Clinical datasets roendocrine tumours in general [54]. Interestingly, exploration within TANRIC showed a strong positive correlation in expres- We used five clinical cohorts from 1) WCM [ 5]; 2) GenomeDx Bio- sion with SSTR5 to SSTR5-AS1 in 14 of 20 cancer types (P < 2.18 sciences (GX) Inc. (MCI and MCII); 3) JHSM; and 4) VPC, cumu- −15 × 10 ; Supplementary Table S29). Furthermore, SSTR5 mRNA latively totalling 927 samples. For the VPC, 80 specimens were is detectable in the blood of neuroendocrine tumours of the lung obtained from patients undergoing RP and snap frozen follow- [156] and could be a valuable non-invasive diagnostic marker for ing a protocol approved by the Clinical Research Ethics Board of NEPC. In fact, clinicians utilize this biological feature in other the University of British Columbia, the BC Cancer Agency, and neuroendocrine tumours (NETs) using Octreoscans to determine Vancouver General Hospital pathology (depending on the sam- tumour stage and/or identification of sites of metastasis. Oc- ple source). All patients signed a formal consent form approved treoscans, when compared with positron emission tomography by the ethics board. A subset of the GX Decipher GRID database (PET) scans (another commonly used approach for this), appear of clinical specimens was selected, totalling 777 patient PCa more sensitive in the detection of well-differentiated NETs [157]. expression profiles (all from formalin-fixed parafin embedded In addition to this, therapeutically, somatostatin analogues are tissue) and were obtained from two RP Mayo Clinic cohorts emerging as a promising treatment option for inoperable or that have been previously described (MCI [116]and MCII [117]). metastatic NETs [158]. However, specifically in NEPC, targeting JHSM samples, totalling 33 samples, were retrieved from surgical SSTR5 and/or SSTR5-AS1 for diagnostic or therapeutic purposes pathology and consultation files of Johns Hopkins Hospital (John is in its infancy. Interestingly, SSTR5 (C terminal) is required for Hopkins Registry) from 1999 to 2013, as previously described Rb induction and G1 cell cycle arrest [159], resulting in anti- [164]. The 33 samples were annotated as 6 morphologically di- proliferative effects. However, without Rb (known to be lost in agnosed pure SCPC samples, 12 high-risk (Gleason 9–10) AD, 10 NEPC), this function is negated. Alternatively, the interaction ev- SCPC (SC-mixed), and 5 AD (AD-mixed) from mixed histology tu- idence for SSTR5-AS1 and KDM4B (Supplementary Fig. S13) pro- mours containing separate adenocarcinoma and small cell com- vides another strong connection to NEPC biology. KDM4B is a ponents. For this cohort, samples were dicotimized into either histone demethylase and a key molecule in AR signaling and AD (AD and AD-mixed samples) or NEPC (SCPC and SCPC-mixed turnover [160]. In NEPC with the absence of the AR, KDM4B could sampels) for the purposes of validating the 122 NEPC lncRNA pa- interact with N-Myc instead, where it has been shown to regu- tient signature. We also explored an externally processed cohort late and epigenetically activate this oncogene in neuroblastoma comprising 114 metastatic CRPC specimens, of which 44 were [161]. N-Myc has been seen to drive the progression of NEPC [5, NEPC [32] and used in this study. Referred to in the text as the 26, 27] and recently through EZH2-mediated transcription [27]. extNEPC cohort, we accessed and visualized this data through However, another mechanism of activation could be facilitated cBioPortal [165, 166] Version 1.9.0 [167]. OncoPrint schematics through SSTR5-AS1 regulation. However, both of these putative were generated for displaying multiple genomic alterations by Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 15 heatmap for the lncRNAs. The extNEPC study samples were clas- then at 70C for 15 minutes. Prior to use in qRT-PCR, products sified using a pathologic classification system [ 6] that included were diluted 10-fold with water. FastStart Essential Green Mas- five catagories: “A,” usual prostate adenocarcinoma without ter kit from Roche (Catalogue #06 402 712 001) was used as de- neuroendocrine differentiation; “B,” usual prostate adenocarci- scribed from their protocol for qRT-PCR reactions. In brief, 2 μl noma with neuroendocrine differentiation >20%; “C,” small-cell of water, 3 μl of a mixture of forward and reverse primers (each carcinoma; “D,” large-cell neuroendocrine carcinoma; and “E,”, at a concentration of 10 μM),and 10 μl of the Roche Master Mix mixed small-cell carcinoma–adenocarcinoma. were aliquoted into each well of a 96-well plate. A mixture of 4 μl of water plus 1 μl of the diluted cDNA was then added to the appropriate wells. Expression was then quantified (as mea- Material collection and processing (VPC Cohort) sured by Ct) through the Roche Light Cycler 96 machine. Each Hematoxylin and eosin (H&E) stained, formalin-fixed paraffin- lncRNA/sample pair was quantified with technical replicates in embedded, and fresh frozen sections were reviewed by a pathol- triplicate. Average and standard deviation of Ct were calculated ogist to identify blocks with highest tumour content. For each across these triplicates, and Ct calculated relative to house frozen block used, a 5-μmslide was first taken for H&E staining; keeper gene PSMB4 (most consistent and highly expressed gene then 4 × 100-μmsections were taken for DNA and RNA isolation vs REEP5 and SNRPD3). Delta Cts were calculated relative to before a second 5-μm slide was taken for H&E staining. Each control samples, and fold changes were plotted using Prisms H&E slide was required to have tumour content >50% for a tu- GraphPad software (Supplemental Table 12). mour to proceed for sequencing. RNA from 100-μmsections of snap frozen tissue was isolated using the mirVana Isolation Kit RNA sequence analysis pipeline from Ambion (AM 1560). RNA sequencing was performed on Il- lumina HiSeq 2000 at BC Cancer Agency Michael Smith Genome We implemented an lncRNA sequence analysis pipeline that in- Sciences Centre according to standard protocols. cludes algorithms catered to the detection of known and novel transcripts (Supplementary Figs S1 and S2). Implemented in- house, this pipeline is modified and extended from the tuxedo Material collection and processing (GRID and JHSM) suite of sequence analysis algorithms [45]. Once received from For GRID (MCI and MCII) and JHSM cohorts, specimen selection, the sequencing centre in bam format, all sequenced model sys- RNA extraction, and microarray hybridization were performed tems and patient samples were de-aligned into raw fastq for- in a Clinical Laboratory Improvement Amendments-certified mat (including flagged reads) using bam2fastq and put through laboratory facility (GenomeDx Biosciences, San Diego, CA, USA) the following pipeline. To ensure high-quality sequence reads, as described previously [116, 117]. Total RNA extraction, purifica- libraries were trimmed using a Sickle, a windowed-adaptive ap- tion, RNA amplification, and labelling were done using the Ova- proach [171]. For each read pair processed together, the algo- tion WTA FFPE system (NuGen, San Carlos, CA, USA). RNA was rithm determines the most optimal inner read sequence by hybridized to Human Exon 1.0 ST GeneChips (Affymetrix, Santa trimming both 3’ and 5’ prime ends based on quality and Clara, CA, USA). After microarray profiling, quality control was length thresholds (for full description, see [172]). Bases with a preformed using the Affymetrix Power Tools package, and probe quality score of <99.0% base call accuracy (corresponding to set normalization was performed using the Single Channel Ar- a Phred quality score of 20) were removed. Reads less than ray Normalization algorithm [168]. approximately two-thirds read length (30 nt in WCM and 60 nt in VPC) post-trimming were discarded. Highly repetitive se- quences (>2% of library) were also discarded post-trimming us- Quantitative Real-Time Polymerase Chain Reaction ing the cutadapt tool. All quality control metrics were gener- (qRT-PCR) ated and quanitified (pre- and post-trimming) using the FASTX- Primers were designed using Primer3 and checked with in sil- Toolkit and the FastQC Windows software. Reads were aligned ico PCR in UCSC Genome Browser (s Supplementary Table S26 to the Hg19 human genome build using an unspliced aligner for forward and reverse primer sequences). Housekeeping genes for handling exonic reads (Bowtie–v2.2.3) in conjunction with PSMB4, REEP5, and SNRPD3 were selected on the basis of high, a spliced aligner to handle reads spanning exon-exon junc- consistent expression levels across many cell and tissue types tions (Tophat 2.0.12). Transcriptome reconstruction using En- and were used in the MiTranscriptome lncRNA study [169, 170]. sembl GRCh37.75 gene tracks for each library was performed us- Two lncRNAs from each NEtD class and three NEPC lncRNAs ing a quasi de novo (genome-guided) approach (Cufflinks v2.2.1), from among the top candidates (Supplementary Fig. S12 and where reads were assembled and abundances estimated us- Supplementary Table S26) were selected (n = 11) for qRT-PCR ing an overlap graph producing a minimal spanning network validation. The cDNA from the PDX LTL331 models three time of transcripts. This version of Ensembl contained 38 transcript points (AD, postTX, and NEPC) were used to validate the NEtD classes grouped by 4 core biotypes. At this stage, transcripts lncRNAs and a subset of the VPC clinical samples for the NEPC were also multi-read and fragment bias corrected. Transcripts lncRNAs. With the rarity of clinical NEPC samples, tumour tissue with highly abundant expression were masked (e.g., rRNAs) and subsequent RNA were extremely limited. Due to this, only from downstream steps to increase transcript quantification ac- three NEPC (V73, V90, and V91) and one AD (V60) clinical sample curacy. Sample transcriptomes, the reference genome, and the were included in this validation. For each lncRNA and sample transcript annotation were then meta-assembled (Cuffmerge) to tested, the following experimental protocol was carried out: 1 produce a single annotation transcriptome model. Based on this μg of total RNA for each sample was diluted to 18 μl with water model, transcript quantification (Cuffquant) and normalization and 1 μl of random hexamers (50μM; Thermo Fisher). The mix- (Cuffnorm). Geometric and FPKM normalization performed in- ture was heated to 65C for 5 minutes and chilled. Afterwards, 5 dependently (Cuffnorm) corrected for uneven library sequencing μl of 5x reverse transcriptase buffer, 1 μlof10mMdNTPs,and depths between samples and variable transcript lengths within 1 μl of Superscript II reverse transcriptase (Thermo Fisher) were samples. Transcript expression displaying computational arti- added. Each sample was then incubated at 42C for 1 hour and facts (expression values <0.1 known to occur with Cufflinks) Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 16 NEPC lncRNAs were converted to zero values. This generated transcript expres- (GREF [includes the androgen receptor and the closely related sion where only lncRNAs (Ensembl and ENCODE-based) were glucocorticoid, mineralocorticoid, and progesterone receptors], extracted and used for all downstream analysis. All algorithms NRSF [REST], SOX, HOX, STEM, E2FF, ETSF, and ETVI1) motifs denoted in brackets are referenced and described in the Trap- included in this study are described in Supplementary Table nell et al. Nature protocol [45]. Each cohort (VPC and WCM) was S25. Core and matrix similarities were again set to 1 and op- processed independently by this pipeline, and then transcrip- timized, respectively. (3) Multiple lncRNA enrichment analysis tome annotations were merged. This was accomplished using was performed using the Overrepresented TFBS algorithm. En- Ensembl transcript IDs combined with transcript lengths to pro- richment of matrix/matrix family was determined by Genomatix duce unique transcript identifiers for each lncRNA across co- calculated z-scores (>2or <-2), which is based on the distance horts. from the population mean (genome or promoter sequence back- ground) in units of the population standard deviation for query sequence/promoter. Genomatix calculates z-scores with a conti- RNA-RNA interaction analysis nuity correction using the formula z = (x-E-0.5)/S, where x is the A genome-wide analysis for SSTR5-AS1 and LINC00514 lncRNA number of found matches in the input data, E is the expected interactions was performed using a multistep systemic ap- value, and S is the standard deviation. This formula is also de- proach [119]. This tool is available publicly within an online scribed in the oPOSSUM algorithm [180]. A z-score < -2 or >2 database [173] hosted by the Computational Biological Research can be considered statistically significant and corresponds to a Center at the National Institute of Advanced Industrial Science P-value of approximately 0.05. and Technology in Japan. The interaction search space included all hg19 annotated lncRNA and mRNA transcripts. This gener- Microarray to sequencing platform lift over/mapping ated top-ranking interaction partners (n = 100), based on local interaction minimum free energy (Supplementary Tables S27– Affymetrix Human Exon 1.0 ST GeneChip probes were mapped S28). R-chie [174, 175] was used to visualize the top-ranking pre- to hg19 coordinates using SMALT v0.76 [181]. Probe set genomic dictions KDM4B and TADA3 for SSTR5-AS1 and LINC00514, re- regions (PSRs) were redefined accordingly. Exons within each spectively, using the double structure feature (Supplementary lncRNA from sequencing cohorts (VPC and WCM) were inte- Figs S13 and S14). All bases that were not within the interac- grated with PSRs to build an overlap table to determine ab- tion site were predicted to form RNA secondary structure by sence/presence of lncRNA transcripts on the affymetrix microar- RNAfold [176–178] selecting to enforce constrained pairing pat- ray. R function iRanges v2.9.18 was used with method findOver- terns for the interacting bases. Minimum free-energy structures lap to build the described table above. Microarray PSRs were re- were predicted by RNAfold on the 300nt sequences upstream quired to be entirely within sequenced exon genome regions; and downstream of the interaction site. otherwise they were excluded. Applying this methodology, 106 of 122 NEPC lncRNAs (87%) and 81 of 100 NEtD lncRNAs (81%) mapped to microarray PSRs for clinicopathological analysis on TFBS identification and enrichment analysis GRID cohorts MCI and MCII. All TFBS analysis was performed using Genomatix software, databases, and algorithms [179]. Three types of TF analysis were Statistical analysis carried out in this study: (1) single lncRNA motif characteriza- tion, (2) multiple lncRNA analysis for select TFs, and (3) multiple For all cohorts, the programming language R v3.0 was used lncRNA enrichment analysis. Prior to any of the above, lncRNA for statistical analysis. For VPC and WCM cohorts, unsuper- transcript(s) were submitted to the Gene2Promoter algorithm for vised hierarchical clustering was performed with the h.clust retrieval of promoter sequences. Databases used with this algo- package with Pearson correlation for distance and average link- rithm included ElDorado 12–2013 and NCBI build 37 (for multiple age used. Only transcripts within the top fifth percentile based lncRNA analysis where genomic background needed to match on their standard deviations were selected. The clustering and sequencing data) or the most recent databases ElDorado 12–2016 heatmaps generated were built using the heatmap.2 function. and GRCh38 (for single lncRNA analysis where genomic back- Similar clustering analysis was performed for GRID cohorts ground was not relevant). Transcripts with alternative isoforms except with Euclidian distance, the ward method for linkage, were required to have gold level (experimentally verified 5’ com- and the use of the heatmap.3 function due to its advanced plete transcript), silver level (transcript with 5’ end confirmed by row/column labelling features. For all cohorts before clustering, PromoterInspector prediction), or bronze level (annotated tran- normalized log2 expression values were standardized/scaled us- script, no confirmation for 5’ completeness) quality for their al- ing a z-score that ranged from -2 to 2. For principal component ternative isoforms. (1) Single lncRNA motif characterization was analysis, the R package prcomp was used to calculate variance performed using the MatInspector algorithm [97–99]withpa- among transcript and sample subsets for the calculation of tran- rameters “core similarity” (degree of similarity for highest con- script weights and principle components. The top three compo- served bases of motif) set to 1 and “matrix similarity” (degree of nents were used for visual inspection. For all clinical group-wise similarity between motif and query sequence) set to optimized comparisons, a standard Student’s t-test was applied to iden- as recommended by Genomatix and as described in MatInspec- tify significantly differentially expressed transcripts between tor referenced papers above. MatInspector uses the best in field groups/phenotypes. Significance thresholds were implemented MatBase database for TFBS motif/matrix annotation, where Ma- by enforcing a strict P-value cut-off of <0.05. Multiple test correc- trix Family Library Version 10.0 was used. (2) Multiple lncRNA tion was applied to P-values using the Bonferroni and Hochberg analysis for select TFs was performed using MatInspector and method to mathematically minimize false discovery rate (FDR) select TF motifs (“matrix”) applied accordingly. All matrix an- and with a cut-off of P-value < 0.05. See Supplementary Ta- notation, descriptions, and matrix family definitions are listed bleS4for theseresults.To biologically minimize FDR, mathemati- in Supplementary Table S25. Select TF matrices (BRN2, STAT3, cal FDR correction was removed and instead followed the filter- NKX3, NMYC, SOX2, and SOX11) and select TF matrix families down workflows in Figs 3A–C and Fig. 1C. Despite mathematical Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 17 FDR being removed, statistical significance of P-value < 0.05 was Additional file still maintained during the filter-down approach using the de- Table S1: Discovery cohorts clinical and sequencing informa- scribed method in each step. See Supplementary Table S5 for tion. these results. For ROC curves and AUC calculations, the R pack- Table S2: Differentially expressed AR- and AR+ CRPC lncRNA. age “pROC” was used. Kaplan-Meier analysis was performed for Table S3: NEPC lncRNA expression signature labels and order for determining survival outcome using the R package “survfit” with Fig. 2. transcripts displaying below background (<0.1) expression being Table S4: Clinical cohort group-wise comparisons. removed from this analysis. Table S5: NEtD and NEPC lncRNA transcript annotation Table S6: TFBS-Enrichment by Matrix for NEtD Class I Table S7: TFBS-Enrichment by Matrix for NEtD Class II Transcript ranking Table S8: TFBS-Enrichment by Matrix for NEtD Class III Table S9: TFBS-Enrichment by Matrix for NEtD Class IV NEtD lncRNA transcripts were ranked based on fold changes ob- Table S10: TFBS-Enrichment by Matrix for NEtD Class III-IV served in the clinical group-wise comparisons. For NEtD Class Table S11: TFBS-Enrichment by Matrix for all NEtD lncRNAs I Deactivated, the three group-wise comparison fold changes Table S12: TFBS-Enrichment by Matrix Family for NEtD Class I were used (NEPC vs AD, CRPC vs AD, and NHT vs AD), where Table S13: TFBS-Enrichment by Matrix Family for NEtD Class II the minimum fold change observed between the three com- Table S14: TFBS-Enrichment by Matrix Family for NEtD Class III parisons was selected and then ranked in decreasing order. For Table S15: TFBS-Enrichment by Matrix Family for NEtD Class IV Class II Activated and Class III Persistent transcripts, NEPC vs Table S16: TFBS-Enrichment by Matrix Family for NEtD Class III- AD fold changes were calculated and ranked in increasing or- IV der for both VPC and WCM cohorts, where the maximum fold Table S17: TFBS-Enrichment by Matrix Family for all NEtD lncR- change between VPC and WCM was selected. For Class IV Tran- NAs sient transcripts, absolute fold changes for AD vs NHT and NHT Table S18: TFBS Overlap Table by Matrix vs NEPC were calculated and ranked in increasing order with the Table S19: TFBS Overlap Table by Matrix Family maximum fold change from either group selected. Similar rank- Table S20: Unique and common TFBS in NEtD lncRNA ing was performed for NEPC lncRNA transcripts, which were or- Table S21: Select TF Identification for NEtD Class I dered by increasing (up-regulated transcripts in NEPC vs. AD) Table S22: Select TF Identification for NEtD Class II-III or decreasing (down-regulated transcripts in NEPC vs. AD) or- Table S23: Select TF Identification for NEtD Class III-IV der to determine the highest- and lowest-expressed transcripts Table S24: Select TF Identification for all NEtD lncRNAs in NEPC vs AD, respectively. Concordantly expressed transcripts were required between VPC and WCM cohorts. The top 20 lncR- Table S25: Genomatix Matrix and Matrix Family definitions Table S26: Top-ranking NEPC and NEtD lncRNAs. NAs (based on fold changes from clinical samples defined above) Table S27: SSTR5-AS1-predicted RNA (mRNA or lncRNA) interac- were taken from each group. This produced 20 × 5 group (n = tions with associated binding energies, predicted transcript En- 100) isoforms representing 76 unique lncRNA transcripts. These sembl ID, name, interaction position, and ranking. represent the top NEtD/NEPC lncRNA candidates from this study Table S28: LINC00514-predicted RNA (mRNA or lncRNA) interac- (Supplementary Table S26). No pseudogenes were included in tions with associated binding energies, predicted transcript En- these rankings. sembl ID, name, interaction position, and ranking. Table S29: TANRIC results for lncRNAs FENDRR and SSTR5-AS1. Spearman rank correlation for protein coding genes that are Availability of supporting data the sense forms to the above antisense transcripts. Numbers A subset of the sequenced samples (n = 70) used in this study in brackets denote P-values. NSC, no significant correlation; NA, was from previous studies with all raw sequencing data rean- mRNA data was not available for this tumour type, therefor the alyzed here using the described pipeline above. These 70 sam- analysis was not applicable. ples have been previously submitted to the European Nucleotide Figure S1: The next-generation sequence analysis pipeline im- Archive (ENA) or NCBIs Gene Expression Omnibus. This in- plemented for the detection and quantification of lncRNAs in cludes the 6 NEPC PDX model samples [24] (ENA accession num- this study. A) The nine-step lncRNA next-generation sequenc- ber PRJEB9660 and GEO accession number GSE59986), 2 CRPC ing analysis pipeline with core algorithms (Bowtie2, Tophat2, PDXmodelsamples[162] (ENA accession number PRJEB19256), Cufflinks2, Cuffmerge, Cuffquant, and Cuffnorm) implemented 4 NEPC (VPC) samples [31, 56], 23 AD (VPC) samples [56](ENA from the Tuxedo suite of analysis tools. Sequencing quality con- accession number PRJEB6530), 30 AD (WCM) samples [5], and trol metrics before and after trimming of data for sample V60 7 NEPC (WCM) samples [5]. The remaining unpublished se- is outlined in B–F. This includes pre-trimming (A) phred quality quenced samples (n = 55) have been submitted to the ENA under scores and (B) percentage of each base type across read library accession number PRJEB21092. Please see Supplementary Table at each base pair position. After quality control corrections were S1 for a summary of sequencing and clinical information on applied, V60 read library had acceptable (D) Phred quality scores these 125 samples. All microarray samples from GX cohorts, in- (∼30 Phred Score) and (E) expected base type percentages (∼25%) cluding 545 AD (MCI [116]) samples and 232 AD (MCII [117]) sam- for T, C, A, and G. F) All over-represented sequences that were ples, are accessible through Gene Expression Omnibus acces- >2% of library were removed from the V60 read library. See Meth- sion numbers GSE46691 and GSE62116, respectively. Additional ods for complete listing and version numbers for all algorithms supporting data and custom code from the sequencing pipeline and tools used in the sequence analysis pipeline. described above are also available from the GigaScience GigaDB Figure S2: Average Phred quality scores for all VPC and WCM database [182]. samples before and after quality control corrections were ap- plied. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 18 NEPC lncRNAs Figure S3: Unsupervised hierarchal clustering (A–D) and princi- Funding ple component analysis (E–H) on the four major Ensembl tran- This work was supported by the Mitacs Accelerate PhD Fellow- script classes detected within the VPC cohort. Samples are la- ship Program (IT04310 to V.R.R.) in collaboration with GenomeDx belled as adenocarcinomas (blue) and neuroendocrine tumours Biosciences, Terry Fox Foundation (201012TFF to C.C.C.), and (gold). Prostate Cancer Canada Team Grant (T2013–01 to C.C.C.). Figure S4: Unsupervised hierarchal clustering (A–D) and princi- ple component analysis (E–H) on the four major Ensembl tran- script classes detected within the WCM cohort. Samples are la- Competing interests belled as adenocarcinomas (blue) and neuroendocrine tumours The authors declare that they have no competing interests. (gold). Figure S5: Detected and differentially expressed lncRNAs among −/+ (A) AR CRPC xenograft models (B) and matched clinical sam- Aknowledgements ples. P, patient; X, xenograft. We are grateful to the following GenomeDx bioinformaticians: Figure S6: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC Mandeep Takhar for her help with GRID statistical analy- Cohort – LEGEND. Clinical, transcriptome, and genomic annota- sis/primer code; Hussam Al-Deen Ashab for his help with IA tion for samples plotted in Supplementary Figs S7–S11. All an- analysis (ultimately was not included in the paper, but effort notations were generated from cBioportal with the exception and knowledge gained from the results were insightful for this of the NEPC molecular subtype, which was assigned from this work); Nicholas Erho for his effort in mapping our sequenced study. Transcripts denoted with 1 in superscript within Supple- data to the GRID microarray, mentoring, and constant support; mentary Figs S8–S11 are lncRNAs that overlap a NEPC lncRNA and Mohammed Alshalalfa for his guidance and supervision on with a NEtD lncRNA. For example, H19 appears in Supplemen- all NEPC/NEtD lncRNA clinicopathological GenomeDx analysis. tary Fig. S8 (NEPC lncRNA) and supplementary 9 (NEtD lncRNA We would also like to deeply thank Daniel Lai and Alex Gawron- Class II). ski for their advice with RNA-RNA visualization and interaction Figure S7: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC analysis algorithms. We would like to thank Faraz Hach for his Cohort - Known NEPC genes and TFs. Select NEPC onco- manuscript insights and advice. Lastly, we are extremely grate- genes, tumour suppressor, and transcription factors that have ful to Stephanie Giles Ramnarine for her manuscript comments, been reported previously or within this study for transcrip- advice, and support. tomic/genomic context with Supplementary Figs S8–S11. Please see Supplementary Fig. S6 for figure legend. Figure S8: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC References Cohort - NEPC lncRNA. Transcripts from the NEPC lncRNA ex- pression signature that are up-regulated (74 of 122), testable (58 1. Torre LA, Bray F, Siegel RL et al. Global cancer statistics, of 74), and altered (25 of 58) in the extNEPC cohort. Please see 2012. CA Cancer J Clin 2015;65:87–108. Supplementary Fig. S6 for figure legend. 2. Karantanos T, Evans CP, Tombal B et al. Understanding the Figure S9: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC mechanisms of androgen deprivation resistance in prostate Cohort - NEtD lncRNA Class II. Transcripts from NEtD lncRNA cancer at the molecular level. Eur Urol 2015;67:470–9. Class II (222), testable (128 of 222), and altered (26 of 128) in the 3. Grasso CS, Wu YM, Robinson DR et al. The mutational land- extNEPC cohort. Please see Supplementary Fig. S6 for figure leg- scape of lethal castration-resistant prostate cancer. Nature end. 2012;487:239–43. Figure S10: NEPC and NEtD lncRNA Oncoprint Plot in the 4. Vlachostergios PJ, Puca L, Beltran H. Emerging variants extNEPC Cohort - NEtD lncRNA Class III. Transcripts from NEtD of castration-resistant prostate cancer. Curr Oncol Rep lncRNA Class III (84), testable (79 of 84), and altered (29 of 84) in 2017;19:32. the extNEPC cohort. Please see Supplementary Fig. S6 for figure 5. Beltran H, Rickman DS, Park K et al. Molecular characteri- legend. zation of neuroendocrine prostate cancer and identification Figure S11: NEPC and NEtD lncRNA Oncoprint Plot in the of new drug targets. Cancer Discov 2011;1:487–95. extNEPC Cohort - NEtD lncRNA Class IV. Transcripts from NEtD 6. Epstein JI, Amin MB, Beltran H et al. Proposed morphologic lncRNA Class IV (45), testable (36 of 45), and altered (11 of 31) in classification of prostate cancer with neuroendocrine dif- the extNEPC cohort. Please see Supplementary Fig. S6 for figure ferentiation. Am J Surg Pathol 2014;38:756–67. legend. 7. Lin D, Wyatt AW, Xue H et al. High fidelity patient-derived Figure S12: Quantitative Real-Time Polymerase Chain Reac- xenografts for accelerating prostate cancer discovery and tionon select NEPC and NEtD lncRNAs. drug development. Cancer Res 2014;74:1272–83. Figure S13: Hypothetical RNA-RNA folding structure for exon 4 8. Aparicio AM, Harzstark AL, Corn PG et al. Platinum-based of SSTR5-AS1 (top) and the 3’UTR of KDM4B (bottom). Predicted chemotherapy for variant castrate-resistant prostate can- base pair binding (green arcs) along the sequence (black arrow) cer. Clin Cancer Res 2013;19:3621–30. are displayed, included predicted interaction site (orange bars). 9. Wang HT, Yao YH, Li BG et al. Neuroendocrine Prostate Can- Figure S14: Hypothetical RNA-RNA folding structure for exon 4 cer (NEPC) progressing from conventional prostatic adeno- of LINC00514 (top) and the 3’UTR of TADA3 (bottom). Predicted carcinoma: factors associated with time to development of base pair binding (green arcs) along the sequence (black arrow) NEPC and survival from NEPC diagnosis-a systematic re- are displayed, included predicted interaction site (orange bars). view and pooled analysis. J Clin Oncol 2014;32:3383–90. 10. Terry S, Beltran H. The many faces of neuroendocrine dif- ferentiation in prostate cancer progression. Front Onco 2014;4:60. 11. Lee JK, , Bangayan NJ, Chai T et al. Systemic surfaceome pro- filing identifies target antigens for immune-based therapy Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 19 in subtypes of advanced prostate cancer. PNAS 2018; 115, Biochim Biophys Acta 2017;1860:1002–12. E4473–E4482. 31. Lapuk AV, Wu C, Wyatt AW et al. From sequence 12. Shen R, Dorai T, Szaboles M et al. Transdifferentiation of to molecular pathology, and a mechanism driving the cultured human prostate cancer cells to a neuroendocrine neuroendocrine phenotype in prostate cancer. J Pathol cell phenotype in a hormone-depleted medium. Urol Oncol 2012;227:286–97. 1997;3:67–75. 32. Beltran H, Prandi D, Mosquera JM et al. Divergent clonal evo- 13. Palmgren JS, Karavadia SS, Wakefield MR. Unusual and un- lution of castration-resistant neuroendocrine prostate can- derappreciated: small cell carcinoma of the prostate. Semin cer. Nat Med 2016;22:298–305. Oncol 2007;34:22–29. 33. Viswanathan VS, Ryan MJ, Dhruv HD et al. Dependency of a 14. Tan HL, Sood A, Rahimi HA et al. Rb loss is characteristic of therapy-resistant state of cancer cells on a lipid peroxidase prostatic small cell neuroendocrine carcinoma. Clin Cancer pathway. Nature 2017;547:453–7. Res 2014;20:890–903. 34. Gibb EA, Brown CJ, Lam WL. The functional role of long non- 15. Chen H, Sun Y, Wu C et al. Pathogenesis of prostatic small coding RNA in human carcinomas. Mol Cancer 2011;10:38. cell carcinoma involves the inactivation of the P53 pathway. 35. Cheetham SW, Gruhl F, Mattick JS et al. Long noncoding Endocr Relat Cancer 2012;19:321–31. RNAs and the genetics of cancer. Br J Cancer 2013;108:2419– 16. Ku SY, Rosario S, Wang Y et al. Rb1 and Trp53 cooperate to 25. suppress prostate cancer lineage plasticity, metastasis, and 36. Marchese FP, Raimondi I, Huarte M. The multidimensional antiandrogen resistance. Science 2017;355:78–83. mechanisms of long noncoding RNA function. Genome Biol 17. Ham WS, Cho NH, Kim WT et al. Pathological effects of 2017;18:206. prostate cancer correlate with neuroendocrine differentia- 37. Sun W, Yang Y, Xu C et al. Regulatory mechanisms of long tion and PTEN expression after bicalutamide monotherapy. noncoding RNAs on gene expression in cancers. Cancer J Urol 2009;182:1378–84. Genet 2017;216–217:105–10. 18. Zou M, Toivanen R, Mitrofanova A et al. Transdifferentia- 38. Gutschner T, Diederichs S. The hallmarks of cancer: a long tion as a mechanism of treatment resistance in a mouse non-coding RNA point of view. RNA Biol 2012;9:703–19. model of castration-resistant prostate cancer. Cancer Dis- 39. Kondo Y, Shinjo K, Katsushima K. Long non-coding RNAs cov 2017;7:736–49. as an epigenetic regulator in human cancers. Cancer Sci 19. Li Y, Donmez N, Sahinalp C et al. SRRM4 drives neuroen- 2017;108:1927–33. docrine transdifferentiation of prostate adenocarcinoma 40. Sahu A, Singhal U, Chinnaiyan AM. Long noncoding RNAs under androgen receptor pathway inhibition. Eur Urol, in cancer: from function to translation. Trends Cancer 2016. 2015;1:93–109. 20. Li Y, Chen R, Bowden M et al. Establishment of a neuroen- 41. Bhan A, Soleimani M, Mandal SS. Long noncoding RNA and docrine prostate cancer model driven by the RNA splicing cancer: a new paradigm. Cancer Res 2017;77:3965–81. factor SRRM4. Oncotarget 2017;8:66878–88. 42. Cheng W, Zhang Z, Wang J. Long noncoding RNAs: new 21. Zhang X, Coleman IM, Brown LG et al. SRRM4 expression players in prostate cancer. Cancer Lett 2013;339:8–14. and the loss of REST activity may promote the emergence 43. Lin D, Dong X, Wang K et al. Identification of DEK as a poten- of the neuroendocrine phenotype in castration-resistant tial therapeutic target for neuroendocrine prostate cancer. prostate cancer. Clin Cancer Res 2015;21:4698–708. Oncotarget 2015;6:1806–20. 22. Bishop JL, Thaper D, Vahid S et al. The master neural tran- 44. Clermont PL, Lin D, Crea F et al. Polycomb-mediated si- scription factor BRN2 is an androgen receptor suppressed lencing in neuroendocrine prostate cancer. Clin Epigenetics driver of neuroendocrine differentiation in prostate cancer. 2015;7:40. Cancer Discov 2016, 7, 54–71;7:54-71. 45. Trapnell C, Roberts A, Goff L et al. Differential gene and 23. Kim J, Jin H, Zhao JC et al. FOXA1 inhibits prostate cancer transcript expression analysis of RNA-seq experiments neuroendocrine differentiation. Oncogene 2017;36:4072–80. with TopHat and Cufflinks. Nat Protoc 2012; 7:562–78. 24. Akamatsu S, Wyatt AW, Lin D et al. The placental gene 46. Gutschner T, Hammerle M, Diederichs S. MALAT1 – a PEG10 promotes progression of neuroendocrine prostate paradigm for long noncoding RNA function in cancer. J Mol cancer. Cell Rep 2015;12:922–36. Med (Berl) 2013;91:791–801. 25. Ci X, Hao J, Dong X et al. Heterochromatin protein 1al- 47. Cross DS, Burmester JK. Functional characterization of the pha mediates development and aggressiveness of neuroen- GDEP promoter and three enhancer elements in retinoblas- docrine prostate cancer. Cancer Res 2018. toma and prostate cell lines. Med Oncol 2008;25:40–49. 26. Lee JK, Phillips JW, Smith BA et al. N-Myc drives neuroen- 48. Reding DJ, Zhang KQ, Salzman SA et al. Identification of docrine prostate cancer initiated from human prostate ep- a gene frequently mutated in prostate tumors. Med Oncol ithelial cells. Cancer Cell 2016;29:536–47. 2001;18:179–87. 27. Dardenne E, Beltran H, Benelli M et al. N-Myc induces an 49. Niknafs YS, Han S, Ma T et al. The lncRNA landscape of EZH2-Mediated transcriptional program driving neuroen- breast cancer reveals a role for DSCAM-AS1 in breast cancer docrine prostate cancer. Cancer Cell 2016;30:563–77. progression. Nat Commun 2016;7:12791. 28. Mu P, Zhang Z, Benelli M et al. SOX2 promotes lineage 50. Wang O, Yang F, Liu Y et al. C-MYC-induced upregulation of plasticity and antiandrogen resistance in TP53- and RB1- lncRNA SNHG12 regulates cell proliferation, apoptosis and deficient prostate cancer. Science 2017; 355:84–88. migration in triple-negative breast cancer. Am J Transl Res 29. Maina PK, Shao P, Liu Q et al. c-MYC drives histone 2017;9:533–45. demethylase PHF8 during neuroendocrine differentiation 51. Chen T, Yang P, He ZY. Long non-coding RNA H19 can pre- and in castration-resistant prostate cancer. Oncotarget dict a poor prognosis and lymph node metastasis: a meta- 2016;7:75585–602. analysis in human cancer. Minerva Med 2016;107:251–8. 30. Maina PK, Shao P, Jia X et al. Histone demethylase PHF8 reg- 52. Raveh E, Matouk IJ, Gilon M et al. The H19 long non-coding ulates hypoxia signaling through HIF1alpha and H3K4me3. RNA in cancer initiation, progression and metastasis - a Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 20 NEPC lncRNAs proposed unifying theory. Mol Cancer 2015;14:184. dogenous miR-342-3p in gallbladder cancer. J Exp Clin Can- 53. Li H, Zhu L, Xu L et al. Long noncoding RNA linc00617 cer Res 2016;35:160. exhibits oncogenic activity in breast cancer. Mol Carcinog 71. Xu TP, Huang MD, Xia R et al. Decreased expression of the 2017;56:3–17. long non-coding RNA FENDRR is associated with poor prog- 54. Pisarek H, Pawlikowski M, Kunert-Radek J et al. SSTR1 and nosis in gastric cancer and FENDRR regulates gastric can- SSTR5 subtypes are the dominant forms of somatostatin re- cer cell metastasis by affecting fibronectin1 expression. J ceptor in neuroendocrine tumors. Folia Histochem Cytobiol Hematol Oncol 2014;7:63. 2010;48:142–7. 72. Fernando TR, Contreras JR, Zampini M et al. The lncRNA 55. Childs A, Vesely C, Ensell L et al. Expression of somatostatin CASC15 regulates SOX4 expression in RUNX1-rearranged receptors 2 and 5 in circulating tumour cells from patients acute leukemia. Mol Cancer 2017;16:126. with neuroendocrine tumours. Br J Cancer 2016;115:1540–7. 73. Diskin SJ, Capasso M, Schnepp RW et al. Common variation 56. Wyatt AW, Mo F, Wang K et al. Heterogeneity in the inter- at 6q16 within HACE1 and LIN28B influences susceptibility tumor transcriptome of high risk prostate cancer. Genome to neuroblastoma. Nat Genet 2012;44:1126–30. Biol 2014;15:426. 74. Russell MR, Penikis A, Oldridge DA et al. CASC15-S is a tu- 57. Ahlgren G, Pedersen K, Lundberg S et al. Regressive changes mor suppressor lncRNA at the 6p22 neuroblastoma suscep- and neuroendocrine differentiation in prostate cancer after tibility locus. Cancer Res 2015;75:3155–66. neoadjuvant hormonal treatment. Prostate 2000;42:274–9. 75. Nakakura EK, Watkins DN, Schuebel KE et al. Mammalian 58. Wolf DA, Herzinger T, Hermeking H et al. Transcriptional Scratch: a neural-specific snail family transcriptional re- and posttranscriptional regulation of human androgen re- pressor. PNAS 2001;98:4010–5. ceptor expression by androgen. Mol Endocrinol 1993;7:924– 76. Ishii J, Sato H, Yazawa T et al. Class III/IV POU transcrip- 36. tion factors expressed in small cell lung cancer cells are in- 59. Cai C, He HH, Chen S et al. Androgen receptor gene ex- volved in proneural/neuroendocrine differentiation. Pathol pression in prostate cancer is directly suppressed by the Int 2014;64:415–22. androgen receptor through recruitment of lysine-specific 77. Beltran H, Wyatt AW, Chedgy EC et al. Impact of therapy on demethylase. Cancer Cell 2011;20:457–71. genomics and transcriptomics in high-risk prostate cancer 60. Knuuttila M, Yatkin E, Kallio J et al. Castration in- treated with neoadjuvant docetaxel and androgen depriva- duces up-regulation of intratumoral androgen biosynthe- tion therapy. Clin Cancer Res 2017;, 23, 6802–6811. sis and androgen receptor expression in an orthotopic 78. Manohar CF, Furtado MR, Salwen HR et al. Hox gene ex- VCaP human prostate cancer xenograft model. Am J Pathol pression in differentiating human neuroblastoma cells. 2014;184:2163–73. Biochem Mol Biol Int 1993;30:733–41. 61. Wang Q, Zhang J, Liu Y et al. A novel cell cycle-associated 79. Manohar CF, Salwen HR, Furtado MR et al. Up-regulation of lncRNA, HOXA11-AS, is transcribed from the 5-prime end HOXC6, HOXD1, and HOXD8 homeobox gene expression in of the HOXA transcript and is a biomarker of progression in human neuroblastoma cells following chemical induction glioma. Cancer Lett 2016;373:251–9. of differentiation. Tumour Biol 1996;17:34–47. 62. Sun M, Nie F, Wang Y et al. LncRNA HOXA11-AS promotes 80. Hessenkemper W, Baniahmad A. Targeting heat shock pro- proliferation and invasion of gastric cancer by scaffold- teins in prostate cancer. Curr Med Chem 2013;20:2731–40. ing the chromatin modification factors PRC2, LSD1, and 81. Azad AA, Zoubeidi A, Gleave ME et al. Targeting heat shock DNMT1. Cancer Res 2016;76:6299–310. proteins in metastatic castration-resistant prostate cancer. 63. Misawa A, Takayama K, Urano T et al. Androgen-induced Nature Rev Urol 2015;12:26–36. long noncoding RNA (lncRNA) SOCS2-AS1 promotes cell 82. Wang B, Lee CW, Witt A et al. Heat shock factor 1 induces growth and inhibits apoptosis in prostate cancer cells. J Biol cancer stem cell phenotype in breast cancer cell lines. Chem 2016;291:17861–80. Breast Cancer Res Treat 2015;153:57–66. 64. Zhao W, Luo J, Jiao S. Comprehensive characterization of 83. Chen WM, Huang MD, Sun DP et al. Long intergenic non- cancer subtype associated long non-coding RNAs and their coding RNA 00152 promotes tumor cell cycle progression clinical implications. Sci Rep 2014;4:6591. by binding to EZH2 and repressing p15 and p21 in gastric 65. Li Z, Yu X, Shen J. ANRIL: a pivotal tumor suppressor cancer. Oncotarget 2016;7:9773–87. long non-coding RNA in human cancers. Tumour Biol 84. Chen QN, Chen X, Chen Z-Y et al. Long intergenic non- 2016;37:5657–61. coding RNA 00152 promotes lung adenocarcinoma prolif- 66. Aguilo F, Zhou MM, Walsh MJ. Long noncoding RNA, poly- eration via interacting with EZH2 and repressing IL24 ex- comb, and the ghosts haunting INK4b-ARF-INK4a expres- pression. Mol Cancer 2017;16:17. sion. Cancer Res 2011;71:5365–9. 85. Jin HJ, Zhao JC, Wu L et al. Cooperativity and equilibrium 67. Kotake Y, Nakagawa T, Kitagawa K et al. Long non-coding with FOXA1 define the androgen receptor transcriptional RNA ANRIL is required for the PRC2 recruitment to and si- program. Nat Commun 2014;5:3972. lencing of p15(INK4B) tumor suppressor gene. Oncogene 86. Jin HJ, Zhao JC, Ogden I et al. Androgen receptor- 2011;30:1956–62. independent function of FoxA1 in prostate cancer metas- 68. Yap KL, Li S, Munoz-Ca ˜ bello AM et al. Molecular interplay tasis. Cancer Res 2013;73:3725–36. of the noncoding RNA ANRIL and methylated histone H3 87. Yang J, Mani SA, Donaher JL et al. Twist, a master regulator lysine 27 by polycomb CBX7 in transcriptional silencing of of morphogenesis, plays an essential role in tumor metas- INK4a. Mol Cell 2010;38:662–74. tasis. Cell 2004;117:927–39. 69. Yu W, Gius D, Onyango P et al. Epigenetic silencing of tu- 88. Wang J, Nikhil K, Viccaro K et al. The Aurora-A-Twist1 axis mour suppressor gene p15 by its antisense RNA. Nature promotes highly aggressive phenotypes in pancreatic car- 2008;451:202–6. cinoma. J Cell Sci 2017;130:1078–93. 70. Wang SH, Ma F, Tang ZH et al. Long non-coding RNA H19 89. Galvan JA, Astudillo A, Vallina A et al. Epithelial- regulates FOXM1 expression by competitively binding en- mesenchymal transition markers in the differential Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 21 diagnosis of gastroenteropancreatic neuroendocrine prostate cancer. EMBO J 2013;32:1665–80. tumors. Am J Clin Pathol 2013;140:61–72. 108. Seifert A, Werheid DF, Knapp SM et al. Role of Hox genes in 90. Fendrich V, Maschuw K, Waldmann J et al. Epithelial- stem cell differentiation. World J Stem Cells 2015;7:583–95. mesenchymal transition is a critical step in tumorgene- 109. Cho KH, Jeong KJ, Shin SC et al. STAT3 mediates TGF-beta1- sis of pancreatic neuroendocrine tumors. Cancers (Basel) induced TWIST1 expression and prostate cancer invasion. 2012;4:281–94. Cancer Lett 2013;336:167–73. 91. Eide T, Ramberg H, Glackin C et al. TWIST1, a novel 110. Yang X, Duan B, Zhou X. Long non-coding RNA FOXD2-AS1 androgen-regulated gene, is a target for NKX3-1 in prostate functions as a tumor promoter in colorectal cancer by regu- cancer cells. Cancer Cell Int 2013;13:4. lating EMT and Notch signaling pathway. Eur Rev Med Phar- 92. Ainechi S, Mann SA, Lin J et al. Paired Box 5 (PAX5) expres- macol Sci 2017;21:3586–91. sion in poorly differentiated neuroendocrine carcinoma of 111. Zhou W, Ye XL, Xu J et al. The lncRNA H19 mediates breast the gastrointestinal and pancreatobiliary tract: diagnostic cancer cell plasticity during EMT and MET plasticity by and potentially therapeutic implications. Appl Immunohis- differentially sponging miR-200b/c and let-7b. Sci Signal tochem Mol Morphol, 2016. 2017;10. 93. Song J, Li M, Tretiakova M et al. Expression patterns of PAX5, 112. Bauderlique-Le Roy H, Vennin C, Brocqueville G et al. En- c-Met, and paxillin in neuroendocrine tumors of the lung. richment of human Stem-Like prostate cells with s-SHIP Arch Pathol Lab Med 2010;134:1702–5. promoter activity uncovers a role in stemness for the long 94. Czapiewski P, Gorczynski ´ A, Haybaeck J et al. Expression noncoding RNA H19. Stem Cells Dev 2015;24:1252–62. pattern of ISL-1, TTF-1 and PAX5 in olfactory neuroblas- 113. Zhao J, Liu Y, Zhang W et al. Long non-coding RNA toma. Pol J Pathol 2016;67:130–5. Linc00152 is involved in cell cycle arrest, apoptosis, epithe- 95. Kanteti R, Nallasura V, Loganathan S et al. PAX5 is ex- lial to mesenchymal transition, cell migration and invasion pressed in small-cell lung cancer and positively regulates in gastric cancer. Cell Cycle 2015;14:3112–23. c-Met transcription. Lab Invest 2009;89:301–14. 114. Emmrich S, Streltsov A, Schmidt F et al. LincRNAs MONC 96. Walter RF, Mairinger FD, Werner R et al. SOX4, SOX11 and and MIR100HG act as oncogenes in acute megakaryoblastic PAX6 mRNA expression was identified as a (prognostic) leukemia. Mol Cancer 2014;13:171. marker for the aggressiveness of neuroendocrine tumors 115. Yan K, Tian J, Shi W et al. LncRNA SNHG6 is associated with of the lung by using next-generation expression analysis poor prognosis of gastric cancer and promotes cell prolif- (NanoString). Future Oncol 2015;11:1027–36. eration and EMT through epigenetically silencing p27 and 97. Quandt K, Frech K, Karas H et al. MatInd and MatInspec- sponging miR-101-3p. Cell Physiol Biochem 2017;42:999– tor: new fast and versatile tools for detection of consen- 1012. sus matches in nucleotide sequence data. Nucleic Acids Res 116. Erho N, Crisan A, Vergara IA et al. Discovery and valida- 1995;23:4878–84. tion of a prostate cancer genomic classifier that predicts 98. Cartharius K, Frech K, Grote K et al. MatInspector and be- early metastasis following radical prostatectomy. PLoS One yond: promoter analysis based on transcription factor bind- 2013;8:e66855. ing sites. Bioinformatics 2005;21:2933–42. 117. Karnes RJ, Bergstralh EJ, Davicioni E et al. Validation of a 99. Markoff A. Analytical Tools for DNA, Genes and Genomes: genomic classifier that predicts metastasis following rad- Nuts & Bolts, 1st ed. Eagleville, Pa. USA, DNA Press, 2005. ical prostatectomy in an at risk patient population. J Urol 100. Yang L, Lin C, Jin C et al. lncRNA-dependent mechanisms 2013;190:2047–53. of androgen-receptor-regulated gene activation programs. 118. Beltran H, Tagawa ST, Park K et al. Challenges in recog- Nature 2013;500:598–602. nizing treatment-related neuroendocrine prostate cancer. 101. Cui Z, Ren S, Lu J et al. The prostate cancer-up-regulated J Clin Oncol 2012;30:e386–389. long noncoding RNA PlncRNA-1 modulates apoptosis and 119. Terai G, Iwakiri J, Kameda T et al. Comprehensive prediction proliferation through reciprocal regulation of androgen re- of lncRNA-RNA interactions in human transcriptome. BMC ceptor. Urol Oncol 2013;31:1117–23. Genomics 2016;17(Suppl 1):12. 102. Crea F, Watahiki A, Quagliata L et al. Identification of a long 120. Kiryu H, Terai G, Imamura O et al. A detailed investigation non-coding RNA as a novel biomarker and potential ther- of accessibilities around target sites of siRNAs and miRNAs. apeutic target for metastatic prostate cancer. Oncotarget Bioinformatics 2011;27:1788–97. 2014, 5, 764–74; 121. Busch A, Richter AS, Backofen R. IntaRNA: efficient predic- 103. Malik R, Patel L, Prensner JR et al. The lncRNA PCAT29 in- tion of bacterial sRNA targets incorporating target site ac- hibits oncogenic phenotypes in prostate cancer. Mol Cancer cessibility and seed regions. Bioinformatics 2008;24:2849– Res 2014;12(8):1081–7. 56. 104. Wan X, Huang W, Yang S et al. Identification of androgen- 122. Kato Y, Sato K, Hamada M et al. RactIP: fast and accurate responsive lncRNAs as diagnostic and prognostic markers prediction of RNA-RNA interaction using integer program- for prostate cancer. Oncotarget 2016;7(37):60503–18. ming. Bioinformatics 2010;26:i460–466. 105. Lu W, Zhou D, Glusman G et al. KLK31P is a novel androgen 123. Helpap B, Kollermann J, Oehler U. Neuroendocrine differen- regulated and transcribed pseudogene of kallikreins that is tiation in prostatic carcinomas: histogenesis, biology, clin- expressed at lower levels in prostate cancer cells than in ical relevance, and future therapeutical perspectives. Urol normal prostate cells. Prostate 2006;66:936–44. Int 1999;62:133–8. 106. Misawa A, Takayama KI, Fujimura T et al. Androgen- 124. Hirano D, Okada Y, Minei S et al. Neuroendocrine differ- induced lncRNA POTEF-AS1 regulates apoptosis-related entiation in hormone refractory prostate cancer following pathway to facilitate cell survival in prostate cancer cells. androgen deprivation therapy. Eur Urol 2004;45:586–92; dis- Cancer Sci 2017;108:373–9. cussion 592. 107. Takayama K, Horie-Inoue K, Katayama S et al. Androgen- 125. Berruti A, Mosca A, Porpiglia F et al. Chromogranin A ex- responsive long noncoding RNA CTBP1-AS promotes pression in patients with hormone naive prostate cancer Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 22 NEPC lncRNAs predicts the development of hormone refractory disease. J 144. Camacho N, Van Loo P, Edwards S et al. Appraising the rel- Urol 2007;178:838–43; quiz 1129. evance of DNA copy number loss and gain in prostate can- 126. Aggarwal R, Zhang T, Small EJ et al. Neuroendocrine cer using whole genome DNA sequence data. PLos Genet prostate cancer: subtypes, biology, and clinical outcomes. 2017;13:e1007001. J Natl Compr Canc NE 2014;12:719–26. 145. Schuettengruber B, Chourrout D, Vervoort M et al. Genome 127. Yuan TC, Veeramani S, Lin MF. Neuroendocrine-like regulation by polycomb and trithorax proteins. Cell prostate cancer cells: neuroendocrine transdifferentiation 2007;128:735–45. of prostate adenocarcinoma cells. Endocr Relat Cancer 146. Khalil AM, Guttman M, Huarte M et al. Many human 2007;14:531–47. large intergenic noncoding RNAs associate with chromatin- 128. Terry S, Maille´ P, Baaddi H et al. Cross modulation between modifying complexes and affect gene expression. PNAS the androgen receptor axis and protocadherin-PC in medi- 2009;106:11667–72. ating neuroendocrine transdifferentiation and therapeutic 147. Rinn JL, Kertesz M, Wang JK et al. Functional demarcation resistance of prostate cancer. Neoplasia 2013;15:761–72. of active and silent chromatin domains in human HOX loci 129. Huss WJ, Gregory CW, Smith GJ. Neuroendocrine cell differ- by noncoding RNAs. Cell 2007;129:1311–23. entiation in the CWR22 human prostate cancer xenograft: 148. Li J, Han L, Roebuck P et al. TANRIC: an interactive open plat- association with tumor cell proliferation prior to recur- form to explore the function of lncRNAs in cancer. Cancer rence. Prostate 2004;60:91–97. Res 2015;75:3728–37. 130. Vashchenko N, Abrahamsson PA. Neuroendocrine differen- 149. Ren X, Ustiyan V, Pradhan A et al. FOXF1 transcription tiation in prostate cancer: implications for new treatment factor is required for formation of embryonic vasculature modalities. Eur Urol 2005;47:147–55. by regulating VEGF signaling in endothelial cells. Circ Res 131. Aparicio A, Tzelepi V. Neuroendocrine (small-cell) carcino- 2014;115:709–20. mas: why they teach us essential lessons about prostate 150. Tamura M, Sasaki Y, Koyama R et al. Forkhead transcription cancer. Oncology 2014;28:831–8. factor FOXF1 is a novel target gene of the p53 family and 132. Beltran H, Tomlins S, Aparicio A et al. Aggressive vari- regulates cancer cell migration and invasiveness. Oncogene ants of castration-resistant prostate cancer. Clin Cancer Res 2014;33:4837–46. 2014;20:2846–50. 151. Sekaric P, Shamanin VA, Luo J et al. hAda3 regulates 133. Bishop JL, Davies A, Ketola K et al. Regulation of tumor cell p14ARF-induced p53 acetylation and senescence. Onco- plasticity by the androgen receptor in prostate cancer. En- gene 2007;26:6261–8. docr Relat Cancer 2015;22:R165–182. 152. Wang T, Kobayashi T, Takimoto R et al. hADA3 is required 134. Davies AH, Beltran H, Zoubeidi A. Cellular plasticity and the for p53 activity. EMBO J 2001;20:6404–13. neuroendocrine phenotype in prostate cancer. Nat Rev Urol 153. Lin N, Chang KY, Li Z et al. An evolutionarily conserved long 2018;15(5):271–86. noncoding RNA TUNA controls pluripotency and neural lin- 135. Zhang W, Liu B, Wu W et al. Targeting the MYCN-PARP-DNA eage commitment. Mol Cell 2014;53:1005–19. damage response pathway in neuroendocrine prostate can- 154. Lee HJ, Choi NY, Lee SW et al. Epigenetic alteration of im- cer. Clin Cancer Res 2017. printed genes during neural differentiation of germline- 136. Wang C, , Peng G, Huang H et al. Blocking the feed- derived pluripotent stem cells. Epigenetics 2016;11:177–83. back loop between neuroendocrine differentiation and 155. Tsuta K, Wistuba II, Moran CA. Differential expression of macrophages improves the therapeutic effects of enza- somatostatin receptors 1–5 in neuroendocrine carcinoma of lutamide (MDV3100) on prostate cancer. Clin Cancer Res the lung. Pathol Res Pract 2012;208:470–4. 2017;, 24, 708–723. 156. Muscarella LA, D’Alessandro V, la Torre A et al. Gene ex- 137. , Brzezniak C, Oronsky B, Aggarwal RA. Complete metabolic pression of somatostatin receptor subtypes SSTR2a, SSTR3 response of metastatic castration-resistant neuroen- and SSTR5 in peripheral blood of neuroendocrine lung can- docrine carcinoma of the prostate after treatment with cer affected patients. Cell Oncol 2011;34:435–41. RRx-001 and reintroduced platinum doublets. Eur Urol 157. Squires MH, 3rd, Volkan Adsay N, Schuster DM et al. Octre- 2017;, 73, 306–307. oscan versus FDG-PET for neuroendocrine tumor staging: a 138. Thakur MK, Heilbrun L, Dobson K et al. Phase I trial of biological approach. Ann Surg Oncol 2015;22:2295–301. the combination of docetaxel, prednisone, and pasireotide 158. Narayanan S, Kunz PL. Role of somatostatin analogues in in metastatic castrate-resistant prostate cancer. Clin Geni- the treatment of neuroendocrine tumors. Hematol/Oncol tourin Can 2018, Epub ahead of print. Clin North Am 2016;30:163–77. 139. Akamatsu S, Inoue T, Ogawa O et al. Clinical and molecu- 159. Sharma K, Patel YC, Srikant CB. C-terminal region of human lar features of treatment-related neuroendocrine prostate somatostatin receptor 5 is required for induction of Rb and cancer. Int J Urol 2018;, 25, 345–351. G1 cell cycle arrest. Mol Endocrinol 1999;13:82–90. 140. Stone L. Prostate cancer: a novel mechanism of neuroen- 160. Coffey K, Rogerson L, Ryan-Munden C et al. The lysine docrine transdifferentiation. Nat Rev Urol 2018;, 15, 262– demethylase, KDM4B, is a key molecule in androgen recep- 263. tor signalling and turnover. Nucleic Acids Res 2013;41:4433– 141. Huarte M. The emerging role of lncRNAs in cancer. Nat Med 46. 2015;21:1253–61. 161. Yang J, AlTahan AM, Hu D et al. The role of histone 142. Szafranski P, Dharmadhikari AV, Brosens E et al. Small non- demethylase KDM4B in Myc signaling in neuroblastoma. J coding differentially methylated copy-number variants, in- Natl Cancer Inst 2015;107:djv080. cluding lncRNA genes, cause a lethal lung developmental 162. Mo F, Lin D, Takhar M et al. Stromal gene expression is disorder. Genome Res 2013;23:23–33. predictive for metastatic primary prostate cancer. Eur Urol, 143. White NM, Cabanski CR, Silva-Fisher JM et al. Transcrip- 2017;, 73, 524–532. tome sequencing reveals altered long intergenic non- 163. Living Tumor Labratory, www.livingtumorlab.com/index.h coding RNAs in lung cancer. Genome Biol 2014;15:429. tml, Accessed January 2018 Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 23 164. Tsai H, Morais CL, Alshalalfa M et al. Cyclin D1 loss distin- 2018 guishes prostatic small-cell carcinoma from most prostatic 174. Lai D, Proctor JR, Zhu JY et al. R-CHIE: a web server and R adenocarcinomas. Clin Cancer Res 2015;21:5619–29. package for visualizing RNA secondary structures. Nucleic 165. Gao J, Aksoy BA, Dogrusoz U et al. Integrative analysis of Acids Res 2012;40:e95. complex cancer genomics and clinical profiles using the 175. R-chie: A web server and R package for potting arc diagrams cBioPortal. Sci Signal 2013;6:pl1. of RNA secondary structures, http://www.e-rna.org/r-chie/, 166. Cerami E, Gao J, Dogrusoz U et al. The cBio cancer genomics Accessed January 2018 portal: an open platform for exploring multidimensional 176. Mathews DH, Disney MD, Childs JL et al. Incorporating cancer genomics data. Cancer Discov 2012;2:401–4. chemical modification constraints into a dynamic program- 167. cBioPortal for Cancer Genomics, www.cbioportal.org, Ac- ming algorithm for prediction of RNA secondary structure. cessed September 2017. PNAS 2004;101:7287–92. 168. Piccolo SR, Sun Y, Campbell JD et al. A single-sample mi- 177. Gruber AR, Lorenz R, Bernhart SH et al. The vienna RNA croarray normalization method to facilitate personalized- websuite. Nucleic Acids Res 2008;36:W70–74. medicine workflows. Genomics 2012; 100:337–44. 178. RNAfold web server, http://rna.tbi.univie.ac.at/cgi-bin/RNA 169. Iyer MK, Niknafs YS, Malik R et al. The landscape of long WebSuite/RNAfold.cgi, Accessed January 2018 noncoding RNAs in the human transcriptome. Nat Genet 179. Genomatix AG, www.genomatix.de, Accessed August 2017 2015;47:199–208. 180. Ho Sui SJ, Mortimer JR, Arenillas DJ et al. oPOSSUM: identifi- 170. Eisenberg E, Levanon EY. Human housekeeping genes, re- cation of over-represented transcription factor binding sites visited. Trends Genet 2013;29:569–74. in co-expressed genes. Nucleic Acids Res 2005;33:3154–64. 171. Sickle: Windowed Adaptive Trimming for fastq files us- 181. SMALT - Wellcome Sanger Institute Tools and Soft- ing quality, https://github.com/ucdavis-bioinformatics/sick ware, http://www.sanger.ac.uk/science/tools/smalt-0, Ac- le, Accessed January 2018 cessed January 2018 172. UC Davis Bioinformatics Core Software, http://bioinforma 182. Ramnarine VR, Alshalalfa M, Mo F et al. Supporting data for tics.ucdavis.edu/research-computing/software/, Accessed “The Long Noncoding RNA Landscape of Neuroendocrine January 2018 Prostate Cancer and its Clinical Implications” GigaScience 173. Bioinformatics Webserver for RNA-RNA Interaction, http:// Database. 2018. http://doi.org/10.5524/100443. rtools.cbrc.jp/cgi-bin/RNARNA/index.pl, Accessed January Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png GigaScience Oxford University Press

Loading next page...
 
/lp/ou_press/the-long-noncoding-rna-landscape-of-neuroendocrine-prostate-cancer-and-20kpquJDgH
Publisher
Oxford University Press
Copyright
© The Author(s) 2018. Published by Oxford University Press.
eISSN
2047-217X
D.O.I.
10.1093/gigascience/giy050
Publisher site
See Article on Publisher Site

Abstract

Background: Treatment-induced neuroendocrine prostate cancer (tNEPC) is an aggressive variant of late-stage metastatic castrate-resistant prostate cancer that commonly arises through neuroendocrine transdifferentiation (NEtD). Treatment options are limited, ineffective, and, for most patients, result in death in less than a year. We previously developed a first-in-field patient-derived xenograft (PDX) model of NEtD. Longitudinal deep transcriptome profiling of this model enabled monitoring of dynamic transcriptional changes during NEtD and in the context of androgen deprivation. Long non-coding RNA (lncRNA) are implicated in cancer where they can control gene regulation. Until now, the expression of lncRNAs during NEtD and their clinical associations were unexplored. Results: We implemented a next-generation sequence analysis pipeline that can detect transcripts at low expression levels and built a genome-wide catalogue (n = 37,749) of lncRNAs. We applied this pipeline to 927 clinical samples and our high-fidelity NEtD model LTL331 and identified 821 lncRNAs in NEPC. Among these are 122 lncRNAs that robustly distinguish NEPC from prostate adenocarcinoma (AD) patient tumours. The highest expressed lncRNAs within this signature are H19, LINC00617, and SSTR5-AS1. Another 742 are Received: 28 April 2017; Revised: 15 December 2017; Accepted: 1 May 2018 The Author(s) 2018. Published by Oxford University Press. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 2 NEPC lncRNAs associated with the NEtD process and fall into four distinct patterns of expression (NEtD lncRNA Class I, II, III, and IV) in our PDX model and clinical samples. Each class has significant (z-scores >2) and unique enrichment for transcription factor binding site (TFBS) motifs in their sequences. Enriched TFBS include (1) TP53 and BRN1 in Class I, (2) ELF5, SPIC, and HOXD1 in Class II, (3) SPDEF in Class III, (4) HSF1 and FOXA1 in Class IV, and (5) TWIST1 when merging Class III with IV. Common TFBS in all NEtD lncRNA were also identified and include E2F, REST, PAX5, PAX9, and STAF. Interrogation of the top deregulated candidates (n = 100) in radical prostatectomy adenocarcinoma samples with long-term follow-up (median 18 years) revealed significant clinicopathological associations. Specifically, we identified 25 that are associated with rapid metastasis following androgen deprivation therapy (ADT). Two of these lncRNAs (SSTR5-AS1 and LINC00514) stratified patients undergoing ADT based on patient outcome. Discussion: To date, a comprehensive characterization of the dynamic landscape of lncRNAs during the NEtD process has not been performed. A temporal analysis of the PDX-based NEtD model has for the first time provided this dynamic landscape. TFBS analysis identified NEPC-related TF motifs present within the NEtD lncRNA sequences, suggesting functional roles for these lncRNAs in NEPC pathogenesis. Furthermore, select NEtD lncRNAs appear to be associated with metastasis and patients receiving ADT. Treatment-related metastasis is a clinical consequence of NEPC tumours. Top candidate lncRNAs FENDRR, H19, LINC00514, LINC00617, and SSTR5-AS1 identified in this study are implicated in the development of NEPC. We present here for the first time a genome-wide catalogue of NEtD lncRNAs that characterize the transdifferentiation process and a robust NEPC lncRNA patient expression signature. To accomplish this, we carried out the largest integrative study that applied a PDX NEtD model to clinical samples. These NEtD and NEPC lncRNAs are strong candidates for clinical biomarkers and therapeutic targets and warrant further investigation. Keywords: neuroendocrine prostate cancer; transdifferentiation; small cell carcinoma; long non-coding RNA ing data suggest drivers include splice factor SRRM4 [19–21], Introduction master neural transcription factor BRN2 [22], and forkhead box Prostate cancer (PCa) is the most common cancer affecting men A1 (FOXA1) [23]. NEPC tumours have been characterized with (1) and is the third highest cause of cancer death in developed gains in MYCN and AURKA [5]; (2) overexpression of PEG10 [24], countries globally [1]. Advances in detection and treatment for HP1α [25], N-Myc [26, 27], SOX2 [28], and SOX11 [18]; (3) down- PCa have translated to many men being successfully treated regulation of PHF8, KDM3A [29, 30], REST [31], and SPEDF [32]; by surgery and/or radiation. Concomitantly, androgen depriva- and lastly (4) disease dependency on GPX4 [33]. Discoveries such tion therapy (ADT) has resulted in significant survival gains for as these continue to define the protein-coding transcriptome men with metastatic PCa. Commonly administered therapeu- of NEPC. The process of transdifferentiation, however, is highly tics include Enzalutamide, Bicalutamide, and Abiraterone [2]. complex and likely involves multiple layers of genetic and epi- These drugs inhibit the androgen signaling axis, a growth and genetic regulation. differentiation-inducing pathway mediated by the androgen re- Dysregulation of long non-coding RNAs (lncRNAs) could pro- ceptor (AR). Despite these successes, with the steady accumula- vide an additional mechanism for the gene expression alter- tion of facilitating genomic and epigenomic aberrations, a more ations that occur during NEtD. LncRNAs are broadly defined as aggressive tumour capable of growing in castrate levels of testos- large (>200 nucleotides/nt) RNA transcripts, with the most abun- terone can develop [3], termed castration-resistant prostate can- dant subtypes classified as antisense RNAs, pseudogenes, and cer (CRPC). Three main classes of treatment resistance to AR- long intergenic noncoding RNAs (lincRNA) [34]. They are impli- targeted therapies exist, falling into two broad categories asso- cated in a variety of diseases, and their association with cancer ciated with AR signaling [4]. The majority of CRPC reactivate the progression is reported through mechanisms such as remod- AR signaling axis (AR CRPC). However, some tumour cells lever- eling of chromatin, transcriptional co-activation or repression, age their inherent plasticity and progress to an AR-negative state modulation of protein activity, post-transcriptional regulation, − − (AR CRPC), circumventing AR dependence. AR CRPC is highly or as decoy elements [35–37]. LncRNAs form an important regu- heterogeneous, but a major established aggressive subtype is latory layer in global gene expression and, as such, alterations of neuroendocrine prostate cancer (NEPC) [5]. NEPC is pathologi- lncRNAs in cancer are identified as one of the driving forces for cally and clinically similar to small cell carcinoma of the prostate tumorigenesis [38, 39], cancer progression, and metastasis [40, (SCPC), which has been defined as a distinct morphological sub- 41]. More specifically in PCa, lncRNAs have been reported to play type of PCa with neuroendocrine differentiation [6]. Xenograft critical roles at every stage, including the transformation of nor- NEPC models have shown expression of a dominant and irre- mal prostate cells to prostate intraepithelial neoplastic cells, the versible neuronal-like phenotype [7] where conventional CRPC development of localized tumours, and finally progression to ad- therapies are ineffective. Platinum-based chemotherapy is only vanced metastatic disease [42]. These roles in initiation and pro- transiently effective, resulting in poor overall survival [8]with gression are due to aberrant lncRNA expression, which changes most patients surviving ∼7months[9]. Molecular pathology the balance of protein-coding genes involved in processes such markers include expression of chromogranin A (CHGA), synap- as proliferation and apoptosis, thereby facilitating cellular trans- tophysin (SYP), neuron-specific enolase (NSE) [ 10], cell-surface formation. marker CEACAM5 [11], and negative (or low) levels of AR and AR- We recently developed a first-in-field transplantable patient- regulated genes such as PSA [7]. NEPC can arise de novo but much derived xenograft (PDX) model of NEtD: a treatment-na¨ıve ade- more commonly occurs as a consequence of ADT via an adap- nocarcinoma (LTL331) that upon host castration initially re- tive process termed neuroendocrine transdifferentiation (NEtD) gresses (LTL331–8 and 12 week) but then rapidly relapses as ter- [7, 12] and frequently metastasizes to visceral organs [13]. Pre- minally differentiated NEPC (LTL331R) [7]. In our previous study disposing aberrations for NEtD include loss of RB1 [14], TP53 [15], using this model, we demonstrated a lack of evidence for NEPC mutation of Trp53 [16], and/or PTEN inactivation [17, 18]. Emerg- cells before host castration and the conservation of genome Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Table 1: Xenograft model systems used in the study Resistance Molecular Characteristics Model System TMPRSS2- PTEN Name Name Source Phenotype AN TE EZ BI AR PSA SYP SPINK1 ERG ERG PTEN GENE p53 RB LTL313B 313 Primary PCa AD ---- + + - - + + - -/- M WT LTL313BR 313 LTL313B CRPC - ++++ + - - + + - -/- M WT LTL418 418 Primary PCa AD ---- + + - + - - + +/+ WT WT LTL418BR 418 LTL418B CRPC - + - - + + - ————————–Same as parental xenograft above————– LTL331-3 331 Primary PCa AD ---- + + - - + + - -/- M M LTL331-7 331 Primary PCa AD ---- + + - - + + - -/- M M LTL331-5-8 week 331 LTL331-5 AD ---- + - - ————————Same as parental xenograft above————— LTL331-5-12 week 331 LTL331-5 AD ---- + - - LTL331-3-R 331 LTL331-3 NEPC - + - - - - + - - + - -/- M M LTL331-7-R 331 LTL313-3-R NEPC - + - - - - + - - + - -/- M M + − AR and AR CRPC xenograft model samples and their associated molecular characteristics. AD, adenocarcinoma; AN, grows in the absence of androgen; BI, resistant to Bicalutamide; CRPC, castration-resistant prostate cancer; EZ, resistant to Enzalutamide; M, mutation present; NEPC, neuroendocrine prostate cancer; TE, grows in the absence of supplemented testosterone; WT, wild type. Ramnarine et al. 3 characteristics pre- and post-castration, strongly suggesting a phase transition or state change from adenocarcinoma to NEPC [24]. With this model we have identified protein-coding tran- scripts such as PEG10 [24], SRRM4 [19], and HP1α [25]thatare active in the phase transition and validated the discovery of BRN2 [22]. In addition to these, our model has led to the iden- tification of candidate biomarkers and therapeutic targets for NEPC, including the DEK proto-oncogene [43] and epigenetic reg- ulators CBX2 and EZH2 [44] (members of the polycomb group family of transcriptional repressors). We now report the com- prehensive characterization of lncRNAs in our NEtD model. In the current study, we used the longitudinal genomic profiling of our PDX-based NEtD model focusing on lncRNA transcripts. We hypothesized that lncRNA expression across the “time se- ries” would associate with the development of NEPC. Our objec- tive was to comprehensively characterize the dynamic lncRNA landscape of NEtD and NEPC, identify putative functional motifs within lncRNA sequences, determine the clinical relevance of lncRNA expression, and identify associated clinicopathological features. To accomplish this, we implemented a sequence anal- ysis pipeline optimized for the detection of lncRNAs, identified a clinical signature that can robustly distinguish NEPC from ade- nocarcinoma (AD) tumours, and identified four NEtD-associated lncRNA expression profiles. We also identified significant en- richment of well-known transcription factor motifs within the lncRNA sequences. Lastly, we observed that a subset of these lncRNAs is associated with rapid metastasis in treated patients and can stratify tumours based on patient outcome. We present here for the first time a comprehensive landscape of NEPC lncR- NAs and their clinical associations. Results Comprehensive catalogue of lncRNAs in NEPC To identify lncRNAs involved in NEPC, we performed next- generation polyadenylated RNA sequencing on the PDX CRPC models and patient samples. We implemented a sequence anal- ysis pipeline composed primarily of algorithms from the Tuxedo suite of analysis tools [45]. Typically, lncRNAs are expressed at low levels, so the pipeline was augmented to include windowed- adaptive quality control corrections (see Methods and Supple- mentary Figs S1 and S2) that increase the ability to detect low abundance transcripts. We applied this pipeline to all PDX (n = 10) and clinical specimens (n = 117) acquired from the Vancou- ver Prostate Centre (VPC) and Weill Cornell Medicine (WCM) (Ta- bles 1 and 2). Using a quasi de novo mapping strategy combined with amalgamating all sample transcriptome assemblies, we identified 210,999 annotated transcripts spanning 38 Ensembl transcript classes. Defined by Ensembl’s core biotypes, tran- scripts are classified as either protein-coding RNAs, lncRNAs, short ncRNAs, or pseudogenes, which totaled 102,334 (48%), 82,846 (39%), 9,803 (5%),and 16,016 (8%), respectively (Fig. 1A, pie chart 1–2). Within lncRNA, seven subclasses exist: processed transcripts (n = 31,142), retained intron (n = 28,455), lincRNA (n = 12,047), antisense (10,012), sense intronic (n = 821), sense overlapping (n = 340), and three prime overlapping ncRNA (n = 29) (Fig. 1A, pie chart 3; due to their small totals, sense in- tronic, sense overlapping, and three prime overlapping ncRNAs are labeled as “Other”). Despite pseudogenes not being included within Ensembl’s lncRNA classes (listed above), they are by def- inition considered under the umbrella of lncRNA [34]. For each of the eight lncRNA subclasses and their corre- sponding transcripts, we performed unsupervised hierarchi- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 4 NEPC lncRNAs Figure 1: Transcriptome composition and study design. A) Proportions and totals of transcripts detected using our sequence analysis pipeline. Transcripts were separated into protein coding (mRNA) or non-coding RNA (ncRNA) and as defined by Ensembl’s core biotypes as either mRNA, lncRNA, short ncRNA, or pseudo gene. Within lncRNA, there exist seven classes, including processed transcripts, retained intron, lincRNA, antisense, sense intronic, sense overlapping, and three prime overlapping ncRNA (the last three labelled as “other”). Transcript totals are denoted around each pie slice. B) The three transcript classes used in this study due to their ability to separate AD and NEPC tumours, which collectively totalled 37,749 lncRNAs. ∗The pseudogene total was the combination of eight pseudogene subclasses and collectively referred to as pseudogene here. These subclasses include processed pseudogene, unprocessed pseudogene, transcribed unprocessed pseudogene, transcribed processed pseudogene, translated processed pseudogene, polymorphic pseudogene, unitary pseudogene, and pseudogene. These lncRNAs formed the basis for all down-stream analysis and C) the studies project workflow and study design. AUC, area under the curve; GRID, GenomeDx Decipher GRID databa se; JHSM, Johns Hopkins School of Medicine; PDX, patient derived xenograft; ROC, receiver-operating characteristic. See Table 2 for cohort clinical features and compositions. cal clustering (UHC) and principle component analysis (PCA) (n = 12,047),are collectively referred to as lncRNAs here on in (n on the VPC and WCM cohorts (see Methods: Statistical analy- = 37,749 transcripts in total; Fig. 1B). It should be noted that im- sis). Five were incapable of distinctly separating NEPC and AD munoglobulin and T cell receptor genes (n = 326) were removed clinical samples due to insufficient transcript counts, incorrect from the pseudogene transcript total. We explored these lncR- transcript classification, or in general poor transcript annota- NAs in our samples through two analytical workflows (model- tion. The remaining three subclasses were capable of separating based discovery and patient-based discovery), which we later NEPC and AD (Supplementary Figs S3 and S4) and became the fo- merged for clinicopathological analysis. The outline presented cus of all downstream analysis. These three lncRNA subclasses, in this figure represents the study’s overall workflow (Fig. 1C). antisense (n = 10,012), pseudogenes (n = 15,690), and lincRNAs Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Table 2: Clinical samples used in the study Treatment Status Gleason Grade Clinical End Points Cohort +RMET -RMET Institute Name Clinical Group Total Naive NHT ADT RT -6 7 8+ BCR MET PCSM +ADT +ADT VPC VPC AD-NAIVE 56 56 0 0 0 23 0 33 VPCVPC AD-NHT 14 0140 0 0 014 VPC VPC NEPC 50150005 VPC VPC CRPC 53250113 WCM RUBIN NEPC 7 NA NA NA WCM RUBIN AD 30 2235 JHSM LOTAN AD2 17 00 12 JHSM LOTAN NEPC 16 NA NA NA GRID MCI AD 545 0 0 124 54 63 271 211 388 212 132 11 113 GRID MCII AD 232 0 0 77 24 18 117 97 124 75 34 19 24 Total 927 59 17 211 78 107 412 380 512 287 166 30 137 Patient samples and their associated clinical variables including treatment status, Gleason grading, and clinical endpoints. AD, adenocarcinoma; ADT, androgen deprivation therapy; BCR, biochemical recurrence; CRPC, castration- resistant prostate cancer; MET, adenocarcinoma metastasis; NAIVE, adenocarcinoma na¨ıvely treated; NEPC, neuroendocrine prostate cancer; NHT, adenocarcinoma with neoadjuvant treatment; PCSM, prostate cancer specific mortality; +RMET+ADT, ADT-treated rapid metastasis (<36 months) with at least 10 years of clinical follow-up; -RMET+ADT, ADT treated no metastasis with at least 10 years of clinical follow-up. a b Patient overlap exists across different sites of metastasis. Contains a subset of mixed histology tumours (see Methods for breakdown). Cells with a dash are clinical features that are unknown for cohort. Ramnarine et al. 5 A lncRNA expression signature for NEPC − + Recently it has been shown that AR and AR CRPC share sub- stantial genomic overlap yet display significant epigenetic dif- ferences [32]. Here, we hypothesized that the lncRNA transcrip- tome would similarly show unique and common expression al- + − terations between AR and AR CRPC (unexplored to date). To +/− investigate this, we used the AR CRPC xenograft models (Ta- ble 1) to identify changes occurring temporally within the same tumour pre- and post-castration. Once castrated the three AD models (LTL313, LTL418, and LTL331) progress to either AR CRPC (LTL313BR and LTL418BR) or AR CRPC/NEPC (LTL331R). This al- lowed for the identification and quantification of differentially expressed transcripts between pre- and post-CRPC. We inte- grated this data with patient tumour data having matched clini- cal information to ensure the results were clinically relevant and to remove any model-based bias. As we suspected, of all lncR- NAs altered between pre- and post-CRPC (>2 fold, P < 0.05), only 8% (n = 300) were commonly deregulated in both CPRC subtypes. The remaining transcripts (n = 2,669) displayed unique changes + − in the AR or AR CRPC subtype (Supplementary Table S2 and Supplementary Fig. S5). These data support the notion that AR and AR CRPC contain largely distinct lncRNA landscapes. LncRNA expression may be useful as additional biomark- ers beyond those currently used in the diagnosis of NEPC (i.e., CGHA, SYP, and NSE). Moreover, an lncRNA expression signa- ture would strongly support the involvement of lncRNAs in NEPC at a molecular and cellular level. These lncRNAs would be candidates for mechanisms in the activation of a develop- mental pathway and/or plasticity involving previously identified protein-coding genes (PEG10, HP1α, NMYC, SOX2, SRRM4, REST, BRN2, etc.) in NEPC/NEtD. Conversely, since some of these genes (NMYC, SOX2, BRN2, and SRRM4) are well-studied transcription or splicing factors, NEPC lncRNA could be under their regulation. To build an lncRNA expression signature for NEPC, we selected the top fifth percentile of transcripts based on standard devia- tions of expression for the VPC and WCM cohorts independently and performed UHC. All uncharacterized transcripts (i.e., RP##- ######.#, AC######.#, etc.) were removed from the analysis at this point. This produced 265 and 490 NEPC lncRNAs in the VPC and WCM cohorts, respectively. Taking the intersection of these lists and then repeating UHC generated an expression signature of 122 lncRNAs (Supplementary Table S5) that distinctly segre- gated NEPC from AD tumours (Fig. 2A and B). To assess the ro- bustness of this signature, we validated it on an external clin- ical cohort of tumours (n = 33 ; Table 2) from Johns Hopkins School of Medicine (JHSM). These tumours contained 17 AD and 16 NEPC samples and were profiled on the Human Exon array 1.0 ST platform (see Methods) compared with the sequenced dis- covery cohorts. Using the same approach (UHC), a clear sepa- ration of NEPC and AD was observed (Fig. 2C). Observing consis- tent results across different technologies, institutes, and clinical samples further strengthens the robust nature of the NEPC sig- nature. To our knowledge, this is the first report of lncRNAs ex- hibiting a unique, unbiased expression classifier capable of seg- regating NEPC and AD patient samples. Some lncRNA from the patient-derived signature have been previously reported as altered in other cancer types. These lncR- NAs include MALAT1 (alias NEAT2), PCA4 (aliases GDEP, PCAN1, or PCAT4), DSCAM-AS1, and SNHG12. MALAT1 is one of the most well-characterized and studied lncRNAs in cancer and has been identified as a regulator of metastasis and cell migration, a prognostic marker, and a transcriptional regulator of alter- native splicing in lung cancer [46]. PCA4 has been identified Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 6 NEPC lncRNAs Colour Key (A) (B) (C) VPC WCM JHSM −2 −1 0 1 2 Row Z−Score Figure 2: NEPC lncRNA expression signature and clinical classifier. Unsupervised hierarchical clustering of the 122 identified lncRNAs from (A) VPC and (B) WCM cohorts. Validation of this signature was shown in the (C) JHSM cohort. Samples (columns) are labelled as adenocarcinomas (blue) or neuroendocrine (yellow) tumours. See Supplementary Fig. S3 for row/lncRNA labels for each plot. as a prostate and retinal specific transcript [ 47] and frequently uate SSTR-targeted therapy for neuroendocrine tumours in cir- mutated in PCa [48]. DSCAM-AS1 mediates tumour progression culating tumour cells [55], and its use in patient management and tamoxifen resistance in breast cancer through interacting is being tested in a Phase IV clinical trial (NCT02075606). Over- protein hnRNPL [49]. SNHG12 is induced by c-MYC and regu- all, the identification of the NEPC lncRNA expression signature lates cell proliferation, apoptosis, and migration in triple neg- has provided a previously unexplored component of the NEPC ative breast cancer [50]. We were interested in identifying the transcriptome, revealed candidate NEPC biomarkers, and asso- most highly expressed lncRNAs in the signature. Therefore we ciations to NEPC biology. ranked each according to their fold changes when compared with AD samples, required concordance in fold changes across both of the cohorts, and >10-fold change in magnitude. H19, Distinct expression profiles of lncRNAs are associated LINC00617 (alias TUNA/TUNAR), NKX2–1-AS1, and SSTR5-AS1 with NEtD were the only four that fit these thresholds and each with pre- A major goal of this study was to characterize the lncRNA land- vious reports in cancer. Of note, H19 is the most studied among scape during the dynamic phase transition from AD to NEPC the four lncRNAs and is implicated in numerous cancer types using our unique PDX model LTL331 [7](Fig. 3A). To accom- [51]. It is involved in proliferation and both differentiation pro- plish this, we sequenced six samples of our PDX NEtD model cesses related to metastasis, epithelial-to-mesenchymal tran- representing three primary time points along disease progres- sition (EMT), and mesenchymal-to-epithelial transition (MET) sion: two samples from each terminal point (AD and NEPC) and [52]. LINC00617 in breast cancer regulates EMT, cancer progres- two samples post-castration (postTX). Time points 8- and 12- sion, and metastasis through activation of the transcription of week post-castration were selected to represent postTX due to SOX2 [53]. SSTR5-AS1 has not been functionally characterized, tumour volume and serum PSA levels reaching nadir (Table 1 but its sense form SSTR5 has and is a biomarker for neuroen- and Fig. 3A). We identified and quantified all lncRNA transcripts docrine tumours [54]. In fact, recently it has been used to eval- that were altered across the time series and defined four pat- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 AD NEPC Ramnarine et al. 7 Figure 3: Xenograft model of neuroendocrine transdifferentiation, phenotype-driven data integration, and NEtD-associated lncRNAs. A) Schematic depicting the time points at which xenograft tumours were collected along the transdifferentiation of AD to NEPC (adapted from Akamatsu et al., 2015 [24]). B) Phenotypes for clinical samples (coloured circles) that align to various time points from above xenograft model and group-wise comparisons (black connector bars) analyzed for clinical samples. C) Four isolated expression profiles (grey triangles) from select time points in A (light grey circles) with appropriate clinical group-wis e comparisons overlaid and integrated. Unsupervised hierarchical clustering with NEtD lncRNAs (Class I, Deactivated: black bars; Class II, Activated: orange bars; and Class III, Persistent: red bars) identified from integration outlined in (C). Distinct clusters of AD and NEPC clinical samples are observed in the (D) VPC and (E) WCM cohorts. Cla ss IV, Transient lncRNAs were excluded from the clustering due to the lack of clinical samples that would represent this intermediate state. terns of transcript expression: (a) continuous decline in expres- tegration, the following patient group-wise comparisons were sion (Class I, Deactivated, n = 1,613); (b) increasing expression performed: (a) NEPC vs AD, (b) NEPC vs NHT, (c) CRPC vs AD, (d) ¨ ¨ from either AD to postTX or postTX to NEPC (Class II, Activated, NHT vs NAIVE (untreated AD), and (e) NHT vs NAIVE in combi- n = 4,281); (c) continuous increased expression (Class III, Persis- nation with NEPC vs NHT. This produced 1,927;713; 975; 1,045; tent, n = 1,054); and (d) maximum expression at postTX (Class and 117 transcripts, respectively (>2 fold with P < 0.05; total n = IV, Transient, n = 2,668)(total n = 7,627;Fig. 3A). The NEtD model 3,154;Fig. 3B; Supplementary Table S4). Integrating these results and postTX state represents a biological process that currently with the PDX NEtD model transcripts above led to 475, 222, 84, is not characterized as a clinical entity but offers invaluable in- and 45 lncRNAs identified within Class I (Deactivated), Class II sight into the transcriptome of transdifferentiating AD cells. (Activated), Class III (Persistent), and Class IV (Transient), respec- To determine the clinical relevance of Class (I–IV) lncRNAs, tively (total n = 742; Fig. 3C). Unsupervised hierarchical cluster- we integrated patient samples (VPC and WCM, Table 2 - col- ing of Class I-III within WCM (Fig. 3D) and VPC (Fig. 3E) cohorts umn “Clinical group”) with time points in our model (see Fig. exhibited (as expected) a distinct separation of AD and NEPC tu- 3B for alignment of time points to patient groups). Terminal mours and a distinct separation between lncRNAs in Class I-III time points were appropriately aligned to AD and NEPC sam- (rows of heat map). Class IV transcripts were excluded from this ples. However, due to the lack of clinical specimens undergo- illustration due to their lack of altered expression between AD ing NEtD, we hypothesized that neoadjuvant hormone ther- and NEPC clinical samples. Collectively these 742 NEtD lncRNAs apy (NHT) given to AD patients may exhibit characteristics of are associated with the pathogenesis of tNEPC. the postTX state. The transcriptomes from these patients have Prominent examples identified by this biological integration been shown to display the effects of therapy response and more of our NEtD model (Fig. 4A), WCM cohort (Fig. 4B), and VPC co- specifically androgen depletion [ 56]. In fact, neuroendocrine dif- hort (Fig. 4C) illustrate each of these NEtD-defining transcript ferentiation has been shown to increase after only three months classes. PCA3, PCAT1, and PCGEM1 were selected as controls of NHT in a retrospective analysis of 103 radical prostatectomy for this study due to their elevated expression in PCa and high specimens [57]. These early events are the specific alterations we level of characterization. As expected, their expression patterns sought to isolate from the postTX time points of our PDX model. followed the trend in the PCa and NEPC samples (Fig. 4A–C, We also postulated that a subset of Class I (down-regulated in NEtD controls; P < 0.001). SOCS2-AS1 and HOXA11-AS are se- our PDX model) would be up-regulated in the (AR ) CRPC clini- lect examples that characterize the NEtD lncRNA Class I De- cal samples due to reactivation of the AR signalling axis in clas- activated (Fig. 4A–C, Deactivated; P < 0.01). HOXA11-AS, asso- sical CRPC [56, 58–60]. Based on this model-to-patient data in- ciated with the cell cycle through E2F1 [61], has been seen to Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 8 NEPC lncRNAs (A) (B) (C) Figure 4: Select NEtD lncRNAs that exemplify each expression pattern are shown from the (A) NEtD PDX LTL331 model, (B) WCM cohort, and (C) VPC cohort. The expression for NEtD lncRNAs within Class IV, Transient was only identified through the VPC cohort due to the presence of NHT samples, which were not pres ent within the WCM cohort. All boxplots showed significant separation ( P < 0.05) between groups based on a standard Student’s t-test with the exception of ∗ lncRNAs. promote gastric cancer proliferation and invasion (with EZH2) and a focus of our future research and functionalization. Taken and can act as a “molecular sponge” for EZH2 by absorbing (via together, these NEtD lncRNAs (n = 742) characterize the transd- direct interaction) miR-1297 [62]. SOCS2-AS1 is another lncRNA ifferentiation that occurs post-castration and is associated with in this class that has been identified as an AR-regulated tran- tNEPC. script [63] and further supports our hypothesis of AR-regulated lncRNAs within NEtD Class I Deactivated. NKX2–1-AS1 exempli- fies the NEtD Class II Activated (Fig. 4A–C, Activated; P < 0.05) NEtD lncRNAs are enriched with distinct transcription and has been previously seen to characterize lung cancer sub- factor binding motifs types AD and squamous [64]. CDKN2B-AS1 (alias ANRIL) and H19 LncRNAs are not translated and carry out their functions post- are prime illustrations for persistently expressed NEtD Class III transcription in their secondary or tertiary RNA form. This is (Fig. 4A–C, Persistent; P < 0.05). Both of these lncRNAs have been unlike protein-coding transcripts that function in their post- identified across a number of cancer studies (H19 [ 51], ANRIL [65, translational form. Thus, identifying sequence motifs within 66]);however, depending on the cancer type, each has functioned lncRNAs could identify interacting transcripts or proteins that as a tumour suppressor (i.e., ANRIL deactivating tumour sup- provide clues to function. Enrichment of transcription factor pressors CDKN2A/B in cis by three different epigenetic mech- (TF) binding sites (TFBS) was determined by calculating z-scores anisms [67–69]) and as an oncogene (i.e., H19 acts as a sponge for overrepresentation of motifs present in the NEtD lncRNA for FOXM1 by absorbing miR-342–3p [70]). Two demonstrations Classes (I–IV) against their genomic background (Supplementary for transiently expressed NEtD lncRNA Class IV are FENDRR Tables S6–S17 and Methods: Genomatix overrepresented TFBS). and CASC15 (Fig. 4A–C, Transient; P < 0.01). These lncRNAs are We also integrated each of these class-specific enrichment re- well studied in cancer: FENDRR for its prognostic value and its sults to identify unique TFBS for each NEtD Class (Supplemen- involvement in gastric cancer metastasis [71] and CASC15 for tary Tables S18 and S19). All TF and TFBS descriptions, family its regulation of SOX4 in RUNX1-rearranged leukemia [72]and classifications, and annotation is described in Supplementary harboring a risk SNP for susceptibility of neuroblastoma [73]. Table 25. CASC15 has also been identified as a mediator of neural growth In NEtD Class I we identified 33 significant and uniquely en- and differentiation [74], which we believe could be occurring in riched TFBS (Supplementary Tables S6, S18, and S20). Interesting our NEtD model based on the data presented here. Each of these results included binding motifs for TP53, scratch family tran- lncRNAs are among the top candidates identified in this study scriptional repressor 2 (SCRT2), and POU Class III homeobox 3 Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 9 (POU3F3) (z-scores = 4.02, 4.27, and 2.41, respectively). TP53 of- With FOXA1 as one of the characterizing TFBS in Class IV, we ten absent in NEPC, could be an activating TF for many of these sought to explore the persistently expressed (Class III) in con- Deactivated lncRNAs, and suggests an apoptosis or cell cycle ar- junction with the transiently expressed transcripts (Class IV). rest role is present here. Scratch family transcriptional repres- We hypothesized that subsets of these lncRNAs have mechanis- sor 2 has been linked as a neural-specific Snail family transcrip- tic involvement in the transdifferentiation process. To investi- tional repressor and critical for neuronal differentiation [75]. gate this we repeated the TFBS enrichment analysis on Class III Similar to REST, this TF is likely causing the down-regulation of a and IV together and identified six significantly enriched TFBS subset of these lncRNAs. Lastly, POU3F3/BRN1 (a member of the (Supplementary Tables S10, S18, and S20). Confirming our hy- POU family of TFs, as is BRN2) is involved in the development pothesis was the presence of TWIST1 (z-score = 4.01), an es- of the nervous system, expressed in small cell lung cancer cells sential member of the EMT transcriptional reprogramming fac- (which has pathological overlaps with NEPC), and involved in tors [87]. Interestingly, TWIST1 and AURKA have very recently proneural/neuroendocrine differentiation [76]. Considering this been seen to form a feedback loop promoting metastasis and and the significant enrichment of these TF motifs, this suggests highly aggressive phenotypes in pancreatic carcinoma [88], and a role in proliferation and differentiation in NEtD Class I. TWIST1 is a marker for EMT in neuroendocrine tumours [89, Performing TFBS enrichment analysis in NEtD Class II and III 90]. Concerning PCa, it has been identified as AR regulated (and identified 12 and 15 significant and distinct TFBS motifs, respec- repressed via NKX3–1), whereas in the absence of AR it is up- tively (Supplementary Tables S7, S8, S18, and S20). Interestingly, regulated and present in metastatic disease [91]. both classes had significant enrichment for at least one ETS and We further investigated global functional characteristics HOX family member, suggesting overlapping functional roles for across all NEtD lncRNAs. Specifically, we wanted to identify TFBS their respective lncRNAs. For Class II this included ELF5, SPIC, that were significantly enriched and common across all classes. and HOXD1 (z-scores = 2.17, 2.43, and 2.22, respectively) and Due to the high number of lncRNAs (n = 2,147), we decided to for Class III, PDEF (alias SPDEF) and HOX/PBX (z-scores = 3.03 perform this analysis at the TF family level; therefore for each and 2.71, respectively). Members of the ETS family fused to TM- class and the full lncRNA set, we repeated the motif enrichment PRSS2 is the most frequent genomic alteration in PCa; therefore analysis and integrated all of their results (Supplementary Ta- the prevalence of their motifs in these classes is not surprising. bles S12–S17, S19, and S20). We identified 62 significantly com- While the ETS fusion transcript is relatively more specific to PCa mon TFBS families (z-score = >2;Supplementary Table S20). Not vs NEPC, ETS TFs on their own are involved in a wide variety surprising were families involving cell cycle regulation, cyclin of functions, including cellular differentiation and angiogene- B2/CCNB2 and the E2F family (z-scores = 40.66 and 67.36, re- sis. In fact, recently SPDEF was found to be down-regulated in spectively). We also observed both the ETS (z-score = 9.9) and metastatic NEPC due to DNA methylation [32] and was also sig- REST (z-score = 20.52) families of TFs, which reaffirmed our hy- nificantly down-regulated in treated vs untreated high-risk PCa pothesis that these lncRNAs are involved in tumour progres- patients [77]. Conversely, the HOX family has never been linked sion and neuronal pathways. Surprising was the presence of two to PCa or NEPC for that matter, and so this result was unex- PAX families, PAX5 (z-score = 10.15) and PAX9 (z-score = 18.18). pected. In neuroblastoma, however, the HOX genes have been The PAX family is known to regulate lineage specification and linked to differentiating cells [78] and specifically HOXD1 iden- progenitor cell maintenance. In developmental biology, PAX5 is tified here (as well as HOXC6 and HOXD8) are associated with involved in B-cell differentiation and PAX9 in neural crest de- differentiation towards a neuronal phenotype [79]. velopment. PAX5 has been observed as overexpressed in other Performing TFBS enrichment analysis in NEtD Class IV iden- neuroendocrine tumours [92, 93], overexpressed in neuroblas- tified enrichment of 17 distinct TFBS motifs (Supplementary Ta- toma [94], and shown to positively regulate c-Met transcription bles S9, S18, and S20). Class IV transcripts are only expressed in small cell lung cancer [95]. In lung NETs, PAX6 is prognostic during treatment (castration) response. Interestingly, heat shock for aggressiveness [96]. Their role in NEPC is yet to be charac- TFs HSF1 (z-score = 2.06) and HSF2 (z-score = 3.87) were within terized; however, evidence here supports their global involve- these results. Heat shock proteins (HSPs) are expressed at low ment in NEtD and lncRNA function. Lastly, the selenocysteine levels under normal conditions, up-regulated by cellular stress, tRNA activating factor (STAF, z-score = 15.18) was very intrigu- and function as molecular chaperones to control client protein ing to us. A recent Nature study by Schreiber et al. suggested stability and function. Their candidacy as therapeutic targets that treatment resistance in NEtD of PCa depends on a drug- has been well studied in PCa [80]and AR CRPC [81]. In breast gable lipid-peroxidase pathway that protects against ferropto- cancer, HSF1 specifically induces a cancer stem cell phenotype sis (a non-apoptotic form of cell death) [33]. The increased lipid in vitro [82]. In PCa, HSPs bind dihydrotestosterone to the AR metabolism creates a dependency on GPX4, which prevents fer- and enhance AR-mediated transcription. One of the functions roptosis. GPX4 is a selenocysteine-containing enzyme and 1 of of lncRNAs is to facilitate this type of mechanism. For exam- only 25 proteins with this rare amino acid in the entire human ple, LINC00152/CYTOR (identified within this class) binds and re- genome. The data presented here suggest that some of these cruits EZH2 to its target promoters p15 and p21 in gastric cancer lncRNAs may be involved in the selenocysteine pathway via [83] and IL24 in lung cancer [84] and thereby causes repression STAF and in selenoprotein biosynthesis of molecules (i.e., GPX4). of their expression. Considering the transient expression of the Identifying and targeting these lncRNAs could be a path for up- lncRNAs in this class, this data suggest a subset may be stress stream inhibition of GPX4 up-regulation and therefore allow cell response mediators via HSPs. Lastly, FOXA1 showed a significant death in these resistant cells to occur naturally by ferroptosis. enrichment (z-score = 2.8) in this class. Recently, FOXA1 loss was Comprehensive in vitro experimentation would need to be car- identified as a driver of NEtD [ 23], which leads to AR reprogram- ried out to confirm this therapeutic avenue. ming [85] and EMT through direct regulation of SLUG expression [86]. This suggests that some of the lncRNAs in this class could NEtD lncRNAs contain NEPC-related TFBS have a functional role in maintaining cellular identity when un- der the control of FOXA1. It is now well established that complex cellular reprogramming occurs during NEtD, and master regulators such as REST [31], Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 10 NEPC lncRNAs BRN2 [22], SOX2 [28], and SOX11 [18] have been identified as key AS1 regulates EMT and Notch signaling to promote colorectal TFs involved in this process. Identification of well-known TFBS cancer [110]. H19 has been identified as a mediator of breast such as these would test our current hypotheses on the func- cancer plasticity during EMT and its reverse process MET [111], tional involvement of individual lncRNAs in NEtD pathogenesis as well as having a role in stemness in prostate cells [112]. (Supplementary Tables S21–S24). TFBS identification was carried LINC00152 is involved in EMT (combined with cell migration and out using MatInspector [97–99] (see Methods: Genomatix Mat- invasion) in gastric cancer [113]. LINC00478 (alias MONC) inter- base and MatInspector) on each of the NEtD Class for select TFs. feres with hematopoietic lineage decisions and enhances pro- As before, all TF and TFBS descriptions, family classifications, liferation of immature progenitor cells in acute megakaryoblas- and annotation is described in Supplementary Table 25. tic leukemia [114]. Lastly, again in gastric cancer, lncRNA SNHG6 With the dominance of AR-regulated genes in AD, the lack of has been seen to promote cell proliferation and EMT [115]. Based expression that defines Class I NEtD lncRNAs is likely caused by on these data, Class III and IV lncRNAs could have a role in de- the absence of androgen (post-castration), and therefore these veloping a cellular “plastic” state during NEtD. lncRNAs are putatively AR-regulated. To test this, we searched To test the involvement of known NEPC-involved TFs in all for androgen and glucocorticoid response elements (ARE and NEtD lncRNA, we searched for BRN2, ARE/GRE, REST, SOX11, GRE, respectively). The results showed that 107 lncRNAs con- SOX2, NMYC, ETVI, ETS, and NKX3 motifs (Supplementary Ta- tained ARE and/or GRE motifs, of which 16 contained only an ble S24). Since each of the classes had different sizes, this would ARE motif, 49 contained only a GRE motif, and 21 contained both influence the distribution and presence/absence of these mo- ARE and GRE motifs (Supplementary Table S21). To further test tifs, so we extracted the top 25 lncRNA within each class (n = and support our AR-regulated lncRNA hypothesis, we explored 100 NEtD lncRNA), ranked by their magnitude of fold change. all previously reported AR-regulated lncRNAs. Currently, the fol- Observing the distribution of these TFs separated by NEtD class lowing 17 lncRNAs have been identified with experimental evi- revealed an interesting pattern (Fig. 5B). TFs SOX2, SOX11, and dence: PCGEM1 [100], PlncRNA-1/CBR3-AS1 [101], PCAT-18 [102], REST had a relatively more balanced distribution across each PCAT29 [103], SOCS2-AS1 [63], RP1–45I4.2 [104], SUZ12P1 [104], class compared to NKX3, ETSF, ETVI, and NMYC, which showed SNHG5 [104], LINC01138 [104], SNHG1 [104], KLKP1 [104, 105], a preference to binding persistent and transiently expressed LINC00969 [104], LINC-PINT [104], TUG1 [104], MIR17HG [104], lncRNA. Interestingly, more than 50% of ARE/GRE motifs were POTEF-AS1 [106], and CTBP1-AS1 [107]). Of these, 4 (PCAT29, present in transiently expressed lncRNA vs relatively few in SUZ12P1, SNHG1, and CTBP1-AS1) were not within our pipeline Class I Deactivated and Class II Activated. Conversely, BRN2 mo- lncRNA class annotation, and 8 of 13 (61%) were represented in tifs were relatively more present in Class I and II. These patterns NEtD Class I Deactivated lncRNAs (PCGEM1, PlincRNA-1, PCAT- suggest a time-dependent or cellular phase-dependent usage of 18, SOCS-AS1, KLKP1, LINC00969, LINC-PINT, and POTEF-AS1). TFs post-castration and during the NEtD process. Due to our integrative study design (Fig. 3A–C), the remaining 5 did not move forward in the analysis. However, removing the NEPC and NEtD lncRNAs identify putative NEPC integrative steps, down-regulation of these lncRNAs did occur subtypes in either our model or patient samples independently. Overall, of the 13 lncRNA annotated by our pipeline and reported as AR To corroborate the lncRNA expression in an external NEPC regulated, all overlapped in this study. (extNEPC) cohort [32], we visualized NEtD lncRNA Classes II-IV Due to the elevated pattern of expression that defines Class and the up-regulated NEPC lncRNA expression, including ge- II and III NEtD lncRNAs, we hypothesized that a subset of these nomic profiles (copy number and mutation), through an Onco- lncRNAs are constituents of the neuronal phenotype present in Print schematic. The cohort consisted of 44 NEPC specimens NEPC. To test this hypothesis, we analyzed these lncRNAs for (largest published to date) from 30 patients that were classi- the presence/absence of the following select TFs known to ac- fied based on their histomorphology [ 6]. Duetothe exomese- tivate this cell type: APOU Class III homeobox 2 (POU3F2), also quencing performed on this cohort, not all lncRNAs were rep- known as BRN2, and RE1 silencing transcription factor (REST). resented/detectable in this sequencing profiling. We also plot- Activation of BRN2 and deactivation of REST are involved in ted previously reported NEPC predisposing genes, oncogenes, neuronal differentiation and regulation of neurogenesis, respec- drivers, and the TFs we identify above to provide “transcriptome tively. Again, using the MatInspector algorithm, we identified 11, context” for the altered lncRNAs (Supplementary Figs S6 and 22, and 21 lncRNAs in Class II or III with TFBS for BRN2, REST, S7). From the 58 NEPC and 243 NEtD lncRNAs represented in or both, respectively (Supplementary Table S22). Taken together, the extNEPC exome sequencing profiling, 43% (25/58) and 27% this evidence supports involvement for a subset of these lncR- (66/243) showed altered expression in 2–34% of NEPC patients, NAs to neuronal function/pathways in NEtD. respectively (Supplementary Figs S6, S8–S11). To further support the hypothesis of mechanistic involve- Surprisingly, these testable lncRNAs (58 from the NEPC lncR- ment for the NEtD process in Classes III and IV, we expected TFBS NAs and 66 from the NEtD lncRNA) in combination with known related to plasticity and stemness to be present. Therefore, we oncogenes/tumour suppressors/transcription factors (Supple- used MatInspector to identify binding motifs for members of the mentary Fig. S7) resulted in identifying three distinct subsets of following well-studied cellular differentiation TF families: HOX NEPC patients within the extNEPC cohort. Group 1 had relatively [108], SOX, STAT3 [109], and “STEM” (STEM members are defined higher mutation frequencies, higher ploidy, mixed tumour sites, by Matbase and include POU5F1/OCT4, SALL4B, SOX2, NANOG, and mixed pathological classifications. Group 2 had a relatively and TCF7L1). We observed 42, 49, 30, and 33 lncRNAs with TFBS low mutation frequency and low ploidy, derived mostly from for HOX, SOX, STAT3, and STEM genes, respectively (Supplemen- pelvic masses and with pathological classification D (large-cell tary Table S23). In fact, some lncRNA had TFBS within more neuroendocrine carcinoma). Group 3 tumours, however, were than one of these TFs (Fig. 5A). Previous studies have linked 6/7 mostly derived from the prostate with a pathological classifica- of these (highlighted in Fig. 5A) to various components of EMT tion B, and likely primary (de novo) NEPC samples where NEtD and/or cellular plasticity. FENDRR (antisense lncRNA to FOXF1) has not occurred. Of note, copy number loss or mutations in regulates gastric cancer metastasis via fibronectin1 [ 71]. FOXD2- TP53 and RB1 were present in 60% of patients (26/44), spread Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 11 (A) SOX STAT3 HOX STEM FOXD2-AS1 5 1 1 0 7 1 SNHG6 TDRG1 LINC00152 H19 0 2 FENDRR LINC00478 0 10 (B) NKX3 ETSF ETVI Class I NMYC Class II SOX2 SOX11 Class III NRSF Class IV ARE+GRE BRN2 0% 20% 40% 60% 80% 100% Figure 5: TFBS Venn diagram and distribution plots. A) Common and unique TFBS for HOX, SOX, STAT3, and STEM families of transcription factors within Class III and IV of the NEtD lncRNAs. B) Distribution of TFBS for known NEPC-involved TFs within NEtD Class lncRNAs. across the cohort, and did not appear to be associated with a marily with adverse pathology (i.e., high grade/stage) and long- particular group (Supplementary Fig. S7). The three groups could term follow-up for treatment and outcomes (median 18 years). be revealing an lncRNA expression signature that is specific for From these cohorts, a subset (n = 211) received adjuvant ADT tumour site, degree of genomic mutations (SNPs or CNVs), and post-radical prostatectomy (RP). To determine the most clini- pathological classification. However, it is important to note this cally relevant lncRNA transcripts, we first ranked the NEtD/NEPC is an observational result requiring statistical validation in a lncRNAs within their respective classes and selected the top larger cohort. The specificity of these genomic and lncRNA tran- deregulated from each. The ranking was performed based on scriptome profiles would need to be explored across a variety fold changes observed within the clinical groups (see Methods). of metastatic sites and NEPC pathologies to validate these three This produced 100 top-ranking NEtD/NEPC lncRNAs that we in- novel NEPC molecular subtypes. vestigated within the GRID cohorts (Fig. 1C and Supplementary Table S26). We validated 11 of these (2 from each NEtD Class and 3 from the NEPC lncRNA signature) by quantitative real-time PCR to confirm expression changes identified in the model and clin- NEPC and NEtD lncRNAs are associated with ical samples (Supplementary Fig. S12). Due to the difference in treatment-related metastasis profiling platforms between GRID (Affymetrix microarray) and Prognostic and predictive biomarkers for NEtD and NEPC are in VPC/WCM cohorts (Illumina Sequencing), it was necessary to dire need since ADT is not effective for a cancer that has un- remap the GRID microarray probes (see Methods) that aligned dergone NEtD and thus circumvents the AR signalling axis. We within NEtD/NEPC lncRNA sequenced regions. This resulted in examined if the NEtD (n = 742) and NEPC (n = 122) lncRNAs are 81/100 being present and quantifiable on the microarray plat- associated with NEPC related clinical outcomes in patients with form. primary prostatic adenocarcinoma. To accomplish this, we ex- A characteristic of NEPC patients in the clinic is the occur- plored the candidates in two cohorts from the Mayo Clinic (MCI rence of rapid metastasis following treatment [118], and so we [116]and MCII [117]) from the Decipher GRID database (GRID) first tested the lncRNAs’ ability to predict rapid metastasis post- (n = 777, Table 2). We could not perform this analysis within ADT. We performed receiver-operating characteristic (ROC) anal- VPC/WCM cohorts due to their small sample sizes and short- ysis to compare the sensitivity and specificity of predicting rapid term clinical follow-up. The GRID cohorts represent tumours pri- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 12 NEPC lncRNAs (A) 0 (events = 15) 0 (events = 13) 1 (events = 22) 1 (events = 20) p−value: 0.005 p−value: 0.010 0 20 40 60 80 100 120 140 020 40 60 80 100 120 140 Time(months) RP−Met Time(months) RP−Met Group=0 21 20 15 15 11 9988 2 Group=0 21 21 16 16 12 10 9 7 6 Group=1 22 18 11 86 5 3 1 Group=1 22 17 10 7 5 4 3 222 SSTR5-AS1 LINC00514 (B) 0 (events = 20) 0 (events = 19) 1 (events = 20) 1 (events = 21) p−value: 0.900 p−value: 0.832 0 20 40 60 80 100 120 140 020 40 60 80 100 120 140 Time(months) RP−Met Time(months) RP−Met Group=0 28 23 20 19 14 11 10 986 Group=0 28 23 20 19 14 11 10 9 9 4 Group=1 28 24 23 17 15 12 10 983 Group=1 28 24 23 17 15 12 10 9 7 5 VPC WCMC (C) SSTR5-AS1 LINC00514 SSTR5-AS1 LINC00514 Figure 6: Kaplan-Meier estimates and expression for SSTR5-AS1 and LINC00514. Kaplan-Meier estimates for metastasis-free survival in the MCII cohort comparing low (blue lines) and high (yellow lines) expression (split by median) in (A) treated patients that received post-prostatectomy adjuvant ADT for SSTR5-AS1 (left) and LINC00514 (right) and (B) patients not receiving ADT treatment. C) Box plot expression for the top two NEPC lncRNA candidates (SSTR5-AS1 and LINC00514) within the VPC and WCM cohorts. metastasis (within 36 months) for each lncRNA. We then cal- The expression of two NEtD/NEPC lncRNA transcripts (SSTR5- culated the area under the curve (AUC) for each lncRNA ROC AS1 and LINC00514) was able to separate patients more likely in both cohorts using probe set region expression summarized to develop metastatic disease from those that did not (P = 0.005 across the full lncRNA transcript (Supplementary Table S26). and P = 0.010, respectively; Fig. 6A). To increase our confidence This identified eight lncRNAs with the highest scores: NR2F1- that the results were associated with treatment status, we gen- AS1, LINC00654, FENDRR, PCAT2, and NKX2–1-AS1 in MCI (AUC erated Kaplan-Meier estimates for these transcripts in untreated >0.70) and LINC00478, LINC00173, and LINC00514 in MCII (AUC patients from the same cohort; neither showed significant sep- > 0.70). These lncRNAs serve as candidates for predicting rapid aration in their performance (P = 0.905 and P = 0.832, respec- metastasis in patients receiving ADT. Selecting all NEtD/NEPC tively; Fig. 6B). Expression for SSTR5-AS1 and LINC00514 in the lncRNAs with AUC >0.65 (n = 25), we performed survival analy- VPC and WCM cohorts illustrate their distinct and elevated ex- sis to ascertain their ability to separate patients for metastasis as pression in NEPC vs. AD patient samples (Fig. 6C). These results an outcome and end-point. Specifically, we calculated Kaplan- suggest a strong association between treatment status and in- Meier estimates for metastatic disease progression stratified by creased probability of metastatic disease in patients with dif- median expression in ADT-treated samples of the MCII cohort. ferential expression of these lncRNAs. This, together with re- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 AD NEPC EXPRESSION (LOG2) Probability of Mets Free Survival Probability of Mets Free Survival 02468 10 12 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 Probability of Mets Free Survival Probability of Mets Free Survival 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 Ramnarine et al. 13 sults from the NEtD model and NEPC clinical samples, implicate In this study, we characterized the unexplored global lncRNA SSTR5-AS1 and LINC00514 in NEPC and these lncRNAs serve as landscape during NEtD to provide insights into the NEPC non- strong candidates as predictive biomarkers for metastatic dis- coding milieu of this lethal and treatment-induced process. This ease post-RP following ADT. required the implementation of a sequence analysis pipeline One of the mechanisms observed with lncRNAs is direct with increased sensitivity towards lower expressed transcripts, RNA-RNA interaction with mRNA, resulting in regulation of their characteristic of lncRNAs. The pipeline was able to detect 37,749 expression (activation or repression). This type of investigation lncRNA transcripts (subclassified as either lincRNA, antisense, is computationally intensive, and there are limited algorithms or pseudogene) and quantify them in the two clinical cohorts available to identify putative mRNA targets genome-wide. How- (VPC and WCM). The novelty of this study lies in the use of pa- ever, a method was recently published to predict lncRNA-mRNA tient samples integrated with the NEtD PDX model to detect clin- interactions genome-wide [119], and so we sought to identify ically relevant lncRNAs involved in the NEtD/phase transition candidate mRNA transcripts interacting with SSTR5-AS1 and process. In this study, we identified 742 lncRNAs associated with LINC00514. The pipeline’s three core algorithms include Rac- NEtD and identified a robust 122 NEPC lncRNA patient signature cess [120] for the identification of accessible regions within the capable of classifying NEPC from AD patient samples. The motif lncRNA, IntaRNA [121] to calculate nucleotide interaction en- analysis identified significantly enriched TFBS that were unique ergies, and RactIP [122] to predict joint secondary structures. to NEtD Classes I Deactivated (TP53 and BRN1), II Activated (ELF5, Applying this methodology to SSTR5-AS1 and LINC00514 pro- SPIC, and HOXD1), III Persistent (SPDEF and HOX), IV Transient duced a list of predicted interacting partners for these lncRNAs (TP53, HSF1, HSF2, and FOXA1), and III Persistent combined with (Supplementary Tables S27–S28). The top-ranked mRNAs were IV Transient (TWIST1). Through similar analysis, we also iden- KDM4B and TADA3 that are predicted to hybridize and form joint tified common TFBS (CCNB2, E2F, ETS, REST, PAX5, PAX9, and structures independently with SSTR5-AS1 and LINC00514, re- STAF) enriched across all of the NEtD lncRNAs. From among the spectively (Supplementary Figs S12 and S13). In the clinical co- 100 top-ranking lncRNA, we observe that a subset have strong horts, TADA3 is down-regulated in NEPC vs AD (>2-fold), while clinical associations with metastatic PCa patients after receiv- KDM4B is up-regulated (>5 fold); however, only the deregulation ing ADT. In previous lncRNA studies in cancer, several have been of TADA3 is statistically significant (VPC P = 0.003 and WCM P = linked to malignant transformation with key roles affecting vari- 0.017). Both genes have NEPC associations (see Discussion), and ous aspects of cellular homeostasis, including proliferation, sur- our data suggest they are being regulated by these lncNRAs. vival, migration, and genomic instability [141]. Similarly, lncR- NAs identified in this study, including SSTR5-AS1 and LINC00514 with their association with poor outcome, FENDRR for its asso- ciation with rapid metastasis, and H19 and LINC00617 for their Discussion concordantly high expression across both of the discovery co- Primary NEPC arises de novo in 0.5% to 2% of all prostate can- horts, could be the missing links in the mechanisms causing cer patients [123]. However, tNEPC can develop in 20–30% of NEtD. These five represent the top candidates discovered in this metastatic castrate-resistant prostrate cancer tumours [124]and study due to this evidence but also for their characterization in increases with disease progression [125]. The real incidence of other cancer types. tNEPC may be higher because of under-recognition due to tu- FENDRR is a top deregulated lncRNA in NEtD Class IV Tran- mour heterogeneity, the limited number of metastatic tumour sient and may have a role in the NEtD process. It is implicated biopsies performed, lack of uniform consensus definition based in a lethal lung development disorder [142], lung cancer [143], within a mutational hotspot that is copy number lost in PCa on histology or biomarker expression, and frequent misclassifi- cation as high-grade PCa (most notable in tumours with mixed [144], and can bind to PRC2 [145, 146]. PRC2 plays a significant histologies) [126]. NEPC can be induced in vitro in AR LNCaP role in tumour progression through binding of HOTAIR (a very cells in androgen-depleted culture conditions [127, 128], sim- well-studied lncRNA). Together, HOTAIR and PRC2 are involved ilarly in vivo [7, 129], and in patient tumours long-term ADT in the control of chromatin structure and associated gene activ- has increased neuroendocrine differentiation [118, 124, 130]. It ity [147]. FENDRR may be involved in tumorigenesis like HOTAIR is now common to observe treatment-resistant tumours with due to its known interaction with PRC2. A recent study showed neuroendocrine features upon metastatic biopsy, and the pre- down-regulation of FENDRR is associated with poor prognosis vailing consensus is that epithelial plasticity enables tumour in gastric cancer and regulates cancer cell metastasis through adaptation in response to AR-targeted therapies [7, 9, 118, 126, fibronectin [ 71]. Functionally, this could be occurring in NEPC as 131–134]. This evidence supports the notion that tNEPC inci- well due to FENDRR’s transient expression in the NEtD model dence through NEtD will increase as new powerful ADTs en- and its association to rapid metastasis in ADT-treated PCa pa- ter the clinic. There is an urgency for therapeutic strategies and tients from the GRID (MCI) cohort. Another putative function of this transcript is through upregulating FOXF1, which is a clinical biomarkers defining NEtD/NEPC. Currently, the only op- tion for patients is the short-lived effects of platinum-based protein-coding gene and the sense form for the antisense tran- chemotherapy. Optimism is on the rise, as there is an AURKA in- script FENDRR. Antisense transcripts are known to regulate their hibitor (MLN8237) in a Phase 2 clinical trial (NCT01799278), com- sense forms (positively or negatively). Using TANRIC, an interac- binatorial approaches using AURKA with PARP inhibitors un- tive resource for the exploration of lncRNAs in large patient co- der investigation [135], indirect methods that resensitize the tu- horts within 20 TCGA cancer types [148], we see that FENDRR ex- mour to Enzalutamide [136] or platinum-based chemotherapy pression is positively correlated to FOXF1 in 16 of 20 cancer types −9 [137] (Phase 2 clinical trial NCT02489903 with a Phase 3 clinical (P < 3.71 × 10 ; Supplementary Table S29). In fact, FOXF1 dele- trial being planned), a SSTR4/5 analogue (Pasireotide/SOM230) in tion has been seen to significantly reduce FENDRR in endothe- four independent clinical trials at various Phases (NCT01646684, lial cells [149]. FOXF1 is also a target gene of p53 and is seen to NCT01313559, NCT01468532, and NCT01794793) with one al- regulate cancer cell migration and invasiveness [150]. Together these transcripts may play a transient coordinated role in NEtD ready reporting promising clinical efficacy [ 138], and increased study of NEPC/NEtD in general [134, 139, 140]. through PRC2 or fibronectin. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 14 NEPC lncRNAs LINC00514 is amongst the highest expressed lncRNAs in functions (SSTR5-AS1: SSTR5 or SSTR5-AS1: KDM4B: N-Myc) re- NEtD Class III Persistent. It has not been characterized. It is pre- quire thorough in vitro and in vivo exploration to ascertain their dicted to bind to TADA3 (Supplementary Fig. S14), potentially validity. causing a reduction of its activity. This is intriguing because Although multiple layers of genetic and epigenetic deregula- TADA3 is involved in the stabilization and activation of p53 [151, tion likely cooperate to facilitate NEtD, understanding the non- 152], and this putative interaction (LINC00514: TADA3) could be coding contribution to this multifarious process is necessary to an alternative mechanism for loss of p53 activity, already known design effective novel therapeutics. Using the five independent to be frequently lost in NEPC [15]. H19 and LINC00617 were two patient cohorts and our proven NEtD PDX LTL331 model, lncR- of the four highest (>10-fold) NEPC-expressed lncRNAs in this NAs such as FENDRR, LINC00514, LINC00617, H19, SSTR5-AS1, study and fortunately (unlike most lncRNAs) have both been and others identified in this study may provide more in-depth thoroughly characterized functionally. LINC00617 is highly con- insights into NEtD and NEPC. Research identifying the relation- served across vertebrate genomes and is required for mainte- ship of these lncRNAs to other known drivers, oncogenes, and nance of pluripotency and neural differentiation in embryonic Activated pathways in NEtD is now required. This study is the stem cells [153]. It controls this lineage commitment through first to report the lncRNA landscape of NEtD, a robust NEPC RNA-binding proteins PTBP1, hnRNP-K, and Nucleolin. These lncRNA expression clinical classifier, and provides numerous RNA-binding protein complexes have been detected at promot- candidate biomarkers and therapeutic targets. ers of NANOG, SOX2 (promoter of lineage plasticity in NEPC [28]), and FGF4 [153]. H19 has also been identified in neural differen- Methods tiation of pluripotent stem cells [154] but with unknown mech- anisms. With such an elevated level of expression in the clini- PDXs cal cohorts (∼30- to 40-fold and ∼20- to 30-fold in VPC/WCM for Animal ethics, care, experiments, xenograft generation, and all LINC00617 and H19, respectively), these lncRNA could be respon- protocols were carried out in accordance with the guidelines of sible for maintaining the neuronal component of NEPC through the Canadian Council of Animal care as previously described [7]. epigenetic regulation. Specific xenograft models used in this study have been previ- SSTR5-AS1 is the highest expressed lncRNA in the NEPC clin- ously published (protein-coding transcriptomes) by Akamatsu ical samples when requiring expression concordance in VPC et al. [24] and Mo et al. [162]. In brief, six LTL331, two LTL313, and WCM cohorts. It is an antisense transcript of SSTR5, which and two LTL418 PDXs were raised in NOD-SCID mice (NOD.CB17- is a member of the superfamily of somatostatin receptors. So- Prkdcscid/J) at the Living Tumor Laboratory [163]. Xenograft tis- matostatins are peptide hormones that regulate diverse cellu- sue was harvested after fixed lengths of time post host castra- lar functions such as neurotransmission, cell proliferation, and tion, tissue was measured, fixed for histopathological analysis, endocrine signalling, as well as inhibiting the release of many and processed for RNA analysis. hormones and other secretory proteins. The SSTR family (1–5) are markers for neuroendocrine tumours of the lung [155], with SSTR1 and SSTR5 the most dominant forms of SSTR in neu- Clinical datasets roendocrine tumours in general [54]. Interestingly, exploration within TANRIC showed a strong positive correlation in expres- We used five clinical cohorts from 1) WCM [ 5]; 2) GenomeDx Bio- sion with SSTR5 to SSTR5-AS1 in 14 of 20 cancer types (P < 2.18 sciences (GX) Inc. (MCI and MCII); 3) JHSM; and 4) VPC, cumu- −15 × 10 ; Supplementary Table S29). Furthermore, SSTR5 mRNA latively totalling 927 samples. For the VPC, 80 specimens were is detectable in the blood of neuroendocrine tumours of the lung obtained from patients undergoing RP and snap frozen follow- [156] and could be a valuable non-invasive diagnostic marker for ing a protocol approved by the Clinical Research Ethics Board of NEPC. In fact, clinicians utilize this biological feature in other the University of British Columbia, the BC Cancer Agency, and neuroendocrine tumours (NETs) using Octreoscans to determine Vancouver General Hospital pathology (depending on the sam- tumour stage and/or identification of sites of metastasis. Oc- ple source). All patients signed a formal consent form approved treoscans, when compared with positron emission tomography by the ethics board. A subset of the GX Decipher GRID database (PET) scans (another commonly used approach for this), appear of clinical specimens was selected, totalling 777 patient PCa more sensitive in the detection of well-differentiated NETs [157]. expression profiles (all from formalin-fixed parafin embedded In addition to this, therapeutically, somatostatin analogues are tissue) and were obtained from two RP Mayo Clinic cohorts emerging as a promising treatment option for inoperable or that have been previously described (MCI [116]and MCII [117]). metastatic NETs [158]. However, specifically in NEPC, targeting JHSM samples, totalling 33 samples, were retrieved from surgical SSTR5 and/or SSTR5-AS1 for diagnostic or therapeutic purposes pathology and consultation files of Johns Hopkins Hospital (John is in its infancy. Interestingly, SSTR5 (C terminal) is required for Hopkins Registry) from 1999 to 2013, as previously described Rb induction and G1 cell cycle arrest [159], resulting in anti- [164]. The 33 samples were annotated as 6 morphologically di- proliferative effects. However, without Rb (known to be lost in agnosed pure SCPC samples, 12 high-risk (Gleason 9–10) AD, 10 NEPC), this function is negated. Alternatively, the interaction ev- SCPC (SC-mixed), and 5 AD (AD-mixed) from mixed histology tu- idence for SSTR5-AS1 and KDM4B (Supplementary Fig. S13) pro- mours containing separate adenocarcinoma and small cell com- vides another strong connection to NEPC biology. KDM4B is a ponents. For this cohort, samples were dicotimized into either histone demethylase and a key molecule in AR signaling and AD (AD and AD-mixed samples) or NEPC (SCPC and SCPC-mixed turnover [160]. In NEPC with the absence of the AR, KDM4B could sampels) for the purposes of validating the 122 NEPC lncRNA pa- interact with N-Myc instead, where it has been shown to regu- tient signature. We also explored an externally processed cohort late and epigenetically activate this oncogene in neuroblastoma comprising 114 metastatic CRPC specimens, of which 44 were [161]. N-Myc has been seen to drive the progression of NEPC [5, NEPC [32] and used in this study. Referred to in the text as the 26, 27] and recently through EZH2-mediated transcription [27]. extNEPC cohort, we accessed and visualized this data through However, another mechanism of activation could be facilitated cBioPortal [165, 166] Version 1.9.0 [167]. OncoPrint schematics through SSTR5-AS1 regulation. However, both of these putative were generated for displaying multiple genomic alterations by Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 15 heatmap for the lncRNAs. The extNEPC study samples were clas- then at 70C for 15 minutes. Prior to use in qRT-PCR, products sified using a pathologic classification system [ 6] that included were diluted 10-fold with water. FastStart Essential Green Mas- five catagories: “A,” usual prostate adenocarcinoma without ter kit from Roche (Catalogue #06 402 712 001) was used as de- neuroendocrine differentiation; “B,” usual prostate adenocarci- scribed from their protocol for qRT-PCR reactions. In brief, 2 μl noma with neuroendocrine differentiation >20%; “C,” small-cell of water, 3 μl of a mixture of forward and reverse primers (each carcinoma; “D,” large-cell neuroendocrine carcinoma; and “E,”, at a concentration of 10 μM),and 10 μl of the Roche Master Mix mixed small-cell carcinoma–adenocarcinoma. were aliquoted into each well of a 96-well plate. A mixture of 4 μl of water plus 1 μl of the diluted cDNA was then added to the appropriate wells. Expression was then quantified (as mea- Material collection and processing (VPC Cohort) sured by Ct) through the Roche Light Cycler 96 machine. Each Hematoxylin and eosin (H&E) stained, formalin-fixed paraffin- lncRNA/sample pair was quantified with technical replicates in embedded, and fresh frozen sections were reviewed by a pathol- triplicate. Average and standard deviation of Ct were calculated ogist to identify blocks with highest tumour content. For each across these triplicates, and Ct calculated relative to house frozen block used, a 5-μmslide was first taken for H&E staining; keeper gene PSMB4 (most consistent and highly expressed gene then 4 × 100-μmsections were taken for DNA and RNA isolation vs REEP5 and SNRPD3). Delta Cts were calculated relative to before a second 5-μm slide was taken for H&E staining. Each control samples, and fold changes were plotted using Prisms H&E slide was required to have tumour content >50% for a tu- GraphPad software (Supplemental Table 12). mour to proceed for sequencing. RNA from 100-μmsections of snap frozen tissue was isolated using the mirVana Isolation Kit RNA sequence analysis pipeline from Ambion (AM 1560). RNA sequencing was performed on Il- lumina HiSeq 2000 at BC Cancer Agency Michael Smith Genome We implemented an lncRNA sequence analysis pipeline that in- Sciences Centre according to standard protocols. cludes algorithms catered to the detection of known and novel transcripts (Supplementary Figs S1 and S2). Implemented in- house, this pipeline is modified and extended from the tuxedo Material collection and processing (GRID and JHSM) suite of sequence analysis algorithms [45]. Once received from For GRID (MCI and MCII) and JHSM cohorts, specimen selection, the sequencing centre in bam format, all sequenced model sys- RNA extraction, and microarray hybridization were performed tems and patient samples were de-aligned into raw fastq for- in a Clinical Laboratory Improvement Amendments-certified mat (including flagged reads) using bam2fastq and put through laboratory facility (GenomeDx Biosciences, San Diego, CA, USA) the following pipeline. To ensure high-quality sequence reads, as described previously [116, 117]. Total RNA extraction, purifica- libraries were trimmed using a Sickle, a windowed-adaptive ap- tion, RNA amplification, and labelling were done using the Ova- proach [171]. For each read pair processed together, the algo- tion WTA FFPE system (NuGen, San Carlos, CA, USA). RNA was rithm determines the most optimal inner read sequence by hybridized to Human Exon 1.0 ST GeneChips (Affymetrix, Santa trimming both 3’ and 5’ prime ends based on quality and Clara, CA, USA). After microarray profiling, quality control was length thresholds (for full description, see [172]). Bases with a preformed using the Affymetrix Power Tools package, and probe quality score of <99.0% base call accuracy (corresponding to set normalization was performed using the Single Channel Ar- a Phred quality score of 20) were removed. Reads less than ray Normalization algorithm [168]. approximately two-thirds read length (30 nt in WCM and 60 nt in VPC) post-trimming were discarded. Highly repetitive se- quences (>2% of library) were also discarded post-trimming us- Quantitative Real-Time Polymerase Chain Reaction ing the cutadapt tool. All quality control metrics were gener- (qRT-PCR) ated and quanitified (pre- and post-trimming) using the FASTX- Primers were designed using Primer3 and checked with in sil- Toolkit and the FastQC Windows software. Reads were aligned ico PCR in UCSC Genome Browser (s Supplementary Table S26 to the Hg19 human genome build using an unspliced aligner for forward and reverse primer sequences). Housekeeping genes for handling exonic reads (Bowtie–v2.2.3) in conjunction with PSMB4, REEP5, and SNRPD3 were selected on the basis of high, a spliced aligner to handle reads spanning exon-exon junc- consistent expression levels across many cell and tissue types tions (Tophat 2.0.12). Transcriptome reconstruction using En- and were used in the MiTranscriptome lncRNA study [169, 170]. sembl GRCh37.75 gene tracks for each library was performed us- Two lncRNAs from each NEtD class and three NEPC lncRNAs ing a quasi de novo (genome-guided) approach (Cufflinks v2.2.1), from among the top candidates (Supplementary Fig. S12 and where reads were assembled and abundances estimated us- Supplementary Table S26) were selected (n = 11) for qRT-PCR ing an overlap graph producing a minimal spanning network validation. The cDNA from the PDX LTL331 models three time of transcripts. This version of Ensembl contained 38 transcript points (AD, postTX, and NEPC) were used to validate the NEtD classes grouped by 4 core biotypes. At this stage, transcripts lncRNAs and a subset of the VPC clinical samples for the NEPC were also multi-read and fragment bias corrected. Transcripts lncRNAs. With the rarity of clinical NEPC samples, tumour tissue with highly abundant expression were masked (e.g., rRNAs) and subsequent RNA were extremely limited. Due to this, only from downstream steps to increase transcript quantification ac- three NEPC (V73, V90, and V91) and one AD (V60) clinical sample curacy. Sample transcriptomes, the reference genome, and the were included in this validation. For each lncRNA and sample transcript annotation were then meta-assembled (Cuffmerge) to tested, the following experimental protocol was carried out: 1 produce a single annotation transcriptome model. Based on this μg of total RNA for each sample was diluted to 18 μl with water model, transcript quantification (Cuffquant) and normalization and 1 μl of random hexamers (50μM; Thermo Fisher). The mix- (Cuffnorm). Geometric and FPKM normalization performed in- ture was heated to 65C for 5 minutes and chilled. Afterwards, 5 dependently (Cuffnorm) corrected for uneven library sequencing μl of 5x reverse transcriptase buffer, 1 μlof10mMdNTPs,and depths between samples and variable transcript lengths within 1 μl of Superscript II reverse transcriptase (Thermo Fisher) were samples. Transcript expression displaying computational arti- added. Each sample was then incubated at 42C for 1 hour and facts (expression values <0.1 known to occur with Cufflinks) Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 16 NEPC lncRNAs were converted to zero values. This generated transcript expres- (GREF [includes the androgen receptor and the closely related sion where only lncRNAs (Ensembl and ENCODE-based) were glucocorticoid, mineralocorticoid, and progesterone receptors], extracted and used for all downstream analysis. All algorithms NRSF [REST], SOX, HOX, STEM, E2FF, ETSF, and ETVI1) motifs denoted in brackets are referenced and described in the Trap- included in this study are described in Supplementary Table nell et al. Nature protocol [45]. Each cohort (VPC and WCM) was S25. Core and matrix similarities were again set to 1 and op- processed independently by this pipeline, and then transcrip- timized, respectively. (3) Multiple lncRNA enrichment analysis tome annotations were merged. This was accomplished using was performed using the Overrepresented TFBS algorithm. En- Ensembl transcript IDs combined with transcript lengths to pro- richment of matrix/matrix family was determined by Genomatix duce unique transcript identifiers for each lncRNA across co- calculated z-scores (>2or <-2), which is based on the distance horts. from the population mean (genome or promoter sequence back- ground) in units of the population standard deviation for query sequence/promoter. Genomatix calculates z-scores with a conti- RNA-RNA interaction analysis nuity correction using the formula z = (x-E-0.5)/S, where x is the A genome-wide analysis for SSTR5-AS1 and LINC00514 lncRNA number of found matches in the input data, E is the expected interactions was performed using a multistep systemic ap- value, and S is the standard deviation. This formula is also de- proach [119]. This tool is available publicly within an online scribed in the oPOSSUM algorithm [180]. A z-score < -2 or >2 database [173] hosted by the Computational Biological Research can be considered statistically significant and corresponds to a Center at the National Institute of Advanced Industrial Science P-value of approximately 0.05. and Technology in Japan. The interaction search space included all hg19 annotated lncRNA and mRNA transcripts. This gener- Microarray to sequencing platform lift over/mapping ated top-ranking interaction partners (n = 100), based on local interaction minimum free energy (Supplementary Tables S27– Affymetrix Human Exon 1.0 ST GeneChip probes were mapped S28). R-chie [174, 175] was used to visualize the top-ranking pre- to hg19 coordinates using SMALT v0.76 [181]. Probe set genomic dictions KDM4B and TADA3 for SSTR5-AS1 and LINC00514, re- regions (PSRs) were redefined accordingly. Exons within each spectively, using the double structure feature (Supplementary lncRNA from sequencing cohorts (VPC and WCM) were inte- Figs S13 and S14). All bases that were not within the interac- grated with PSRs to build an overlap table to determine ab- tion site were predicted to form RNA secondary structure by sence/presence of lncRNA transcripts on the affymetrix microar- RNAfold [176–178] selecting to enforce constrained pairing pat- ray. R function iRanges v2.9.18 was used with method findOver- terns for the interacting bases. Minimum free-energy structures lap to build the described table above. Microarray PSRs were re- were predicted by RNAfold on the 300nt sequences upstream quired to be entirely within sequenced exon genome regions; and downstream of the interaction site. otherwise they were excluded. Applying this methodology, 106 of 122 NEPC lncRNAs (87%) and 81 of 100 NEtD lncRNAs (81%) mapped to microarray PSRs for clinicopathological analysis on TFBS identification and enrichment analysis GRID cohorts MCI and MCII. All TFBS analysis was performed using Genomatix software, databases, and algorithms [179]. Three types of TF analysis were Statistical analysis carried out in this study: (1) single lncRNA motif characteriza- tion, (2) multiple lncRNA analysis for select TFs, and (3) multiple For all cohorts, the programming language R v3.0 was used lncRNA enrichment analysis. Prior to any of the above, lncRNA for statistical analysis. For VPC and WCM cohorts, unsuper- transcript(s) were submitted to the Gene2Promoter algorithm for vised hierarchical clustering was performed with the h.clust retrieval of promoter sequences. Databases used with this algo- package with Pearson correlation for distance and average link- rithm included ElDorado 12–2013 and NCBI build 37 (for multiple age used. Only transcripts within the top fifth percentile based lncRNA analysis where genomic background needed to match on their standard deviations were selected. The clustering and sequencing data) or the most recent databases ElDorado 12–2016 heatmaps generated were built using the heatmap.2 function. and GRCh38 (for single lncRNA analysis where genomic back- Similar clustering analysis was performed for GRID cohorts ground was not relevant). Transcripts with alternative isoforms except with Euclidian distance, the ward method for linkage, were required to have gold level (experimentally verified 5’ com- and the use of the heatmap.3 function due to its advanced plete transcript), silver level (transcript with 5’ end confirmed by row/column labelling features. For all cohorts before clustering, PromoterInspector prediction), or bronze level (annotated tran- normalized log2 expression values were standardized/scaled us- script, no confirmation for 5’ completeness) quality for their al- ing a z-score that ranged from -2 to 2. For principal component ternative isoforms. (1) Single lncRNA motif characterization was analysis, the R package prcomp was used to calculate variance performed using the MatInspector algorithm [97–99]withpa- among transcript and sample subsets for the calculation of tran- rameters “core similarity” (degree of similarity for highest con- script weights and principle components. The top three compo- served bases of motif) set to 1 and “matrix similarity” (degree of nents were used for visual inspection. For all clinical group-wise similarity between motif and query sequence) set to optimized comparisons, a standard Student’s t-test was applied to iden- as recommended by Genomatix and as described in MatInspec- tify significantly differentially expressed transcripts between tor referenced papers above. MatInspector uses the best in field groups/phenotypes. Significance thresholds were implemented MatBase database for TFBS motif/matrix annotation, where Ma- by enforcing a strict P-value cut-off of <0.05. Multiple test correc- trix Family Library Version 10.0 was used. (2) Multiple lncRNA tion was applied to P-values using the Bonferroni and Hochberg analysis for select TFs was performed using MatInspector and method to mathematically minimize false discovery rate (FDR) select TF motifs (“matrix”) applied accordingly. All matrix an- and with a cut-off of P-value < 0.05. See Supplementary Ta- notation, descriptions, and matrix family definitions are listed bleS4for theseresults.To biologically minimize FDR, mathemati- in Supplementary Table S25. Select TF matrices (BRN2, STAT3, cal FDR correction was removed and instead followed the filter- NKX3, NMYC, SOX2, and SOX11) and select TF matrix families down workflows in Figs 3A–C and Fig. 1C. Despite mathematical Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 17 FDR being removed, statistical significance of P-value < 0.05 was Additional file still maintained during the filter-down approach using the de- Table S1: Discovery cohorts clinical and sequencing informa- scribed method in each step. See Supplementary Table S5 for tion. these results. For ROC curves and AUC calculations, the R pack- Table S2: Differentially expressed AR- and AR+ CRPC lncRNA. age “pROC” was used. Kaplan-Meier analysis was performed for Table S3: NEPC lncRNA expression signature labels and order for determining survival outcome using the R package “survfit” with Fig. 2. transcripts displaying below background (<0.1) expression being Table S4: Clinical cohort group-wise comparisons. removed from this analysis. Table S5: NEtD and NEPC lncRNA transcript annotation Table S6: TFBS-Enrichment by Matrix for NEtD Class I Table S7: TFBS-Enrichment by Matrix for NEtD Class II Transcript ranking Table S8: TFBS-Enrichment by Matrix for NEtD Class III Table S9: TFBS-Enrichment by Matrix for NEtD Class IV NEtD lncRNA transcripts were ranked based on fold changes ob- Table S10: TFBS-Enrichment by Matrix for NEtD Class III-IV served in the clinical group-wise comparisons. For NEtD Class Table S11: TFBS-Enrichment by Matrix for all NEtD lncRNAs I Deactivated, the three group-wise comparison fold changes Table S12: TFBS-Enrichment by Matrix Family for NEtD Class I were used (NEPC vs AD, CRPC vs AD, and NHT vs AD), where Table S13: TFBS-Enrichment by Matrix Family for NEtD Class II the minimum fold change observed between the three com- Table S14: TFBS-Enrichment by Matrix Family for NEtD Class III parisons was selected and then ranked in decreasing order. For Table S15: TFBS-Enrichment by Matrix Family for NEtD Class IV Class II Activated and Class III Persistent transcripts, NEPC vs Table S16: TFBS-Enrichment by Matrix Family for NEtD Class III- AD fold changes were calculated and ranked in increasing or- IV der for both VPC and WCM cohorts, where the maximum fold Table S17: TFBS-Enrichment by Matrix Family for all NEtD lncR- change between VPC and WCM was selected. For Class IV Tran- NAs sient transcripts, absolute fold changes for AD vs NHT and NHT Table S18: TFBS Overlap Table by Matrix vs NEPC were calculated and ranked in increasing order with the Table S19: TFBS Overlap Table by Matrix Family maximum fold change from either group selected. Similar rank- Table S20: Unique and common TFBS in NEtD lncRNA ing was performed for NEPC lncRNA transcripts, which were or- Table S21: Select TF Identification for NEtD Class I dered by increasing (up-regulated transcripts in NEPC vs. AD) Table S22: Select TF Identification for NEtD Class II-III or decreasing (down-regulated transcripts in NEPC vs. AD) or- Table S23: Select TF Identification for NEtD Class III-IV der to determine the highest- and lowest-expressed transcripts Table S24: Select TF Identification for all NEtD lncRNAs in NEPC vs AD, respectively. Concordantly expressed transcripts were required between VPC and WCM cohorts. The top 20 lncR- Table S25: Genomatix Matrix and Matrix Family definitions Table S26: Top-ranking NEPC and NEtD lncRNAs. NAs (based on fold changes from clinical samples defined above) Table S27: SSTR5-AS1-predicted RNA (mRNA or lncRNA) interac- were taken from each group. This produced 20 × 5 group (n = tions with associated binding energies, predicted transcript En- 100) isoforms representing 76 unique lncRNA transcripts. These sembl ID, name, interaction position, and ranking. represent the top NEtD/NEPC lncRNA candidates from this study Table S28: LINC00514-predicted RNA (mRNA or lncRNA) interac- (Supplementary Table S26). No pseudogenes were included in tions with associated binding energies, predicted transcript En- these rankings. sembl ID, name, interaction position, and ranking. Table S29: TANRIC results for lncRNAs FENDRR and SSTR5-AS1. Spearman rank correlation for protein coding genes that are Availability of supporting data the sense forms to the above antisense transcripts. Numbers A subset of the sequenced samples (n = 70) used in this study in brackets denote P-values. NSC, no significant correlation; NA, was from previous studies with all raw sequencing data rean- mRNA data was not available for this tumour type, therefor the alyzed here using the described pipeline above. These 70 sam- analysis was not applicable. ples have been previously submitted to the European Nucleotide Figure S1: The next-generation sequence analysis pipeline im- Archive (ENA) or NCBIs Gene Expression Omnibus. This in- plemented for the detection and quantification of lncRNAs in cludes the 6 NEPC PDX model samples [24] (ENA accession num- this study. A) The nine-step lncRNA next-generation sequenc- ber PRJEB9660 and GEO accession number GSE59986), 2 CRPC ing analysis pipeline with core algorithms (Bowtie2, Tophat2, PDXmodelsamples[162] (ENA accession number PRJEB19256), Cufflinks2, Cuffmerge, Cuffquant, and Cuffnorm) implemented 4 NEPC (VPC) samples [31, 56], 23 AD (VPC) samples [56](ENA from the Tuxedo suite of analysis tools. Sequencing quality con- accession number PRJEB6530), 30 AD (WCM) samples [5], and trol metrics before and after trimming of data for sample V60 7 NEPC (WCM) samples [5]. The remaining unpublished se- is outlined in B–F. This includes pre-trimming (A) phred quality quenced samples (n = 55) have been submitted to the ENA under scores and (B) percentage of each base type across read library accession number PRJEB21092. Please see Supplementary Table at each base pair position. After quality control corrections were S1 for a summary of sequencing and clinical information on applied, V60 read library had acceptable (D) Phred quality scores these 125 samples. All microarray samples from GX cohorts, in- (∼30 Phred Score) and (E) expected base type percentages (∼25%) cluding 545 AD (MCI [116]) samples and 232 AD (MCII [117]) sam- for T, C, A, and G. F) All over-represented sequences that were ples, are accessible through Gene Expression Omnibus acces- >2% of library were removed from the V60 read library. See Meth- sion numbers GSE46691 and GSE62116, respectively. Additional ods for complete listing and version numbers for all algorithms supporting data and custom code from the sequencing pipeline and tools used in the sequence analysis pipeline. described above are also available from the GigaScience GigaDB Figure S2: Average Phred quality scores for all VPC and WCM database [182]. samples before and after quality control corrections were ap- plied. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 18 NEPC lncRNAs Figure S3: Unsupervised hierarchal clustering (A–D) and princi- Funding ple component analysis (E–H) on the four major Ensembl tran- This work was supported by the Mitacs Accelerate PhD Fellow- script classes detected within the VPC cohort. Samples are la- ship Program (IT04310 to V.R.R.) in collaboration with GenomeDx belled as adenocarcinomas (blue) and neuroendocrine tumours Biosciences, Terry Fox Foundation (201012TFF to C.C.C.), and (gold). Prostate Cancer Canada Team Grant (T2013–01 to C.C.C.). Figure S4: Unsupervised hierarchal clustering (A–D) and princi- ple component analysis (E–H) on the four major Ensembl tran- script classes detected within the WCM cohort. Samples are la- Competing interests belled as adenocarcinomas (blue) and neuroendocrine tumours The authors declare that they have no competing interests. (gold). Figure S5: Detected and differentially expressed lncRNAs among −/+ (A) AR CRPC xenograft models (B) and matched clinical sam- Aknowledgements ples. P, patient; X, xenograft. We are grateful to the following GenomeDx bioinformaticians: Figure S6: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC Mandeep Takhar for her help with GRID statistical analy- Cohort – LEGEND. Clinical, transcriptome, and genomic annota- sis/primer code; Hussam Al-Deen Ashab for his help with IA tion for samples plotted in Supplementary Figs S7–S11. All an- analysis (ultimately was not included in the paper, but effort notations were generated from cBioportal with the exception and knowledge gained from the results were insightful for this of the NEPC molecular subtype, which was assigned from this work); Nicholas Erho for his effort in mapping our sequenced study. Transcripts denoted with 1 in superscript within Supple- data to the GRID microarray, mentoring, and constant support; mentary Figs S8–S11 are lncRNAs that overlap a NEPC lncRNA and Mohammed Alshalalfa for his guidance and supervision on with a NEtD lncRNA. For example, H19 appears in Supplemen- all NEPC/NEtD lncRNA clinicopathological GenomeDx analysis. tary Fig. S8 (NEPC lncRNA) and supplementary 9 (NEtD lncRNA We would also like to deeply thank Daniel Lai and Alex Gawron- Class II). ski for their advice with RNA-RNA visualization and interaction Figure S7: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC analysis algorithms. We would like to thank Faraz Hach for his Cohort - Known NEPC genes and TFs. Select NEPC onco- manuscript insights and advice. Lastly, we are extremely grate- genes, tumour suppressor, and transcription factors that have ful to Stephanie Giles Ramnarine for her manuscript comments, been reported previously or within this study for transcrip- advice, and support. tomic/genomic context with Supplementary Figs S8–S11. Please see Supplementary Fig. S6 for figure legend. Figure S8: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC References Cohort - NEPC lncRNA. Transcripts from the NEPC lncRNA ex- pression signature that are up-regulated (74 of 122), testable (58 1. Torre LA, Bray F, Siegel RL et al. Global cancer statistics, of 74), and altered (25 of 58) in the extNEPC cohort. Please see 2012. CA Cancer J Clin 2015;65:87–108. Supplementary Fig. S6 for figure legend. 2. Karantanos T, Evans CP, Tombal B et al. Understanding the Figure S9: NEPC and NEtD lncRNA Oncoprint Plot in the extNEPC mechanisms of androgen deprivation resistance in prostate Cohort - NEtD lncRNA Class II. Transcripts from NEtD lncRNA cancer at the molecular level. Eur Urol 2015;67:470–9. Class II (222), testable (128 of 222), and altered (26 of 128) in the 3. Grasso CS, Wu YM, Robinson DR et al. The mutational land- extNEPC cohort. Please see Supplementary Fig. S6 for figure leg- scape of lethal castration-resistant prostate cancer. Nature end. 2012;487:239–43. Figure S10: NEPC and NEtD lncRNA Oncoprint Plot in the 4. Vlachostergios PJ, Puca L, Beltran H. Emerging variants extNEPC Cohort - NEtD lncRNA Class III. Transcripts from NEtD of castration-resistant prostate cancer. Curr Oncol Rep lncRNA Class III (84), testable (79 of 84), and altered (29 of 84) in 2017;19:32. the extNEPC cohort. Please see Supplementary Fig. S6 for figure 5. Beltran H, Rickman DS, Park K et al. Molecular characteri- legend. zation of neuroendocrine prostate cancer and identification Figure S11: NEPC and NEtD lncRNA Oncoprint Plot in the of new drug targets. Cancer Discov 2011;1:487–95. extNEPC Cohort - NEtD lncRNA Class IV. Transcripts from NEtD 6. Epstein JI, Amin MB, Beltran H et al. Proposed morphologic lncRNA Class IV (45), testable (36 of 45), and altered (11 of 31) in classification of prostate cancer with neuroendocrine dif- the extNEPC cohort. Please see Supplementary Fig. S6 for figure ferentiation. Am J Surg Pathol 2014;38:756–67. legend. 7. Lin D, Wyatt AW, Xue H et al. High fidelity patient-derived Figure S12: Quantitative Real-Time Polymerase Chain Reac- xenografts for accelerating prostate cancer discovery and tionon select NEPC and NEtD lncRNAs. drug development. Cancer Res 2014;74:1272–83. Figure S13: Hypothetical RNA-RNA folding structure for exon 4 8. Aparicio AM, Harzstark AL, Corn PG et al. Platinum-based of SSTR5-AS1 (top) and the 3’UTR of KDM4B (bottom). Predicted chemotherapy for variant castrate-resistant prostate can- base pair binding (green arcs) along the sequence (black arrow) cer. Clin Cancer Res 2013;19:3621–30. are displayed, included predicted interaction site (orange bars). 9. Wang HT, Yao YH, Li BG et al. Neuroendocrine Prostate Can- Figure S14: Hypothetical RNA-RNA folding structure for exon 4 cer (NEPC) progressing from conventional prostatic adeno- of LINC00514 (top) and the 3’UTR of TADA3 (bottom). Predicted carcinoma: factors associated with time to development of base pair binding (green arcs) along the sequence (black arrow) NEPC and survival from NEPC diagnosis-a systematic re- are displayed, included predicted interaction site (orange bars). view and pooled analysis. J Clin Oncol 2014;32:3383–90. 10. Terry S, Beltran H. The many faces of neuroendocrine dif- ferentiation in prostate cancer progression. Front Onco 2014;4:60. 11. Lee JK, , Bangayan NJ, Chai T et al. Systemic surfaceome pro- filing identifies target antigens for immune-based therapy Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 19 in subtypes of advanced prostate cancer. PNAS 2018; 115, Biochim Biophys Acta 2017;1860:1002–12. E4473–E4482. 31. Lapuk AV, Wu C, Wyatt AW et al. From sequence 12. Shen R, Dorai T, Szaboles M et al. Transdifferentiation of to molecular pathology, and a mechanism driving the cultured human prostate cancer cells to a neuroendocrine neuroendocrine phenotype in prostate cancer. J Pathol cell phenotype in a hormone-depleted medium. Urol Oncol 2012;227:286–97. 1997;3:67–75. 32. Beltran H, Prandi D, Mosquera JM et al. Divergent clonal evo- 13. Palmgren JS, Karavadia SS, Wakefield MR. Unusual and un- lution of castration-resistant neuroendocrine prostate can- derappreciated: small cell carcinoma of the prostate. Semin cer. Nat Med 2016;22:298–305. Oncol 2007;34:22–29. 33. Viswanathan VS, Ryan MJ, Dhruv HD et al. Dependency of a 14. Tan HL, Sood A, Rahimi HA et al. Rb loss is characteristic of therapy-resistant state of cancer cells on a lipid peroxidase prostatic small cell neuroendocrine carcinoma. Clin Cancer pathway. Nature 2017;547:453–7. Res 2014;20:890–903. 34. Gibb EA, Brown CJ, Lam WL. The functional role of long non- 15. Chen H, Sun Y, Wu C et al. Pathogenesis of prostatic small coding RNA in human carcinomas. Mol Cancer 2011;10:38. cell carcinoma involves the inactivation of the P53 pathway. 35. Cheetham SW, Gruhl F, Mattick JS et al. Long noncoding Endocr Relat Cancer 2012;19:321–31. RNAs and the genetics of cancer. Br J Cancer 2013;108:2419– 16. Ku SY, Rosario S, Wang Y et al. Rb1 and Trp53 cooperate to 25. suppress prostate cancer lineage plasticity, metastasis, and 36. Marchese FP, Raimondi I, Huarte M. The multidimensional antiandrogen resistance. Science 2017;355:78–83. mechanisms of long noncoding RNA function. Genome Biol 17. Ham WS, Cho NH, Kim WT et al. Pathological effects of 2017;18:206. prostate cancer correlate with neuroendocrine differentia- 37. Sun W, Yang Y, Xu C et al. Regulatory mechanisms of long tion and PTEN expression after bicalutamide monotherapy. noncoding RNAs on gene expression in cancers. Cancer J Urol 2009;182:1378–84. Genet 2017;216–217:105–10. 18. Zou M, Toivanen R, Mitrofanova A et al. Transdifferentia- 38. Gutschner T, Diederichs S. The hallmarks of cancer: a long tion as a mechanism of treatment resistance in a mouse non-coding RNA point of view. RNA Biol 2012;9:703–19. model of castration-resistant prostate cancer. Cancer Dis- 39. Kondo Y, Shinjo K, Katsushima K. Long non-coding RNAs cov 2017;7:736–49. as an epigenetic regulator in human cancers. Cancer Sci 19. Li Y, Donmez N, Sahinalp C et al. SRRM4 drives neuroen- 2017;108:1927–33. docrine transdifferentiation of prostate adenocarcinoma 40. Sahu A, Singhal U, Chinnaiyan AM. Long noncoding RNAs under androgen receptor pathway inhibition. Eur Urol, in cancer: from function to translation. Trends Cancer 2016. 2015;1:93–109. 20. Li Y, Chen R, Bowden M et al. Establishment of a neuroen- 41. Bhan A, Soleimani M, Mandal SS. Long noncoding RNA and docrine prostate cancer model driven by the RNA splicing cancer: a new paradigm. Cancer Res 2017;77:3965–81. factor SRRM4. Oncotarget 2017;8:66878–88. 42. Cheng W, Zhang Z, Wang J. Long noncoding RNAs: new 21. Zhang X, Coleman IM, Brown LG et al. SRRM4 expression players in prostate cancer. Cancer Lett 2013;339:8–14. and the loss of REST activity may promote the emergence 43. Lin D, Dong X, Wang K et al. Identification of DEK as a poten- of the neuroendocrine phenotype in castration-resistant tial therapeutic target for neuroendocrine prostate cancer. prostate cancer. Clin Cancer Res 2015;21:4698–708. Oncotarget 2015;6:1806–20. 22. Bishop JL, Thaper D, Vahid S et al. The master neural tran- 44. Clermont PL, Lin D, Crea F et al. Polycomb-mediated si- scription factor BRN2 is an androgen receptor suppressed lencing in neuroendocrine prostate cancer. Clin Epigenetics driver of neuroendocrine differentiation in prostate cancer. 2015;7:40. Cancer Discov 2016, 7, 54–71;7:54-71. 45. Trapnell C, Roberts A, Goff L et al. Differential gene and 23. Kim J, Jin H, Zhao JC et al. FOXA1 inhibits prostate cancer transcript expression analysis of RNA-seq experiments neuroendocrine differentiation. Oncogene 2017;36:4072–80. with TopHat and Cufflinks. Nat Protoc 2012; 7:562–78. 24. Akamatsu S, Wyatt AW, Lin D et al. The placental gene 46. Gutschner T, Hammerle M, Diederichs S. MALAT1 – a PEG10 promotes progression of neuroendocrine prostate paradigm for long noncoding RNA function in cancer. J Mol cancer. Cell Rep 2015;12:922–36. Med (Berl) 2013;91:791–801. 25. Ci X, Hao J, Dong X et al. Heterochromatin protein 1al- 47. Cross DS, Burmester JK. Functional characterization of the pha mediates development and aggressiveness of neuroen- GDEP promoter and three enhancer elements in retinoblas- docrine prostate cancer. Cancer Res 2018. toma and prostate cell lines. Med Oncol 2008;25:40–49. 26. Lee JK, Phillips JW, Smith BA et al. N-Myc drives neuroen- 48. Reding DJ, Zhang KQ, Salzman SA et al. Identification of docrine prostate cancer initiated from human prostate ep- a gene frequently mutated in prostate tumors. Med Oncol ithelial cells. Cancer Cell 2016;29:536–47. 2001;18:179–87. 27. Dardenne E, Beltran H, Benelli M et al. N-Myc induces an 49. Niknafs YS, Han S, Ma T et al. The lncRNA landscape of EZH2-Mediated transcriptional program driving neuroen- breast cancer reveals a role for DSCAM-AS1 in breast cancer docrine prostate cancer. Cancer Cell 2016;30:563–77. progression. Nat Commun 2016;7:12791. 28. Mu P, Zhang Z, Benelli M et al. SOX2 promotes lineage 50. Wang O, Yang F, Liu Y et al. C-MYC-induced upregulation of plasticity and antiandrogen resistance in TP53- and RB1- lncRNA SNHG12 regulates cell proliferation, apoptosis and deficient prostate cancer. Science 2017; 355:84–88. migration in triple-negative breast cancer. Am J Transl Res 29. Maina PK, Shao P, Liu Q et al. c-MYC drives histone 2017;9:533–45. demethylase PHF8 during neuroendocrine differentiation 51. Chen T, Yang P, He ZY. Long non-coding RNA H19 can pre- and in castration-resistant prostate cancer. Oncotarget dict a poor prognosis and lymph node metastasis: a meta- 2016;7:75585–602. analysis in human cancer. Minerva Med 2016;107:251–8. 30. Maina PK, Shao P, Jia X et al. Histone demethylase PHF8 reg- 52. Raveh E, Matouk IJ, Gilon M et al. The H19 long non-coding ulates hypoxia signaling through HIF1alpha and H3K4me3. RNA in cancer initiation, progression and metastasis - a Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 20 NEPC lncRNAs proposed unifying theory. Mol Cancer 2015;14:184. dogenous miR-342-3p in gallbladder cancer. J Exp Clin Can- 53. Li H, Zhu L, Xu L et al. Long noncoding RNA linc00617 cer Res 2016;35:160. exhibits oncogenic activity in breast cancer. Mol Carcinog 71. Xu TP, Huang MD, Xia R et al. Decreased expression of the 2017;56:3–17. long non-coding RNA FENDRR is associated with poor prog- 54. Pisarek H, Pawlikowski M, Kunert-Radek J et al. SSTR1 and nosis in gastric cancer and FENDRR regulates gastric can- SSTR5 subtypes are the dominant forms of somatostatin re- cer cell metastasis by affecting fibronectin1 expression. J ceptor in neuroendocrine tumors. Folia Histochem Cytobiol Hematol Oncol 2014;7:63. 2010;48:142–7. 72. Fernando TR, Contreras JR, Zampini M et al. The lncRNA 55. Childs A, Vesely C, Ensell L et al. Expression of somatostatin CASC15 regulates SOX4 expression in RUNX1-rearranged receptors 2 and 5 in circulating tumour cells from patients acute leukemia. Mol Cancer 2017;16:126. with neuroendocrine tumours. Br J Cancer 2016;115:1540–7. 73. Diskin SJ, Capasso M, Schnepp RW et al. Common variation 56. Wyatt AW, Mo F, Wang K et al. Heterogeneity in the inter- at 6q16 within HACE1 and LIN28B influences susceptibility tumor transcriptome of high risk prostate cancer. Genome to neuroblastoma. Nat Genet 2012;44:1126–30. Biol 2014;15:426. 74. Russell MR, Penikis A, Oldridge DA et al. CASC15-S is a tu- 57. Ahlgren G, Pedersen K, Lundberg S et al. Regressive changes mor suppressor lncRNA at the 6p22 neuroblastoma suscep- and neuroendocrine differentiation in prostate cancer after tibility locus. Cancer Res 2015;75:3155–66. neoadjuvant hormonal treatment. Prostate 2000;42:274–9. 75. Nakakura EK, Watkins DN, Schuebel KE et al. Mammalian 58. Wolf DA, Herzinger T, Hermeking H et al. Transcriptional Scratch: a neural-specific snail family transcriptional re- and posttranscriptional regulation of human androgen re- pressor. PNAS 2001;98:4010–5. ceptor expression by androgen. Mol Endocrinol 1993;7:924– 76. Ishii J, Sato H, Yazawa T et al. Class III/IV POU transcrip- 36. tion factors expressed in small cell lung cancer cells are in- 59. Cai C, He HH, Chen S et al. Androgen receptor gene ex- volved in proneural/neuroendocrine differentiation. Pathol pression in prostate cancer is directly suppressed by the Int 2014;64:415–22. androgen receptor through recruitment of lysine-specific 77. Beltran H, Wyatt AW, Chedgy EC et al. Impact of therapy on demethylase. Cancer Cell 2011;20:457–71. genomics and transcriptomics in high-risk prostate cancer 60. Knuuttila M, Yatkin E, Kallio J et al. Castration in- treated with neoadjuvant docetaxel and androgen depriva- duces up-regulation of intratumoral androgen biosynthe- tion therapy. Clin Cancer Res 2017;, 23, 6802–6811. sis and androgen receptor expression in an orthotopic 78. Manohar CF, Furtado MR, Salwen HR et al. Hox gene ex- VCaP human prostate cancer xenograft model. Am J Pathol pression in differentiating human neuroblastoma cells. 2014;184:2163–73. Biochem Mol Biol Int 1993;30:733–41. 61. Wang Q, Zhang J, Liu Y et al. A novel cell cycle-associated 79. Manohar CF, Salwen HR, Furtado MR et al. Up-regulation of lncRNA, HOXA11-AS, is transcribed from the 5-prime end HOXC6, HOXD1, and HOXD8 homeobox gene expression in of the HOXA transcript and is a biomarker of progression in human neuroblastoma cells following chemical induction glioma. Cancer Lett 2016;373:251–9. of differentiation. Tumour Biol 1996;17:34–47. 62. Sun M, Nie F, Wang Y et al. LncRNA HOXA11-AS promotes 80. Hessenkemper W, Baniahmad A. Targeting heat shock pro- proliferation and invasion of gastric cancer by scaffold- teins in prostate cancer. Curr Med Chem 2013;20:2731–40. ing the chromatin modification factors PRC2, LSD1, and 81. Azad AA, Zoubeidi A, Gleave ME et al. Targeting heat shock DNMT1. Cancer Res 2016;76:6299–310. proteins in metastatic castration-resistant prostate cancer. 63. Misawa A, Takayama K, Urano T et al. Androgen-induced Nature Rev Urol 2015;12:26–36. long noncoding RNA (lncRNA) SOCS2-AS1 promotes cell 82. Wang B, Lee CW, Witt A et al. Heat shock factor 1 induces growth and inhibits apoptosis in prostate cancer cells. J Biol cancer stem cell phenotype in breast cancer cell lines. Chem 2016;291:17861–80. Breast Cancer Res Treat 2015;153:57–66. 64. Zhao W, Luo J, Jiao S. Comprehensive characterization of 83. Chen WM, Huang MD, Sun DP et al. Long intergenic non- cancer subtype associated long non-coding RNAs and their coding RNA 00152 promotes tumor cell cycle progression clinical implications. Sci Rep 2014;4:6591. by binding to EZH2 and repressing p15 and p21 in gastric 65. Li Z, Yu X, Shen J. ANRIL: a pivotal tumor suppressor cancer. Oncotarget 2016;7:9773–87. long non-coding RNA in human cancers. Tumour Biol 84. Chen QN, Chen X, Chen Z-Y et al. Long intergenic non- 2016;37:5657–61. coding RNA 00152 promotes lung adenocarcinoma prolif- 66. Aguilo F, Zhou MM, Walsh MJ. Long noncoding RNA, poly- eration via interacting with EZH2 and repressing IL24 ex- comb, and the ghosts haunting INK4b-ARF-INK4a expres- pression. Mol Cancer 2017;16:17. sion. Cancer Res 2011;71:5365–9. 85. Jin HJ, Zhao JC, Wu L et al. Cooperativity and equilibrium 67. Kotake Y, Nakagawa T, Kitagawa K et al. Long non-coding with FOXA1 define the androgen receptor transcriptional RNA ANRIL is required for the PRC2 recruitment to and si- program. Nat Commun 2014;5:3972. lencing of p15(INK4B) tumor suppressor gene. Oncogene 86. Jin HJ, Zhao JC, Ogden I et al. Androgen receptor- 2011;30:1956–62. independent function of FoxA1 in prostate cancer metas- 68. Yap KL, Li S, Munoz-Ca ˜ bello AM et al. Molecular interplay tasis. Cancer Res 2013;73:3725–36. of the noncoding RNA ANRIL and methylated histone H3 87. Yang J, Mani SA, Donaher JL et al. Twist, a master regulator lysine 27 by polycomb CBX7 in transcriptional silencing of of morphogenesis, plays an essential role in tumor metas- INK4a. Mol Cell 2010;38:662–74. tasis. Cell 2004;117:927–39. 69. Yu W, Gius D, Onyango P et al. Epigenetic silencing of tu- 88. Wang J, Nikhil K, Viccaro K et al. The Aurora-A-Twist1 axis mour suppressor gene p15 by its antisense RNA. Nature promotes highly aggressive phenotypes in pancreatic car- 2008;451:202–6. cinoma. J Cell Sci 2017;130:1078–93. 70. Wang SH, Ma F, Tang ZH et al. Long non-coding RNA H19 89. Galvan JA, Astudillo A, Vallina A et al. Epithelial- regulates FOXM1 expression by competitively binding en- mesenchymal transition markers in the differential Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 21 diagnosis of gastroenteropancreatic neuroendocrine prostate cancer. EMBO J 2013;32:1665–80. tumors. Am J Clin Pathol 2013;140:61–72. 108. Seifert A, Werheid DF, Knapp SM et al. Role of Hox genes in 90. Fendrich V, Maschuw K, Waldmann J et al. Epithelial- stem cell differentiation. World J Stem Cells 2015;7:583–95. mesenchymal transition is a critical step in tumorgene- 109. Cho KH, Jeong KJ, Shin SC et al. STAT3 mediates TGF-beta1- sis of pancreatic neuroendocrine tumors. Cancers (Basel) induced TWIST1 expression and prostate cancer invasion. 2012;4:281–94. Cancer Lett 2013;336:167–73. 91. Eide T, Ramberg H, Glackin C et al. TWIST1, a novel 110. Yang X, Duan B, Zhou X. Long non-coding RNA FOXD2-AS1 androgen-regulated gene, is a target for NKX3-1 in prostate functions as a tumor promoter in colorectal cancer by regu- cancer cells. Cancer Cell Int 2013;13:4. lating EMT and Notch signaling pathway. Eur Rev Med Phar- 92. Ainechi S, Mann SA, Lin J et al. Paired Box 5 (PAX5) expres- macol Sci 2017;21:3586–91. sion in poorly differentiated neuroendocrine carcinoma of 111. Zhou W, Ye XL, Xu J et al. The lncRNA H19 mediates breast the gastrointestinal and pancreatobiliary tract: diagnostic cancer cell plasticity during EMT and MET plasticity by and potentially therapeutic implications. Appl Immunohis- differentially sponging miR-200b/c and let-7b. Sci Signal tochem Mol Morphol, 2016. 2017;10. 93. Song J, Li M, Tretiakova M et al. Expression patterns of PAX5, 112. Bauderlique-Le Roy H, Vennin C, Brocqueville G et al. En- c-Met, and paxillin in neuroendocrine tumors of the lung. richment of human Stem-Like prostate cells with s-SHIP Arch Pathol Lab Med 2010;134:1702–5. promoter activity uncovers a role in stemness for the long 94. Czapiewski P, Gorczynski ´ A, Haybaeck J et al. Expression noncoding RNA H19. Stem Cells Dev 2015;24:1252–62. pattern of ISL-1, TTF-1 and PAX5 in olfactory neuroblas- 113. Zhao J, Liu Y, Zhang W et al. Long non-coding RNA toma. Pol J Pathol 2016;67:130–5. Linc00152 is involved in cell cycle arrest, apoptosis, epithe- 95. Kanteti R, Nallasura V, Loganathan S et al. PAX5 is ex- lial to mesenchymal transition, cell migration and invasion pressed in small-cell lung cancer and positively regulates in gastric cancer. Cell Cycle 2015;14:3112–23. c-Met transcription. Lab Invest 2009;89:301–14. 114. Emmrich S, Streltsov A, Schmidt F et al. LincRNAs MONC 96. Walter RF, Mairinger FD, Werner R et al. SOX4, SOX11 and and MIR100HG act as oncogenes in acute megakaryoblastic PAX6 mRNA expression was identified as a (prognostic) leukemia. Mol Cancer 2014;13:171. marker for the aggressiveness of neuroendocrine tumors 115. Yan K, Tian J, Shi W et al. LncRNA SNHG6 is associated with of the lung by using next-generation expression analysis poor prognosis of gastric cancer and promotes cell prolif- (NanoString). Future Oncol 2015;11:1027–36. eration and EMT through epigenetically silencing p27 and 97. Quandt K, Frech K, Karas H et al. MatInd and MatInspec- sponging miR-101-3p. Cell Physiol Biochem 2017;42:999– tor: new fast and versatile tools for detection of consen- 1012. sus matches in nucleotide sequence data. Nucleic Acids Res 116. Erho N, Crisan A, Vergara IA et al. Discovery and valida- 1995;23:4878–84. tion of a prostate cancer genomic classifier that predicts 98. Cartharius K, Frech K, Grote K et al. MatInspector and be- early metastasis following radical prostatectomy. PLoS One yond: promoter analysis based on transcription factor bind- 2013;8:e66855. ing sites. Bioinformatics 2005;21:2933–42. 117. Karnes RJ, Bergstralh EJ, Davicioni E et al. Validation of a 99. Markoff A. Analytical Tools for DNA, Genes and Genomes: genomic classifier that predicts metastasis following rad- Nuts & Bolts, 1st ed. Eagleville, Pa. USA, DNA Press, 2005. ical prostatectomy in an at risk patient population. J Urol 100. Yang L, Lin C, Jin C et al. lncRNA-dependent mechanisms 2013;190:2047–53. of androgen-receptor-regulated gene activation programs. 118. Beltran H, Tagawa ST, Park K et al. Challenges in recog- Nature 2013;500:598–602. nizing treatment-related neuroendocrine prostate cancer. 101. Cui Z, Ren S, Lu J et al. The prostate cancer-up-regulated J Clin Oncol 2012;30:e386–389. long noncoding RNA PlncRNA-1 modulates apoptosis and 119. Terai G, Iwakiri J, Kameda T et al. Comprehensive prediction proliferation through reciprocal regulation of androgen re- of lncRNA-RNA interactions in human transcriptome. BMC ceptor. Urol Oncol 2013;31:1117–23. Genomics 2016;17(Suppl 1):12. 102. Crea F, Watahiki A, Quagliata L et al. Identification of a long 120. Kiryu H, Terai G, Imamura O et al. A detailed investigation non-coding RNA as a novel biomarker and potential ther- of accessibilities around target sites of siRNAs and miRNAs. apeutic target for metastatic prostate cancer. Oncotarget Bioinformatics 2011;27:1788–97. 2014, 5, 764–74; 121. Busch A, Richter AS, Backofen R. IntaRNA: efficient predic- 103. Malik R, Patel L, Prensner JR et al. The lncRNA PCAT29 in- tion of bacterial sRNA targets incorporating target site ac- hibits oncogenic phenotypes in prostate cancer. Mol Cancer cessibility and seed regions. Bioinformatics 2008;24:2849– Res 2014;12(8):1081–7. 56. 104. Wan X, Huang W, Yang S et al. Identification of androgen- 122. Kato Y, Sato K, Hamada M et al. RactIP: fast and accurate responsive lncRNAs as diagnostic and prognostic markers prediction of RNA-RNA interaction using integer program- for prostate cancer. Oncotarget 2016;7(37):60503–18. ming. Bioinformatics 2010;26:i460–466. 105. Lu W, Zhou D, Glusman G et al. KLK31P is a novel androgen 123. Helpap B, Kollermann J, Oehler U. Neuroendocrine differen- regulated and transcribed pseudogene of kallikreins that is tiation in prostatic carcinomas: histogenesis, biology, clin- expressed at lower levels in prostate cancer cells than in ical relevance, and future therapeutical perspectives. Urol normal prostate cells. Prostate 2006;66:936–44. Int 1999;62:133–8. 106. Misawa A, Takayama KI, Fujimura T et al. Androgen- 124. Hirano D, Okada Y, Minei S et al. Neuroendocrine differ- induced lncRNA POTEF-AS1 regulates apoptosis-related entiation in hormone refractory prostate cancer following pathway to facilitate cell survival in prostate cancer cells. androgen deprivation therapy. Eur Urol 2004;45:586–92; dis- Cancer Sci 2017;108:373–9. cussion 592. 107. Takayama K, Horie-Inoue K, Katayama S et al. Androgen- 125. Berruti A, Mosca A, Porpiglia F et al. Chromogranin A ex- responsive long noncoding RNA CTBP1-AS promotes pression in patients with hormone naive prostate cancer Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 22 NEPC lncRNAs predicts the development of hormone refractory disease. J 144. Camacho N, Van Loo P, Edwards S et al. Appraising the rel- Urol 2007;178:838–43; quiz 1129. evance of DNA copy number loss and gain in prostate can- 126. Aggarwal R, Zhang T, Small EJ et al. Neuroendocrine cer using whole genome DNA sequence data. PLos Genet prostate cancer: subtypes, biology, and clinical outcomes. 2017;13:e1007001. J Natl Compr Canc NE 2014;12:719–26. 145. Schuettengruber B, Chourrout D, Vervoort M et al. Genome 127. Yuan TC, Veeramani S, Lin MF. Neuroendocrine-like regulation by polycomb and trithorax proteins. Cell prostate cancer cells: neuroendocrine transdifferentiation 2007;128:735–45. of prostate adenocarcinoma cells. Endocr Relat Cancer 146. Khalil AM, Guttman M, Huarte M et al. Many human 2007;14:531–47. large intergenic noncoding RNAs associate with chromatin- 128. Terry S, Maille´ P, Baaddi H et al. Cross modulation between modifying complexes and affect gene expression. PNAS the androgen receptor axis and protocadherin-PC in medi- 2009;106:11667–72. ating neuroendocrine transdifferentiation and therapeutic 147. Rinn JL, Kertesz M, Wang JK et al. Functional demarcation resistance of prostate cancer. Neoplasia 2013;15:761–72. of active and silent chromatin domains in human HOX loci 129. Huss WJ, Gregory CW, Smith GJ. Neuroendocrine cell differ- by noncoding RNAs. Cell 2007;129:1311–23. entiation in the CWR22 human prostate cancer xenograft: 148. Li J, Han L, Roebuck P et al. TANRIC: an interactive open plat- association with tumor cell proliferation prior to recur- form to explore the function of lncRNAs in cancer. Cancer rence. Prostate 2004;60:91–97. Res 2015;75:3728–37. 130. Vashchenko N, Abrahamsson PA. Neuroendocrine differen- 149. Ren X, Ustiyan V, Pradhan A et al. FOXF1 transcription tiation in prostate cancer: implications for new treatment factor is required for formation of embryonic vasculature modalities. Eur Urol 2005;47:147–55. by regulating VEGF signaling in endothelial cells. Circ Res 131. Aparicio A, Tzelepi V. Neuroendocrine (small-cell) carcino- 2014;115:709–20. mas: why they teach us essential lessons about prostate 150. Tamura M, Sasaki Y, Koyama R et al. Forkhead transcription cancer. Oncology 2014;28:831–8. factor FOXF1 is a novel target gene of the p53 family and 132. Beltran H, Tomlins S, Aparicio A et al. Aggressive vari- regulates cancer cell migration and invasiveness. Oncogene ants of castration-resistant prostate cancer. Clin Cancer Res 2014;33:4837–46. 2014;20:2846–50. 151. Sekaric P, Shamanin VA, Luo J et al. hAda3 regulates 133. Bishop JL, Davies A, Ketola K et al. Regulation of tumor cell p14ARF-induced p53 acetylation and senescence. Onco- plasticity by the androgen receptor in prostate cancer. En- gene 2007;26:6261–8. docr Relat Cancer 2015;22:R165–182. 152. Wang T, Kobayashi T, Takimoto R et al. hADA3 is required 134. Davies AH, Beltran H, Zoubeidi A. Cellular plasticity and the for p53 activity. EMBO J 2001;20:6404–13. neuroendocrine phenotype in prostate cancer. Nat Rev Urol 153. Lin N, Chang KY, Li Z et al. An evolutionarily conserved long 2018;15(5):271–86. noncoding RNA TUNA controls pluripotency and neural lin- 135. Zhang W, Liu B, Wu W et al. Targeting the MYCN-PARP-DNA eage commitment. Mol Cell 2014;53:1005–19. damage response pathway in neuroendocrine prostate can- 154. Lee HJ, Choi NY, Lee SW et al. Epigenetic alteration of im- cer. Clin Cancer Res 2017. printed genes during neural differentiation of germline- 136. Wang C, , Peng G, Huang H et al. Blocking the feed- derived pluripotent stem cells. Epigenetics 2016;11:177–83. back loop between neuroendocrine differentiation and 155. Tsuta K, Wistuba II, Moran CA. Differential expression of macrophages improves the therapeutic effects of enza- somatostatin receptors 1–5 in neuroendocrine carcinoma of lutamide (MDV3100) on prostate cancer. Clin Cancer Res the lung. Pathol Res Pract 2012;208:470–4. 2017;, 24, 708–723. 156. Muscarella LA, D’Alessandro V, la Torre A et al. Gene ex- 137. , Brzezniak C, Oronsky B, Aggarwal RA. Complete metabolic pression of somatostatin receptor subtypes SSTR2a, SSTR3 response of metastatic castration-resistant neuroen- and SSTR5 in peripheral blood of neuroendocrine lung can- docrine carcinoma of the prostate after treatment with cer affected patients. Cell Oncol 2011;34:435–41. RRx-001 and reintroduced platinum doublets. Eur Urol 157. Squires MH, 3rd, Volkan Adsay N, Schuster DM et al. Octre- 2017;, 73, 306–307. oscan versus FDG-PET for neuroendocrine tumor staging: a 138. Thakur MK, Heilbrun L, Dobson K et al. Phase I trial of biological approach. Ann Surg Oncol 2015;22:2295–301. the combination of docetaxel, prednisone, and pasireotide 158. Narayanan S, Kunz PL. Role of somatostatin analogues in in metastatic castrate-resistant prostate cancer. Clin Geni- the treatment of neuroendocrine tumors. Hematol/Oncol tourin Can 2018, Epub ahead of print. Clin North Am 2016;30:163–77. 139. Akamatsu S, Inoue T, Ogawa O et al. Clinical and molecu- 159. Sharma K, Patel YC, Srikant CB. C-terminal region of human lar features of treatment-related neuroendocrine prostate somatostatin receptor 5 is required for induction of Rb and cancer. Int J Urol 2018;, 25, 345–351. G1 cell cycle arrest. Mol Endocrinol 1999;13:82–90. 140. Stone L. Prostate cancer: a novel mechanism of neuroen- 160. Coffey K, Rogerson L, Ryan-Munden C et al. The lysine docrine transdifferentiation. Nat Rev Urol 2018;, 15, 262– demethylase, KDM4B, is a key molecule in androgen recep- 263. tor signalling and turnover. Nucleic Acids Res 2013;41:4433– 141. Huarte M. The emerging role of lncRNAs in cancer. Nat Med 46. 2015;21:1253–61. 161. Yang J, AlTahan AM, Hu D et al. The role of histone 142. Szafranski P, Dharmadhikari AV, Brosens E et al. Small non- demethylase KDM4B in Myc signaling in neuroblastoma. J coding differentially methylated copy-number variants, in- Natl Cancer Inst 2015;107:djv080. cluding lncRNA genes, cause a lethal lung developmental 162. Mo F, Lin D, Takhar M et al. Stromal gene expression is disorder. Genome Res 2013;23:23–33. predictive for metastatic primary prostate cancer. Eur Urol, 143. White NM, Cabanski CR, Silva-Fisher JM et al. Transcrip- 2017;, 73, 524–532. tome sequencing reveals altered long intergenic non- 163. Living Tumor Labratory, www.livingtumorlab.com/index.h coding RNAs in lung cancer. Genome Biol 2014;15:429. tml, Accessed January 2018 Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018 Ramnarine et al. 23 164. Tsai H, Morais CL, Alshalalfa M et al. Cyclin D1 loss distin- 2018 guishes prostatic small-cell carcinoma from most prostatic 174. Lai D, Proctor JR, Zhu JY et al. R-CHIE: a web server and R adenocarcinomas. Clin Cancer Res 2015;21:5619–29. package for visualizing RNA secondary structures. Nucleic 165. Gao J, Aksoy BA, Dogrusoz U et al. Integrative analysis of Acids Res 2012;40:e95. complex cancer genomics and clinical profiles using the 175. R-chie: A web server and R package for potting arc diagrams cBioPortal. Sci Signal 2013;6:pl1. of RNA secondary structures, http://www.e-rna.org/r-chie/, 166. Cerami E, Gao J, Dogrusoz U et al. The cBio cancer genomics Accessed January 2018 portal: an open platform for exploring multidimensional 176. Mathews DH, Disney MD, Childs JL et al. Incorporating cancer genomics data. Cancer Discov 2012;2:401–4. chemical modification constraints into a dynamic program- 167. cBioPortal for Cancer Genomics, www.cbioportal.org, Ac- ming algorithm for prediction of RNA secondary structure. cessed September 2017. PNAS 2004;101:7287–92. 168. Piccolo SR, Sun Y, Campbell JD et al. A single-sample mi- 177. Gruber AR, Lorenz R, Bernhart SH et al. The vienna RNA croarray normalization method to facilitate personalized- websuite. Nucleic Acids Res 2008;36:W70–74. medicine workflows. Genomics 2012; 100:337–44. 178. RNAfold web server, http://rna.tbi.univie.ac.at/cgi-bin/RNA 169. Iyer MK, Niknafs YS, Malik R et al. The landscape of long WebSuite/RNAfold.cgi, Accessed January 2018 noncoding RNAs in the human transcriptome. Nat Genet 179. Genomatix AG, www.genomatix.de, Accessed August 2017 2015;47:199–208. 180. Ho Sui SJ, Mortimer JR, Arenillas DJ et al. oPOSSUM: identifi- 170. Eisenberg E, Levanon EY. Human housekeeping genes, re- cation of over-represented transcription factor binding sites visited. Trends Genet 2013;29:569–74. in co-expressed genes. Nucleic Acids Res 2005;33:3154–64. 171. Sickle: Windowed Adaptive Trimming for fastq files us- 181. SMALT - Wellcome Sanger Institute Tools and Soft- ing quality, https://github.com/ucdavis-bioinformatics/sick ware, http://www.sanger.ac.uk/science/tools/smalt-0, Ac- le, Accessed January 2018 cessed January 2018 172. UC Davis Bioinformatics Core Software, http://bioinforma 182. Ramnarine VR, Alshalalfa M, Mo F et al. Supporting data for tics.ucdavis.edu/research-computing/software/, Accessed “The Long Noncoding RNA Landscape of Neuroendocrine January 2018 Prostate Cancer and its Clinical Implications” GigaScience 173. Bioinformatics Webserver for RNA-RNA Interaction, http:// Database. 2018. http://doi.org/10.5524/100443. rtools.cbrc.jp/cgi-bin/RNARNA/index.pl, Accessed January Downloaded from https://academic.oup.com/gigascience/article-abstract/7/6/giy050/4994835 by Ed 'DeepDyve' Gillespie user on 21 June 2018

Journal

GigaScienceOxford University Press

Published: May 10, 2018

There are no references for this article.

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off