Main

Pancreatic cancer is the third-highest cause of cancer-related mortality and has the lowest 5-year overall survival (OS) rate of all cancer types1. For most patients with pancreatic adenocarcinoma (PDAC), the most common pancreatic cancer histology, chemotherapy remains the mainstay of therapy2,3,4,5. More than 90% of PDACs exhibit activating mutations in KRAS hotspot residues, a majority of which have eluded targeted therapeutic approaches until recently6. KRAS wild-type (WT) tumors are enriched for actionable alterations in mitogen-activated protein kinase (MAPK) pathway genes such as BRAF, NTRK1, NTRK3 and NRG1 (refs. 7,8,9,10,11,12,13,14).

Superior outcomes are observed in molecularly selected cohorts, such as patients with pathogenic variants in BRCA1, BRCA2 and PALB2, and with tumors harboring targetable alterations2,3,15,16,17,18. Molecular profiling studies have also identified several prognostic features associated with poor outcomes in PDAC19,20. In addition, gene expression-based stratification has identified two main subtypes of PDAC, classical and basal-like, that differ in molecular pathology, therapeutic vulnerabilities and outcomes21,22,23,24,25,26,27.

More recently, allelic imbalance (unequal number of maternal and paternal copies) at the KRAS locus has been associated with poor outcomes in PDAC28. Moreover, KRAS mutant allele gains promote aggressive phenotypes in mouse models of pancreatic cancer29. Whole-genome doubling (WGD), a hallmark of advanced cancer and a negative prognostic factor for OS, is a key driver promoting allelic imbalances leading to KRAS mutant dosage gains30,31. However, the prevalence of KRAS mutant allele dosage gains in PDAC and its association with disease progression remains poorly understood.

In this study, we leverage germline and somatic profiling of n = 2,336 patients to study the genomic and clinical correlates of outcomes for patients with PDAC. We incorporate clinical histories of n = 1,480 patients with long-term follow-up. We demonstrate that KRAS mutant dosage gains are a hallmark of disease progression and are prognostic of poor outcomes across all stages of PDAC.

Results

MSK-IMPACT PDAC study cohort

This study included n = 2,336 patients with PDAC whose tumors were prospectively sequenced as part of standard care at Memorial Sloan Kettering Cancer Center (MSK; Methods; Supplementary Tables 1 and 2). At diagnosis, 31% (n = 731) of patients had resectable tumors, 25% (n = 581) had borderline resectable/locally advanced tumors and 44% (n = 1,024) had metastatic disease (Table 1). Detailed clinical information including lines of treatments, time on treatment and best overall responses were manually curated for 63% (n = 1,480) of patients. The majority (61%) of the sequenced specimens were from primary PDAC, and 39% were from distant metastases. Median age at diagnosis was 67 years. Tumor specimens were sequenced to median depth of 606× using the FDA-authorized MSK-IMPACT clinical sequencing assay that encompasses up to 505 cancer genes32. Somatic substitutions, insertions, deletions, focal copy number amplifications, homozygous deletions and fusions in select genes were identified using a clinically validated pipeline and annotated using the FDA-recognized precision oncology knowledge base, OncoKB33,34.

Table 1 MSK-IMPACT PDAC study cohort characteristics
Full size table

Genomic characteristics of PDACs

Overall, 95% (n = 2,209) of tumors harbored oncogenic alterations in KRAS (Fig. 1a). Notably, this includes 1% (n = 22) of tumors in which KRAS mutations were identified with sequencing read evidence below thresholds for clinical reporting (Methods). Expectedly, tumors WT for KRAS (KRASWT) were significantly enriched for oncogenic alterations in other MAPK pathway genes such as BRAF, NRAS, NF1, NTRK1, NTRK3, FGFR2, ERBB2, MAP2K1, ROS1, MET and RAF1 (collectively, 60% in KRASWT versus 7% in KRASMUT, P = 1.6 × 10−47; Fig. 1a). Seven of 26 oncogenic BRAF alterations were in-frame deletions between amino acids N486 and P490, which were nearly absent in BRAFMUT tumors of melanoma (2 of 749) and thyroid cancer (0 of 473)14,35 (Extended Data Fig. 1a). Moreover, gain-of-function oncogenic fusions involving MAPK genes were nearly exclusive to KRASWT tumors (3.2% in KRASWT versus 0.04% in KRASMUT, P = 2 × 10−35). Collectively, these other MAPK pathway-altered KRASWT (referred to as other-MAPKMUT) tumors comprise 3% (n = 76) of tumors. Finally, 2% (n = 51) of tumors were WT for any MAPK pathway alteration (referred to as MAPKWT). Hypothesizing that these MAPKWT tumors may harbor occult MAPK alterations that eluded detection by genomic sequencing, we performed transcriptome sequencing on 11 tumors where sufficient quality material was available. We identified activating fusions involving MAPK pathway genes BRAF and NRG1 in two of 11 tumors (18%), suggesting that a substantial fraction of these tumors may be driven by alterations in non-MAPK pathway genes.

Fig. 1: Somatic alteration landscape in PDAC.
figure 1

a, Oncoprint of somatic oncogenic alterations in selected genes (Methods) across the following three genomic groups: KRASMUT, other-MAPKMUT and MAPKWT. Tile plot on the left indicates gene-level alteration enrichment in other-MAPKMUT and MAPKWT subtypes compared to KRASMUT using two-sided Fisher exact test. An asterisk indicates that NRG1 was not included in the enrichment analysis as it was profiled only in a subset of samples (Methods). Other MAPK pathway genes include FGFR1, ERBB3, FGFR4, EGFR, RASA1, CBL, MAPK1, ALK, MAP2K2, ERRFI1, FLT3, JAK2, KIT, PDGFRA, RAC1, RET, RRAS2, SOS1 and SPRED1. b, TMB-H and MSI-H prevalence by genomic subtype among tumor samples with ≥30% purity (n = 1,126). c, Age at diagnosis by genomic subtype with statistical comparison by two-sided Wilcoxon rank sum test (n = 2,336). Boxes represent the 25th, 50th (median) and 75th percentiles. Whiskers represent the minimum and maximum values, no further than 1.5× the interquartile range from the respective quartiles, with points beyond this range plotted individually. d, Kaplan–Meier curves showing OS for three genomic subtypes (Methods). e, Forest plot of multivariable Cox proportional hazards model of OS in overall cohort (top, corresponding to d; n = 2,270) and among patients who did not receive targeted therapies (bottom; n = 2,187). Models were stratified by stage at diagnosis and adjusted for sex, age, ancestry, disease status, resection and interval between diagnosis and sample collection (full model shown in Extended Data Table 1). f, Distribution of stage at diagnosis by genomic subtypes (n = 2,336). g, Prevalence of oncogenic alterations in GNAS and CDKN2A/CDKN2B by stage at diagnosis among tumor samples with purity >30% (n = 1076). h, Genetic ancestry by stage at diagnosis (n = 2,336). i, Tumor ___location (body/tail versus head) by stage (n = 1,462). Statistical significance is displayed as nominal P value for significant results after multiple test correction by FDR by Wilcoxon rank sum test for c, two-sided chi-squared test for h and two-sided Fisher’s exact test for b, g and i. Error bars represent 95th percentile binomial CI around the mean for b, g and i, and 95th percentile confidence intervals of the HR estimate (colored squares) for e. mo, month.

Full size image

Several significant differences in oncogenic alterations were identified among these three genomic subtypes of PDAC (KRASMUT, other-MAPKMUT and MAPKWT; Fig. 1a and Extended Data Fig. 1b). TP53 mutations were significantly more frequent in KRASMUT (78%) compared to other-MAPKMUT (38%, P = 1.4 × 10−9) and MAPKWT (45%, P = 1.7 × 10−5) tumors. Interestingly, the TP53 alteration rate among BRAF-mutated tumors was indistinguishable from that of KRASMUT tumors, supporting prior observations that BRAF mutations phenocopy KRAS mutations in pancreas cancers (78% in KRASMUT versus 73% in other-MAPKMUT/BRAFMUT, P = 0.6)36. In contrast, GNAS, SMARCB1 and PIK3CA alterations were significantly enriched in MAPKWT tumors compared to KRASMUT (16% versus 2%, P = 1.5 × 10−5; 8% versus <1%, P = 2.4 × 10−6; 10% versus 2%, P = 0.002, respectively). SMARCB1 loss in KRASWT tumors has previously been associated with the monomorphic anaplastic subtype of undifferentiated carcinomas with rhabdoid features37. Here, three of four SMARCB1-altered tumors in MAPKWT subtype presented with aggressive histologic features of either poor differentiation or high grade, although none exhibited undifferentiated rhabdoid features. FOXP1 and CREBBP alterations were nearly exclusive to other-MAPKMUT (P = 5.3 × 10−5 and P = 1.3 × 10−5, respectively, compared to KRASMUT; 0% prevalence in MAPKWT). ARID1A alterations were enriched in other-MAPKMUT compared to KRASMUT tumors (21% versus 8% prevalence, P = 2.4 × 10−4). High tumor mutational burden (10 or more nonsynonymous mutations per megabase, TMB-H) and microsatellite instability-high (MSI-H) tumors were infrequent (1.5% and 0.4%, respectively) and were enriched in KRASWT tumors (TMB-H—6.1% versus 1.8% in KRASMUT tumors, P = 0.03; MSI-H—3.7% versus 0.3%, P = 0.006; Fig. 1b)11.

We next explored the association between various clinical characteristics and genomic subtypes of PDAC. Age at diagnosis significantly varied by KRAS alteration status. Compared to KRASMUT tumors (median age = 67 years), patients in other-MAPKMUT (64 years, P = 0.03) and MAPKWT (58 years, P = 2 × 10−4) subtypes presented with a significantly earlier age at diagnosis (Fig. 1c and Extended Data Fig. 1c). Patients with MAPKWT tumors also had a significantly different ancestry composition compared to patients with KRASMUT (two-sided chi-squared test, P = 0.004) or other-MAPKMUT tumors (P = 0.04), marked by elevated rates of East Asian patients in other-MAPKMUT and patients with African ancestry in MAPKWT subtypes, indicating the possibility of unrecognized driver alterations in underrepresented patient populations (Fig. 1c). No significant differences were observed in gene-level alteration frequencies across sex or ancestry group (Extended Data Fig. 1d,e).

To evaluate differences in OS across the three genomic subtypes, we used a multivariable Cox proportional hazards model stratified by clinical stage accounting for age, resection status, disease status, sex and genetic ancestry (Methods). We observed that patients with other-MAPKMUT and MAPKWT tumors had significantly longer OS compared to patients with KRASMUT tumors (KRASMUT versus MAPKWT—adjusted hazard ratio (HRadj) = 0.69, CI = 0.48–0.98, P = 0.041; KRASMUT versus other-MAPKMUT—HRadj = 0.69, CI = 0.51–0.93, P = 0.014; Fig. 1d,e, top, and Extended Data Table 1). Improved outcomes among KRASWT patients have been attributed to higher prevalence of targetable alterations in this group14. After excluding patients who received targeted therapies (n = 83), we noted that OS was indistinguishable between patients with KRASMUT and other-MAPKMUT tumors (HRadj = 0.95, P = 0.7), whereas patients with MAPKWT tumors had significantly better OS (HRadj = 0.68, CI = 0.47–0.97, P = 0.035; Fig. 1e, bottom, and Extended Data Table 1). This suggests that the well-recognized OS advantage among patients with KRASWT tumors extends beyond those with targetable alterations.

Somatic and clinical characteristics across clinical stages

The stage at diagnosis was similar across genomic subtypes (Fig. 1f). Among other genes, GNAS, which is associated with IPMN precursor lesions with improved outcomes over PanIN-derived PDACs, was significantly enriched in resectable tumors (7.9%) while CDKN2A/CDKN2B alteration rate increased with disease progression ranging from 44% in resectable disease to 60% in metastatic tumors (P = 6.9 × 10−5; Fig. 1g and Extended Data Fig. 2a)38,39. No pathway-level differences in alteration rates were observed across disease stages (Extended Data Fig. 2b). Patients with African ancestry presented with more advanced disease (87%) compared to other ancestries (67% in European and 63% in East Asian; Fig. 1h). The ___location of the tumor in the pancreas was also strongly associated with stage—resectable tumors most commonly arose in the head (71%), while metastatic tumors frequently arose in the body/tail (57%; Fig. 1i and Extended Data Fig. 2d).

Germline mutations and concomitant somatic alterations

Ten percent of all patients harbored germline pathogenic mutations in high- and moderate-penetrance genes, including BRCA2 (n = 86, 3.7%), BRCA1 (n = 41, 1.8%), ATM (n = 41, 1.8%) and PALB2 (n = 11, 0.5%; Fig. 2a). Lynch syndrome with germline mutations in mismatch repair genes MLH1, MSH2, MSH6 and PMS2 was identified in 17 patients (0.7%), of which six (35%) presented with MSI inferred from sequencing40,41. Although pathogenic germline variants in BRCA1, BRCA2 and CHEK2 were more frequent among those with Ashkenazi Jewish (ASJ) ancestry, ATM and PALB2 germline variants were more common among those without ASJ ancestry (Extended Data Fig. 3a). Patients with metastatic disease presented with slightly elevated rates of germline pathogenic variants compared to those with earlier stage disease (12.3% versus 9.6%, P = 0.01), driven primarily by an enrichment in BRCA1/BRCA2 alterations in tumors from patients with metastatic disease (Extended Data Fig. 3b). Germline pathogenic variants were more frequent among MAPKWT (25%) compared to other-MAPKMUT (17%, P = 0.05) and KRASMUT (10%, P = 0.001). This elevated rate was underpinned by higher prevalence of pathogenic germline mutations in ATM (gATM) in MAPKWT patients compared to KRASMUT (18% versus 1.4%, P = 2 × 10−6; Fig. 2a). Germline ATM mutations also co-occurred with somatic GNAS mutations, possibly attributed to IPMN-derived PDAC42,43. Notably, this increased gATM burden in MAPKWT tumors was specific to mutations of germline origin as the rate of somatic ATM mutations (sATM) in MAPKWT tumors was not different from the other genomic subtypes (P > 0.05; Fig. 2a)43.

Fig. 2: Germline alteration landscape in PDAC.
figure 2

a, Oncoprint of pathogenic germline variants and somatic oncogenic alterations by genomic subtype. Displayed sBRCA2, sBRCA1 and sATM alterations are exclusively in sporadic tumors without a pathogenic germline variant in these genes. An asterisk indicates that MMR includes MSH2, MSH6, MLH1 and PMS2. ‘Other’ includes all other pathogenic variants in high- and moderate-penetrance genes—BRIP1, CDKN2A, CHEK2, FLCN, HOXB13, MITF, NBN, NF1, RAD51D, SDHA, SMARCA4, STK11, TP53 and TSC1. Tiles at left show gene-level enrichment in other-MAPKMUT and MAPKWT subtypes compared to KRASMUT by two-sided Fisher’s exact test. b, Rates of loss of heterozygosity (Methods) at BRCA2, BRCA1, ATM, PALB2, MMR and other loci (as in a) in n = 1,946 patients with germline pathogenic variants, somatic mutations in sporadic cancers and patients WT for any alteration in corresponding genes (comparisons by two-sided Fisher's exact test). Error bars represent 95th percentile binomial confidence intervals around the mean. c, Pattern of germline and somatic ATM and TP53 alterations, with monoallelic or biallelic zygosity status indicated. d, Pattern of germline and somatic BRCA1 and TP53 alterations, as in c.

Full size image

Strong selection for biallelic inactivation through copy number loss of heterozygosity (LOH) was observed in gBRCA1, gBRCA2, gPALB2 and gATM carriers (Fig. 2b; Methods). The rate of LOH was higher in sporadic tumors with somatic mutations in BRCA2 and ATM compared to tumors without mutations in these genes, but lower compared to LOH rates observed in germline carriers18,44,45,46. As previously reported, gATM and sTP53 mutations were mutually exclusive (P = 4 × 10−11), and sTP53 alterations and gBRCA1 mutations co-occurred (P = 4 × 10−7)44,47. Interestingly, taking zygosity into account, irrespective of origin of mutation (somatic or germline), tumors with biallelic loss of ATM were notably depleted for sTP53 mutations (Fig. 2c)44. In contrast, both monoallelic and biallelic BRCA1 mutations showed equal co-occurrence with sTP53 mutations (Fig. 2d).

Mutant allelic imbalance at KRAS locus

We next evaluated the extent to which KRASMUT tumors harbor copy number allelic imbalance (unequal number of maternal and paternal alleles) at KRAS locus and mechanism by which allelic imbalance occurs28 (Fig. 3a). To increase sensitivity and specificity to infer allelic state, we restricted our analysis to n = 1,157 KRASMUT tumors with sufficiently high-quality copy number fits (Methods). Of these 1,157 tumors, 42% presented with allelic imbalance at KRAS locus including focal or arm-level amplifications (4%, n = 48), shallow gains (16%, n = 186), copy-neutral LOH (CNLOH; 5%, n = 62), LOH (11%, n = 129) and losses after WGD (5%, n = 56; Fig. 3b). Selection for KRAS mutant allele was widespread with 93% of all imbalance events preferentially gaining or retaining the mutant allele. One in five (19%) KRASMUT tumors harbored WGD, and the rate of allelic imbalance was substantially higher in WGD tumors (75%) versus non-WGD tumors (30%; Fig. 3c). Consequently, KRASMUT allele gain was more than twice as common in WGD tumors (43% with three or more mutant copies) compared to non-WGD tumors (20% with two or more mutant copies; Fig. 3d). Moreover, the magnitude of mutant-allele gains also varied by WGD status. Tumors with two or more mutant-allele gains (that is, ≥4 mutant-alleles in WGD and ≥3 mutant-alleles in non-WGD tumors) were substantially more frequent in tumors with WGD (20% versus 6% in non-WGD tumors).

Fig. 3: KRAS mutant allele dosage gains and their prognostic implications.
figure 3

a, Schematic of copy number states as related to mutant copy gain and retention of the WT allele in cases of allelic imbalance. The ‘Loss after WGD’ state indicates any copy number losses of the minor allele following WGD but excludes complete losses of the minor allele which are considered as CNLOH (Methods). b, Overall prevalence of copy number states in KRASMUT tumors. c, Prevalence of copy number state with allele selection by WGD status. d, Estimated number of gained mutant KRAS copies by WGD status. Tumors in which the WT allele was gained, or the mutant allele was lost are not shown here (n = 32, ‘WT selection’ in c). e, Kaplan–Meier curves of OS stratified by the number of mutant KRAS copies in diploid (non-WGD) KRASMUT tumors excluding tumors with gain of WT allele (n = 865; Extended Data Table 2). f, Prevalence of KRAS copy number states (left) and mutant copy gain (two or more mutant copies; right) by clinical stage at diagnosis (n = 874). Statistical comparisons show pairwise two-sided Fisher's exact tests. g, Kaplan–Meier curves of OS stratified by the number of mutant KRAS copies as in e, within each clinical stage at diagnosis (Extended Data Table 2). h, Kaplan–Meier curves of OS stratified by copy number state in diploid (non-WGD) KRASMUT tumors excluding tumors with gain of WT allele (n = 865; Extended Data Table 3). i, Forest plots of multivariable Cox proportional hazards model of OS by KRAS copy number state as in h, within each clinical stage at diagnosis (n = 865; Extended Data Table 3). Error bars represent 95th percentile binomial CI in f, and 95th percentile CI of the HR in i. Displayed P values in e, g–i are two-sided nominal P values from multivariable Cox proportional hazards models that include age, sex, ancestry and time from diagnosis to sample collection as covariates. Models for e and h are stratified by clinical stage at diagnosis (Extended Data Tables 3 and 4). R, resectable; M, metastatic.

Full size image

KRAS mutant dosage and prognostic effect

KRASMUT PDAC tumors with a gain of mutant-allele present with aggressive phenotypes in mice and have worse OS in patients28,29. Herein we sought to evaluate the prognostic effect of dosage gains of KRASMUT on OS across disease stages. To mitigate the confounding effects of WGD, a notable negative predictor of OS, we limited our analysis to non-WGD tumors (n = 934; Methods)30. In a multivariable Cox model of OS stratified by stage at diagnosis and adjusted for sex, age at diagnosis, genetic ancestry, and time from diagnosis to sample collection, we observed that patients with tumors with any gain of KRASMUT allele had significantly lower OS compared to those with tumors with one mutant allele (P = 3.5 × 10−7, HRadj = 1.7, CI = 1.4–2.0; Fig. 3e and Extended Data Table 2). Although no significant difference in OS was observed between patients with tumors harboring two mutant copies (n = 130) and those with 3 or more mutant copies (n = 49) when stratifying by stage, a larger sample size may be required to identify incremental effects of additional dosage gains on prognosis.

KRASMUT allele dosage gains were more frequent in metastatic tumors (29%) compared to tumors from patients with locally advanced (14%) or resectable (8%) disease and were strongly correlated with advanced disease (P = 5.7 × 10−11, chi-squared trend test; Fig. 3f). This supports prior observations in mouse models with KRAS G12D dosage gains, which demonstrated amplified RAS transcriptional programs and exhibited rapid disease progression29. Notably, we observe the poor prognostic effect of mutant allele dosage across disease stages (Fig. 3g). Among patients with resectable disease, the median OS was significantly lower in patients with tumors harboring dosage gains of the mutant-allele compared to those whose tumors had one mutant copy (23 months versus 32 months; HRadj = 2.16, CI = 1.1–4.3, P = 0.03). Similarly, among patients with metastatic disease, KRASMUT dosage gains were associated with 5 months shorter OS compared to patients whose tumors had a single KRASMUT copy (8.5 months versus 13 months; HRadj = 1.63, CI = 1.3–2.1, P = 4.9 × 10−5; Extended Data Table 2).

The WT KRAS allele in KRASMUT tumors has previously been shown to act as a tumor suppressor14,48,49. Tumors adapt during disease progression by losing the WT allele or acquiring additional copies of the mutant allele50,51. We therefore hypothesized that patients with KRASMUT tumors with either gain of mutant or with loss of WT should have substantially worse prognosis compared to patients with balanced KRASMUT tumors. However, we observed that only the KRASMUT allele dosage is a notable predictor of poor prognosis, independent of loss or retention of the WT allele. Among patients with KRASMUT tumors that have a single copy of the mutant allele, OS was indistinguishable between patients with tumors that retained the WT allele and those that lost the WT allele (Fig. 3h and Extended Data Table 3). Both groups of patients also demonstrated significantly improved OS compared to patients with KRASMUT tumors with mutant-allele gains, irrespective of WT allele status. However, interestingly, among patients with tumors with KRASMUT dosage gains, patients with tumors with loss of WT had significantly worse OS compared to patients with tumors with WT retained (HRadj = 1.64, 95% confidence interval (CI) = 1.1–2.5, P = 0.016 for CNLOH with gain-of-mutant as reference). This effect was most prominent among metastatic tumors (Fig. 3i and Extended Data Table 3). Taken together, our findings indicate a synergistic effect of losing the WT allele in tumors with KRASMUT dosage gains in promoting worse disease outcomes.

KRAS mutant allele-specific differences

Over 98% of all codon substitutions in KRAS were at G12 (91%) and Q61 (7%) residues (Fig. 4a). Among KRASMUT tumors, G12D (41%) was the most abundant hotspot mutation followed by G12V (32%), and G12R (16%). Fourteen tumors (0.6%) harbored multiple KRAS hotspot mutations (Extended Data Fig. 4a,b). No significant differences were observed in KRAS variant prevalence by stage at diagnosis, sex, genetic ancestry or age at diagnosis (Fig. 4a and Extended Data Fig. 4c). Patients with G12R-mutant tumors compared to those with G12D-mutant tumors were less likely to have a smoking history (40% versus 55%, P = 2 × 10−4; Extended Data Fig. 4c).

Fig. 4: Differential genomic and prognostic features of KRAS variants.
figure 4

a, Top, prevalence of the most common KRAS variants among all KRASMUT tumors. Bottom, disease stage composition, distributions of allelic imbalances and KRASMUT allele dosage gains across patients with different KRAS variants. Allelic imbalances and dosage gains are shown only for non-WGD tumors as indicated. Displayed nominal P value denotes statistical comparison of distribution of gene-level copy number state by two-sided chi-squared test. Tumors with multiple driver KRAS mutations are not shown here (n = 14; Extended Data Fig. 4a,b). Error bars for KRASMUT CN represent 95th percentile binomial CIs around the mean. b, Kaplan–Meier curve of OS among KRAS G12D, G12V and G12R variants. Displayed P values are nominal two-sided P values from a multivariable Cox proportional hazards model stratified by stage at diagnosis and accounting for KRAS mutant allele copy number, sex, age, ancestry, disease status and interval from diagnosis to sample collection (Extended Data Table 4). c, Oncoprint showing the prevalence of co-occurring TP53, CDKN2A/CDKN2B, SMAD4, ARID1A, AKT2, and RB1 mutations and TGFβ, SWI/SNF, PI3K and RTK–Ras signaling pathway alterations among KRAS G12D, G12R and G12V tumors. Tiles on the left indicate pairwise enrichment testing of co-occurring alterations between KRAS G12D versus G12R and G12V mutant tumors by two-sided Fisher's exact test.

Full size image

In a multivariable Cox model stratified by stage at diagnosis and accounting for mutant KRAS gain, we evaluated the differences in OS among patients with tumors harboring the most common KRAS variants (Extended Data Table 4). Patients with KRAS G12R-mutant tumors had significantly better OS compared with G12D-driven cancers (HRadj = 0.78, CI = 0.67–0.92, P = 0.003; Fig. 4b). Among patients with de novo metastatic disease, those with KRAS G12D, G12V or G12R mutations did not show significant differences in progression-free survival between those who received first-line FOLFIRINOX compared to gemcitabine-based therapy (Extended Data Fig. 4e,f). Despite the noted differences in functional effects of different KRAS alleles, the rates of genome doubling and mutant-allele dosage gains were similar across alleles (Fig. 4a and Extended Data Fig. 4d)52. However, in comparing G12D and G12R tumors, SMAD4 alterations were significantly more common in G12R tumors (30% in G12R versus 21% in G12D, P = 0.001), while ARID1A alterations were more frequent in G12D tumors (10% in G12D versus 5% in G12R, P = 0.002; Fig. 4c)53.

Clinical actionability in PDAC

Approximately 10% of patients with PDAC harbored standard-of-care biomarkers of response to targeted therapies as defined by OncoKB level 1 (for example, MSI-H, TMB-H and oncogenic MAPK pathway alterations in NTRK1, NTRK3 and RET) and level 2 (KRAS G12C, BRCA1/BRCA2, PALB2; OncoKB version 4.12, December 2023)34. Most strikingly, an additional 78% of patients harbored biomarkers with compelling clinical evidence of response to specific drugs (OncoKB level 3A), nearly all of which (98%) were attributed to KRAS G12D/V/R/A/S mutations for which favorable responses to RAS inhibitors have been observed in recent phase I trials (Fig. 5a)6,34. The clinical actionability varied across the genomic subtypes with other-MAPKMUT patients, expectedly, harboring the highest fraction (18%) of OncoKB level 1 alterations (Extended Data Fig. 5). For patients with metastatic PDAC who were fit for systemic therapy, chemotherapy with either a 5-fluorouracil (5-FU) or gemcitabine-based backbone was administered (Table 1 and Fig. 5b)5. Treatment selection reflected changes in standard-of-care guidelines over the course of the study54.

Fig. 5: Clinically actionable alterations and treatment landscape.
figure 5

a, Highest OncoKB level of evidence by patient, with actionable (levels 1, 2 and 3A) alterations labeled. b, Alluvial diagram of treatment sequence for all 1,480 patients with treatment annotation by clinical stage at diagnosis showing prevalence of neoadjuvant therapy, resection surgery, adjuvant therapy and up to 5 rounds of systemic therapy. Different lines of treatment are shown on the x axis. c, Clinicogenomic characterization of n = 29 patients who were metastatic at presentation and received PARP-inhibitor therapy as part of systemic treatment. Cumulative times on PARP-inhibitor therapy and platinum therapy are shown as bar charts, along with demographic information (age at diagnosis, sex and genetic ancestry), germline and somatic alterations in HRD genes with associated zygosity, and somatic alterations in other commonly altered genes or genes of interest. The patient marked with a plus (+) received PARP-inhibitor therapy to target a germline RAD50 mutation (not shown). d, OS for n = 304 patients who were metastatic at presentation and received either 5-FU or gemcitabine-based first-line treatments (top left). Top-right, bottom-left and bottom-right, Kaplan–Meier curves of OS by alteration status of indicated genes. P values are nominal two-sided P values from a multivariable Cox model that accounts for sex, age, ancestry and interval from diagnosis to sample collection for each gene (Methods; Extended Data Table 5).

Full size image

We next evaluated molecular correlates of response to poly-ADP ribose polymerase inhibitors (PARPi). In our clinically curated cohort, n = 29 patients with stage IV disease at diagnosis received PARPi therapy at any time during their treatment. Of these, 25 patients had germline (n = 23) or somatic (n = 2) mutations in BRCA2 (n = 18), BRCA1 (n = 6) and PALB2 (n = 1; Fig. 5c). Twenty-eight of 29 patients had KRASMUT tumors, 4 of which also had mutant allele dosage gains. Overall, 38% (11 of 29) of patients received PARPi therapy for more than six months (median 17 months; range = 8–56). Ten of 11 (91%) patients had BRCA2MUT; the remaining patient had a BRCA1WT/BRCA2WT and MAPKWT tumor with a BAP1 loss-of-function fusion. While all ten BRCA2-mutant tumors deriving benefit from PARPi (defined as >6 months on treatment) had biallelic inactivation, no benefit was observed among six other BRCA2-mutant tumors with biallelic losses, suggesting that biallelic loss of BRCA1/BRCA2 is an important but insufficient biomarker alone of response to PARPi55.

We next sought to identify the genomic correlates of OS in KRASMUT patients with de novo metastatic disease who received first-line standard-of-care chemotherapy (n = 304, median OS = 10.5 months; Methods; Fig. 5d). We evaluated associations between genes altered in at least 3% of patients (n = 13 genes) and OS using a multivariable Cox regression model accounting for sex, age, ancestry and interval from diagnosis to biospecimen collection (Methods; Extended Data Table 5). Only patients with BRCA2 alterations (including germline carriers) showed modestly improved outcomes (HRadj = 0.66, 95% CI = 0.44–0.98, P = 0.038), in line with the established response to platinum therapy16,18,56 (Fig. 5d). After correcting for multiple hypothesis testing, alterations in two genes remained significant. Both were associated with a significantly shorter OS—RNF43 (HRadj = 2.79, 95% CI = 1.40–5.56, Padj = 0.047; median OS = 6.6 months RNF43MUT versus 10.8 months RNF43WT) and AKT2 amplifications (HRadj = 2.03, 95% CI = 1.26, 3.29, Padj = 0.048; median OS = 8.6 months in AKT2MUT versus 10.8 months in AKT2WT). Although not commonly altered in PDAC (<5% in this cohort), RNF43 loss-of-function activity has been shown to be KRAS dependent57. While KRAS mutant gains were more frequent in RNF43MUT tumors, this did not reach statistical significance in our cohort (67% of RNF43MUT tumors had KRAS mutant gains versus 30% in RNF43WT, P = 0.2). Together, these results suggest that prognosis in patients with PDAC is most linked to KRAS mutation and dosage in the majority of patients, and that somatic alterations in other genes do not explain the wide variation in response to standard therapies.

Discussion

Herein we study the detailed clinicogenomic profiles of n = 2,336 patients to characterize prognostic biomarkers guiding clinical outcomes in PDAC. Compared to prior studies, our cohort has several unique strengths including many molecularly profiled metastatic tumors, detailed clinical histories for two-thirds of patients, and matched tumor and normal sequencing that allowed for incorporation of germline results and robust inference of allele-specific copy number states at the mutated loci. Collectively, these strengths enabled improved classification of PDAC into genomic subtypes and the demonstration of prognostic potential of mutant allele dosage gains in KRASMUT tumors across all disease stages.

Our findings have several clinical implications. The predictive role of KRAS mutant dosage gains in impacting antitumor response and therapeutic outcomes remains unknown, although we expect it to be important. As many RAS-targeted therapies are advancing in clinical trials, it is imperative to conduct analyses to evaluate whether KRAS mutant dosage has a role in mediating efficacy as well as depth and duration of response to treatment. Findings from such studies may provide impetus for future trials to stratify patients by KRAS mutant dosage in evaluating therapeutic outcomes. For example, randomized clinical trials that evaluate patients with KRAS mutant dosage gains in both the intervention as well as control arms are key to evaluate response to RAS inhibitors as well as current standard-of-care therapies. Second, shallow gains of KRAS mutant alleles, including single-copy gains, which are associated with poor outcomes and can potentially guide care, are currently not reported by any clinical sequencing assays. There is an imminent need for evaluation, standardization and incorporation of these low-level gains into reports for clinical interpretation. Third, acknowledging the challenge of small biopsy specimens in the diagnosis and management of patients with PDAC, it is essential to understand whether circulating tumor DNA can provide a noninvasive means of assessing KRAS mutant dosage and whether sequential analyses could provide a dynamic and real-time assessment to inform response and resistance.

Although PDACs are genomically grouped by KRAS status, our data argue that the KRASWT tumors comprise two distinct molecular groups stratified by whether or not they harbor oncogenic alterations in other MAPK pathway genes. The other-MAPKMUT and MAPKWT subtypes comprise 3% and 2% of all PDACs, respectively, and reveal distinct etiologic and molecular factors underpinning tumorigenesis. However, further molecular and biochemical studies that measure the activity of Ras–ERK signaling as well as the upregulation of downstream components including MYC or YAP are needed to establish MAPK independence among these MAPKWT tumors58.

In this large clinicogenomically characterized cohort of patients with PDAC, several important limitations need highlighting. Patients were evaluated and treated at a major referral center for PDAC with inherent biases of fitness, age and other characteristics. With current trends of increasing incidence of early-onset PDAC, especially in younger women, these characteristics are important to note59,60. We report a high proportion of patients undergoing resection of primary pancreatic cancer, 41%, which is higher than the broader population of all patients with PDAC (Table 1). Relatedly, only 15% of our cohort was non-European, highlighting the need for these results to be evaluated in subgroups with greater ethnic and racially diverse populations. This study primarily focuses on genomic features derived from MSK-IMPACT. Additional transcriptomic assessments are needed to fully characterize pathway-level activity and to identify the relationship between genomic biomarkers described herein and previously described transcriptomic subtypes (classical and basal) that have been associated with prognosis and for which predictive implications are under evaluation in clinical trials (NCT04469556)25. Given the prognostic potential of both KRAS dosage and the basal-like transcriptomic subtype, a comprehensive analysis of subtyping is needed to determine whether these are indicative of the same or independent pathways of tumor progression.

Comprehensive analysis of germline and somatic alterations is critically important in PDAC and both have therapeutic implications. Additionally, integration with demographic and clinical information is required to link understanding of underlying tumor biology with treatment implications and outcomes. We identify distinct molecular and clinical features of KRASMUT, other-MAPKMUT and MAPKWT genomic subtypes and further stratify KRASMUT tumors by allele-specific copy number. Collectively, these findings demonstrate that increased KRAS mutant allele dosage is an important negative prognostic feature that will need to be integrated into clinical practice as we move into the era of RAS-directed therapeutics for PDAC.

Methods

The research was reviewed and overseen by the Institutional Review Board at Memorial Sloan Kettering Cancer Center.

Study cohort

Between January 2014 and September 2021, 2,671 tumor and matched normal DNA samples from 2,566 patients with PDAC were subjected to molecular testing using MSK-IMPACT, an FDA-authorized clinical next-generation sequencing panel. Patients provided consent for tumor profiling under an institutional Memorial Sloan Kettering Cancer Center (MSK) IRB-approved research prospective protocol, ‘Tumor Genomic Profiling in Patients Evaluated for Targeted Cancer Therapy’ (NCT01775072). Patients were not compensated financially for participation in the study. Tumor samples with low sequence coverage (<100×) or no detectable somatic alterations likely due to low tumor content quality were excluded (112 samples from 99 patients). Upon further review, 62 patients were excluded for having a diagnosis other than PDAC, and an additional 69 patients were excluded for having incomplete medical records. From the remaining 2,336 patients, one representative sample from each patient was selected for further analysis based on several criteria including higher tumor purity, higher sequence coverage, sample passing of allele-specific copy number (FACETS) quality control criteria (https://github.com/taylor-lab/facets-preview/) and the specific IMPACT panel size on which the specimen was sequenced. Of these 2,336 patients, 1,480 patients who had at least one year of clinical follow-up at our center between January 2014 and March 2021 underwent comprehensive clinical annotation after manual curation of medical health records. Genetic ancestry was inferred from IMPACT as described previously60. Comprehensive demographic and clinical information are included in Table 1. Self-reported sex information was collected from patient medical records. No sex-based criteria were used to select patients in this study. Gender information was not considered in the study.

Tumor sequencing and mutation assessment

Tumor samples were sequenced using the MSK-IMPACT targeted sequencing panel as described previously33. Briefly, tumor type and purity were assessed by pathology from H&E-stained slides of tumor samples. Genomic DNA from formalin-fixed, paraffin-embedded (FFPE) tumor and matched normal (peripheral blood) samples was extracted, and targeted sequencing was performed using custom DNA probes against all exons and selected introns of a given panel of cancer genes. Tumors were sequenced using four different generations of MSK-IMPACT panels containing 341 genes (n = 17 samples), 410 genes (n = 438), 468 genes (n = 1,536) and 505 genes (n = 345). The median sequencing depth of tumors was 606× (25th percentile: 469×, 75th percentile: 749×). The median purity of the tumors, estimated by FACETS, was 31%. For tumors without FACETS estimated purity, the pathologist estimated tumor purity was used to determine if the sample meets sufficient quality criteria to be included in the cohort (see above). Somatic mutations, copy number alterations and structural rearrangements in select genes were identified using a previously described pipeline validated for use in a Clinical Laboratory Improvement Amendments (CLIA)-compliant laboratory33. All somatic alterations were annotated for clinical actionability using OncoKB version. 4.12, December 202334. All somatic alterations that were identified as ‘oncogenic’ or ‘likely oncogenic’ by OncoKB were considered as drivers and used exclusively in all analyses (annotated as [Gene]MUT, regardless of alteration type) unless otherwise noted. Genes were included in the oncoprint in Fig. 1a if they met any of the following criteria: mutated in >5% of patients, MAPK pathway, MYC pathway, clinically actionable genes mutated in ≥1% of patients, or significantly enriched in either other-MAPKMUT or MAPKWT tumors. Pathway gene lists were defined by the oncogenic signaling pathways described in ref. 61. Germline variant discovery and pathogenicity assessment was performed as previously described62,63. MSI was determined by MSIsensor, with tumors identified with an MSIsensor score of 10 or higher classified as microsatellite instable (MSI-H, microsatellite instability high)64,65. A select subset of the cohort (n = 90 patients) received additional clinical testing for fusions using the custom RNA-seq panel (MSK-Fusion) that utilizes Archer Multiplex PCR technology66. All tumors wild-type for KRAS mutations were evaluated for read evidence at subdetection thresholds. In these tumors, we genotyped known KRAS hotspot mutations using a custom tool (https://github.com/mskcc/GetBaseCountsMultiSample) and identified reads with mapping quality scores of at least 20. All mapping reads were manually reviewed using IGV (http://software.broadinstitute.org/software/igv/) to identify high-quality reads. A read is considered high quality if there are no mismatches and if the evidence supporting the mutated base is not in the leading or trailing ten bases of the sequencing read.

RNA sequencing for fusion detection

We performed RNA sequencing on 11 tumor specimens derived from FFPE tissue blocks. Fusions were called using FusionCatcher v1.20 (https://github.com/ndaniel/fusioncatcher) and Arriba v2.1.0 (https://github.com/suhrig/arriba). Results from each caller were merged into a single output file according to matching gene IDs and breakpoints. Noncoding fusions from FusionCatcher were excluded from the final output of the pipeline. All fusion calls were annotated by OncoKB. Two tumors harbored activating fusions involving MAPK pathway genes. In one sample, the ATP1B1–NRG1 rearrangement fused the third exon of ATP1B1 to the second exon of NRG1. In another sample, the fourteenth exon of GIT2 was fused to the fourth exon of BRAF. In both events, the entire kinase ___domain of the 3′ gene was retained.

Allele-specific copy number analysis

Allele-specific copy number was inferred using the FACETS algorithm (v0.5.14)67. Briefly, FACETS was run on each tumor sample in a two-step mode as previously described28. The first step aims to identify the normalized tumor to normal sequencing coverage ratio corresponding to the diploid state. Using this, the second step aims to identify the focal gains and losses. FACETS-derived purity and ploidy estimate for the tumor sample and the allele-specific integer copy number for each locus was inferred. To identify and exclude tumors with poor quality copy number fits, we applied a series of quality control criteria (https://github.com/taylor-lab/facets-preview/) that include degree of evidence supporting the diploid state, the fraction of the genome with homozygous deletions, fraction of the genome that is estimated to be subclonal, hypersegmentation, concordance between integer copy number estimate and the allelic configuration and whether FACETS was able to estimate purity (C.B., D.M., Ino de Bruijn, Mingxuan Zhang, Michael V. Gormally et al., in revision). Overall, 1,555 of 2,322 tumor samples for which we were able to generate FACETS profiles had sufficiently high-quality copy number fits evaluable for further analysis. WGD status is inferred as previously described (https://github.com/mskcc/facets-suite)30. Briefly, a tumor is deemed to have undergone genome doubling if more than 50% of the autosomal genome has a major copy number of 2 or higher. In our cohort, 339 of 1,555 tumors were WGD-positive. Clonality of somatic mutations was determined as described previously68. Briefly, we first infer the cancer cell fraction (CCF) of mutation using the variant allele frequency of the somatic mutation, the integer copy number at the locus, the read coverage and the FACETS estimated tumor purity. A mutation is deemed to be clonal if either the CCF estimate is 80% or higher, or if the CCF estimate is 70% or higher and the upper bound of the 95th percentile CI of the CCF estimate is greater than 90% (https://github.com/mskcc/facets-suite). LOH at specific loci harboring germline/somatic mutations or WT alleles was assessed using allele-specific copy number inference from FACETS. A locus is considered to have LOH if the lcn is 0.

KRAS allelic imbalance and mutant allele dosage gains

Of the 1,555 tumors with sufficiently high-quality copy number profiles, 1,157 were KRASMUT with determinable integer copy number. We then assigned each of these tumors to six different allelic copy number states based on the total copy number (tcn) and lower copy number (lcn) at the KRAS locus as well as the WGD status (Fig. 3a). In non-WGD tumors, diploid tumors with tcn:lcn of 2:1 were considered as ‘Balanced’ (heterozygous). Tumors with complete loss of one inherited allele (lcn of 0) were considered as ‘LOH’ if the tcn is 1, and ‘CNLOH’ if tcn is >1. Tumors were considered to harbor KRAS ‘amplifications’ if our clinical pipeline identified focal KRAS amplifications or the tcn in FACETS was ≥5 in the absence of WGD, or ≥6 with WGD. The tumors that retained both inherited maternal and paternal alleles but have gained additional copies of one or both inherited alleles were considered as ‘Gains’. Among WGD-positive tumors, we accounted for the higher balanced state (tcn:lcn of 4:2) in ascribing allelic copy number states associated with losses and gains. WGD-positive tumors that acquired single-copy losses (irrespective of whether the loss targeted the WT or the mutant allele) were considered as ‘Loss after WGD’. To be considered as a ‘Gain’ in WGD-positive tumors, the tcn is required to be 5 or higher.

The expected number of KRAS mutant alleles is estimated using the observed variant allele frequency, tumor purity and total copy number at the KRAS locus as previously described28,69. We considered tumors with a higher number of mutant copies than WT copies to have undergone mutant allele selection. These include tumors identified as having ‘LOH’ (tcn:lcn of 1:0 for non-WGD tumors and 2:0 for WGD-positive tumors), ‘CNLOH’ (complete loss of wild-type with either ≥2 mutant copies in non-WGD tumors, or, ≥3 mutant copies in WGD-positive tumors), ‘Gain of mutant’ (with higher number of mutant copies than wild-type copies), and finally, among WGD-positive tumors, ‘WGD+Loss of WT’ tumors in which a copy number losses after WGD event targeted the wild-type allele (Fig. 3c). Tumors in which selection following copy number gains/losses resulted in preferential retention of the KRAS WT allele were rare (n = 32, ‘WT selection’, Fig. 3c) and were not considered further.

Statistics and reproducibility

Enrichment testing was performed using two-sided Fisher’s exact tests with multiple test correction by the FDR method where appropriate. For gene-level enrichment analyses, we evaluated all cancer genes that are captured on all versions of the MSK-IMPACT panels (n = 341 genes) and were altered in at least 3 samples in the dataset. Unless otherwise noted, enrichment tests were limited to alterations detected by MSK-IMPACT and annotated as drivers (oncogenic or likely oncogenic) by OncoKB. Error bars on all figures showing counts/prevalence data represent 95% CIs of the binomial probability. Two-sided chi-squared tests were used for testing differences in prevalence of nonbinary variables such as ancestry and clinical stage. Two-sided Wilcoxon rank sum tests were used for comparison of numeric variables such as age. In boxplots, the center line represents the median, lower and upper hinges show the 25th and 75th percentiles, and lower and upper whiskers show the minimum and maximum values no further than 1.5× the interquartile range, respectively, with outlier points beyond this range plotted individually. Analyses were conducted using R version 4.3.0 with the rstatix (v0.7.2), stats (v4.3.0) and binom (v1.1.1.1) packages and visualized using ggplot2 (v3.4.4). Oncoprints were generated using cBioPortal for Cancer Genomics70,71,72.

Survival

OS was calculated as the interval from date of IMPACT sample collection to date of death or last follow-up. Patients with no follow-up after IMPACT sample collection were censored on day 1 (n = 9). Patients with missing date of IMPACT specimen collection (n = 2) or those with IMPACT specimen collection before January 2014 (n = 64) were excluded from survival analyses. All P values and adjusted HRs associated with survival analyses are from multivariable Cox proportional hazards models accounting for covariates including age, sex, ancestry and interval from diagnosis to sample collection. Full models are shown in Extended Data Tables. Covariates were selected based on statistically significant univariate associations with OS. MSI type and TMB were not found to have a statistically significant association with OS in a univariate model and therefore were not included as covariates in multivariable models. Models were either stratified by clinical stage at diagnosis or evaluated within each stage to account for known differences in OS by stage. Median OS values were calculated using the univariate Kaplan–Meier estimator. Analyses were conducted using R version 4.3.0 with the tidyverse (v2.0.0), survival (v3.5-5) and gtsummary (v1.7.2) packages.

For stage and survival analyses associated with KRAS copy number (Fig. 3), resectable patients who had experienced a recurrence by the time of sample collection were excluded (n = 95) to limit variation in disease progression within the stage. For OS comparison by binary gene alteration status (Fig. 5d), genes were considered if they had a somatic alteration prevalence of at least 3% of metastatic patients who had received standard-of-care chemotherapy with either 5-FU to gemcitabine backbones who could be evaluated for OS (n = 473 patients, 13 genes). Samples were then limited to those with high purity (30% or higher) to reduce the bias of false negatives in the WT group (n = 304 patients). Chemotherapy backbone (5-FU versus gemcitabine/nab-paclitaxel) was not significantly associated with OS and was therefore not included as a covariate (P > 0.05). Multivariable Cox proportional hazards models were constructed for each gene (Extended Data Table 5), and P values of the genomic alterations in each model were adjusted by FDR with n = 13.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.