Genome-wide association study identifies multiple susceptibility loci for pancreatic cancer.

We performed a multistage genome-wide association study including 7,683 individuals with pancreatic cancer and 14,397 controls of European descent. Four new loci reached genome-wide significance: rs6971499 at 7q32.3 (LINC-PINT, per-allele odds ratio (OR) = 0.79, 95% confidence interval (CI) 0.74-0.84, P = 3.0 × 10(-12)), rs7190458 at 16q23.1 (BCAR1/CTRB1/CTRB2, OR = 1.46, 95% CI 1.30-1.65, P = 1.1 × 10(-10)), rs9581943 at 13q12.2 (PDX1, OR = 1.15, 95% CI 1.10-1.20, P = 2.4 × 10(-9)) and rs16986825 at 22q12.1 (ZNRF3, OR = 1.18, 95% CI 1.12-1.25, P = 1.2 × 10(-8)). We identified an independent signal in exon 2 of TERT at the established region 5p15.33 (rs2736098, OR = 0.80, 95% CI 0.76-0.85, P = 9.8 × 10(-14)). We also identified a locus at 8q24.21 (rs1561927, P = 1.3 × 10(-7)) that approached genome-wide significance located 455 kb telomeric of PVT1. Our study identified multiple new susceptibility alleles for pancreatic cancer that are worthy of follow-up studies.

Pancreatic cancer is the fourth leading cause of cancer-related death in the United States and the fifth leading cause in the European Union[1,2]. Over 80% of patients have incurable disease at the time of diagnosis, and the majority live for less than 12 months[3]. Rare, moderately- to highly-penetrant mutations account for a small fraction of the familial aggregation of pancreatic cancer[4]. In two previous GWAS called PanScan I[5] and PanScan II[6], we identified common variants at four loci associated with risk of sporadic pancreatic cancer in European populations. Subsequent GWAS demonstrated five distinct susceptibility loci among individuals of Chinese descent[7] and three suggestive loci among individuals of Japanese descent[8].

In the current study (designated PanScan III), we performed a multistage GWAS of 7,683 individuals diagnosed with pancreatic cancer and 14,397 control individuals of European descent (Online Methods, Table 1, Supplementary Table 1 and Supplementary Fig. 1). In stage 1, we newly genotyped 1,582 cases from 13 prospective cohort studies, 2 case series, and 1 case-control study using Illumina OmniExpress Beadchip array. The control population included 5,203 cancer-free individuals previously genotyped using second-generation Illumina SNP microarrays (e.g. OmniExpress, Omni 1M or Omni 2.5M) and drawn from PanScan III prospective cohorts and a Spanish case-control study of bladder cancer. Of newly genotyped cases, 94% passed quality control criteria (Online Methods, Supplementary Tables 2 and 3), and 712,704 SNPs were included with a minimum call rate of 94%. In stage 2, we used the primary whole genome scan data from the reported PanScan I[5] (1,757 cases and 1,801 controls from 12 cohort studies and 1 case-control study typed on Illumina HumanHap550 array) and PanScan II[6] (1,768 cases and 1,841 controls from 8 case-control studies typed with Illumina Human 610-Quad array) studies. To address differences in typed SNPs across the arrays, we utilized the DCEG Imputation Reference Set[9] to fill in missing genotypes (Online Methods).

Table 1

Subject numbers and characteristics of pancreatic cancer cases and controls

	CasesNo. (%)	ControlsNo. (%)
No. of subjects
Stage 1	1,582	5,203
Stage 2	3,525	3,642
Replication	2,576	5,552
Full study population	7,683	14,397
Geographic region
United States	4,387 (57.1)	7,962 (55.3)
Central/Northern Europe	2,264 (29.5)	3,853 (26.8)
Southern Europe	1,032 (13.4)	2,582 (17.9)
Sex
Male	4,107 (53.5)	8,841 (61.4)
Female	3,576 (46.5)	5,556 (38.6)
Age, years
≤ 60	1,972 (25.7)	4,577 (31.8)
61 – 70	2,688 (35.0)	5,906 (41.0)
> 70	3,023 (39.3)	3,914 (27.2)
Smoking status*
Current / past	2,634 (51.6)	4,541 (51.3)
Never	1,642 (32.2)	3,186 (36.0)
Unknown	831 (16.3)	1,118 (12.6)

Smoking status was available for subjects in Stages 1 and 2.

In a meta-analysis of stages 1 and 2, we observed robust associations for the four previously identified loci in individuals of European descent: rs687289 at 9q34.2 (ABO, OR = 1.27; 95% CI = 1.20–1.35; P = 1.6×10−16), rs9543325 at 13q22.1 (KLF5/KLF12, OR = 1.23; 95% CI = 1.18–1.30; P = 4.3×10−14), rs10919791 at 1q32.1 (NR5A2, OR = 0.79; 95% CI = 0.75–0.85; P = 1.4×10−11), and rs31490 at 5p15.33 (CLPTM1L, OR = 1.20; 95% CI = 1.14–1.27 P = 2.0×10−11).

We observed two new SNPs below genome-wide significance (P<5×10−8) in the meta-analysis of stages 1 and 2, plus 11 additional promising SNPs (P<5×10−5) from distinct regions (Supplementary Table 4). These 13 SNPs were carried forward for replication (stage 3) in 2,576 cases and 5,552 controls, drawn from: (a) cases in stage 1 with DNA quantity insufficient for full GWAS, (b) cases and controls from the PANDoRA consortium[10], and (c) cases enrolled to CALGB 80303, a U.S. cooperative group clinical trial[11] (Supplementary Table 5). Additional control subjects were selected from cancer-free individuals previously genotyped using Illumina HumanHap550 array (Online Methods). Of 13 SNPs advanced to replication, nine SNPs were associated with pancreatic cancer risk (P <0.05) in the replication stage (Supplementary Table 6).

For the complete study of 7,683 cases and 14,397 controls, we applied a fixed-effect meta-analysis to the results from the three stages. Overall, six SNPs had P-values below genome-wide significance: rs2736098 at 5p15.33 (a second signal in TERT, P=9.8×10−14), rs6971499 at 7q32.3 (LINC-PINT, P=3.0×10−12), rs7190458 at 16q23.1 (BCAR1/CTRB1/CTRB2, P=1.1×10−10), rs9581943 at 13q12.2 (PDX1, P=2.4×10−9), rs16986825 at 22q12.1 (ZNRF3, P=1.2×10−8) (Table 2 and Figure 1), and rs4962153 at 9q34.2 (ADAMTS13, P=1.5×10−8). In a subsequent conditional analysis described below, rs4962153 in ADAMTS13 marked the same signal as rs687289 in ABO identified in PanScan I and II. An additional locus at 8q24.21 was close to genome-wide significance (rs1561927, P=1.3×10−7) and located in a region previously associated with multiple cancers (Table 2 and Figure 1).

Table 2

Association results for five new pancreatic cancer susceptibility loci and one suggestive locus

Chr	Nearest gene(s)a	SNP	Positionb	Minorallelec	Majorallelec	Stage	Allelic OR(95% CI)	Minor allelefrequency		Pf
								Controls	Cases
5p15.33	TERT, MIR4457, CLPTM1L	rs2736098	1,294,086	T	C	Stage 1	0.76 (0.68–0.86)	0.268	0.216	8.22 × 10−6
						Stage 2	0.82 (0.74–0.90)	0.284	0.259	2.63 × 10−5
						Replicationd	0.81 (0.74–0.89)			1.36 × 10−5
						Combinede	0.80 (0.76–0.85)			9.78 × 10−14
7q32.3	LINC-PINT	rs6971499	130,680,521	C	T	Stage 1	0.79 (0.68–0.90)	0.155	0.127	6.58 × 10−4
						Stage 2	0.81 (0.74–0.90)	0.147	0.124	4.69 × 10−5
						Replicationd	0.77 (0.69–0.86)			4.37 × 10−6
						Combinede	0.79 (0.74–0.84)			2.98 × 10−12
16q23.1	BCAR1, CTRB1, CTRB2	rs7190458	75,263,661	A	G	Stage 1	1.61 (1.32–1.96)	0.042	0.065	4.14 × 10−6
						Stage 2	1.47 (1.20–1.82)	0.039	0.049	2.17 × 10−4
						Replicationd	1.33 (1.10–1.61)			3.14 × 10−3
						Combinede	1.46 (1.30–1.65)			1.13 × 10−10
13q12.2	PDX1	rs9581943	28,493,997	A	G	Stage 1	1.23 (1.12–1.35)	0.397	0.441	1.34 × 10−5
						Stage 2	1.12 (1.05–1.20)	0.406	0.434	5.51 × 10−4
						Replicationd	1.11 (1.03–1.20)			4.80 × 10−3
						Combinede	1.15 (1.10–1.20)			2.35 × 10−9
22q12.1	ZNRF3	rs16986825	29,300,306	T	C	Stage 1	1.25 (1.10–1.42)	0.150	0.184	4.96 × 10−4
						Stage 2	1.15 (1.05–1.26)	0.149	0.168	2.11 × 10−3
						Replicationd	1.18 (1.08–1.30)			5.13 × 10−4
						Combinede	1.18 (1.12–1.25)			1.18 × 10−8
8q24.21	MIR1208, PVT1	rs1561927	129,568,078	C	T	Stage 1	0.88 (0.78–0.97)	0.269	0.251	1.59 × 10−2
						Stage 2	0.86 (0.80–0.93)	0.279	0.250	1.11 × 10−4
						Replicationd	0.89 (0.82–0.97)			6.44 × 10−3
						Combinede	0.87 (0.83–0.92)			1.30 × 10−7

Results from unconditional logistic regression of the genotypes generated in Stage 1, Stage 2, and replication (total of 7,683 individuals diagnosed with pancreatic cancer and 14,397 controls).

Closest RefSeq gene(s). Genes located within 25 kb of SNP are listed in black in order of closest gene to those further away; closest genes outside this 50 kb window are listed in grey.

Position of SNP in NCBI genome build 37 (Hg19).

Minor and major alleles.

The replication is a meta-analysis of three groups and thus, minor allele frequency (MAF) is not listed.

Number of case and control subjects in joint analysis of Stage 1, Stage 2 and replication: rs2736098 (7,199/13,121), rs6971499 (7,435/13,289), rs7190458 (7,412/13,291), rs9581943 (7,415/13,286), rs16986825 (7,413/13,196), rs1561927 (7,486/13,274).

1 d.f. score test;

Chr: chromosome and band; OR, per-allele OR for the minor allele adjusted for age, sex, geographic region and significant principal components for Stage 1; per-allele OR adjusted for age, sex, study, arm and significant principal components for Stage 2; per-allele OR adjusted for age, sex and study for Replication. Text in bold indicates the combined meta-analysis results.

Figure 1

Association results, recombination hotspots and LD plots for new pancreatic cancer susceptibility regions (a–e) and one suggestive region (f)

Top, association results of GWAS data from the stage 1 (gray diamonds), stage 2 (purple diamonds), replication (blue diamonds) and the combined data from stages 1–3 (red diamonds) plotted against −log10 P values (left y axis). Overlaid are likelihood ratio statistics (right y axis) estimating putative recombination hotspots across the region on the basis of five unique sets of 100 randomly selected control samples. Bottom, LD heat map based on r2 values from the total control populations for all SNPs included in the GWAS. The data are based on a total number of 7,683 individuals with pancreatic cancer and 14,397 controls of European descent. Shown are results for 5p15.33 (a), 7q32.3 (b), 16q23.1 (c), 13q12.2 (d), 22q12.1 (e), and 8q24.1 (f).

The SNP, rs6971499, at 7q32.3 maps to an intron in LINC-PINT, which is a p53-induced long intergenic non-protein coding RNA located in a 375 kb region between Muskelin 1 (MKLN1) and KLF14 (Supplementary Table 7 and Figure 1). Muskelin is an intracellular protein that mediates cell response to the extracellular matrix, particularly influencing cell adhesion and cytoskeleton organization[12]. KLF14 is a member of the Kruppel-like family of transcription factors, which have been implicated as tumor suppressors, including in mutant KRAS-driven tumors[13]. KLF14 has also been identified as a regulator of several metabolic phenotypes, including type 2 diabetes[14]. Notably, the previously established susceptibility locus at 13q22.1 is located in an intergenic region between KLF5 and KLF12, two other members of the Kruppel-like family of transcription factors.

The SNP, rs7190458, at 16q23.1 is a synonymous SNP residing in the last exon of BCAR1 (also known as p130Cas) and close to two chymotrypsinogen genes, CTRB1 (5 kb) and CTRB2 (23kb) (Supplementary Table 7 and Figure 1). Aberrant expression of BCAR1 has been linked with transformation and progression of multiple cancer types, and BCAR1 functions as an adaptor protein that coordinates cell cycle control, cytoskeleton organization, and cell migration[15,16]. The chymotrypsinogens are members of a family of serine proteases that are secreted by the pancreas into the gastrointestinal tract[17]. Mutations in the related genes PRSS1 (trypsin 1) and CTRC have been associated with hereditary pancreatitis[18], a known risk factor for pancreatic cancer[19]. In addition, a susceptibility locus for types 1 and 2 diabetes[20,21] is located 16 kb centromeric to rs7190458 (rs7202877, r2=0.32 in 1000G CEU data). Functional analyses indicate that this variant (rs7202877) leads to impaired pancreatic beta-cell function[22] and influences expression of CTRB1 and CTRB2 in pancreatic tissue[23].

At chromosome 13q12.2, the newly identified SNP, rs9581943, is approximately 200bp upstream of PDX1 (pancreatic and duodenal homeobox1 protein 1) and intronic to PDX1-AS1 (PDX1 antisense RNA 1), a recently identified noncoding RNA (Supplementary Table 7 and Figure 1). PDX1 is critical for early pancreatic development, plays a role in differentiation of exocrine pancreas, and regulates beta-cell function in the mature pancreas[24,25]. Mutations in PDX1 have been linked to agenesis of the pancreas[24] and maturity onset diabetes of the young (MODY)[26], a dominantly inherited disorder of non-autoimmune diabetes. Furthermore, PDX1 has been implicated in glucose-dependent regulation of insulin gene transcription[27], and GWAS have identified a SNP (rs2293941, r2=0.20 in 1000G CEU data) at the PDX1 locus associated with fasting glucose levels[28].

The signal at 22q12.1, rs16986825, maps to an intron in ZNRF3 (zinc and ring finger 3) (Supplementary Table 7 and Figure 1), encoding a cell surface transmembrane E3 ubiquitin protein ligase that is a negative regulator of the Wnt signaling pathway[29]. Additionally, CHEK2 is located 162 kb centromeric to the marker SNP and encodes a cell-cycle checkpoint kinase that cooperates with p53, BRCA1 and ATM in response to DNA damage[30]. Alterations in CHEK2 have been implicated in susceptibility to several cancer types[31].

We performed conditional analyses to assess whether the newly identified SNPs at 5p15.33 (CLPTM1L/TERT) and 9q34.2 (ABO/ADAMTS13) were independent from those identified previously. After conditioning on the reported SNP within intron 13 of CLPTM1L, the newly identified synonymous SNP within the second exon of TERT (rs2736098) remained statistically significant (P=2.4×10−3) (Table 3). Two strong recombination hotspots lie between the established and new SNPs in 1000G CEU data (likelihood ratios, LR of 27.1 and 261.0)[32], and the two SNPs are in modest linkage disequilibrium (LD; r2=0.22 in 1000G CEU data) (Figure 1).

Table 3

Conditional analyses of SNPs at chromosomes 5p15.33 and 9q34.2

New SNP	Chr	Positiona	Gene	ORb	Pb	Conditional ORc	Conditional Pc	Established SNP	r2,d	ORe	Pe	Conditonal ORf	Conditional Pf
rs2736098	5p15.33	1,294,086	TERT, CLPTM1L	0.80 (0.74–0.86)	1.42 × 10−9	0.88 (0.80–0.95)	2.44 × 10−3	rs401681 / rs31490	0.22	0.83 (0.79–0.88)	1.97 × 10−11	0.86 (0.81–0.92)	7.55 × 10−6
rs4962153	9q34.2	136,323,754	ADAMTS13, ABO	1.20 (1.10–1.30)	1.97 × 10−5	1.05 (0.96–1.16)	0.28	rs505902 / rs687289	0.17	1.27 (1.20–1.35)	1.64 × 10−16	1.25 (1.18–1.33)	6.28 × 10−12

Position of SNP in NCBI genome build 37 (Hg19).

Per-allele ORs for the minor allele and P values for the new SNP from the unconditional meta-analysis of Stage 1 and Stage 2.

Per-allele ORs for the minor allele and P values for the new SNP from the conditional meta-analysis.

r2 LD values between the new SNP and the established SNP at the locus based on 1000 Genomes Project CEU data.

Per-allele ORs for the minor allele and P values for the established SNP from the unconditional meta-analysis of Stage 1 and Stage 2.

Per-allele ORs for the minor allele and P values for the established SNP from the conditional meta-analysis.

Chr: chromosome.

TERT encodes the catalytic subunit of telomerase reverse transcriptase, a component of the ribonucleoprotein complex that maintains integrity of chromosome ends. Inherited mutations affecting TERT underlie cases of dyskeratosis congenita, aplastic anemia, acute myeloid leukemia, familial melanoma, and pulmonary fibrosis[33,34]. CLPTM1L encodes the cleft lip and palate associated transmembrane 1 like protein involved in mediating apoptosis, aneuploidy, cisplatin resistance, and RAS-mediated malignant transformation[35,36]. Variants across the TERT/CLPTM1L region have previously been associated with risk of multiple cancers. Furthermore, independent signals within TERT and CLPTM1L have been identified for bladder cancer[37], CLL[38], and lung cancer[37,39], and fine-mapping studies have identified at least four independent signals across the TERT/CLPTM1L region associated with cancer[40,41]. The new SNP identified in PanScan III (rs2736098) is located in a region of LD spanning ~4kb from the promoter region to exon 2 of TERT. This SNP and several correlated SNPs have been associated with telomere length in white blood cells and TERT promoter activity[37,40,41]. The minor allele of rs2736098 that is associated with a lower risk of pancreatic cancer in PanScan was associated with longer telomeres and lower risk of breast cancer[40]. Although further characterization of this region will be necessary, the new SNP in exon 2 of TERT appears to mark an independent risk locus for pancreatic cancer.

After conditioning on the established SNP at 9q34.2 in the first intron of ABO, the SNP rs4962153 in ADAMTS13 identified in PanScan III was not statistically significant (P=0.28), indicating that these two SNPs point to the same susceptibility haplotype (Table 3).

A promising risk locus was identified at 8q24.21 (rs1561927; P = 1.3×10−7) in a nongenic region between PVT1 and LINC00977 (Supplementary Table 7 and Figure 1). 8q24.21 is known to contain multiple cancer susceptibility loci that span over 2Mb[42,43]. The promising pancreatic cancer SNP is in LD with a SNP associated with ovarian cancer risk (rs10088218, r2=0.37 in 1000G CEU data, 24kb upstream)[44], and the closest genes are centromeric to rs1561927: MIR1208 (406 kb), PVT1 (455 kb), and MYC (814 kb). Several 8q24.21 risk loci have been shown to interact with MYC or PVT1 promoters through long range chromosomal interaction, and allele-specific effects on the expression of both genes have been reported[42,45,46]. However, these loci are located more than 1 Mb upstream of rs1561927 on 8q24.21 (r2<0.03 in 1000G CEU data).

In stratified analyses, no statistically significant heterogeneity was noted by geographic region or smoking status (Supplementary Tables 8 and 9). In a preliminary analysis that included 173 cases and 430 controls of Asian ancestry (Supplementary Table 10), we examined the susceptibility loci identified in individuals of European descent[5,6] (Table 2). We also assessed previously published pancreatic cancer risk loci from individuals of Chinese[7] and Japanese[8] ancestry, noting no loci and one locus, respectively, as nominally statistically significant in PanScan (Supplementary Table 11).

To pursue the first steps towards understanding the functional underpinnings of the newly identified risk alleles, we conducted bioinformatic analyses using HaploReg[47] (Supplementary Table 12). We also evaluated expression quantitative trait locus (eQTL) effects[48-50] (Supplementary Table 12). Cis-eQTLs were noted on chr16q23.1 in peripheral blood (CFDP1), chr13q12.2 in skin and liver (POMP), chr22q12.1 in liver and peripheral blood (CCDC117) and peripheral blood (XBP1), and chr8q24.21 in adipose tissue (PVT1). XBP1 at chr22q12.1 regulates pancreatic beta-cell function with effects on systemic glucose control[51] and modulates acinar cell homeostasis[25]. In gene set enrichment analysis[52] of genes within 100 kb of the 10 index SNPs identified in PanScan, the only statistically significant pathway was maturity onset diabetes of the young (P=3.3×10−4). Understanding the functional consequences of pancreatic cancer susceptibility variants will require further laboratory investigation.

In a linear-mixed model analysis[53] (Online Methods), we estimated that the heritability for pancreatic cancer due to common SNPs present on GWAS arrays was 13% (95% CI, 4–22%). Furthermore, we estimated that the nine loci identified in individuals of European ancestry account for approximately 9% of total heritability tagged by common SNPs. We also evaluated the cumulative association with pancreatic cancer of risk alleles at susceptibility loci identified in individuals of European descent. Compared to individuals with the most prevalent number of risk alleles in controls (n=10), those with ≤6 risk alleles had an OR of 0.55 (95% CI, 0.44–0.68) and those with ≥14 risk alleles had an OR of 2.24 (95% CI, 1.80–2.80) for pancreatic cancer (Supplementary Figure 2).

In conclusion, our multistage GWAS revealed new loci associated with pancreatic cancer risk, as well as promising loci that merit follow-up. Several of the new loci harbor plausible candidate genes implicated in pancreas development, pancreatic beta-cell function, and predisposition to diabetes. Further investigation is warranted to understand the biological underpinnings of these common pancreatic cancer susceptibility alleles.

ONLINE METHODS

Stage 1: GWAS for PanScan III

We conducted a GWAS of pancreatic cancer using case and control subjects from 17 studies (Supplementary Table 1). Pancreatic cancer cases included individuals newly identified from nine cohort studies that participated in PanScan I[5], as well as those from five new cohort studies, two new case series, and one new case-control study. The new cohort studies included the Agricultural Health Study (AHS)[54], Melbourne Collaborative Cohort Study (MCCS)[55], Multiethnic Cohort Study (MEC)[56], Selenium and Vitamin E Cancer Prevention Trial (SELECT)[57], and Vitamins and Lifestyle Study (VITAL)[58]. The new case-based studies were the Gastrointestinal Cancer Clinic of Dana-Farber Cancer Institute (DFCI-GCC), Spanish Pancreatic Cancer Study PANKRAS-II[59], and PANDoRA-Heidelberg pancreatic cancer case-control study[10]. Cases were defined as those individuals having primary adenocarcinoma of the exocrine pancreas (ICD-O-3 code C250–C259). Those with non-exocrine pancreatic tumors (histology types 8150, 8151, 8153, 8155 and 8240) were excluded. Each participating study obtained informed consent from study participants, approval from its institutional review board (IRB) for this study, and IRB certification permitting data sharing in accordance with the NIH Policy for Sharing of Data Obtained in NIH-Supported or NIH-Conducted Genome-Wide Association Studies.

All samples from pancreatic cancer cases with sufficient DNA (n=1,894) were genotyped on the Illumina OmniExpress chip at the NCI Cancer Genomic Research Laboratory (CGR) (Supplementary Table 2). Genotypes were called using the Illumina GenomeStudio software. Genotype clusters for new cases were estimated using samples with a completion rate of 98% to optimize accuracy. Genotypes for all samples, including those initially excluded, were subsequently called based on the optimized cluster file. Extensive quality-control metrics were applied to the data: SNPs with a call rate <94% or Hardy-Weinberg Proportion p value <1×10−7 were excluded (n=18,765); samples with a call rate <94% (n=78), mean heterozygosity <26% or >33% (n=2) based on autosomal SNPs or gender discordance (>5% heterozygosity based on the X chromosome SNPs for males or <20% heterozygosity based on the X chromosome SNPs for females, n=5) were excluded. Unexpected duplicates (>99.9% concordance, n=3) and first-degree relatives (n=2, on the basis of identity-by-descent sharing with Pi-hat >0.40) were removed. Quality-control duplicate samples in PanScan III (n=38 pairs) showed >99.9% genotype concordance. Duplicates with PanScan I or II were removed (>99.9% concordance, n=21). Ancestry was assessed using the Genotyping Library and Utilities (GLU) struct.admix module. Participants with <80% European ancestry (n=199) were excluded for the primary analysis of individuals of European ancestry (Supplementary Fig. 3). After exclusions, 1,582 cases of European ancestry were available for analysis (Supplementary Tables 2 and 3).

Controls of ≥80% European ancestry were drawn from ten of the studies included in PanScan III (Alpha-Tocopherol, Beta-Carotene Lung Cancer Prevention Study (ATBC), American Cancer Society Cancer Prevention Study-II Cohort (CPS-II), European Prospective Investigation into Cancer and Nutrition (EPIC), Health Professionals Follow-Up Study (HPFS), Melbourne Collaborative Cohort Study (MCCS), Multiethnic Cohort Study (MEC), Nurses' Health Study (NHS), Prostate, Lung, Colorectal and Ovarian Cancer Screening Trial (PLCO), Spanish Pancreatic Cancer Study (SPCS, PANKRAS-II and Spanish Bladder Cancer SBC/EPICURO studies)[60] and Women's Health Initiative (WHI)). These controls had no history of cancer, were not included in PanScan I or PanScan II, and had been previously genotyped at CGR on the Illumina OmniExpress, Omni1M or Omni2.5M arrays, with extensive quality-control as previously described[9,61-64]. In total, 5,203 controls were included in the analysis (Supplementary Tables 2 and 3). A total of 608,202 SNPs with overall completion rate >80% in both cases and controls were advanced to the association analysis. To evaluate population substructure, a principal components analysis was performed using the struct.pca module of GLU, version 1.0, which is similar to EIGENSTRAT[65]. Plots of the first six principal components are shown in Supplementary Figure 4. The estimated inflation of the test statistic, λ, was 1.02[66]; a Quantile-quantile (QQ) plot is shown in Supplementary Figure 5.

Association analysis was performed assuming a log-additive genetic model and adjusting for age, sex, geographic region and 6 significant eigenvectors (i.e. EVs that were nominally significant in a baseline risk model adjusting for age, sex, and geographic region). Geographic region was defined as REGION_US (United States): AgHealth, CPS-II, DFCI, HPFS, MEC, NHS, NYU-WHS, PHS, PLCO, SELECT, VITAL, WHI; REGION_CNE (Central and Northern Europe): ATBC, EPIC, PANDoRA-Heidelberg, MCCS (Melbourne); and REGION_SE (Southern Europe): SBCS (Spain controls), PANKRAS-II (cases). All data analyses and management were conducted using GLU.

Stage 2: PanScan I and II data

The second stage involved the primary whole genome scan data from the previously reported PanScan I[5] and PanScan II[6] studies. PanScan I and PanScan II were genotyped on the Illumina HumanHap550 Infinium II and the Human 610-Quad chips, respectively, whereas PanScan III was genotyped on the OmniExpress chip. As the number of overlapping SNPs between the three chips is moderate (~300K), imputation of missing genotypes was performed using phased haplotypes from the DCEG Reference Set and IMPUTE2[9,67]. The DCEG reference set is well-designed for “filling in” missing genotypes across chip designs in PanScan since it is based on several of the same studies included in PanScan and the imputation accuracy is improved over 1,000 Genomes and HapMap data[9]. Imputed SNPs with low minor allele frequencies (MAF <0.01) or low-quality scores (IMPUTE2 information score <0.3) were removed prior to association analysis. The same quality thresholds as described above for stage 1 were applied for stage 2. Final numbers of cases and controls included in stage 2 were 1,757 cases and 1,801 controls from PanScan I and 1,768 cases and 1,841 controls from PanScan II.

To combine data from PanScan I, II, and III, meta-analyses were performed using the fixed-effects inverse-variance method based on the β estimates and standard errors. No heterogeneity was observed across stages 1 and 2 for the SNPs identified as GWAS significant or suggestive in the full study (P heterogeneity ≥0.11; Supplementary Table 4). Manhattan plot for the results of the meta-analysis of stage 1 and stage 2 is shown in Supplementary Figure 6.

Association analysis was also performed in 173 cases and 430 controls of Asian ancestry from the Shanghai Men’s and Women’s Health Study (SMWHS) (Supplementary Table 10). This analysis included case and control subjects from stages 1 and 2 of PanScan III and previously genotyped control subjects from SMWHS[68]. Quality control and association analysis were performed as described above for European ancestry subjects.

Stage 3: Replication studies

Thirteen SNPs (P-value threshold of <5×10−5) were taken forward for de novo replication in an additional 2,576 cases and 5,552 controls. The replication samples were analyzed individually as three groups: (A) CGR: pancreatic cancer case and control subjects from CARET[69] plus samples from cases that did not have sufficient DNA for full GWAS and control subjects previously genotyped at CGR; (B) PANDoRA: case and control subjects from the PANDoRA pancreatic cancer case-control consortium[10] (no overlap with the PANDoRA-Heidelberg cases genotyped in stage 1); and (C) CALBG/Alliance 80303: cases from a randomized clinical trial of gemcitabine plus placebo versus gemcitabine plus bevacizumab[11], and control subjects previously genotyped at CGR (Supplementary Table 5).

Genotyping for cases in group A was performed using custom TaqMan genotyping assays (Applied Biosystems) at CGR. Genotyping for cases and controls from PANDoRA (group B) was performed in the same manner but at the German Cancer Research Center (DKFZ) in Heidelberg, Germany. Quality-control duplicate samples in the replication (CGR, n=20 pairs; PANDoRA, n=512 pairs) showed >99.9% genotype concordance. Patients enrolled on CALGB/Alliance 80303 (group C) were previously genotyped using the Illumina HumanHap550v3 Genotyping BeadChip array[11]. Control subjects from PLCO previously genotyped at CGR using the Illumina HumanHap550v3 Genotyping BeadChip array[70] were used for groups (A) and (C) (Supplementary Table 5), and did not overlap with control subjects included in PanScan I, II or III. CALBG/Alliance 80303 and control genotypes were imputed to OmniExpress SNP content in the same manner as described above for stage 1. Quality control thresholds and exclusions for sample and loci in the replication are listed in Supplementary Table 5B. Association results for the replication studies were adjusted for age, sex and study, and a meta-analysis of the three replication groups was performed using the fixed-effects inverse-variance method based on the β estimates and standard errors (Supplementary Table 6). This was followed by a meta-analysis of stages 1, 2 and replication for the 13 SNPs using the same fixed-effects inverse-variance method.

Technical validation

A comparison of the genotyping calls from the imputation of PanScan I and II into OmniExpress array contents and confirmatory TaqMan assays (n= 511 samples from PanScan I and II) yielded an r2 of 0.74, 0.96, 0.56, 0.99, 0.98 and 1.00 for rs2736098, rs6971499, rs7190458, rs9581943, rs16986825 and rs1561927, respectively.

Estimate of recombination hotspots

To identify recombination hotspots, we used SequenceLDhot[32], a program that uses the approximate marginal likelihood method[71] and calculates likelihood-ratio statistics at a set of possible hotspots. We tested five unique sets of 100 control samples. The PHASE v2.1 program was used to calculate background recombination rates[72,73], and LD heat maps were visualized using the snp.plotter program[74]. For estimation of recombination hotspots between loci in TERT and CLPTM1L on chr5p15.33, we used the 1000G (version 3) CEU data.

Heritability analysis

To estimate heritability explained by common SNPs present on GWAS arrays on the liability scale (lifetime disease risk of 0.015), we used GCTA[53,75] on a set of LD-pruned SNPs (r2<0.5) that passed the following stringent quality control thresholds: MAF>1%, SNP missing rates <5%, subject missing rate <1%, and HWE P-values >10−4. Non-autosomal SNPs and pairs of subjects with genetic relatedness >5% were removed. These analyses were run separately in PanScan I, II and III, adjusting for age, sex, study (or geographic region in PanScan III) and the significant principal components in each study. PanScan III analyses were restricted to participating studies that contributed both cases and controls. PanScan I, II, and III results were combined via meta-analysis. We repeated the analyses restricted to the genome-wide significant SNPs in individuals of European ancestry to estimate the proportion of heritability tagged by these nine SNPs.

Further follow-up analyses

We constructed a genetic risk score for pancreatic cancer, incorporating the susceptibility loci identified in PanScan I, II, and III. For this analysis, subjects could possess zero to 20 risk alleles, based on their genotypes at the 10 identified loci. Odds ratios were calculated using multivariable-adjusted unconditional logistic regression with meta-analysis to combine data from stages 1 and 2, as done in the analyses of individual SNPs. Replication samples were not genotyped for the four susceptibility loci identified in PanScan I and II, and therefore, these subjects could not be included in the risk score analysis. Subjects with missing genotypes for one or more of the 10 SNPs (n=898) were assigned the most common genotype at that SNP among cases or controls. In sensitivity analyses, results were unchanged if these subjects were excluded. Using 1000 Genomes CEU data, we identified SNPs with r2 >0.7 with our lead SNP. We used HaploReg v2[47], a tool for exploring noncoding functional annotation using ENCODE data, to evaluate the genome surrounding our SNPs (Supplementary Table 12). In addition, we evaluated cis associations between all new and promising SNPs discovered in this study and the expression of nearby genes in skin biopsies, adipose biopsies and non-transformed peripheral blood samples from subjects of European descent from publically available data sets[48,50] (Supplementary Table 12). Gene set enrichment analysis was also performed for genes in pancreatic cancer risk loci identified in subjects of European descent (in a window of 100 kb centered on the most significant SNP in each locus) based on KEGG (Kyoto Encyclopedia of Genes and Genomes) annotations using GeneCodis3 with reporting of the corrected hypergeometric P-value[52].