Prostate cancer remains the second leading cause of cancer death in American men and there is an unmet need for biomarkers to identify patients with aggressive disease. In an effort to identify biomarkers of recurrence, we performed global RNA sequencing on 106 formalin-fixed, paraffin-embedded prostatectomy samples from 100 patients at three independent sites, defining a 24-gene signature panel. The 24 genes in this panel function in cell-cycle progression, angiogenesis, hypoxia, apoptosis, PI3K signaling, steroid metabolism, translation, chromatin modification, and transcription. Sixteen genes have been associated with cancer, with five specifically associated with prostate cancer (BTG2, IGFBP3, SIRT1, MXI1, and FDPS). Validation was performed on an independent publicly available dataset of 140 patients, where the new signature panel outperformed markers published previously in terms of predicting biochemical recurrence. Our work also identified differences in gene expression between Gleason pattern 4 + 3 and 3 + 4 tumors, including several genes involved in the epithelial-to-mesenchymal transition and developmental pathways. Overall, this study defines a novel biomarker panel that has the potential to improve the clinical management of prostate cancer. Cancer Res; 74(12); 3228–37. ©2014 AACR.

Prostate cancer remains the most common cancer diagnosed for U.S. males, and ranks second among tumor site–specific mortality with over 28,000 deaths per year (1). Currently, both clinicians and patients are faced with difficult treatment decisions when a diagnosis is made of prostate cancer because it is very hard to predict whether a patient is likely to progress to aggressive, metastatic disease. The most common clinical parameters that are used to make prognoses for patients with prostate cancer include age, stage, prostate-specific antigen (PSA) level, positive margins, perineural invasion, and Gleason score. Gleason score is generated by summing the grades of the two most prominent patterns, and somewhat takes into account the heterogeneity of prostate cancers (2, 3). Patients with a Gleason score of eight or higher generally have a poor prognosis, and those with Gleason scores of six or lower generally have favorable prognoses, but it is much more difficult to predict outcomes for the very large number of patients with a Gleason score of seven. Several studies have supported the concept that the Gleason 4 + 3 pattern tumors have worse outcomes than Gleason 3 + 4 tumors (4–6), but clinical parameters alone are not sufficient for accurate prediction of recurrence. Thus, biomarkers that can predict the likely clinical outcome and recurrence for patients after surgical therapy are urgently needed to aid clinicians in deciding on the appropriate course of treatment and in accurate prognostication for patients.

We recently published a panel of 10 protein-coding genes and two miRNA genes that could be used to separate patients with and without biochemical recurrence (BCR) after prostatectomy (7). Postoperative BCR is defined as two detectable PSA readings (>0.2 ng/mL), and generally precedes clinical recurrence with metastases, although BCR does not always predict clinical recurrence, and patients with BCR do not always have clinical metastases. Importantly, this panel of 12 biomarkers could significantly predict clinical recurrence for patients with a Gleason score of seven, but this analysis was limited to a selected set of approximately 522 prostate cancer–relevant genes and 1,146 human miRNAs. To take a more global approach to identify biomarkers of recurrence, we collected formalin-fixed paraffin-embedded (FFPE) prostatectomy samples from three independent sites (Atlanta Veterans Administration Medical Center, AVAMC; Sunnybrook Health Sciences Centre at the University of Toronto; and Moffitt Cancer Center, MCC), and prepared libraries for RNA sequencing. We performed RNA sequencing on 106 prostatectomy samples from 100 patients (61 from AVAMC, 29 from Toronto, and 10 from Moffitt) using the Illumina HiSeq2000 instrument, and identified a new set of biomarkers of BCR composed of a 24-gene panel. This panel showed significant improvement of prediction of BCR over clinical parameters alone and over previously described biomarker panels based on cellular proliferation (8).

Sample selection

We identified 150 cases at the AVAMC hospital between 1990 and 2000 that could potentially be used for this project. We were able to locate slides and FFPE blocks for 100 of those cases, and identified regions of cancer and benign tissue in slides for each of them. These samples were then submitted for processing to obtain 1-mm tissue cores. We prepared RNA from 99 AVAMC samples, of which 79 passed our RNA quality control (QC) analysis for genomic profiling. We also obtained 81 samples from the Sunnybrook Health Science Center (Toronto, ON) and 99 RNA samples from the MCC in Tampa, FL. Approximately 80% of samples passed QC analysis.

Those samples that were included had sufficient RNA yield, met specific inclusion criteria, had available tissue specimens, documented long-term follow-up, and consented to participate or were included by Institutional Review Board consent waiver. MCC subjects were patients with prostate cancer who underwent radical prostatectomy at the MCC between 1987 and 2003. Prostate cancer cases were men 21 years and older who had surgery for their disease at the MCC and had pathologically confirmed primary prostate cancer. AVAMC cases were patients with prostate cancer who underwent radical prostatectomy between 1990 and 2000. University of Toronto (UT) cases were patients with prostate cancer who underwent radical prostatectomy at the Sunnybrook Health Science Center between 1998 and 2006. Exclusion criterion included recurrent patients who underwent hormonal or radiation treatment before radical prostatectomy. The cases were assigned prostate ID numbers to protect their identities. These patients did not receive neoadjuvant or concomitant hormonal therapy. Their demographic, treatment, and long-term clinical outcome data have been collected and recorded in an electronic database. Clinical data recorded include PSA measurements, radiological studies and findings, clinical findings, tissue biopsies, and additional therapies that the subjects may have received. In total, 106 samples from 100 patients were sequenced by next-generation sequencing: 61 from the AVAMC, 35 from Toronto, and 10 from MCC. Of these cases, 49 had BCR and 51 had no BCR. Of those 100 patients, 97 had complete PSA data and were included in construction of biomarker prediction models. We used only 10 of the 99 MCC cases for several reasons including limited RNA yield, the fact that only three cases had BCR, and funding limitations to conduct sequencing analyses.

RNA preparation

Tissue cores (1 mm) were used for RNA preparation (two cores per case) rather than sections because of the heterogeneity of samples and the opportunity for obtaining cores with very high percentage tumor content, except for the 10 samples from MCC. MCC supplied total FFPE RNA that was prepared from 5-μm unstained sections. Hematoxylin and eosin-stained slides were reviewed by a board-certified anatomic pathologist with expertise in genitourinary pathology (A.O. Osunkoya for AVAMC and L. Sugar for UT) to identify regions of cancer to select corresponding areas for cutting of cores from paraffin blocks as previously described (9). Total RNA was prepared at the Winship Cancer Genomics Shared Resource from FFPE cores as previously described (10), using the Omega Biotek FFPE RNA methodology in 96-well format on a MagMax 96 Liquid Handler Robot (Life Technologies). FFPE RNA was quantitated using a Nanodrop spectrophotometer, and tested for RNA integrity and quality by Taqman analysis of the RPL13a ribosomal protein on a HT7900 real-time PCR instrument (Applied Biosystems). Samples with sufficient yield (>500 ng), A260/A280 ratio > 1.8, and RPL13a CT values less than 31 cycles were used for preparation of RNAseq libraries.

RNA sequencing

We prepared RNA sequencing libraries using 106 prostatectomy RNA samples from 100 patients and performed 50-bp paired-end sequencing using the Illumina HiSeq platform. RNA sequencing libraries were prepared from 2 μg of total RNA using the TruSeq Kit (Illumina, Inc.) with the following modification. Instead of purifying poly-A RNA using poly-dT primer beads, we removed ribosomal RNA using the Ribominus Kit (Ambion). All other steps were performed according to the manufacturer's protocols. RNAseq libraries were analyzed for QC and average size of inserts was approximately 200 to 300 bp. Samples were multiplexed into three samples per lane on the Illumina Version 3 Flowcells on the HiSeq2000 platform. Samples were sequenced at the Emory GRA Genomics Core Facility, Hudson Alpha Institute, and at the Southern California Genotyping Consortium. The complete dataset can be accessed in the Gene Expression Omnibus database with accession number GSE54460.

Bioinformatics analysis

FASTQ files generated from the Illumina HiSeq2000 were analyzed for quality using FASTQC, mapped to the human genome using TopHat software (version 2.0.8) and Bowtie (version 2.1.0; ref. 11), and Cufflinks software (12) was used to generate fragment per kilobase of transcript per million reads (FPKM) values. Genes were filtered to determine if they were detected (defined as FPKM > 1) in each sample. Genes that were detected in 80% of BCR samples or 80% of non-BCR samples were retained, leaving a set of 5,265 genes for further analysis. Duplicates from six patient samples were removed, and three patients had incomplete clinical data on PSA and pathologic stage, leaving 97 samples for biomarker analysis. Pathway and upstream analyses were performed using the Ingenuity Pathway Analysis (IPA) Knowledgebase (13). Transcriptome coverage was computed based on a transcriptome size of 81 Mb from the University of California at Santa Cruz (UCSC) known genes set (build hg19). We performed differential gene expression analysis between Gleason pattern 4 + 3 and Gleason pattern 3 + 4 tumors using DESeq (14) in R Bioconductor.

TaqMan validation

Individual TaqMan assays were obtained for TMSB10 (Hs1005565), PTN (Hs383235), and SYNM (Hs2326749) and FFPE RNA from 22 samples were tested for correlation with FPKM values. Direct Ct values were used for correlation analysis without a reference gene, because substantial variability was observed for both beta actin (ACTB) and 18S with constant RNA inputs.

Statistical analysis

PSA data were skewed with extreme outliers and were log-transformed, log(PSA + 1), in all analyses. To identify important RNA biomarkers and build and evaluate prediction models for prostate cancer recurrence (i.e., BRC), we adopted the following strategy, similar to our earlier work (7). In the training step, the prediction models were built for time to BRC using our RNAseq data. Specifically, we first fit a univariate Cox proportional hazard (PH) model for each individual RNA biomarker, and a set of important RNA biomarkers were then preselected based on a false discovery rate (FDR) threshold of 0.10. Next, we fit a lasso Cox PH model (15, 16) using the preselected RNA biomarkers to identify the optimal panel of RNA biomakers and the resulting prediction score, where the tuning parameter for lasso was selected using a leave-one-out cross-validation technique. On the basis of the RNA biomarker panel, a final prediction model for recurrence was built to also include relevant clinical biomarkers, including pathologic stage, PSA, Gleason score, surgical margin status, and age, through fitting a Cox PH model. For comparison, we also built a prediction model using only clinical information, including pathologic stage, PSA, Gleason score, and age, through fitting a Cox PH model. In the validation step, we first used our samples to evaluate each prediction model; specifically, we calculated the predictive scores for subjects based on each prediction model, divided subjects into high (poor score) and low (good score) risk groups based on the median predictive score, and performed the log-rank test to compare BCR between the two risk groups. Subsequently, we repeated the validation procedure using an independent gene expression microarray study with data from 140 patients with prostate cancer (17). To compare performance of different prediction models, we used three metrics for evaluting predictive accuracy, integrated discrimination improvement (IDI), net reclassification improvement (NRI), and median improvement in risk score for censored survival outcomes as described (18).

Genome-wide RNA sequencing analysis

We identified prostatectomy cases and FFPE blocks at the AVAMC hospital between 1990 and 2000, and also obtained prostatectomy FFPE cores from the Sunnybrook Health Science Center (Toronto, ON) and prostatectomy FFPE RNA samples from the MCC. We prepared RNA from two FFPE 1-mm cores per block, and following QC analysis, prepared a total of 106 RNAseq libraries from 100 patients: 61 from the AVAMC, 35 libraries from 29 patients from Toronto, and 10 from MCC. Of these cases, 49 had BCR and 51 had no BCR. Of these 100 cases, 97 had complete associated clinical data and were used for predictive model building. A summary of the clinical characteristics of these samples is provided in Supplementary Table S1. Detailed clinical data are provided in Supplementary Table S2.

Because RNA from FFPE samples is typically fragmented, conventional methods for sequencing library preparation that use poly-dT hybridization to capture fully processed mRNA are not optimal. We tested alternative library preparation methods and developed a protocol for robust RNAseq analysis of FFPE samples using the Ribominus Kit (Ambion) to remove ribosomal RNA, followed by library preparation using the Illumina TruSeq Kit. We determined that multiplexing three samples per lane gave us adequate coverage for RNAseq analysis. In total, we generated approximately 490 billion base pairs (Gbp) of sequence, of which 294 Gbp mapped uniquely to the human genome (build 19, hg19). In total, we obtained 5.874 billion mapped reads. The average number of mapped reads was 55.4M reads per sample, providing an average coverage of 34.2× for the human transcriptome (UCSC known genes hg19). A summary of the sequencing output is provided in Supplementary Table S3.

Possibly due to the fact that we did not perform a poly-A selection step, we observed a significant number of reads that mapped to gene introns, likely from partially processed mRNAs. We did not observe large numbers of reads mapping to intergenic regions, indicating that there was no DNA contamination in our samples. A representative screenshot for the mapped reads is given in Fig. 1. Moreover, the level of intronic reads was similar to that observed in RNAseq data derived from fresh-frozen samples analyzed by The Cancer Genome Atlas project (Supplementary Fig. S1). To validate and verify the accuracy of our RNAseq data, we performed TaqMan analyses on a few select genes and observed high correlation of TaqMan and RNAseq data for these genes (r2 = 0.80–0.97; Supplementary Fig. S2). In addition, we analyzed expression of genes that are typically involved in chromosomal translocations such as ERG, ETV1, ETV4, ETV5, and SPINK1. We observed mutually exclusive high levels of expression of ERG (45% of samples), ETV1 (6%), ETV4 (4%), and SPINK1 (10%) in proportions similar to previous studies (19–22; Fig. 2).

Figure 1.

Representative tracks of coverage of RNAseq reads from the IGFBP3 gene from the Integrated Genomics Viewer (IGV). Coverage is primarily focused in exons, but some intronic coverage is evident. Little, if any, coverage is observed in intergenic regions, indicating that intronic reads are likely from partially processed RNAs due to the lack of a poly-A selection step.

Figure 1.

Representative tracks of coverage of RNAseq reads from the IGFBP3 gene from the Integrated Genomics Viewer (IGV). Coverage is primarily focused in exons, but some intronic coverage is evident. Little, if any, coverage is observed in intergenic regions, indicating that intronic reads are likely from partially processed RNAs due to the lack of a poly-A selection step.

Close modal
Figure 2.

FPKM expression values for the ERG (blue), ETV1 (red), ETV4 (green), ETV5 (magenta), and SPINK1 (cyan) genes that have been previously shown to exhibit translocations in prostate cancers. Outlier samples are clearly evident, and increased expression is mutually exclusive as expected. Proportions of numbers of cases with increased expression of each gene also are in accordance with previous reports, with high expression of ERG in 45%, ETV1 in 6%, ETV4 in 4%, and SPINK1 (10%) of samples.

Figure 2.

FPKM expression values for the ERG (blue), ETV1 (red), ETV4 (green), ETV5 (magenta), and SPINK1 (cyan) genes that have been previously shown to exhibit translocations in prostate cancers. Outlier samples are clearly evident, and increased expression is mutually exclusive as expected. Proportions of numbers of cases with increased expression of each gene also are in accordance with previous reports, with high expression of ERG in 45%, ETV1 in 6%, ETV4 in 4%, and SPINK1 (10%) of samples.

Close modal

We also prepared RNA from separate cores from the same 6 patients, and prepared separate sequencing libraries on different days. Analysis of the fragment per kilobase of transcript per million mapped reads (FPKM) values from these replicate sequence analyses indicated very strong correlation (r2 = 0.70–0.96) for the 5,265 genes that were robustly detected in at least 80% of samples and used in our biomarker analyses (Supplementary Fig. S3). The pair of samples with the lowest correlation (UTPC034) had the greatest difference in number of mapped reads (18M vs. 112M), whereas the paired samples with the highest correlation (UTPC004) both had very deep coverage (94M and 101M mapped reads each). Differential gene expression analysis using DESeq indicated very few differentially expressed genes between replicate sequencing libraries (14 genes on average). For our biomarker analysis, in each case, we used the library with the higher number of mapped reads.

Biomarker analysis

FASTQ files generated from the Illumina HiSeq2000 were mapped to the human genome using TopHat software (version 2.0.8) and Bowtie (version 2.1.0; ref. 11), and Cufflinks software (12) was used to generate FPKM values. Genes were filtered to determine if they were detected (defined as FPKM > 1) in each sample. First, we analyzed the distribution of detected genes, and observed a peak of genes that were detected in 80% or more of the samples (Supplementary Fig. S4). To avoid excluding those genes that were not expressed in one of the two groups, we elected to retain genes detected in either 80% of BCR or 80% of non-BCR samples. Genes that were detected in 80% of BCR samples or 80% of non-BCR samples included a set of 5,265 protein-coding or noncoding genes for further analysis. Duplicates from six patient samples were removed, and three patients had incomplete clinical data on PSA and pathologic stage, leaving 97 samples for biomarker analysis.

Using the set of 5,265 genes from 97 samples, we built a 24-gene prediction model using a preselection step and a lasso Cox PH model (Table 1), and the final prediction model was built to include the predictive score based on this panel of 24 markers as well as the relevant clinical biomarkers, including pathologic stage, PSA, Gleason score, surgical margin status, and age (Table 2). For comparison, we also built a prediction model using only clinical information, namely, pathologic stage, PSA, Gleason score, surgical margin status, and age, through fitting a Cox PH model (Supplementary Table S4). We then performed log-rank tests to compare BCR between the low-risk (good score) and high-risk (poor score) groups with clinical variables only (Fig. 3A and B), and with both clinical parameters and RNAseq data (Fig. 3C and D). The Kaplan–Meier analysis (Fig. 3A and C) demonstrated that these markers could significantly discriminate patients at higher and lower risk of recurrence by the log-rank test (P = 1.45e−21) in our training data, more significant than using clinical variables alone (P = 5.39e−8). The improvements of the full model in Table 2 over the model using only clinical parameters in prediction as measured by IDI, NRI, and median improvement in risk score for censored survival outcomes (18) were all statistically significant (Table 3).

Figure 3.

Kaplan–Meier survival curves for 97 cases in the training set analyzed by RNAseq (A, C, and E) and for 140 cases in the validation set from Taylor and colleagues (17; B, D, and F). Kaplan–Meier curves using clinical parameters alone are shown in A and B. Kaplan–Meier curves using clinical parameters combined with the 24 RNA biomarker genes as described in Table 2 are shown in C and D. Kaplan–Meier curves using 31 biomarkers from Myriad Genetics are shown in E and F.

Figure 3.

Kaplan–Meier survival curves for 97 cases in the training set analyzed by RNAseq (A, C, and E) and for 140 cases in the validation set from Taylor and colleagues (17; B, D, and F). Kaplan–Meier curves using clinical parameters alone are shown in A and B. Kaplan–Meier curves using clinical parameters combined with the 24 RNA biomarker genes as described in Table 2 are shown in C and D. Kaplan–Meier curves using 31 biomarkers from Myriad Genetics are shown in E and F.

Close modal
Table 1.

RNA biomarkers of BCR following prostatectomy identified by RNAseq analysis. Selected genes and estimated coefficients using all expressed genomic markers

SymbolGene nameFunctionCancersCoefficient
BTG2 BTG family, member 2 Cell cycle, apoptosis Prostate, breast, lung −0.08831 
CDC37L1 Cell division cycle 37 homolog (S. cerevisiae)-like 1 Cell cycle Head and neck −0.13272 
COL15A1 Collagen, type XV, alpha 1 Angiongenesis, extracellular matrix Kidney 0.007957 
COL3A1 Collagen, type III, alpha 1 Extracellular matrix Breast, ovarian, lung, brain, gastric 0.27692 
EIF2D Eukaryotic translation initiation factor 2D Translation na 0.403303 
FDPS Farnesyl diphosphate synthase Sterol biosynthesis Breast, ovarian, head and neck 0.019765 
HIST1H1C Histone cluster 1, H1c Chromatin na 0.088821 
HIST1H2BG Histone cluster 1, H2bg Chromatin na 0.067189 
IFT57 Intraflagellar transport 57 homolog (Chlamydomonas) Apoptosis GBM 0.274924 
IGFBP3 Insulin-like growth factor binding protein 3 Insulin signaling PI3K Prostate, breast, colon, many 0.076922 
ITPR1 Inositol 1,4,5-triphosphate receptor, type 1 PIP3, apoptosis Prostate, breast, leukemia −0.27849 
LBH Limb bud and heart development homolog (mouse) Transcription Breast 0.098016 
LOC284801 LOC284801 Unknown na 0.022612 
MARCH5 Membrane-associated ring finger (C3HC4) 5 E3 ubiquitin ligase of the mitochondria na −0.25456 
MED4 Mediator complex subunit 4 Transcription, chromatin Breast, cervical −0.14 
MEMO1 Mediator of cell motility 1; similar to mediator of cell motility 1 Insulin signaling, EMT, estrogen signaling Breast, ALL 0.116525 
MXI1 MAX interactor 1 Transcription, MYC activity Prostate, melanoma, kidney, bladder −0.03462 
PTN Pleiotrophin Angiongenesis Prostate, colorectal, bladder, breast −0.16271 
RPL23AP53 Ribosomal protein L23a pseudogene 53 Ribosomes na −0.12988 
SACM1L SAC1 suppressor of actin mutations 1-like (yeast) Mitotic spindles, golgi ALL −0.02007 
SIRT1 Sirtuin (silent mating type information regulation 2 homolog) 1 (S. cerevisiae) Transcription, chromatin Prostate, breast, many others −0.08227 
SNORA20 Small nucleolar RNA, H/ACA box 20 Ribosomes na −0.10131 
SRSF3 Serine/arginine-rich splicing factor 3 Splicing, mRNA transport Colon, ovarian −0.02622 
SYNM Synemin, intermediate filament protein Intermediate filament Breast, GBM, liver −0.0574 
SymbolGene nameFunctionCancersCoefficient
BTG2 BTG family, member 2 Cell cycle, apoptosis Prostate, breast, lung −0.08831 
CDC37L1 Cell division cycle 37 homolog (S. cerevisiae)-like 1 Cell cycle Head and neck −0.13272 
COL15A1 Collagen, type XV, alpha 1 Angiongenesis, extracellular matrix Kidney 0.007957 
COL3A1 Collagen, type III, alpha 1 Extracellular matrix Breast, ovarian, lung, brain, gastric 0.27692 
EIF2D Eukaryotic translation initiation factor 2D Translation na 0.403303 
FDPS Farnesyl diphosphate synthase Sterol biosynthesis Breast, ovarian, head and neck 0.019765 
HIST1H1C Histone cluster 1, H1c Chromatin na 0.088821 
HIST1H2BG Histone cluster 1, H2bg Chromatin na 0.067189 
IFT57 Intraflagellar transport 57 homolog (Chlamydomonas) Apoptosis GBM 0.274924 
IGFBP3 Insulin-like growth factor binding protein 3 Insulin signaling PI3K Prostate, breast, colon, many 0.076922 
ITPR1 Inositol 1,4,5-triphosphate receptor, type 1 PIP3, apoptosis Prostate, breast, leukemia −0.27849 
LBH Limb bud and heart development homolog (mouse) Transcription Breast 0.098016 
LOC284801 LOC284801 Unknown na 0.022612 
MARCH5 Membrane-associated ring finger (C3HC4) 5 E3 ubiquitin ligase of the mitochondria na −0.25456 
MED4 Mediator complex subunit 4 Transcription, chromatin Breast, cervical −0.14 
MEMO1 Mediator of cell motility 1; similar to mediator of cell motility 1 Insulin signaling, EMT, estrogen signaling Breast, ALL 0.116525 
MXI1 MAX interactor 1 Transcription, MYC activity Prostate, melanoma, kidney, bladder −0.03462 
PTN Pleiotrophin Angiongenesis Prostate, colorectal, bladder, breast −0.16271 
RPL23AP53 Ribosomal protein L23a pseudogene 53 Ribosomes na −0.12988 
SACM1L SAC1 suppressor of actin mutations 1-like (yeast) Mitotic spindles, golgi ALL −0.02007 
SIRT1 Sirtuin (silent mating type information regulation 2 homolog) 1 (S. cerevisiae) Transcription, chromatin Prostate, breast, many others −0.08227 
SNORA20 Small nucleolar RNA, H/ACA box 20 Ribosomes na −0.10131 
SRSF3 Serine/arginine-rich splicing factor 3 Splicing, mRNA transport Colon, ovarian −0.02622 
SYNM Synemin, intermediate filament protein Intermediate filament Breast, GBM, liver −0.0574 

Abbreviations: na, not applicable; GBM, glioblastoma.

Table 2.

The prediction model using RNA biomarkers and clinical variables is shown, in which the predictive score is calculated from Table 1 

VariablesCoefficient
RNA predictive score 2.35147 
Gleason 0.123248 
Log(prePSA+1) 0.492755 
Age 0.012269 
SMS 0.651851 
pStage > 2 0.679813 
VariablesCoefficient
RNA predictive score 2.35147 
Gleason 0.123248 
Log(prePSA+1) 0.492755 
Age 0.012269 
SMS 0.651851 
pStage > 2 0.679813 

Abbreviation: SMS, surgical margin status.

Table 3.

Comparison of prediction performance of the full model, including the RNA biomarkers and clinical parameters with other models that include clinical variables alone, RNA biomarkers alone, or the MYRIAD model in terms of IDI, NRI, and median improvement in risk score for censored survival outcomes (18)

IDIP valueNRIP valueMedian improvement in risk scoreP value
Training set 
 Full model vs. clinical variables only 0.469 0.024 0.875 0.022 0.394 0.004 
 Full model vs. RNA biomarkers only 0.218 0.036 0.718 0.043 0.2 0.005 
 Full model vs. MYRIAD model 0.676 0.03 0.875 0.045 0.634 <0.001 
Validation set 
 Full model vs. clinical variables only 0.699 0.019 0.669 0.042 0.678 <0.001 
 Full model vs. RNA biomarkers only 0.005 0.643 −0.439 0.246 −0.002 0.986 
 Full model vs. MYRIAD model 0.601 0.027 0.625 0.051 0.652 <0.001 
IDIP valueNRIP valueMedian improvement in risk scoreP value
Training set 
 Full model vs. clinical variables only 0.469 0.024 0.875 0.022 0.394 0.004 
 Full model vs. RNA biomarkers only 0.218 0.036 0.718 0.043 0.2 0.005 
 Full model vs. MYRIAD model 0.676 0.03 0.875 0.045 0.634 <0.001 
Validation set 
 Full model vs. clinical variables only 0.699 0.019 0.669 0.042 0.678 <0.001 
 Full model vs. RNA biomarkers only 0.005 0.643 −0.439 0.246 −0.002 0.986 
 Full model vs. MYRIAD model 0.601 0.027 0.625 0.051 0.652 <0.001 

NOTE: A positive value in all three metrics, including IDI and NRI, indicates an improvement over the second model. Significant P values (bold) indicate a statistically significant improvement of the full model over other models in prediction of BCR.

To validate this panel of biomarkers, we identified an independent gene expression microarray study with data from 140 patients with prostate cancer (17). Using the data from Taylor and colleagues (17), we evaluated the final prediction models obtained from the training phase. Each prediction model from the training phase was used to generate a predictive score for each subject in the testing data set, and subjects were subsequently divided into the high- and low-risk groups using the median predictive score. Log-rank tests were performed to compare time to BCR between the high (poor score) and low (good score) risk groups. We first tested our prediction model based on clinical variables alone (Fig. 3B), showing significant discriminative performance in the validation data as well (P = 2.85e−3). In addition, we observed that the full panel including RNA biomarkers and clinical variables was very significantly prognostic in this independent validation set (P = 7.87e−5, Fig. 3D).

Comparison with existing biomarkers

In addition, we performed a direct comparison of our 24 biomarker genes with a set of 31 cell-cycle progression genes developed by Myriad Genetics (8). The Myriad panel performed less well in our dataset (P = 4.94e−8, Fig. 3E) than our biomarker panel at a level similar to clinical parameters alone. In the validation set from Taylor and colleagues (17), our panel of biomarker genes outperformed the 31 cell-cycle progression genes, which had less significant P value for the log-rank test (P = 1.4e−4, Fig. 3F). Furthermore, analysis for IDI, NRI, and median improvement in risk score demonstrated that our panel was statistically significantly better in prediction of recurrence than the Myriad panel (Table 3).

Pathway analysis

To evaluate the pathways that might be associated with BCR, we selected a list of 894 genes with a univariate Cox PH q value < 0.10 (Supplementary Table S5) and subjected it to pathway enrichment and upstream analyses. The most highly enriched biologic function pathways in this gene set (Supplementary Table S6) included RNA processing (P = 5.1e−13), RNA expression (P = 9.4e−13), apoptosis (P = 1.3e−8), cellular proliferation (P = 5.0e−8), and cell autophagy (P = 1.4e−6). The most significantly enriched canonical pathways from this analysis (Supplementary Table S7) included mTOR signaling (P = 6e−6), p70S6K signaling (P = 5.3e−5), PI3K/AKT signaling (P = 6.0e−5), actin regulation by actin-related protein-Wiskott-Aldrich Syndrome Protein (ARP-WASP) (P = 1.0e−4), integrin-linked kinase (ILK) signaling (P = 1.1e−4), embryonic stem cell pluripotency (P = 6.5e−4), and integrin signaling (P = 6.6e−4). In addition, an upstream analysis of the compounds and transcription factors that affect this set of genes identified a number of compounds and pathways that are associated with prostate cancer progression (Supplementary Table S8), including sirolimus, topotecan, HIF1A, dihydrotestosterone, and TP53. Of the top 12 upstream regulators, 10 included the androgen receptor in the mechanistic network, and all 12 included TP53. Finally, we performed IPA analysis on the final set of 24 biomarker genes, and found that it was significantly associated with malignant neoplasm of the abdomen (P = 7.59e−05) with 15 genes annotated for this disease (BTG2, COL15A1, COL3A1, FDPS, HIST1H1C, IFT57, IGFBP3, ITPR1, MARCH5, MED4, MEMO1, MXI1, PTN, SIRT1, and SYNM), as well as for pelvic cancer (P = 2.98e−06) with a subset of 10 genes having this annotation.

Gleason score analysis

The grading system for prostate cancer is unique in that the final pathologic grade is a Gleason sum obtained by assigning a single Gleason grade to the most prevalent pattern, known as the primary patterns and then adding this to another single Gleason grade assigned to the next most prevalent pattern, known as the secondary pattern, to obtain a sum known as the Gleason score. It has been suggested that primary Gleason 4 pattern and Gleason 3 pattern tumors represent different disease states (23), and several studies have supported the concept that the primary Gleason pattern of patients with a Gleason score of seven is predictive of outcome (4–6). To investigate the differences in gene expression, we performed DEseq differential gene expression analysis comparing 43 samples with Gleason 3 + 4 (primary pattern 3) with 22 samples with Gleason 4 + 3 (primary pattern 4). We identified 304 genes differentially expressed between these two patient groups (Supplementary Table S12). There were several genes differentially expressed relevant to prostate cancer, including upregulated (miR10A, Twist, HOXC6, and AR) and downregulated (ERG, EGF, SOX9, WIF1, WNT5A, and SHH) genes in Gleason 4 + 3 compared with 3 + 4 cases. IPA pathway analysis also determined that the top biologic functions associated with the 304 differentially expressed genes were abnormal bone morphology (P = 7.72e−8), cell differentiation (P = 4.05e−7), genital tumors (P = 1.66e−6), and prostatic bud formation (P = 7.48e−6).

There is a great need for robust biomarkers to predict which tumors are more likely to result in different clinical outcomes to optimize treatment decisions. Through RNAseq analysis of 100 FFPE prostatectomy samples, we have identified a new set of 24 biomarker genes that are highly predictive of BCR in patients with prostate cancer. Biomarker studies in prostate cancer have faced a challenge in which the number of available frozen specimens with very long-term follow-up (> 7 years) is fairly limited. In this study, we have demonstrated that RNAseq analysis of FFPE prostate tissues is feasible and can provide important insight into prostate cancer progression.

Of the 24 biomarker genes, at least 16 (BTG2, CDC27L1, COL15A1, COL3A1, FDPS, IFT57, IGFBP3, ITPR1, LBH, MED4, MEMO1, MXI1, PTN, SACM1L, SIRT1, SRSF3, and SYNM) have been previously associated with some type of cancer, and five (BTG2, ref. 24; IGFBP3, refs. 25–27; MXI1, ref. 28; FDPS, ref. 29; and SIRT1, refs. 30, 31) have been previously associated with prostate cancer. There are several biologic functions that can be used to functionally categorize these biomarker genes. First, there is a set of genes associated with the cytoskeleton, extracellular matrix, angiogenesis, and hypoxia, including SYNM (32), COL15A1, COL3A1, PTN (33), IFT57 (34), and MXI1 (35). Next, there is a set of cell-cycle–related genes, including CDC37L1, BTG2 (24), and SACM1L. In addition, there are several genes associated with apoptosis, PI3K, insulin, metabolism, and steroid hormone signaling, which includes IFT57 (34), ITPR1 (36), IGFBP3 (37), MEMO1 (38), FDPS (29), SRSF3 (39), MARCH5 (40), and MED4 (41–43). Another functional category contains genes associated with ribosomes and translation, including RPL23AP53, SNORA20, and EIF2D. Finally, there is a set of genes relevant to transcription and chromatin modification, which includes HIST1H1C, HIST1H2BG, MED4 (41–43), LBH (44), SIRT1 (30, 31), and MXI1 (35).

FFPE specimens with long-term follow-up clinical data are much more abundant than frozen tissues, and successful sequence analysis of this resource opens up this large resource for further analysis. To date, very few studies have reported RNAseq analysis of FFPE tissues with the notable exception of breast cancer (45). Development of biomarkers from FFPE tissues also increases the likelihood that they can be translated into clinically useful biomarkers that can be used in traditional clinical practice that may not have access to frozen specimens. Nevertheless, challenges remain in translating an RNAseq biomarker test to the clinic without additional validation using alternative platforms. Thus, further validation using TaqMan, NanoString, Digital Droplet PCR, BioMark PCR, or gene-focused RNAseq methods such as PCR-amplicon or capture sequencing in a Clinical Laboratory Improvement Amendments laboratory environment would be important before developing this biomarker panel into a clinical test. This would be important both to simplify the assay and to provide independent validation that the data from RNAseq of FFPE samples are sufficiently robust across technologies. Future experiments are currently planned to determine the optimal approach for developing this biomarker panel.

RNAseq may be too complex and expensive for clinical translation, and thus, development of alternative assays such as TaqMan, digital PCR, or Nanostring type assays may be more realistic. In addition, the fact that approximately 20% of FFPE RNA samples failed to meet our QC criteria is another challenge. It remains to be seen whether alternative assays can accommodate more challenging specimens with more degraded RNA. Furthermore, additional validation of this panel in another large, independent set of patient samples will be important to increase the robustness of this panel. Validation studies on separate groups of patients will speed translation of these biomarkers into a clinical lab test that could be translated to widespread clinical application that would give physicians an idea about what is the best course of treatment for patients with prostate cancer and help avoid unnecessary treatments. Future studies that apply these biomarkers to biopsy or biofluid samples from patients who undergo radiation or active surveillance could determine whether they are also useful for discriminating aggressive from indolent disease. In the long run, this will result in better patient outcomes and reduced healthcare costs and treatment side effects.

The comparison of Gleason 3 + 4 with 4 + 3 cases identified several interesting genes that suggest that increasing miR10A, Twist, HOXC6, and AR expression and decreasing SOX9, WIF1, and WNT5A are associated with increasing tumor grade. Moreover, the fact that the most significant biologic annotation was abnormal bone morphology suggests that higher Gleason primary pattern tumors intrinsically express genes that may facilitate metastasis to the bone. Future studies will be needed to determine the functional role of these genes associated with bone morphology in prostate cancer progression and metastasis.

No potential conflicts of interest were disclosed.

Conception and design: Q. Long, A.O. Osunkoya, T. Gillespie, A.K. Seth, J.A. Petros, C.S. Moreno

Development of methodology: Q. Long, A.O. Osunkoya, S. Sannigrahi, C.S. Moreno

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): A.O. Osunkoya, W. Zhou, T. Gillespie, J.Y. Park, R.K. Nam, A. Stanimirovic, J.A. Petros, C.S. Moreno

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): Q. Long, J. Xu, A.O. Osunkoya, B.A. Johnson, R.K. Nam, C.S. Moreno

Writing, review, and/or revision of the manuscript: Q. Long, J. Xu, A.O. Osunkoya, B.A. Johnson, W. Zhou, T. Gillespie, J.Y. Park, R.K. Nam, L. Sugar, A.K. Seth, J.A. Petros, C.S. Moreno

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): Q. Long, A.O. Osunkoya, S. Sannigrahi, T. Gillespie, L. Sugar, A.K. Seth, C.S. Moreno

Study supervision: A.O. Osunkoya, J.A. Petros, C.S. Moreno

The authors thank Selvi Ramalingam for technical assistance.

This research was supported by Department of Defense (DOD) grant W81XWH-010-1-0090, NIH R01CA106826, NIH U01 CA168449, NIH R03CA173770, NIH R03CA183006, R01CA128813, Canadian Cancer Society Research Institute grant 019038, as well as a National Cancer Institute Cancer Center Support grant P30CA138292 for the Winship Cancer Genomics Shared Resource and the Emory Integrated Genomics Core.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

1.
American Cancer Society [Internet]
. 
Cancer facts and figures 2013: American Cancer Society
;
[updated 2013; cited 2013]. Available from
: http://www.cancer.org/research/cancerfactsfigures/cancerfactsfigures/cancer-facts-figures-2013.
2.
Gleason
DF
. 
Histologic grading of prostate cancer: a perspective
.
Hum Pathol
1992
;
23
:
273
9
.
3.
Gleason
DF
,
Mellinger
GT
. 
Prediction of prognosis for prostatic adenocarcinoma by combined histological grading and clinical staging
.
J Urol
1974
;
111
:
58
64
.
4.
Berg
KD
,
Roder
MA
,
Brasso
K
,
Vainer
B
,
Iversen
P
. 
Primary Gleason pattern in biopsy Gleason score 7 is predictive of adverse histopathological features and biochemical failure following radical prostatectomy
.
Scand J Urol
2014
;
48
:
168
76
.
5.
Ro
YK
,
Lee
S
,
Jeong
CW
,
Hong
SK
,
Byun
SS
,
Lee
SE
. 
Biochemical recurrence in Gleason score 7 prostate cancer in Korean men: significance of the primary Gleason grade
.
Korean J Urol
2012
;
53
:
826
9
.
6.
Herman
CM
,
Kattan
MW
,
Ohori
M
,
Scardino
PT
,
Wheeler
TM
. 
Primary Gleason pattern as a predictor of disease progression in Gleason score 7 prostate cancer: a multivariate analysis of 823 men treated with radical prostatectomy
.
Am J Surg Pathol
2001
;
25
:
657
60
.
7.
Long
Q
,
Johnson
BA
,
Osunkoya
AO
,
Lai
YH
,
Zhou
W
,
Abramovitz
M
, et al
Protein-coding and microRNA biomarkers of recurrence of prostate cancer following radical prostatectomy
.
Am J Pathol
2011
;
179
:
46
54
.
8.
Cuzick
J
,
Swanson
GP
,
Fisher
G
,
Brothman
AR
,
Berney
DM
,
Reid
JE
, et al
Prognostic value of an RNA expression signature derived from cell cycle proliferation genes in patients with prostate cancer: a retrospective study
.
Lancet Oncol
2011
;
12
:
245
55
.
9.
Nam
RK
,
Sugar
L
,
Yang
W
,
Srivastava
S
,
Klotz
LH
,
Yang
LY
, et al
Expression of the TMPRSS2:ERG fusion gene predicts cancer recurrence after surgery for localised prostate cancer
.
Br J Cancer
2007
;
97
:
1690
5
.
10.
Abramovitz
M
,
Ordanic-Kodani
M
,
Wang
Y
,
Li
Z
,
Catzavelos
C
,
Bouzyk
M
, et al
Optimization of RNA extraction from FFPE tissues for expression profiling in the DASL assay
.
Biotechniques
2008
;
44
:
417
23
.
11.
Trapnell
C
,
Pachter
L
,
Salzberg
SL
. 
TopHat: discovering splice junctions with RNA-Seq
.
Bioinformatics
2009
;
25
:
1105
11
.
12.
Trapnell
C
,
Williams
BA
,
Pertea
G
,
Mortazavi
A
,
Kwan
G
,
van Baren
MJ
, et al
Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation
.
Nat Biotechnol
2010
;
28
:
511
5
.
13.
www.ingenuity.com
[Internet]. Ingenuity Knowledgebase [updated 2013; cited 2013]. Available from
: http://www.ingenuity.com/.
14.
Anders
S
,
Huber
W
. 
Differential expression analysis for sequence count data
.
Genome Biol
2010
;
11
:
R106
.
15.
Tibshirani
R
. 
The lasso method for variable selection in the Cox model
.
Stat Med
1997
;
16
:
385
95
.
16.
Goeman
JJ
. 
L1 penalized estimation in the Cox proportional hazards model
.
Biom J
2010
;
52
:
70
84
.
17.
Taylor
BS
,
Schultz
N
,
Hieronymus
H
,
Gopalan
A
,
Xiao
Y
,
Carver
BS
, et al
Integrative genomic profiling of human prostate cancer
.
Cancer Cell
2010
;
18
:
11
22
.
18.
Uno
H
,
Tian
L
,
Cai
T
,
Kohane
IS
,
Wei
LJ
. 
A unified inference procedure for a class of measures to assess improvement in risk prediction systems with survival data
.
Stat Med
2013
;
32
:
2430
42
.
19.
Tomlins
SA
,
Rhodes
DR
,
Perner
S
,
Dhanasekaran
SM
,
Mehra
R
,
Sun
XW
, et al
Recurrent fusion of TMPRSS2 and ETS transcription factor genes in prostate cancer
.
Science
2005
;
310
:
644
8
.
20.
Tomlins
SA
,
Laxman
B
,
Dhanasekaran
SM
,
Helgeson
BE
,
Cao
X
,
Morris
DS
, et al
Distinct classes of chromosomal rearrangements create oncogenic ETS gene fusions in prostate cancer
.
Nature
2007
;
448
:
595
9
.
21.
Tomlins
SA
,
Rhodes
DR
,
Yu
J
,
Varambally
S
,
Mehra
R
,
Perner
S
, et al
The role of SPINK1 in ETS rearrangement-negative prostate cancers
.
Cancer Cell
2008
;
13
:
519
28
.
22.
Helgeson
BE
,
Tomlins
SA
,
Shah
N
,
Laxman
B
,
Cao
Q
,
Prensner
JR
, et al
Characterization of TMPRSS2:ETV5 and SLC45A3:ETV5 gene fusions in prostate cancer
.
Cancer Res
2008
;
68
:
73
80
.
23.
Lavery
HJ
,
Droller
MJ
. 
Do Gleason patterns 3 and 4 prostate cancer represent separate disease states?
J Urol
2012
;
188
:
1667
75
.
24.
Ficazzola
MA
,
Fraiman
M
,
Gitlin
J
,
Woo
K
,
Melamed
J
,
Rubin
MA
, et al
Antiproliferative B cell translocation gene 2 protein is down-regulated post-transcriptionally as an early event in prostate carcinogenesis
.
Carcinogenesis
2001
;
22
:
1271
9
.
25.
Sarma
AV
,
Dunn
RL
,
Lange
LA
,
Ray
A
,
Wang
Y
,
Lange
EM
, et al
Genetic polymorphisms in CYP17, CYP3A4, CYP19A1, SRD5A2, IGF-1, and IGFBP-3 and prostate cancer risk in African-American men: The Flint Men's Health Study
.
Prostate
2008
;
68
:
296
305
.
26.
Ramachandran
S
,
Liu
P
,
Young
AN
,
Yin-Goen
Q
,
Lim
SD
,
Laycock
N
, et al
Loss of HOXC6 expression induces apoptosis in prostate cancer cells
.
Oncogene
2005
;
24
:
188
98
.
27.
McCabe
CD
,
Spyropoulos
DD
,
Martin
D
,
Moreno
CS
. 
Genome-wide analysis of the homeobox C6 transcriptional network in prostate cancer
.
Cancer Res
2008
;
68
:
1988
96
.
28.
Prochownik
EV
,
Eagle Grove
L
,
Deubler
D
,
Zhu
XL
,
Stephenson
RA
,
Rohr
LR
, et al
Commonly occurring loss and mutation of the MXI1 gene in prostate cancer
.
Genes Chromosomes Cancer
1998
;
22
:
295
304
.
29.
Todenhofer
T
,
Hennenlotter
J
,
Kuhs
U
,
Gerber
V
,
Gakis
G
,
Vogel
U
, et al
Altered expression of farnesyl pyrophosphate synthase in prostate cancer: evidence for a role of the mevalonate pathway in disease progression?
World J Urol
2013
;
31
:
345
50
.
30.
Byles
V
,
Zhu
L
,
Lovaas
JD
,
Chmilewski
LK
,
Wang
J
,
Faller
DV
, et al
SIRT1 induces EMT by cooperating with EMT transcription factors and enhances prostate cancer cell migration and metastasis
.
Oncogene
2012
;
31
:
4619
29
.
31.
Powell
MJ
,
Casimiro
MC
,
Cordon-Cardo
C
,
He
X
,
Yeow
WS
,
Wang
C
, et al
Disruption of a Sirt1-dependent autophagy checkpoint in the prostate results in prostatic intraepithelial neoplasia lesion formation
.
Cancer Res
2011
;
71
:
964
75
.
32.
Noetzel
E
,
Rose
M
,
Sevinc
E
,
Hilgers
RD
,
Hartmann
A
,
Naami
A
, et al
Intermediate filament dynamics and breast cancer: aberrant promoter methylation of the Synemin gene is associated with early tumor relapse
.
Oncogene
2010
;
29
:
4814
25
.
33.
Papadimitriou
E
,
Mikelis
C
,
Lampropoulou
E
,
Koutsioumpa
M
,
Theochari
K
,
Tsirmoula
S
, et al
Roles of pleiotrophin in tumor growth and angiogenesis
.
Eur Cytokine Netw
2009
;
20
:
180
90
.
34.
Gervais
FG
,
Singaraja
R
,
Xanthoudakis
S
,
Gutekunst
CA
,
Leavitt
BR
,
Metzler
M
, et al
Recruitment and activation of caspase-8 by the Huntingtin-interacting protein Hip-1 and a novel partner Hippi
.
Nat Cell Biol
2002
;
4
:
95
105
.
35.
Corn
PG
,
Ricci
MS
,
Scata
KA
,
Arsham
AM
,
Simon
MC
,
Dicker
DT
, et al
Mxi1 is induced by hypoxia in a HIF-1-dependent manner and protects cells from c-Myc-induced apoptosis
.
Cancer Biol Ther
2005
;
4
:
1285
94
.
36.
Szlufcik
K
,
Bultynck
G
,
Callewaert
G
,
Missiaen
L
,
Parys
JB
,
De Smedt
H
. 
The suppressor domain of inositol 1,4,5-trisphosphate receptor plays an essential role in the protection against apoptosis
.
Cell Calcium
2006
;
39
:
325
36
.
37.
Cortes-Sempere
M
,
de Miguel
MP
,
Pernia
O
,
Rodriguez
C
,
de Castro Carpeno
J
,
Nistal
M
, et al
IGFBP-3 methylation-derived deficiency mediates the resistance to cisplatin through the activation of the IGFIR/Akt pathway in non-small cell lung cancer
.
Oncogene
2013
;
32
:
1274
83
.
38.
Sorokin
AV
,
Chen
J
. 
MEMO1, a new IRS1-interacting protein, induces epithelial-mesenchymal transition in mammary epithelial cells
.
Oncogene
2013
;
32
:
3130
8
.
39.
Kim
J
,
Park
RY
,
Chen
JK
,
Kim
J
,
Jeong
S
,
Ohn
T
. 
Splicing factor SRSF3 represses the translation of programmed cell death 4 mRNA by associating with the 5′-UTR region
.
Cell Death Differ
2014
;
21
:
481
90
.
40.
Fang
L
,
Hemion
C
,
Goldblum
D
,
Meyer
P
,
Orgul
S
,
Frank
S
, et al
Inactivation of MARCH5 prevents mitochondrial fragmentation and interferes with cell death in a neuronal cell model
.
PLoS One
2012
;
7
:
e52637
.
41.
Lando
M
,
Holden
M
,
Bergersen
LC
,
Svendsrud
DH
,
Stokke
T
,
Sundfor
K
, et al
Gene dosage, expression, and ontology analysis identifies driver genes in the carcinogenesis and chemoradioresistance of cervical cancer
.
PLoS Genet
2009
;
5
:
e1000719
.
42.
Rachez
C
,
Lemon
BD
,
Suldan
Z
,
Bromleigh
V
,
Gamble
M
,
Naar
AM
, et al
Ligand-dependent transcription activation by nuclear receptors requires the DRIP complex
.
Nature
1999
;
398
:
824
8
.
43.
Wang
Q
,
Sharma
D
,
Ren
Y
,
Fondell
JD
. 
A coregulatory role for the TRAP-mediator complex in androgen receptor-mediated gene expression
.
J Biol Chem
2002
;
277
:
42852
8
.
44.
Rieger
ME
,
Sims
AH
,
Coats
ER
,
Clarke
RB
,
Briegel
KJ
. 
The embryonic transcription cofactor LBH is a direct target of the Wnt signaling pathway in epithelial development and in aggressive basal subtype breast cancers
.
Mol Cell Biol
2010
;
30
:
4267
79
.
45.
Prokopec
SD
,
Watson
JD
,
Waggott
DM
,
Smith
AB
,
Wu
AH
,
Okey
AB
, et al
Systematic evaluation of medium-throughput mRNA abundance platforms
.
RNA
2013
;
19
:
51
62
.

Supplementary data