Transcriptome sequencing across a prostate cancer cohort identifies PCAT-1, an unannotated lincRNA implicated in disease progression.

Noncoding RNAs (ncRNAs) are emerging as key molecules in human cancer, with the potential to serve as novel markers of disease and to reveal uncharacterized aspects of tumor biology. Here we discover 121 unannotated prostate cancer-associated ncRNA transcripts (PCATs) by ab initio assembly of high-throughput sequencing of polyA(+) RNA (RNA-Seq) from a cohort of 102 prostate tissues and cells lines. We characterized one ncRNA, PCAT-1, as a prostate-specific regulator of cell proliferation and show that it is a target of the Polycomb Repressive Complex 2 (PRC2). We further found that patterns of PCAT-1 and PRC2 expression stratified patient tissues into molecular subtypes distinguished by expression signatures of PCAT-1-repressed target genes. Taken together, our findings suggest that PCAT-1 is a transcriptional repressor implicated in a subset of prostate cancer patients. These findings establish the utility of RNA-Seq to identify disease-associated ncRNAs that may improve the stratification of cancer subtypes.

Introduction

Recently, next generation transcriptome sequencing (RNA-Seq) has provided a method to delineate the entire set of transcriptional aberrations in a disease, including novel transcripts and non-coding RNAs (ncRNAs) not measured by conventional analyses[1-5]. To facilitate interpretation of sequence read data, existing computational methods typically process individual samples using either short read gapped alignment followed by ab initio reconstruction[2, 3], or de novo assembly of read sequences followed by sequence alignment[4, 5]. These methods provide a powerful framework to uncover uncharacterized RNA species, including antisense transcripts, short RNAs <250 bps, or long ncRNAs (lincRNAs) >250 bps.

While still largely unexplored, ncRNAs, particularly lincRNAs, have emerged as a new aspect of biology, with evidence suggesting that they are frequently cell-type specific, contribute important functions to numerous systems[6, 7], and may interact with known cancer genes such as EZH2[8]. Indeed, several well-described examples, such as HOTAIR[8, 9] and ANRIL[10, 11], indicate that ncRNAs may be essential actors in cancer biology, typically facilitating epigenetic gene repression via chromatin modifying complexes[12, 13]. Moreover, ncRNA expression may confer clinical information about patient outcomes and have utility as diagnostic tests[9, 14]. The characterization of RNA species, their functions, and their clinical applicability is therefore a major area of biological and clinical importance.

Here, we describe a comprehensive analysis of lincRNAs in 102 prostate cancer tissue samples and cell lines by RNA-Seq. We employ ab initio computational approaches to delineate the annotated and unannotated transcripts in this disease, and we find 121 ncRNAs, termed Prostate Cancer Associated Transcripts (PCATs), whose expression patterns distinguish benign, localized cancer, and metastatic cancer samples. Notably, we discover PCAT-1, a novel prostate cancer ncRNA alternately demonstrating either repression by PRC2 or an active role in promoting cell proliferation through transcriptional regulation of target genes. Our findings describe the first comprehensive study of lincRNAs in prostate cancer, provide a computational framework for large-scale RNA-Seq analyses, and describe PCAT-1 as a novel prostate cancer ncRNA functionally implicated in disease progression.

Results

RNA-Seq analysis of the prostate cancer transcriptome

Over two decades of research has generated a genetic model of prostate cancer based on numerous neoplastic events, such as loss of the PTEN[15] tumor suppressor gene and gain of oncogenic ETS transcription factor gene fusions[16-18] in large subsets of prostate cancer patients. We hypothesized that prostate cancer similarly harbored disease-associated ncRNAs in molecular subtypes.

To pursue this hypothesis, we employed transcriptome sequencing on a cohort of 102 prostate tissues and cell lines (20 benign adjacent prostates (benign), 47 localized tumors (PCA), and 14 metastatic tumors (MET) and 21 prostate cell lines). From a total of 1.723 billion sequence fragments from 201 lanes of sequencing (108 paired-end, 93 single read on the Illumina Genome Analyzer and Genome Analyzer II), we performed short read gapped alignment[19] and recovered 1.41 billion mapped reads, with a median of 14.7 million mapped reads per sample ( for sample information). We used the Cufflinks ab initio assembly approach[3] to produce, for each sample, the most probable set of putative transcripts that served as the RNA templates for the sequence fragments in that sample ( and and ).

As expected from a large tumor tissue cohort, individual transcript assemblies may exhibit sources of “noise”, such as artifacts of the sequence alignment process, unspliced intronic pre-mRNA, and genomic DNA contamination. To exclude these from our analyses, we trained a decision tree to classify transcripts as “expressed” versus “background” on the basis of transcript length, number of exons, recurrence in multiple samples, and other structural characteristics ( and ). The classifier demonstrated a sensitivity of 70.8% and specificity of 88.3% when trained using transcripts that overlapped genes in the AceView database[20], including 11.7% of unannotated transcripts that were classified as “expressed” (). We then clustered the “expressed” transcripts into a consensus transcriptome and applied additional heuristic filters to further refine the assembly (). The final ab initio transcriptome assembly yielded 35,415 distinct transcriptional loci ( and ).

Discovery of prostate cancer non-coding RNAs

We compared the assembled prostate cancer transcriptome to the UCSC, Ensembl, Refseq, Vega, and ENCODE gene databases to identify and categorize transcripts (). While the majority of the transcripts (77.3%) corresponded to annotated protein coding genes (72.1%) and non-coding RNAs (5.2%), a significant percentage (19.8%) lacked any overlap and were designated “unannotated” (). These included partially intronic antisense (2.44%), totally intronic (12.1%), and intergenic transcripts (5.25%), consistent with previous reports of unannotated transcription[21, 22, 23]. Due to the added complexity of characterizing antisense or partially intronic transcripts without strand-specific RNA-Seq libraries, we focused on totally intronic and intergenic transcripts.

Global characterization of novel intronic and intergenic transcripts demonstrated that they were more highly expressed (), had greater overlap with expressed sequence tags (ESTs) (), and displayed a clear but subtle increase in conservation over randomly permuted controls (novel intergenic transcripts p = 2.7 × 10-4 ± 0.0002 for 0.4 < ω < 0.8; novel intronic transcripts p = 2.6 × 10-5 ± 0.0017 for 0 < ω < 0.4, Fisher's exact test, ). By contrast, unannotated transcripts scored lower than protein-coding genes for these metrics, which corroborates data in previous reports[2, 24]. Interestingly, a small subset of novel intronic transcripts showed a profound degree of conservation (, insert). Finally, analysis of coding potential revealed that only 5 of 6,144 transcripts harbored a high quality open reading frame (ORF), indicating that the vast majority of these transcripts represent ncRNAs ().

To determine whether our unannotated transcripts were supported by histone modifications defining active transcriptional units, we used published prostate cancer ChIP-Seq data for two prostate cell lines[25], VCaP and LNCaP (). After filtering our dataset for transcribed repetitive elements known to display alternative patterns of histone modifications[26], we observed a strong enrichment for histone modifications characterizing transcriptional start sites (TSSs) and active transcription, including H3K4me2, H3K4me3, Acetyl-H3 and RNA polymerase II () but not H3K4me1, which characterizes enhancer regions[27] ( and ). Interestingly, intergenic ncRNAs showed greater enrichment compared to intronic ncRNAs in these analyses ().

To elucidate global changes in transcript abundance in prostate cancer, we performed a differential expression analysis for all transcripts. We found 836 genes differentially-expressed between benign samples and localized tumors (FDR < 0.01), with annotated protein-coding and ncRNA genes constituting 82.8% and 7.4% of differentially-expressed genes, respectively, including known prostate cancer biomarkers such AMACR[28], HPN[29], and PCA3[14] ( and ). Finally, 9.8% of differentially-expressed genes corresponded to unannotated ncRNAs, including 3.2% within gene introns and 6.6% in intergenic regions.

Characterization of Prostate Cancer Associated Transcripts

As ncRNAs may contribute to human disease[6-9], we identified aberrantly expressed uncharacterized ncRNAs in prostate cancer. We found a total of 1,859 unannotated lincRNAs throughout the human genome. Overall, these intergenic RNAs resided approximately half-way between two protein coding genes (), and over one-third (34.1%) were ≥10kb from the nearest protein-coding gene, which is consistent with previous reports[30] and supports the independence of intergenic ncRNAs genes. For example, visualizing the Chr15q arm using the Circos program (http://mkweb.bcgsc.ca/circos) illustrated genomic positions of eighty-nine novel intergenic transcripts, including one differentially-expressed gene centromeric to TLE3 ().

A focused analysis of the 1,859 unannotated intergenic RNAs yielded 106 that were differentially expressed in localized tumors (FDR < 0.05, ). A cancer outlier expression analysis () similarly nominated numerous unannotated ncRNA outliers () as well as known prostate cancer outliers, such as ERG[18], ETV1[[31] and CRISP3[32]. Merging these results produced a set of 121 unannotated transcripts that accurately discriminated benign, localized tumor, and metastatic prostate samples by unsupervised clustering (). Indeed, clustering analyses using novel ncRNA outliers also suggested disease subtypes (). These 121 unannotated transcripts were ranked and named as Prostate Cancer Associated Transcripts (PCATs) according to their fold change in localized tumor versus benign tissue ( and ).

Validation of novel ncRNAs

To gain confidence in our transcript nominations, we validated multiple unannotated transcripts in vitro by reverse transcription PCR (RT-PCR) and quantitative real-time PCR (qPCR) (). qPCR for four transcripts (PCAT-114, PCAT-14, PCAT-43, PCAT-1) on two independent cohorts of prostate tissues confirmed predicted cancer-specific expression patterns ( and ). Interestingly, all four are prostate-specific, with minimal expression seen by qPCR in breast (n=14) or lung cancer (n=16) cell lines or in 19 normal tissue types (). This is further supported by expression analysis of these transcripts in our RNA-Seq compendium of 13 tumor types, representing 325 samples (). This tissue specificity was not necessarily due to regulation by androgen receptor signaling, as only PCAT-14 expression was induced when androgen responsive VCaP and LNCaP cells were treated with the synthetic androgen R1881, consistent with previous data from this locus[17] (). PCAT-1 and PCAT-14 also showed cancer-specific upregulation when tested on a panel of matched tumor-normal samples ().

Of note, PCAT-114, which ranks as the #5 best outlier, just ahead of ERG ( and ), appears as part of a large, >500 kb locus of expression in a gene desert in Chr2q31. We termed this region Second Chromosome Locus Associated with Prostate-1 (SChLAP1) (). Careful analysis of the SChLAP1 locus revealed both discrete transcripts and intronic transcription, highlighting this region as an intriguing aspect of the prostate cancer transcriptome.

PCAT-1, a novel prostate cancer lincRNA

To explore several transcripts more closely, we performed 5’ and 3’ rapid amplification of cDNA ends (RACE) for PCAT-1 and PCAT-14. Interestingly, the PCAT-14 locus contained components of viral ORFs from the HERV-K endogenous retrovirus family (), whereas PCAT-1 incorporates portions of a mariner family transposase[33, 34], an Alu, and a viral long terminal repeat (LTR) promoter region ( and ). While PCAT-14 was upregulated in localized prostate cancer but largely absent in metastases (), PCAT-1 was strikingly upregulated in a subset of metastatic and high-grade localized (Gleason score ≥7) cancers ( and ). Because of this notable profile, we hypothesized that PCAT-1 may have coordinated expression with the oncoprotein EZH2, a core PRC2 protein that is upregulated in solid tumors and contributes to a metastatic phenotype[35, 36]. Surprisingly, we found that PCAT-1 and EZH2 expression were nearly mutually exclusive (), with only one patient showing outlier expression of both. This suggests that outlier PCAT-1 and EZH2 expression may define two subsets of high-grade disease.

PCAT-1 is located in the chromosome 8q24 gene desert approximately 725 kb upstream of the c-MYC oncogene. To confirm that PCAT-1 is a non-coding gene, we cloned the full-length PCAT-1 transcript and performed in vitro translational assays, which were negative as expected (). Next, since Chr8q24 is known to harbor prostate cancer-associated single nucleotide polymorphisms (SNPs) and to exhibit frequent chromosomal amplification[37-42], we evaluated whether the relationship between EZH2 and PCAT-1 was specific or generalized. To address this, we measured expression levels of c-MYC and NCOA2, two proposed targets of Chr8q amplification[39, 42], by qPCR. Neither c-MYC nor NCOA2 levels showed striking expression relationships to PCAT-1, EZH2, or each other (). Likewise, PCAT-1 outlier expression was not dependent on Chr8q24 amplification, as highly expressing localized tumors often did not have 8q24 amplification and high copy number gain of 8q24 was not sufficient to upregulate PCAT-1 ( and ).

PCAT-1 Function and Regulation

Despite reports showing that upregulation of the ncRNA HOTAIR participates in PRC2 function in breast cancer[9], we do not observe strong expression of this ncRNA in prostate (), suggesting that other ncRNAs may be important in this cancer. To determine the mechanism for the expression profiles of PCAT-1 and EZH2, we inhibited EZH2 activity in VCaP cells, which express low-to-moderate levels of PCAT-1. Knockdown of EZH2 by shRNA or pharmacologic inhibition of EZH2 with the inhibitor 3-deazaneplanocin A (DZNep) caused a dramatic upregulation in PCAT-1 expression levels (), as did treatment of VCaP cells with the demethylating agent 5’deoxyazacytidine, the histone deacetylase inhibitor SAHA, or both (). Chromatin immunoprecipitation (ChIP) assays also demonstrated that SUZ12, a core PRC2 protein, directly binds the PCAT-1 promoter approximately 1kb upstream of the TSS (). Interestingly, RNA immunoprecipitation (RIP) similarly showed binding of PCAT-1 to SUZ12 protein in VCaP cells (). RIP assays followed by RNase A, RNase H, or DNase I treatment either abolished, partially preserved, or totally preserved this interaction, respectively (). This suggests that PCAT-1 exists primarily as a single-stranded RNA and secondarily as a RNA/DNA hybrid.

To explore the functional role of PCAT-1 in prostate cancer, we stably overexpressed full length PCAT-1 or controls in RWPE benign immortalized prostate cells. We observed a modest but consistent increase in cell proliferation when PCAT-1 was overexpressed at physiological levels ( and ). Next, we designed siRNA oligos to PCAT-1 and performed knockdown experiments in LNCaP cells, which express higher levels of PCAT-1 without PRC2-mediated repression (). Supporting our overexpression data, knockdown of PCAT-1 with three independent siRNA oligos resulted in a 25% - 50% decrease in cell proliferation in LNCaP cells (), but not control DU145 cells lacking PCAT-1 expression () or VCaP cells, in which PCAT-1 is expressed but repressed by PRC2 ().

Gene expression profiling of LNCaP knockdown samples on cDNA microarrays indicated that PCAT-1 modulates the transcriptional regulation of 370 genes (255 upregulated, 115 downregulated; FDR ≤ 0.01) ( and ). Gene ontology analysis of the upregulated genes showed preferential enrichment for cellular processes such as mitosis and cell cycle, whereas the downregulated genes had no concepts showing statistical significance ( and ). These results suggest that PCAT-1's function is predominantly repressive in nature, similar to other lincRNAs. We next validated expression changes in three key PCAT-1 target genes (BRCA2, CENPE and CENPF) whose expression is upregulated upon PCAT-1 knockdown () in LNCaP and VCaP cells, the latter of which appear less sensitive to PCAT-1 knockdown likely due to lower overall expression levels of this transcript.

PCAT-1 signatures in prostate cancer

Because of the regulation of PCAT-1 by PRC2 in VCaP cells, we hypothesized that knockdown of EZH2 would also downregulate PCAT-1 targets as a secondary phenomenon due to the subsequent upregulation of PCAT-1. Simultaneous knockdown of PCAT-1 and EZH2 would thus abrogate expression changes in PCAT-1 target genes. Performing this experiment in VCaP cells demonstrated that PCAT-1 target genes were indeed downregulated by EZH2 knockdown, and that this change was either partially or completely reversed using siRNA oligos to PCAT-1 (), lending support to the role of PCAT-1 as a transcriptional repressor. Taken together, these results suggest that PCAT-1 biology may exhibit two distinct modalities: one in which PRC2 represses PCAT-1 and a second in which active PCAT-1 promotes cell proliferation. PCAT-1 and PRC2 may therefore characterize distinct subsets of prostate cancer.

To examine our clinical cohort, we used qPCR to measure expression of BRCA2, CENPE, and CENPF in our tissue samples. Consistent with our model, we found that PCAT-1-expressing samples tended to have low expression of PCAT-1 target genes (). Moreover, comparing EZH2-outlier and PCAT-1-outlier patients (see ), we found that two distinct patient phenotypes emerged: those with high EZH2 tended to have high levels of PCAT-1 target genes; and those with high PCAT-1 expression displayed the opposite expression pattern (). Network analysis of the top 20 upregulated genes following PCAT-1 knockdown with the HefaLMP tool[43] further suggested that these genes form a coordinated network (), corroborating our previous observations. Taken together, these results provide initial data into the composition and function of the prostate cancer ncRNA transcriptome.

Discussion

This study represents the largest RNA-Seq analysis to date and the first to comprehensively analyze a common epithelial cancer from a large cohort of human tissue samples. As such, our study has adapted existing computational tools intended for small-scale use[3] and developed new methods in order to distill large numbers of transcriptome datasets into a single consensus transcriptome assembly that reflects a coherent biological picture.

Among the numerous uncharacterized ncRNA species detected by our study, we have focused on 121 prostate cancer-associated PCATs, which we believe represent a set of uncharacterized ncRNAs that may have important biological functions in this disease. In this regard, these data contribute to a growing body of literature supporting the importance of unannotated ncRNA species in cellular biology and oncogenesis[6-12], and broadly our study confirms the utility of RNA-Seq in defining functionally-important elements of the genome[2-4].

Of particular interest is our discovery of the prostate-specific ncRNA gene PCAT-1, which is markedly overexpressed in a subset of prostate cancers, particularly metastases, and may contribute to cell proliferation in these tumors. It is also notable that PCAT-1 resides in the 8q24 “gene desert” locus, in the vicinity of well-studied prostate cancer risk SNPs and the c-MYC oncogene, suggesting that this locus—and its frequent amplification in cancer—may be linked to additional aspects of cancer biology. In addition, the interplay between PRC2 and PCAT-1 further suggests that this ncRNA may have an important role in prostate cancer progression (). Other ncRNAs identified by this analysis may similarly contribute to prostate cancer as well. Furthermore, recent pre-clinical efforts to detect prostate cancer non-invasively through the collection of patient urine samples have shown promise for several urine-based prostate cancer biomarkers, including the ncRNA PCA3[44, 45]. While additional studies are needed, our identification of ncRNA biomarkers for prostate cancer suggests that urine-based assays for these ncRNAs may also warrant investigation, particularly for those that may stratify patient molecular subtypes.

Taken together, our findings support an important role for tissue-specific ncRNAs in prostate cancer and suggest that cancer-specific functions of these ncRNAs may help to “drive” tumorigenesis. We further speculate that specific ncRNA signatures may occur universally in all disease states and applying these methodologies to other diseases may reveal key aspects of disease biology and clinically important biomarkers.