Abstract
Purpose: To carry out an integrative profile of human pancreatic ductal adenocarcinoma (PDAC) to identify prognosis-significant genes and their related pathways.
Experimental Design: A concordant survival-based whole genome in silico array analysis of DNA copy number, and mRNA and miRNA expression in 25 early-stage PDAC was carried out. A novel composite score simultaneously integrated gene expression with regulatory mechanisms to identify the signature genes with the most levels of prognosis-significant evidence. The predominant signaling pathways were determined via a pathway-based approach. Independent patient cohorts (n = 148 and 42) were then used as in vitro validation of the array findings.
Results: The composite score identified 171 genes in which expressions were able to define two prognosis subgroups (P = 3.8e-5). Eighty-eight percent (151 of 171) of the genes were regulated by prognosis-significant miRNAs. The phosphoinositide 3-kinase/AKT pathway and SRC signaling were densely populated by prognosis-significant genes and driven by genomic amplification of SRC and miRNA regulation of p85α and CBL. On tissue microarray validation (n = 148), p85α protein expression was associated with improved survival for all patients (P = 0.02), and activated P-SRC (Y418) was associated shorter survival for patients with low-grade histology tumors (P = 0.04). Interacting P-SRC and p85α revealed that they define two distinct PDAC patient subgroups (P = 0.0066). Furthering the importance of these pathways, CBL protein expression was associated with improved survival (P = 0.03) on a separate cohort (n = 42).
Conclusions: These pathways and related genes may represent putative clinical biomarkers and possible targets of individualized therapy in the distinct patient subgroups they define. Clin Cancer Res; 18(5); 1352–63. ©2012 AACR.
Despite evidence of the heterogeneous genomic and molecular changes associated with pancreatic ductal adenocarcinoma (PDAC) tumorigenesis, little is known about biologic subsets of tumors to guide patient stratification and individualized therapies. Recent phase III trials with molecular agents have neither resulted in robust survival benefits nor attempted to stratify patients prior to randomization. We use a novel integrative survival–based genomic and molecular array analysis of human PDACs to derive a composite score that ranks genes based on expression and regulatory mechanisms. We find that miRNA regulation plays a critical role in the malignant phenotype of PDAC. This approach shows that dysregulated phosphoinositide 3-kinase/AKT or SRC signaling is significantly associated with distinct patient subgroups with more aggressive disease. These pathway-specific genes and miRNAs represent potentially useful clinical biomarkers and targets of individualized therapy for well-defined patient subgroups.
Introduction
Pancreatic ductal adenocarcinoma (PDAC) is the fourth leading cause of cancer-related deaths in the United States (1) and has an extremely poor overall 5-year survival rate of only 4%. Most patients present with advanced stage disease and have a median survival of less than 1 year (2). Cytotoxic chemotherapy is marginally effective with standard gemcitabine- or 5-fluorouracil–based regimens increasing PDAC median survival by less than 2 months in advanced disease (3–5). Published phase III clinical trials of targeted molecular agents in unselected PDAC populations have also not shown robust survival benefits (6–10). Ultimately, our evolving understanding of significant genomic diversity in PDAC must be used to better inform targeted drug design and delivery.
Recent in-depth exome sequencing showed individual PDAC tumors average more than 60 distinct alterations, the majority of which occur at low frequencies across all tumors. Only a few high prevalence genomic changes were detected, including expected mutations in KRAS and loss or inactivation of known tumor suppressor genes (e.g., TP53 and SMAD4). Despite this genomic heterogeneity, all tumors had genetic alterations that were linked to 12 core signaling pathways (11). Follow-up work comparing patient-matched primary PDAC tumors and subsequent metastases revealed acquisition of further mutations that varied by metastatic site. Strikingly, founder mutations of each metastatic subclone could be traced back to sequenced geographic subregions of the primary tumor (12, 13), providing new insights into the genetic events and timing of PDAC initiation and malignant phenotype. Notably, SMAD4 deletion has been the only genetic alteration from this work that has been linked to patient survival (14).
Others have used gene expression microarray analyses to define molecular signatures associated with PDAC disease progression. Stratford and colleagues (15) identified a 6 gene signature in primary tumors that was associated with metastatic disease and predicted shorter survival in an independent set of 67 patients. Collisson and colleagues (16) analyzed primary PDAC from cell lines and a combination of clinical data sets to classify 3 distinct PDAC molecular subtypes that were able to predict clinical survival, as well as response to therapy in experimental models. Although such molecular profiling has provided valuable information, the remarkable genomic diversity of PDAC and the small size of most patient cohorts has clearly hindered the discovery of additionally biologically important molecular changes.
As a means to effectively study diverse genomic alterations in a small patient data set, we hypothesized that the identification and refinement of prognosis-related genes in PDAC should be improved by increasing the depth of analysis for each tumor using multiple array platforms. The potential for this type of multidimensional analysis was shown in a recent prostate cancer study in which several pathways of known prognostic significance were validated and new ones were additionally implicated (17). For our own survival-based analysis of PDAC, individual gene expression changes associated with survival were matched to potential genomic or epigenetic modes of regulation by integrating microarray results of mRNA expression with DNA copy number variation and miRNA levels. This approach validated pathways implicated in pancreatic tumorigenesis and uncovered previously unrecognized molecular events associated with poor prognosis. The expressions of many identified genes were found to have associated miRNA alterations linked to survival. These genes and their regulatory mechanisms represent promising candidates for future studies addressing their function and evaluating their efficacy as predictive biomarkers and/or targets for molecular-based therapies.
Materials and Methods
Patients and samples
All work was carried out with University of California, Los Angeles (UCLA) Institutional Review Board approval. Three independent, nonoverlapping patient cohorts were used in this study. The initial test cohort of 42 PDAC tumors and 7 nonmalignant pancreas samples snap frozen at the time of surgery were used for microarrays. Of these, only samples with tumor cell content more than 30% were chosen for final multiplatform analysis (n = 25) as determined on representative hematoxylin and eosin (H&E) sections by a practicing gastrointestinal pathologist (DWD). The second patient cohort (n = 42) was tumors isolated from formalin-fixed paraffin-embedded (FFPE) tissue blocks and used as a validation cohort for quantitative PCR (qPCR). The third patient data set (n = 148) was tumors on a TMA used as an immunohistochemistry (IHC) validation cohort. All clinicopathologic and survival information for each patient cohort were extracted from a prospectively maintained UCLA surgical database of pancreatic patients. Disease recurrence was assessed based on biopsy, radiographic evidence, or death. The electronic medical record was used to determine associated clinical and pathologic features, as well as disease-free and disease-specific survival (DSS). Search of the Social Security Death Index was used to determine overall survival. Survival analysis of the tissue microarray (TMA) cohort was limited to overall survival. Disease-free, disease-specific and overall survival times were examined for the microarray and qPCR validation cohorts. Survival intervals were calculated from date of surgery to date of confirmed death or last patient contact.
Gene expression analysis
The gene expression was investigated by Affymetrix HGU133 Plus 2 Array on which multiple probe sets might be used to measure expression of a single gene. Therefore, we analyzed the data based on probe set IDs and used the highest absolute value/score among multiple probe sets for the gene-based interpretation. Details of the array procedure and normalization and filtering are detailed in the Supplementary Information. In brief, GCRMA was used for normalization, and probe sets having presence calls in less than 30% of samples were eliminated before further analyses. Cox scores (18) were used to determine the correlation between individual probe set–based expression and disease-free survival (DFS) time. Probe sets were then sorted based on the absolute values of their Cox scores. Prediction analysis for microarrays (19), implemented by Bioconductor pamr package, was then used to determine (i) if expression of the top-ranked probe sets could be used to predict the sample outcome and (ii) the minimal probe sets to define 2 prognosis groups with the highest statistical significance.
miRNA expression analysis
The Exiqon miRNA arrays (miRCURY LNA microRNA Array v.11.0 -hsa, mmu & rno) were used for measuring genome miRNA expression. Array intensities were adjusted and normalized by variance stabilizing transformation implemented in the Bioconductor vsn packages. The expression profile of each miRNA was represented by the average expression of its multiple spots on the array. Statistical analyses for the association between miRNAs and survival was conducted as detailed for gene expression analysis.
DNA copy number analysis
Affymetrix SNP 6.0 arrays were used to detect copy number aberrations (CNA) in tumor samples. CEL files produced by Affymetrix GeneChip Command Console (AGCC) software were imported into Affymetrix GTC 3.0.1 and analyzed by the Copy Number Analysis workflow with HapMap270 as the reference model. Regions with CNAs were then annotated with gene symbols base on the annotation file from the UCSC genome browser (build hg18). Only those genes in which loci had CNA present in at least 20% of the tumors were included for survival analysis in which Cox proportional hazards model was used to determine whether the group with CNA at a given locus was at higher risk than the group without CNA.
Integrating multiple dimension data to identify signature genes with multiple levels of evidence
In this meta-analysis, a composite score was generated to quantitatively measure gene prognosis significance based on the multiple array platforms and is described in detail in the Supplementary Information. In brief, the composite score was developed based on the assumptions that besides being correlated to survival, expression of the gene of interest is (i) correlated to its copy number, (ii) anticorrelated to its regulating miRNA, and (iii) such changes in copy number and miRNA expression are also associated with prognosis. Mathematically, the operators of the survival-based composite score for each gene included the following: (i) rank of gene expression Cox score; (ii) binary CNA score (0 or 1) requiring both a HR P < 0.2 and concordance with gene expression change; (iii) rank of miRNA Cox score with required anticoncordance between miRNA and corresponding gene expression. The miRNA operator was further strengthened if a prognosis-significant CNA overlapped with the miRNA coding region.
Pathway and gene ontology analysis
Survival signature genes were annotated by databases in the public domain: the KEGG and Molecular Signatures Databases (20). Significant enrichment of pathway/gene set was determined by Fisher exact test.
Validation with TMA
The PDAC TMA has been previously detailed (21) and represents a totally separate, nonoverlapping cohort of patients. Immunohistochemistry was visualized with the VECTASTAIN ABC Elite Kit (Vector Laboratories) following heat-induced antigen retrieval (0.01mol/L citric acid buffer, pH 6.0) and overnight incubation with primary antibodies (1:100 dilution), including P-SRC (Y418; Abcam, ab47411) or p85α (Epitomics, catalog: 1675-1). Three separate 1.0 mm cores for each tumor in the TMA were independently scored by 2 blinded observers using semiquantitative histoscores (range 0–300). Histoscores were the product of staining intensity (0–3) and percentage of tumor cells staining at that intensity (0–100). If any core's histoscore differed by more than 30 points between observers, a revised score was assigned by consensus evaluation. Median histoscore of both observers was used for analysis.
Real-time qPCR of pancreatic cancer resection samples
For qPCR validation, H&E slides from a separate nonoverlapping cohort of 42 additional PDAC patients resected at UCLA between 2002 and 2009 were reviewed by a practicing gastrointestinal pathologist (DWD) to target the extraction of three 3 mm cores from areas of viable tumor in the corresponding FFPE blocks. RNA was isolated with the RecoverAll Total Nucleic Acid Isolation Kit for FFPE (Ambion). Total RNA was reverse-transcribed by random hexamer primers (High Capacity cDNA Reverse Transcription Kit; Ambion). SYBR green real-time qPCR assays were conducted on generated cDNA using a Roche LightCycler 480 Real-Time PCR System. Reaction parameters and primers sequences are available upon request and were optimized for FFPE-derived cDNA based on the use of short (<90 bp), intron-spanning amplicons. Relative gene expression was normalized to ACTB as housekeeping gene.
Statistical analysis
Survival estimate for each subgroup was generated by the Kaplan–Meier method. Log-rank test was used to compare Kaplan–Meier curves. Multivariate Cox proportional hazards models were used to test statistical independence and significance of multiple predictors.
Results
Survival-based integrated genomic and gene expression analyses
The overall strategy of our survival-based study is outlined in Supplementary Fig. S1. Integrative multiplatform array analysis was used to simultaneously examine gene expression and regulatory mechanisms from an initial cohort of 25 patients to identify and refine genes and pathways with biologic significance. Array analysis was conducted on snap-frozen primary patient tissue samples to directly link gene expression to genomic and epigenomic changes in the in situ context of primary tumor and to extract a sufficient source of high quality nucleic acid for simultaneous mRNA, miRNA, and SNP microarray analysis. As our goal was to identify genomic alterations and expression changes in neoplastic ductal epithelial cells, we limited our analysis to primary tumor samples in which estimated tumor cell content exceeded 30% (median 60%, ranging from 35%–90%). This cutoff was chosen after initial unsupervised clustering of microarrays comparing mRNA expression from normal pancreas, chronic pancreatitis, and pancreas tumor samples showed tumors with low tumor cell content (<30%) more frequently clustered with normal and chronic pancreatitis controls (data not shown), consistent with previous work showing a large stromal content can confound array readouts in PDAC (22). Specific molecular alterations identified through this integrated genomic analysis were then independently validated in separate, nonoverlapping patient cohorts by immunohistochemistry (n = 148) or qPCR (n = 42), in which clinicopathologic characteristics and survival outcomes largely overlap (Supplementary Table S1 and Supplementary Fig. S2).
Table 1 provides the clinical and histopathologic parameters for the 25 patients. All patients had early-stage PDAC and received adjuvant chemotherapy after surgery. At the time of analysis, 16 patients had recurrent disease with a median DFS of 13.3 months, and 10 of them had died of disease with a median DSS of 20.6 months. Given the number of deaths and relatively short follow-up time (median 12.4 months for survivors) in our cohort, DFS was chosen as the outcome measure for our analysis. DFS is considered an accurate metric of survival outcome, as most patients with recurrent PDAC will succumb to disease. The clinicopathologic characteristics and survival outcomes of our cohort was similar to other published large cohorts of early-stage PDAC (23).
Factor . | Subcategory . | n (%) . |
---|---|---|
Total samples | 25 | |
Age, y | Median (range) | 67 (49–85) |
<65 | 9 (36%) | |
≥65 | 16 (64%) | |
Gender | Male | 12 (48%) |
Female | 13 (52%) | |
Tumor location | Head | 24 (96%) |
Tail | 1 (4%) | |
Operative technique | Whipple | 22 (92%) |
Distal pancreatectomy | 1 (4%) | |
Total | 1 (4%) | |
Tumor diameter, cm | Median (range) | 3.4 (1.8–5.9) |
<2.5 | 8 (32%) | |
≥2.5 | 17 (68%) | |
T stage | 2 | 2 (8%) |
3 | 23 (92%) | |
Resection margins | Positive | 4 (16%) |
Negative | 21 (84%) | |
Tumor differentiation | Well | 1 (4%) |
Moderate | 13 (52%) | |
Poor | 11 (44%) | |
Lymph nodes | Negative | 8 (32%) |
Positive | 17 (68%) | |
Lymphovascular invasion | Absent | 5 (20%) |
Present | 18 (72%) | |
Perineural invasion | Absent | 2 (8%) |
Present | 21 (84%) | |
AJCC stage | 1 | 2 (8%) |
2 | 23 (92%) | |
3 | 0 (0%) | |
DFS, months | Median | 13.2 |
DSS, months | Median | 20.6 |
Follow-up survivors, mo | Median (range) | 12.4 (3.8–32.5) |
Factor . | Subcategory . | n (%) . |
---|---|---|
Total samples | 25 | |
Age, y | Median (range) | 67 (49–85) |
<65 | 9 (36%) | |
≥65 | 16 (64%) | |
Gender | Male | 12 (48%) |
Female | 13 (52%) | |
Tumor location | Head | 24 (96%) |
Tail | 1 (4%) | |
Operative technique | Whipple | 22 (92%) |
Distal pancreatectomy | 1 (4%) | |
Total | 1 (4%) | |
Tumor diameter, cm | Median (range) | 3.4 (1.8–5.9) |
<2.5 | 8 (32%) | |
≥2.5 | 17 (68%) | |
T stage | 2 | 2 (8%) |
3 | 23 (92%) | |
Resection margins | Positive | 4 (16%) |
Negative | 21 (84%) | |
Tumor differentiation | Well | 1 (4%) |
Moderate | 13 (52%) | |
Poor | 11 (44%) | |
Lymph nodes | Negative | 8 (32%) |
Positive | 17 (68%) | |
Lymphovascular invasion | Absent | 5 (20%) |
Present | 18 (72%) | |
Perineural invasion | Absent | 2 (8%) |
Present | 21 (84%) | |
AJCC stage | 1 | 2 (8%) |
2 | 23 (92%) | |
3 | 0 (0%) | |
DFS, months | Median | 13.2 |
DSS, months | Median | 20.6 |
Follow-up survivors, mo | Median (range) | 12.4 (3.8–32.5) |
Copy number variations and mRNA and miRNA expression alterations as predictors of PDAC survival
We first individually analyzed mRNA expression, miRNA expression, and DNA copy number to determine DFS-associated changes in each microarray platform. For mRNA analysis, we used a nongrouping-based approach addressing survival as a continuous variable. The advantage of this approach over a grouping-based approach is that it does not set an arbitrary survival time threshold for dichotomization and thus better accounts for tumors with a continuum of clinical behavior (18). A semisupervised clustering methodology was used whereby a Cox score was first generated for each gene based on a direct correlation between its expression and DFS. Signature genes were then prioritized and selected based on the absolute value of their Cox scores. Using this method, the most accurate gene set is selected from a series of prediction analyses whereby different cutoff points are iterated to establish an optimal panel of genes accurately separating different prognosis subgroups. The top 500 scored probe sets (Supplementary Table S2) from this approach included 186 upregulated and 314 downregulated transcripts able to segregate 2 highly significant and distinct prognosis groups (median DFS = 7.7 vs. 25.3 months; log-rank test P = 0.000038; Fig. 1A). The nongrouping approach was also used for separate miRNA microarray analysis, which yielded a panel of 31 miRNAs (1.8% of the total) able to robustly segregate 2 prognosis groups (median DFS = 9 vs. >36 months; P = 0.00047; Fig. 1B and Supplementary Table S3).
For survival-based analysis of DNA CNAs, CNAs for each tumor were determined by comparing Affymetrix SNP microarray to the human HapMap reference model using Affymetric GTC software. CNAs were then mapped to specific gene loci using UCSC genome build hg18. Because of the size of our cohort, we focused on CNAs occurring at higher frequency (≥20% of the samples) and relaxed the cutoff to include Cox values approaching significance (P < 0.2). High frequency CNAs were clustered on chromosomes 1, 7, 8, and 20 for amplifications and chromosomes 6, 9, 17, and 18 for deletions (Fig. 1C, top). However, among them, the prognosis-significant CNAs (Fig. 1C, bottom) were located on specific loci, including amplifications on chromosomes 7, 11, 19, 20, and 22, and deletions on chromosomes 6 and 9. We identified CNAs associated with a total of 68 genes that optimally segregated patients into 2 statistically significant prognosis groups (median DFS = 9.7 vs. 25.0 months, log-rank test P = 0.0063, Fig. 1D and Supplementary Table S4).
Integrated molecular analyses further refines genes of prognostic importance in PDAC
Although our initial analysis showed a large number and variety of molecular alterations correlated with DFS in PDAC, it was less clear which alterations mechanistically link to the malignant phenotype or represent clinically useful biomarkers. Although analysis of larger patient data sets can help address these questions, the cost and availability of large PDAC patient data sets are limiting. As an alternative, we hypothesized that relevant associations could be strengthened and refined in our cohort through an integrated assessment of gene expression in conjunction with concordant changes in genetic and epigenetic regulation.
We merged results from all 3 microarray platforms to generate a single composite score that integrates CNA and levels of mRNA and miRNA expression (see Supplementary Fig. S1 and Supplementary Information). To associate genes with miRNAs, a list of 2.1 million potential miRNA–mRNA seed-match pairs were first generated from data in the public domain, including experimental evidence in TarBase (24) and predictive sequence analysis using TargetScanS (25) and miRbase (26). Next, each component in the composite score was given a weight based on its ability to independently identify 2 prognosis groups. Individually, mRNA expression was the most robust predictor of DFS in our analysis, followed by miRNA expression and CNA changes, as determined by the log-rank P values generated from each platform (Fig. 1). The log10 transformations of these P values were used to weight each platform's contribution to survival in the final composite score calculation for each gene, whereby a gene's mRNA expression was weighted most heavily, followed by its linked miRNA alterations and finally CNA (see Supplementary Methods). The composite score also took into account whether there was concordance between CNA gene expression (e.g., amplification associated with higher gene expression) or anticoncordance between miRNA gene expression (e.g., as miRNAs predominantly will repress target gene expression). If these conditions were not met, the component was not included. Based on these criteria, our integrated composite score yielded a refined list of 171 signature genes (represented by 200 probe sets) that accurately segregated patients into 2 statistically significant prognostic groups (median DFS = 8.6 vs. 20.6, P = 0.001, Fig. 2A and B and Supplementary Table S5). From this list, 134 genes had corresponding changes in miRNA expression, 20 genes had corresponding changes in CNA, and 17 genes shared corresponding changes in both miRNA expression and CNA (Fig. 2C). Importantly, these results highlight a strong link between miRNAs regulated mechanism and survival outcome in PDAC.
Pathway-based analysis of signature genes links PDAC survival to PI3K and SRC signaling pathways
We next carried out survival-based pathway and gene ontology analyses of our integrated composite gene signature. Pathways most highly populated with survival correlated genes included ERBB signaling, apoptosis, purine metabolism, focal adhesion and insulin signaling (Supplementary Table S6). Therefore, our data offers survival-based correlations to support an existing literature that links these pathways to pancreatic tumorigenesis (27–31). Our data also offer potential insights into novel mechanisms regulating these pathways and their relationship to clinical disease progression. Several survival-correlated genes in our integrated composite signature uniformly link poor prognosis to changes that will result in upregulation of PI3K/AKT/mTOR and SRC signaling (Fig. 3). For instance, the expression of EGFR, a potential activator of both AKT and SRC signaling, was associated with worse prognosis in our analysis (Cox score +2.77), whereas that of CBL, a ubiquitin ligase able to negatively impact EGFR (32) or SRC (33) expression, correlated with improved survival in our analysis (Cox score −2.4). PIK3R1, which encodes the class Ia PI3K regulatory subunit p85α that antagonizes PI3K/AKT signaling, also correlated with improved survival (Cox score −2.5). Prognosis-linked changes in gene expression were more frequently correlated with prognosis-significant miRNA changes, implicating miRNAs as a critical factor in the malignant phenotype of PDAC (see Supplementary Table S3). For example, our data suggested prognosis-linked expression PIK3R1/p85α mRNA could be mediated by miRNA 519d, which is (i) broadly conserved among vertebrates, (ii) predicted to bind to the 3′ untranslated region of PI3KR1, (iii) inversely correlated with PIK3R1 expression in our analysis, and (iv) an independent predictor shorter DFS (miRNA expression Cox score +2.42).
Highlighting the value of the integrated composite score over gene expression alone, SRC was not identified on independent mRNA array analysis, as it showed only a weak nonsignificant trend toward worse prognosis (Cox expression score +1.22). However, SRC was subsequently captured as a high ranking gene (#21) in our integrated composite signature when its prognosis-significant (Cox HR = 4.1, P = 0.016) genomic amplification was considered with mRNA expression. Prognosis-related dysregulation of SRC was also inferred from our subsequent pathway analysis. SRC activation can occur through integrin-FAK–dependent (34, 35) or FAK-independent (36) mechanisms facilitated by PTPRA, which was associated with worse prognosis in our integrated composite signature (Cox score +2.24). ARFGAP1, which also correlated with worse prognosis in our analysis (Cox Score +1.62), has been shown to potentate SRC downstream signaling via its regulation of the actin cytoskeleton and, leading to enhanced cell motility (37). Interestingly, SRC, PTPRA, and ARFGAP1 are all located on chromosome 20. Although SRC is regulated by amplification, both PTPRA and ARFGAP1 can be regulated by amplification and by miR541, which was inversely correlated with DFS (Cox score −3.12), again revealing the strength of our integrated approach.
Validation of molecular signatures
To validate the predictive value of genes in our integrated composite signature, we chose to focus on SRC and the PI3K/AKT signaling given that (i) they were populated with centrally located prognosis-significant genes, (ii) they are known to be dysregulated in human PDAC (35, 38), and (iii) clinically available small-molecule inhibitors targeting each have shown promise in preclinical models and phase I and II clinical trials (39, 40). For validation, we carried out IHC to determine levels of p85α and activated SRC phosphorylated at tyrosine position (Y419) using a TMA consisting of a separate large cohort (n = 148) of treatment-naïve, resected stages I and II PDACs. This TMA is previously detailed (21) and is similar to other large patient cohorts of early-stage PDAC in which pathologic stage, lymph node status, and histologic tumor grade were each significantly correlated with survival. The overall survival of patients in this cohort is similar to those in the in silico analysis (median 24.2 vs. 19.8 months, respectively), given the difference in follow-up times (Table 1 and Supplementary Table S1). IHC staining for each tumor was determined semiquantitatively by histoscore, with a broad staining distribution seen for each antibody across the TMA (Fig. 4B and C). For each IHC marker, patients were dichotomized on a histoscore cutoff value of 150, chosen to identify groups sufficiently populated for statistical power analysis and that minimized the chance of assigning patients with small differences to opposite groups. Neither P-SRC nor p85α expression was associated with various clinicopathologic factors, with the exception that low (includes well and moderately differentiated) histologic grade tumors have high p85α protein expression (P = 0.00062, χ2 test).
By Kaplan–Meier survival analysis, high p85α protein expression significantly correlated with better survival in the TMA (median survival 28.7 vs. 19.4 months, log-rank test P = 0.02, Fig. 5A), whereas high P-SRC expression did not significantly correlate with prognosis (median survival 21.1 vs. 25.6 months, P = 0.45). We next examined whether P-SRC expression significantly correlated with survival in subgroups of patients stratified first on the strong independent prognostic factors of node status (pN0 or pN1) or tumor grade (low or high). High P-SRC staining significantly associated with worse survival in the subgroup of patients with low-grade tumor histology (median survival 23.2 vs. 42.5, P = 0.04, Fig. 5B), but not in the subgroup with high-grade tumor histology or either subgroup based on node status (data not shown).
We next addressed a potential interaction between p85α and P-SRC status by examining various combinations of dichotomized groups. Tumors with combined low P-SRC and high p85α had significantly better survival relative to any of the other paired combinations when considered individually (Fig. 5C) or in aggregate (median survival = 36.6 vs. 20.4 months, P = 0.0066, Fig. 5D). Multivariate Cox proportional hazards analyses, controlling for LN and grade, showed that high p85α was a near significant independent predictor of improved overall survival (HR 0.690, P = 0.068), while the combination of low P-SRC and high p85α was a significant independent predictor of improved survival in the TMA (HR = 0.53, P = 0.02; Table 2). Of note, further IHC staining of whole tissue sections of recently resected PDACs found that either p85α or PTEN loss appeared to correlate with enhanced PI3K/AKT signaling, as detected by P-AKT and P-S6 (Supplementary Fig. S3). Although in need of prospective validation and further mechanistic evaluation, these results suggest tumors with combined lower signaling activity for both SRC (as detected by the surrogate of reduced P-SRC) and PI3K-AKT (as detected by the surrogate of increased p85α) may define a subset of patients with more favorable clinical outcome. These or other surrogate markers of pathway activation may be especially useful in discriminating patients with more or less aggressive clinical disease and assessing the use of drugs targeting these pathways.
Parameters . | p85α based . | p85α and P-SRC based . | ||
---|---|---|---|---|
. | HRa . | Pb . | HRa . | Pb . |
Positive pN | 1.75 | 0.005 | 1.6 | 0.02 |
High grade | 1.73 | 0.007 | 1.8 | 0.003 |
High p85α | 0.69 | 0.068 | — | — |
High p85α & low P-SRC | — | — | 0.53 | 0.02 |
Parameters . | p85α based . | p85α and P-SRC based . | ||
---|---|---|---|---|
. | HRa . | Pb . | HRa . | Pb . |
Positive pN | 1.75 | 0.005 | 1.6 | 0.02 |
High grade | 1.73 | 0.007 | 1.8 | 0.003 |
High p85α | 0.69 | 0.068 | — | — |
High p85α & low P-SRC | — | — | 0.53 | 0.02 |
aHR less than or more than 1 indicates decreased or increased risk of death for the listed variable, respectively.
bP from the Wald statistic in which null hypothesis is the corresponding coefficient β = 0 (i.e., HR = 1) in Cox proportional hazard model.
We finally sought to validate the finding that CBL was linked to prognosis in our integrated composite signature (Fig. 3). In the absence of a reliable IHC assay, we pursued survival-based quantitative PCR analysis of CBL normalized to ACTB in a separate cohort of 42 PDAC samples; these patients had a similar DSS with the array and TMA cohorts (25.4 vs. 19.8 vs. 24.2 months, respectively), given the differences in follow-up times (Table 1 and Supplementary Table 1). Patients were again dichotomized into groups with low versus high CBL expression, the latter of which significantly correlated with better survival (median survival 44.3 vs. 20.6, P = 0.03, Fig. 5E).
Discussion
The large number and wide range of genetic alterations that characterize pancreatic cancer present both significant challenges and opportunities for improving our understanding and treatment of this highly aggressive and lethal malignancy. Seminal large-scale genomic sequencing studies of PDAC offer a tantalizing spatial and temporal picture of the genomic alterations occurring in both primary tumors and metastatic lesions, but must also now be examined in greater detail to establish their association with the malignant phenotype of PDAC (11–13). We have adopted an integrative approach to identify and prioritize genes of potential importance in PDAC. Our survival-based approach involved multidimensional analysis of gene expression and genomic and epigenomic regulatory mechanisms. This novel strategy allowed us to identify and refine prognosis-significant genes, some of which would not have been identified based on expression alone (e.g., SRC). We also highlight several observations based on our integrated composite gene signature and subsequent pathway and gene ontology–based analysis.
Many of the pathways enriched for survival-correlated genes in our analysis have already been implicated in pancreatic cancer including ERBB, focal adhesion, insulin signaling, and MAPK pathways (Supplementary Table S6). Several of these pathways are linked by EGFR, which itself seems in our integrated composite signature. This is not unexpected as it has been shown that EGFR is overexpressed and associated with disease progression and poor prognosis in PDAC (30, 31). Apart from EGFR, our integrated composite signature contained several additional prognosis-associated genes linked to both SRC signaling and the PI3K/AKT pathway. These included SRC, PIK3R1/p85α, and CBL, important regulatory components that we further linked to PDAC survival in separate validation cohorts by either protein or gene expression.
The PI3K/AKT pathway can promote both PDAC initiation and invasive cancer progression. AKT activity is enhanced in up to 60% of PDAC tissues and cell lines (38). More recently, we showed PI3K pathway activation is critical for the onset and acceleration of tumors in mice with conditional Kras activation and Pten deletion (41). The factors responsible for PI3K/AKT pathway dysregulation in PDAC remain unresolved. Activating mutations of the p110 subunit of PI3K are rare in PDAC (42), as are mutations or deletions of PTEN (43). To address this issue, Ying and colleagues (44) recently found that AKT activation is increased in 68% of PDAC, but only less than half could be explained by genomic variation of AKT or PTEN. This leaves open other potential mechanisms of AKT activation in PDAC. Although its expression may be silenced via DNA methylation (45), PTEN expression did not correlate with survival in our analysis. Instead, our data offer potential alternative mechanisms for PI3K/AKT/mTOR dysregulation in PDAC. These include increased expression of the upstream receptor tyrosine kinase EGFR and downregulation of the p110 regulatory subunit p85α encoded by PI3KR1, possibly through silencing mediated by prognostically-linked miR519d. As validation of these observations, we showed p85α protein expression to be inversely correlated with overall survival in our large PDAC TMA cohort. Likewise, either p85α or PTEN loss was found to further correlate with increased P-AKT and P-S6 status as measured in a small cohort of PDAC tumors by IHC (Supplementary Fig. S3).
Our findings are consistent with previous studies reporting higher levels of SRC expression are associated with worse DFS in PDAC (35, 46). SRC is the signature member of a family of nonreceptor tyrosine kinases able to mediate diverse effects on cellular proliferation, differentiation, survival, motility, and angiogenesis (47). It is activated through multiple mechanisms, including via integrins and membrane-bound receptor tyrosine kinases (e.g., EGFR; ref. 36). SRC overexpression and its activation (as detected by Y419 phosphorylation) is seen in most PDACs (48, 35). Likewise, mice with conditional Kras activation and deletion of the SRC inhibitory kinase Csk develop invasive PDAC more rapidly and at higher prevalence than those with intact Csk (49). Our integrated composite signature and subsequent IHC validation study of P-SRC identify a potentially important association between SRC activation and clinical disease progression in PDAC. Our results here also provide potential mechanisms of SRC dysregulation in PDAC. Genomic amplification of the SRC locus on chromosome 20 was found in 6 of the 25 patients (25%) and was independently correlated with worse prognosis (HR = 4.1, P = 0.016). In addition to genomic amplification or its possible activation by EGFR, SRC signaling could presumably be dysregulated in PDAC via changes in CBL (Fig. 3), which was correlated with better prognosis in both our initial analysis and separate validation cohort. Although our data are correlative at this point, it raises the intriguing possibility that CBL may act as linchpin molecule regulating SRC and/or PI3K/AKT signaling based on its ability to act as ubiquitin ligase targeting both EGFR and SRC (32, 33).
To our knowledge, this represents the first survival-based integrated analysis of molecular changes in PDAC that considers the multiple dimensions of mRNA, miRNA, and CNA. Our approach offers a paradigm for future larger and more complex multidimensional studies seeking to link clinical phenotype with the highly diverse molecular alterations that define PDAC or other cancer types. Although our study is a preliminary and retrospective analysis of PDAC patients with resected disease, it provides several candidate biomarkers with the potential to stratify risk for disease progression or predict response to molecular targeted therapy. Further prospective and mechanistic studies are not only needed to validate the prognostic or predictive value of these markers following surgical resection, but also to establish their potential use in the nonoperative or neoadjuvant setting. Regardless, we have provided multiple lines of correlative data showing that dysregulation of the PI3K/AKT pathway and SRC signaling are linked to PDAC clinical disease progression. These data are strong rationale for future studies seeking to link prognostically significant signature genes mechanistically to PI3K/AKT or SRC dysregulation and explore their utility as predictive biomarkers and targets of molecular therapy in the subsets of PDAC patients they define.
Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.
Acknowledgments
The authors thank the UCLA Institute for Molecular Medicine and Hirshberg Foundation for Pancreatic Cancer Research for support of array analysis and the UCLA Pancreas Tissue Bank; the UCLA Department of Pathology Tissue Array Core Facility for TMA construction; the UCLA Department of Pathology Clinical Microarray Core Facility for array processing; and Howard Reber, Joe Hines, and James Tomlinson for surgical specimens used in the study.
Grant Support
This work was supported by UCLA Scholars in Translational Medicine Program (T.R. Donahue); the UCLA Tumor Biology Program (USHHS Ruth L. Kirschstein Institutional NRSA # T32 CA009056; L.M. Tran and R. Hill), the Damon-Runyon Cancer Foundation (R. Hill), and American Association for Cancer Research/Pancreatic Action Network Career Development Award (D.W. Dawson).
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.