Introduction Warfarin is the most widely prescribed anticoagulant for reducing thromboembolic events that often give rise to stroke, deep vein thrombosis, pulmonary embolism or serious coronary malfunctions [1]. A combination of genetic and non-genetic factors cause Caucasians to exhibit 20-fold interindividual variation in required warfarin dose needed to achieve the usual therapeutic level of anticoagulation as measured by the prothrombin international normalized ratio or INR [2]–[4]. Thus, in the absence of information (genotypic, clinical, etc.) for predicting each patient's required warfarin dose, initial prescribed doses may be too low (risking thrombosis) or too high (risking over-anticoagulation and severe bleeding). Warfarin's risk of serious side effects, narrow therapeutic range, and wide interindividual variation in warfarin dose have focused attention on the need to better predict dose in the initial stage(s) of treatment. We and others have shown that the warfarin drug target VKORC1 (vitamin K epoxide reductase complex, subunit 1) contains common polymorphisms that account for a major portion (∼30%) of the variance in required warfarin dose [5],[6], and we have recently evaluated ∼1500 Swedish patients of the Warfarin Genetics (WARG) cohort in the largest study to date showing likely patient benefit from genetic forecasting of dose [3]. The study confirmed that SNPs in VKORC1 and in the warfarin-metabolizing gene CYP2C9 (cytochrome P450, family 2, subfamily C, polypeptide 9) predict ∼40% of dose variance while non-genetic factors (age, sex, etc.) jointly account for another ∼15%. The robust and now widely replicated associations of warfarin dose with VKORC1 and CYP2C9 have provided one of the most successful applications of pharmacogenetics to date [7] and offer promise for genetic predication of required dose in a clinical setting [3]. Knowledge of major predictors of warfarin dose also impacts the methodology for finding further dose-related genes. In early candidate gene work with a small sample of 201 patients [8], we noted that univariate regression (with tested SNP as the only dose predictor) could statistically detect warfarin association with VKORC1 and with one of two non-synonymous CYP2C9 SNPs (*3) known to influence warfarin dose (Table 1 in [8]). However, a second non-synonymous CYP2C9 SNP (*2) with known but weaker influence on warfarin dose was not detected by univariate regression, but *2 was statistically significant in multivariate regression adjusted for the other known genetic and non-genetic predictors of dose (Table 3 in [8]). These empirical results in a small warfarin sample provided a signpost underscoring the potential importance of multivariate regression for detecting weak effects in studies now searching for additional warfarin genes across the genome. 10.1371/journal.pgen.1000433.t001 Table 1 Association (p-value) of SNPs tested by univariate regression or multiple regression with progressive addition of known dose predictorsa. Predictors in regression analysis Tested SNP VKORC1 CYP2C9*3 CYP2C9*2 CYP4F2 Distribution of all SNPs rs9923231 rs1057910 rs1799853 rs2108622 None 5.4E-78 4.5E-17 8.8E-13 1.6E-05 Figure 1A Age, sex 7.3E-97 1.2E-24 2.4E-14 4.8E-06 – Age, sex, VKORC1 – 3.8E-43 1.0E-15 4.6E-07 – Age, sex, VKORC1, CYP2C9*3 – – 1.4E-26 8.3E-08 – Age, Sex, VKORC1, CYP2C9*3 and *2 – – – 8.3E-10 Figure 1B a Linear regression on warfarin dose was calculated for the 1,053 GWAS subjects. A genome-wide association study (GWAS) enables a systematic search of the entire genome for genetic factors that cause any inherited trait. This method has successfully identified susceptibility loci for common diseases [9], and is beginning to be applied to pharmacogenomics. A recent warfarin GWAS in 181 patients did not detect other genetic factors with major effects on warfarin dose beyond VKORC1 [10] but was underpowered for identifying loci with a moderate contribution. We have now genotyped 325,997 SNPs in 1053 patients of the WARG cohort and here report the first GWAS that is sufficiently powered to detect additional genetic factors that may only modestly influence warfarin dose. Results Figure 1A and the first line of Table 1 summarize results of testing 325,997 GWAS SNPs for association with warfarin dose by univariate regression. The strongest associations were at multiple SNPs in and near VKORC1 (Figure 1A) with the lowest p-value given by rs9923231 (P = 5.4×10−78). In prior fine-mapping of the VKORC1 locus [8], we identified rs9923231 as one of three SNPs located in introns or immediately flanking VKORC1 that exhibit almost perfectly concordant genotypes yielding pairwise linkage disequilibrium (LD) r 2≈1 and which define the warfarin-sensitive A-T-T haplotype at rs9923231-rs9934438-rs2359612 (see also [11]). These highly concordant SNPs were the best predictors of warfarin dose in our previous study and in this GWAS analysis (p T, 0.402)a −0.96 (−1.03, −0.89) 0.283 1.6E-122 −0.99 (−1.09, −0.88) 0.284 5.0E-62 −0.97 (−1.02, −0.91) 0.283 2.7E-181 CYP2C9*3 rs1057910 (Ile359Leu, 0.070)a −1.13 (−1.26, −1.00) 0.075 2.6E-55 −1.08 (−1.27, −0.89) 0.089 2.3E-26 −1.11 (−1.22, −1.00) 0.080 2.6E-79 CYP2C9*2 rs1799853 (Arg144Cys, 0.109)a −0.63 (−0.74, −0.52) 0.048 1.7E-28 −0.40 (−0.55, −0.24) 0.023 5.5E-07 −0.54 (−0.63, −0.45) 0.038 1.1E-31 CYP4F2 rs2108622 (Val433Met, 0.240)a 0.25 (0.17, 0.33) 0.016 8.3E-10 0.16 (0.05, 0.27) 0.005 c0.0029 0.21 (0.14, 0.27) 0.011 3.3E-10 Age −0.04 (−0.04, −0.03) 0.170 1.9E-63 −0.03 (−0.04, −0.03) 0.129 1.7E-31 −0.04 (−0.04, −0.03) 0.155 1.2E-111 Sex (male) 0.35 (0.25, 0.45) 0.017 7.6E-12 0.25 (0.10, 0.40) 0.009 0.001 0.30 (0.22, 0.38) 0.013 1.6E-12 a In parenthesis are major/minor allele, and minor allele frequency. b Effect of individual predictor on dose is indicated by regression coefficient and 95% confidence interval, proportion of explained variance (R 2) and P-value. c Association in same direction as GWAS was assessed by a one-tailed test. To increase the power of our multivariate regression model and possibly detect additional weak effects, we added CYP4F2 (rs2108622) to the model as a predictor and conducted further analyses. First, we retested the GWAS SNPs, but no new SNPs reached genome-wide significance and there was also no apparent excess of SNPs at lower significance thresholds (Figure S1). We also tested warfarin association with haplotypes and with ungenotyped SNPs imputed at 2.2 million HapMap SNPs, but no haplotype or imputed SNP approached genome-wide significance in a genomic region not containing VKORC1, CYP2C9 or CYP4F2. To explore whether copy number variations (CNVs) detectable by the HumanCNV370 array might influence warfarin dose, we used rigorous quality control and retained 879 samples calling 2530 CNVs (see Materials and Methods). None of the CNV loci were significantly associated with dose after correction for multiple testing (lowest CNV p-value was 1.1×10−4 which exceeds 0.05/2530≈2.0×10−5). We note that probe density in many of the detected CNVs is not optimal for conducting association analyses and these results should therefore be viewed as preliminary. Finally, after excluding SNPs near VKORC1, CYP2C9 and CYP4F2, we identified 40 other loci containing one or more GWAS SNPs with p-values below 2.0×10−4 and we genotyped 40 SNPs representing these loci in a follow-up sample of 588 Swedish warfarin patients. However none of the 40 loci replicated for association with warfarin dose, the lowest p-value being 0.04 which is not significant after correction for 40 tests (Table S1). Having not found evidence for any additional genetic modulators of dose, we examined the entire data set (GWAS plus followup samples) for evidence of statistical interaction between pairs of the established dose predictors (VKORC1, CYP2C9, CYP4F2, age, sex). None of the pairs exhibited statistically significant interaction after p-values were corrected for the 15 interaction tests (Table S2). We also performed a GWAS for a secondary trait (“over-anticoagulation”) which we previously found was associated with VKORC1 and CYP2C9 in a candidate gene study [3]. By titrating warfarin dose, physicians attempt to achieve a target level of anticoagulation determined by a reading of 2.0 to 3.0 for the prothrombin international normalized ratio (INR), which is the ratio of time required for a patient's blood to coagulate relative to that of a reference sample. However over-anticoagulation (defined as an INR above 4.0) sometimes occurs and, using Cox regression, our GWAS tested for SNP association with the occurrence of over-anticoagulation in patients during the first 5 weeks of treatment (see Materials and Methods: Association testing of SNPs and haplotypes). We observed genome-wide significant association (p 0.8 or r 2>0.5 for ∼60% or ∼80% respectively of common SNPs [18] and r 2>0.9 for ∼90% of non-synonymous common SNPs [19] in HapMap Caucasians). We therefore conclude that our GWAS probably detected most common SNP variants explaining 1.5% or more of the warfarin dose variance, but may have failed to detect rarer variants that could individually explain up to 5% of dose variance. We further note that the HumanCNV370 array used in this study does not have the required marker complement to undertake a comprehensive GWAS of common CNVs. As noted in the Introduction, the widely replicated warfarin dose associations with VKORC1 and CYP2C9 represent one of the most successful applications of pharmacogenetics to date. Our study together with that of Caldwell et al. [12] now also clearly demonstrates that CYP4F2 (rs2108622) is a third gene that influences warfarin dose, but our GWAS and statistical analysis also implies that additional common SNP variants that influence dose may not exist in Caucasian populations. However, Caucasians might carry common variants with effects smaller than CYP4F2 or rare variants whose effects are substantially larger than the ∼1% of dose variance explained by CYP4F2. Furthermore, other unidentified genes may influence warfarin dose in other ethnicities such as Asians or Africans, and some rare dose-altering variants in known genes such as VKORC1 may exist in only a subset of populations of European descent [20]. Hence, future research could address ethnic differences in the genetic variants that influence warfarin dose as well as subtle intra-ethnic differences and admixture that may exist in European or other populations. In a recent study [3], we highlighted the potential benefit of pre-treatment forecasting of required warfarin dose based on patient genotypes at VKORC1 and CYP2C9 together with non-genetic predictors of dose. Indeed, in August 2007, the US Food and Drug Administration (FDA) updated warfarin labeling to recommend initiating lower warfarin dose in some patients based on VKORC1 and CYP2C9 genotypes. However this recommendation is not a requirement due to a lack of large trials demonstrating warfarin patient benefit from dose forecasting (though two small trials [21],[22] do support such benefit; see also [23]–[27] for reviews and other trials). The results of our GWAS provide further impetus for conducting large-scale dose-forecasting trials by identifying CYP4F2 as a third genetic predictor of dose and also by showing that additional major genetic predictors may not exist in Caucasians or may not emerge in the near-term. Hence, large-scale trials of patient benefit from dose forecasting based on VKORC1 and CYP2C9 (with possible inclusion of CYP4F2 as a minor predictor) are likely to provide state-of-the-art clinical benchmarks for warfarin use during the foreseeable future. Materials and Methods Subjects and Clinical Data The study subjects were 1053 Swedish patients collected for the WARG study [3] (http://www.druggene.org/). This is a multi-centre study of warfarin bleeding complications and response to warfarin treatment [28]. Anticoagulant response is measured by INR, which is the ratio of the time required for a patient's blood to coagulate relative to that of a reference sample. By titrating warfarin dose, physicians aim for a therapeutic INR reading between 2.0 and 3.0; thus the primary quantitative outcome for the GWAS was the mean warfarin dose (mg/week) given to a patient during a minimum series of three consecutive INR measurements between 2 and 3 [3]. As a secondary GWAS outcome, we also catalogued each patient for the occurrence or non-occurrence of “over-anticoagulation” during the first 5 weeks of treatment (defined as an INR reading above 4.0) and tested for genetic association which adjusted for the treatment day (1 to 35) of the over-anticoagulation event (see “Association testing” below). The clinical data collected by the WARG protocol included gender and age since each is a known non-genetic predictor of warfarin dose but did not include bodyweight and dietary information (e.g. vitamin K intake). Regression analysis of prescribed medication which can potentiate or inhibit warfarin action was not a statistically significant predictor of warfarin dose in the 1053 WARG GWAS subjects and hence was not included as a predictor variable in the multivariate regression analyses. The WARG study samples were previously described elsewhere [3],[4],[28],[29] as were the Uppsala followup samples [8]. The WARG and Uppsala studies received ethical approval from the Ethics Committee of the Karolinska Institute and the Research Ethics Committee at Uppsala University, respectively. Genotyping of SNPs and Sample Quality Control From approximately 1500 WARG samples [3] examined for non-degradation and appropriate concentration of DNA (∼50 ng/µl), we randomly selected 1208 subjects for genotyping SNPs and CNV probes using the HumanCNV370 BeadChip array (Illumina). We excluded SNPs with MAF below 1%, call rate below 95%, or if call rate fell below 99% when MAF was below 5%. SNPs that departed from Hardy-Weinberg equilibrium (P 4.0) during treatment days 1–35, we performed Cox proportional hazard regression on survival time (day of over-anticoagulation) using the survival library of R software. The GWAS data set of 1053 WARG subjects contained 215 subjects whose INR exceeded 4.0 during days 1–35 while the entire dataset of 1489 WARG subjects contained 312 such subjects. Association Testing of CNVs For each CNV locus, association was tested with square root of warfarin dose by multivariate regression analysis in which subject copy number intensity was the CNV predictor of dose. This analysis differs from association testing with SNP genotypes since the two CNV alleles on homologous chromosomes generate one copy number intensity rather than a separate allele for each chromosome. As a QC strategy, we determined each subject's rank in the dataset for copy number intensity at each CNV on chromosome 17. This enabled us to differentiate the majority of subjects (whose individual distribution of ranks were approximately random and uniform) from 174 obvious outliers due to poor quality DNA (whose ranking distributions were “U-shaped” since their intensities strongly clustered at both high and low ranks). These 174 subjects were excluded from the primary CNV association analysis (with further confirmation of lower quality DNA for these subjects being their rough correspondence to the subjects with lower (<99%) SNP call rates). However, we also crosschecked the primary CNV analysis by conducting association testing on the dataset without excluding the 174 subjects and found no statistically significant association with warfarin dose at any CNV whether the dataset excluded or included the subjects. Association testing of the CNVs was executed using R software [37]. Replication of CYP4F2 For the replication of CYP4F2 rs2108622, we genotyped a panel of 588 warfarin patients consisting of 410 subjects from the WARG cohort [3] and 178 from the Uppsala cohort [38]. Table 2 shows regression on this pooled sample of 588 subjects. Separate results for each of the two panels are given in Table S4. Follow-Up of Moderately Significant SNPs To possibly identify SNPs with genuine but weak associations to warfarin dose, we excluded VKORC1, CYP2C9, CYP4F2 and identified 40 other GWAS loci for follow-up genotyping exhibiting multivariate regression p-values below 0.0002, and selected 40 SNPs representing these loci for genotyping. Only genotyped (not imputed) SNPs were chosen for follow-up. We genotyped the same 558 patients as in the CYP4F2 replication using the iPLEX MassARRAY. Power Calculation Suppose multiple regression analysis is conducted in N total samples by testing a SNP with coefficient of determination (i.e., explained variance) R 2 test after adjustment for known predictors whose total of coefficient of determination is R 2 knw . The probability (power) to detect the tested SNP at a significance level α equals: (1) where F′(1, N–2, θ 2) is the probability density function for an F distribution with 1, N-2 degrees of freedom and non-centrality parameter θ 2 (Section 28.28 in [39], Example 8.4 in [40]). Here the constant c satisfies the equation: (2) where α is the significance level, and F(1, N–2) is the probability density function for a F-distribution of degree of freedom one and N–2. How Much Does Linkage Disequilbrium (LD) Attenuate Association with a Quantitative Trait? Association with a quantitative trait (QT) becomes weaker for a marker SNP in LD with a SNP that alters the QT, and hence the association becomes more difficult to detect at the marker than at the QT-altering SNP. Here we quantify the LD attenuation for a QT when testing for association by linear regression (which includes the Cochran-Armitage trend test for dichotomous traits), and we obtain a result analogous to the LD attenuation for the Pearson Chi-square test for allelic association to dichotomous traits as in cases and controls [17]. If a causative QT-altering SNP has a coefficient of determination (i.e., explained variance) and is in pairwise LD of r 2 with a marker SNP, then the coefficient of determination for the marker SNP ( ) is approximated by: (3) In other words, when testing a marker, the proportion of explained variance decreases by a factor of r 2. To begin the proof of Equation 3, let the QT be represented by the random variable “q”, and let “m” and “x” be SNP genotypes (coded 0, 1, or 2) representing the marker and causative (QT-altering) SNP, respectively. The coefficients of determination are equal to the square of two correlation coefficients (denoted by “Corr”) measuring the correlation of m or x with q: (4) (5) Also note that correlation between genotypes at the marker and causative SNP is given by another correlation coefficient: (6) It is well known that the partial correlation coefficient of m and q conditioned on x is (equation 16.20, p. 649 in [41]): (7) However, conditional on genotype at the causative SNP, marker m and the QT q would be uncorrelated (assuming m is not in LD with a second causative polymorphism) and thus the numerator of Equation 7 would be zero implying that: (8) Based on prior work [42]–[44], we show in Text S1 that the squares of the genotypic correlation coefficient and LD correlation coefficient r 2 are approximately equal if the population is in Hardy-Weinberg equilibrium. Therefore, substituting r 2 for in Equation 8 gives Equation 3. Supporting Information Figure S1 QQ plot for association of each GWAS SNP with warfarin dose. SNPs were tested for association with warfarin by regression analysis that adjusted for age, sex, and genotype at VKORC1, CYP2C9*2 and *3, CYP4F2. The QQ plot omits SNPs in loci already known to be associated with warfarin dose (VKORC1, CYP2C9, CYP4F2). The excess of SNPs with small p-values is minor: whereas 65.4 SNPs with p<0.0002 are expected, 70 were observed (1.069 times inflated). (0.20 MB PDF) Click here for additional data file. Figure S2 Manhattan plot of GWAS results of testing for association with warfarin-induced over-anticoagulation. Horizontal axis is the genomic position, and vertical axis is minus log of p-value. Red dots above the gray line indicate association of genome-wide significance (p<1.5×10−7) at SNPs in the VKORC1 locus such as rs9923231 (p = 8.9×10−9). However, no other loci achieved genome-wide significance. See main text for more details. (0.07 MB PDF) Click here for additional data file. Table S1 Multivariate regression results for 40 SNPs followed-up after GWAS. Unlike rs2108622 of CYP4F2, none of these 40 SNPs exhibited statistical significance after correction for multiple testing. See main text for more details. (0.03 MB XLS) Click here for additional data file. Table S2 Testing for statistical interaction between predictors of warfarin dose. After correcting for the 15 interaction tests, no pair of predictors exhibited statistically significant interaction. Data is for the combined panel of subjects (N = 1641). See main text for more details. (0.02 MB XLS) Click here for additional data file. Table S3 Survival time analysis for incidence of over-anticoagulation within the first 5 weeks of treatment in the whole WARG cohort (N = 1489). See main text for more details. (0.02 MB XLS) Click here for additional data file. Table S4 Multiple regression analysis of warfarin dose in WARG replication samples, WARG GWAS plus replication, or Uppsala replication samples. This table shows the same results displayed in Table 2 of the main text except that WARG and Uppsala subjects are separated into different subsets. (0.02 MB XLS) Click here for additional data file. Text S1 For populations in Hardy-Weinberg Equilibrium, Linkage Disequilibrium r2 and Genotypic R2 are approximately equal. (0.03 MB DOC) Click here for additional data file.
