Migraine polygenic risk score associates with efficacy of migraine-specific drugs

Study population—the target sample

The study population consisted of 2,591 patients with migraine who were recruited at the Danish Headache Center in 1999–2002, 2005–2006, and 2010–2011.^13,14 All patients with migraine were interviewed face to face or by telephone by a trained physician or trained senior medical student using a semistructured interview. The interview was designed by head of classification committee Prof. Jes Olesen to phenotype and classify migraine diagnosis according to the International Classification of Headache Disorders, second edition.¹⁵

Migraine drug response

The semistructured interview included questions covering the necessary clinical data for migraine diagnoses and information on the effect of migraine treatment. Acute treatment effect was considered to be positive in cases where the patient reported at least 50% pain reduction within 2 hours after taking medication. Prophylactic treatment effect was considered to be positive in cases where the patient reported a reduction of over 50% in migraine attacks. For acute treatment, the patient was asked about efficiency of (1) triptans and (2) ergotamine, which are both migraine-specific drugs, and (3) weak analgesics, which is nonspecific for migraine treatment. For prophylactic treatment, the patient was asked about the efficiency of (1) β-blockers, (2) Ca²⁺ antagonists, (3) angiotensin II receptor blockers, (4) angiotensin-converting enzyme (ACE) inhibitors, (5) anticonvulsants, (6) antidepressants, and (7) hormone treatment. Both generic and commercial names were mentioned, where the interviewer used pro.medicin.dk as reference. The questioned drug needed to be taken specifically for treatment of migraine. For all questions, the answer “Do not know” was considered as missing data.

Genotyping

All patients with migraine were genotyped on the Illumina HumanOmniExpress 12v1 (n = 2,152) or Illumina HumanOmniExpress 24v1 (n = 439) chip. For each data set, the quality control of genotypes was performed using PLINK 1.9.¹⁶ We used genotypes for 2,766 ethnicity-sensitive single nucleotide polymorphisms (SNPs) common to all Illumina SNP arrays to estimate European, Asian, and African ancestry probabilities with STRUCTURE¹⁷ and excluded individuals with <90% European ancestry. SNPs with <0.95 genotyping rate, <0.01 minor allele frequency, or p < 1 × 10⁻⁶ for Hardy-Weinberg Equilibrium were excluded, and individuals with <0.98 genotype rate were removed. Next, we created a subset of markers independent of each other with respect to linkage disequilibrium (LD) using a window size of 100 markers shifting by 25 markers at a time and removed 1 half of every SNP pair with genotypic r² > 0.1. This was performed to avoid overestimating the effect by including mutually dependent SNPs, i.e., SNPs in LD. Using this subset of markers, we calculated heterozygosity (HET) and sex and removed (1) all individuals with outlying HET values (>5 SD from the median of the whole sample) and (2) all individuals where sex determined from genotype did not match reported sex. We then removed all A/T and C/G markers to avoid strand issues. Related individuals were detected based on their genotype data, and 1 random individual per related couple (Pihat > 0.10) was removed. After filtering and quality control, 542,168 SNPs and 2,219 individuals were retained for analyses. The total data set of 2,219 patients with migraine consisted of 1,201 patients who had migraine without aura (MO) and 1,018 patients who had migraine with aura (MA). Cases with probable migraine (with or without aura) were included in the analysis, and in case the individual had both MA and MO, they were assigned to the MA subgroup.

PRS calculation

The PRS was calculated using LDpred, which adjusts for LD between markers and further rescales allelic effects based on the likelihood of each marker belonging to the fraction of markers assumed to be causal.¹⁸ We calculated PRSs using the default models for causal variant fraction in LDpred (i.e., 1, 0.3, 0.1, 0.03, 0.01, 0.003, and 0.001). The LD information was retrieved from the subjects with migraine and 929 unrelated Danish controls, who were genotyped on the same genotype chip (Illumina HumanOmniExpress 12v1). To calculate PRSs for migraine, we used p values and log₁₀ odds ratios (ORs) from a subset of the International Headache Genetics Consortium (IHGC) migraine GWAS meta-analysis (n_case = 59,674; n_control = 316,078)¹² from which all individuals of Danish origin (1,771 cases and 1,000 controls) had been removed to avoid overlap between the discovery and target sample and a resulting overestimation of allelic effects.

To investigate which fraction of causal variants gives the best prediction of migraine, we compared the PRSs of our migraine sample (n = 2,219) with the 929 Danish controls. Migraine was most significantly predicted with a model assuming the fraction of causal variants to be 0.03 (p = 6.91 × 10⁻²⁷). The PRS generated under this model predicted both MO and MA significantly (p = 3.69 × 10⁻²⁶ and 4.98 × 10⁻¹⁷). Next, we investigated whether the PRSs of the migraine subtypes could predict the respective subtype better than the PRS of migraine, using the GWAS on the clinical subset (5,557 MA and 7,352 MO) of the IHGC meta-analysis.¹³ This was not the case, likely because of the limited sample size of the discovery cohort; the PRS of MO was predicted with a p value of 2.70 × 10⁻¹² and the PRS of MA predicted MA with a p value of 1.61 × 10⁻³. Therefore, all analyses were conducted using the PRS of migraine. We then rescaled the migraine PRS to a mean of zero and a unit corresponding to a twofold genetic increased risk for migraine in the target population; this was done by first subtracting the mean PRS from each subject's PRS and then multiplying it by log(OR)/log(2), where the OR was extracted from the model predicting migraine using the 2,219 cases and 929 controls.

Statistical analysis

The rescaled PRS for migraine was tested for its association with drug response using a logistic regression model including age, sex, genotype chip, and the first 10 principal components (PCs) of the genotypes as covariates. The PCs were calculated in PLINK¹⁶ and included in the model to correct for population stratification. As triptan and ergotamine are both migraine-specific drugs used for acute treatment and act through the same serotonin receptors (5-HT_1B and 5-HT_1D), they were analyzed together to increase the statistical power. The mode of action of prophylactic treatments is unknown; therefore, they were analyzed together. The association of migraine-specific acute, migraine nonspecific acute, and prophylactic treatment with the PRS was corrected for multiple testing (n = 3) using Bonferroni correction resulting in adjusted p values (p_adj). As prophylactic treatment is potentially confounded by comorbid hypertension or epilepsy, we tested whether these comorbidities had a significant effect on treatment response. In case they were statistically significantly associated with treatment response, they were included as covariates.

All analyses were performed for the complete set of patients with migraine. Subsequently, it was tested whether there was a statistically significant difference between the migraine subtypes by including an interaction term between the PRS and migraine subtype. We presented the area under the curve (AUC), representing the prediction accuracy and ORs, using the partial Receiver Operating Caracteristic R-package.¹⁹ The AUC was calculated for both the model including only the covariates and the full model including the PRS and the covariates. The difference between the 2 AUCs was tested using the DeLong test in the partial Receiver Operating Caracteristic R-package. ORs were presented with 95% confidence intervals (CIs). All analyses were performed in R (version 3.4.3.).²⁰

Replication cohort

The Danish Blood Donor Study (DBDS) genomic cohort (n = 79,595) was used as the replication cohort (see detailed description elsewhere).²¹ The PRS for migraine was calculated as done for the clinical cohort. Using a subpopulation of the DBDS genomic cohort (n = 17,222) with information on self-reported migraine (n_migraine = 3,906), we estimated the OR for migraine within the DBDS genomic cohort (OR = 3.98, 95% CI = 3.18–4.98). The OR for migraine was used for subsequent normalization of the PRS as done for the clinical cohort. Using the prescription register of the 79,595 participants, we identified 5,616 users of migraine-specific treatment (1,372 males and 4,244 females). Positive triptan responders were defined as having 10 or more purchases of triptans, as previously suggested to be a reliable indicator of positive treatment response.²² This resulted in 1,246 triptan responders (213 males and 1,033 females). In the regression model, age, sex, and the 10 first PCs were included as covariates.

Standard protocol approval, registrations, and patient consents

Written informed consent was obtained from all participants. The study was approved by the Danish Ethical Standards Committee (H-2-2010-122) and the Danish Data Protection Agency (01080/GLO-2010-10).

Data availability

Summary statistics of the GWAS are available upon agreement with the IHGC due to embargo with 23andMe. Genotype data of our clinical cohort are available upon agreement with the senior author and upon material transfer agreement.