Precision Target Discovery for Migraine: An Integrated GWAS-eQTL-PheWAS Pipeline

SMR analysis results using genetic tool data, including target gene expression in whole blood from the eQTLGen Consortium, prefrontal cortex cis-eQTL data from the PsychENCODE project, and data from 13 brain regions (amygdala, anterior cingulate cortex, caudate nucleus, cerebellar hemisphere, cerebellum, cerebral cortex, frontal cortex, hippocampus, hypothalamus, nucleus accumbens, putamen, spinal cord (as part of the central nervous system), and substantia nigra) [30,31], and whole blood of the GTEx V8 project, are presented in Supplementary Tables S1–S3.

Using FDR < 0.05 and HEIDI > 0.05 as threshold criteria, 60 genes (after deduplication) were identified as significantly associated with migraine. Among these, 31 were expressed in whole blood, 20 in brain regions, and 9 were shared between whole blood and brain regions, as detailed in Table 1.

In order to validate the robustness and reproducibility of the primary MR results, we used two external datasets (BrainMeta v1 and BrainMeta v2) to replicate the analysis. The significant results of the main analysis, cross-validated by ≥2 databases, are presented in. Supplementary Table S4a,b

To clarify whether there are potential repurposable drugs, we extracted the data on druggable genes from the studies of Wei-Ming Su (2023) [32] and Kang-Fu Yin (2024) [33], and screened out eight blood genes (AURKC, CHD4, GBA2, HVCN1, NCOR2, NR1D1, TGFB3, THRA) and five brain region genes (BICD1, IPO8, BACE2, REV3L, CHRNB1) (as shown in Figure 1). The SMR results of the 13 druggable genes are presented in Table 2. Among them, three genes (AURKC, BACE2, IPO8) were validated through independent validation in at least two databases (Supplementary Table S4b).

To assess the consistency of causal estimates, five MR methods were used to validate the significant genes identified by SMR analysis, including inverse variance weighted, weighted mode, simple mode, weighted median, and MR Egger (Supplementary Table S5). Data analysis was performed using the "TwoSampleMR", "MR-PRESSO", and "MendelianRandomization" packages. Data cleaning and statistical analyses were conducted using R version 4.2.2 (https://www.r-project.org/↗, accessed on 4 May 2025). As shown in Figure 2, the results demonstrate that 13 druggable genes that were significantly associated with migraine were simultaneously validated by multiple MR methods. For NCOR2, which did not pass the Egger intercept or global test, the "MR-PRESSO outlier test" and "MR-PRESSO distortion test" were used to verify the consistency of the results after removing outlier SNPs (as shown in Figure 3).

To further identify the pleiotropic association between druggable genes and migraine, we performed a colocalization analysis to confirm that the eQTL and GWAS signals likely originate from the same causal variant. Colocalization analysis was conducted on 13 significant druggable genes, and we found that 10 genes passed the colocalization analysis. The posterior probability results show that the causal variation probability shared by the three genes with migraine is relatively low (GBA2 (PP.H4.abf = 0.53), HVCN1 (PP.H4.abf = 0.63), REV3L (PP.H4.abf = 0.44)), which is insufficient to support the causal effect of this gene on migraine (as shown in Figure 4).

We obtained 74 novel migraine-related targets from Pharos v3.19 (https://pharos.nih.gov/↗), 42 of which are already in clinical development (Supplementary Table S6). For instance, targets such as TRPM8, TNF, and ESR1 have been used in clinical drug development. To further explore the pathological mechanisms of migraine and the druggability of key target genes identified by multi-tissue SMR analysis, we used the STRING database for PPI network analysis. The organism was set to Homo sapiens, and the minimum interaction score was defined as medium confidence (0.40). The results show (as shown in Figure 5) that seven blood genes (THRA, NCOR2, NR1D1, CHD4, TGFB3, PABPC4, PPP1CC) and four brain genes (BACE2, CHRNB1, ID4, MICU1) are moderately associated with ten novel migraine targets (interaction score > 0.4), among which seven are druggable genes (THRA, NCOR2, NR1D1, CHD4, TGFB3, BACE2, CHRNB1). The results of evaluating the drug development potential of the remaining four genes (MICU1, ID4, PPP1CC, PABPC4) are shown in Figure 6, with PPP1CC having the highest standardized score. MICU1 had the highest PubMed score (354.26), followed by ID4. The target development level of both was "Tbio", indicating no related drug development has been conducted yet (specific scores are provided in Supplementary Table S7).

The results of the Gene Ontology enrichment (GO enrichment) analysis show that the 10 migraine-related druggable genes are primarily enriched in gland development, intracellular receptor signaling pathways, and responses to steroid hormones, as well as biological processes, such as steroid hormone biosynthesis (Figure 7a). These findings suggest that the pathogenesis of migraine may be closely linked to endocrine regulation, cellular signal transduction, and biological rhythms, providing critical insights for further elucidating the molecular mechanisms of migraine and developing novel therapeutic strategies.

The GO enrichment analysis of the drug potential genes (REV3L, GBA2, HVCN1) revealed that they were significantly enriched in multiple metabolic processes (e.g., glycosylceramide metabolism, glycolipid catabolic metabolism, lipid glycosylation) and cellular responses (e.g., cellular response to pH, acrosome reaction, and regulation of superoxide anion generation) (Figure 7b). Although the colocalization analysis of these genes with migraine GWASs did not reach statistical significance (PP.H4.abf < 0.8), their enrichment in specific biological processes suggests potential involvement in migraine pathogenesis through indirect mechanisms, thereby providing novel potential therapeutic targets for migraine treatment.

PheWAS analysis further confirmed the non-horizontal pleiotropy of these genes. As shown in Figure 8, when we take NCOR2 as an example, both continuous and binary traits fall within the threshold range (−log₁₀(p-value) < 8).

First-line drug target information for the prevention and acute treatment of migraine was obtained from the DrugBank database (v5.1.7) (encompassing 53 targets; see Supplementary Table S8), and the STRING database (v11.5) was used for network construction. Interaction confidence levels were classified as medium (0.40) or high (0.70). A gene interaction network was constructed between the 13 potential drug targets and 53 known migraine therapeutic targets. As shown in Figure 9, among the 13 identified potential drug targets, 8 (TGFB3, CHRNB1, BACE2, THRA, NCOR2, NR1D1, CHD4, REV3L) interacted with the targets of two acute migraine treatments (aspirin, ibuprofen) and two preventive drugs (topiramate, valproic acid). Among these, five (THRA, NCOR2, NR1D1, CHD4, REV3L) exhibited high-confidence associations with first-line migraine therapeutic targets.

Existing drugs (including marketed and unapproved drugs) related to the potential drug targets were retrieved via the DGIdb (https://dgidb.org↗) and interaction scores were calculated to evaluate the association strength between the drugs and targets (larger nodes in the network indicate higher interaction scores) (specific data are provided in Supplementary Table S9). The results show (as shown in Figure 10) that these 8 druggable genes are associated with 41 existing drugs, revealing that drugs such as 4 TGFB3 inhibitors and 7 THRA agonists may have off-label use for treating or preventing migraine. Furthermore, THRA exhibits a significant interaction (interaction score of 0.168) with aspirin and is the only gene associated with marketed drugs for migraine.

Through SMR analysis, we found that migraine-associated genes are primarily concentrated on Chromosome 6 (UFL1, LY6G5C, REV3L, etc.) and Chromosome 12 (PPP1CC, IPO8, BICD1, NCOR2, CHD4, etc.). Brain region localization analysis further revealed that the cerebellum and cerebellar hemispheres are key regions associated with migraine, enriched for genes such as UFL1 and LY6G5C. This finding aligns with previous neuroimaging studies reporting abnormalities in these brain regions among migraine patients [34], suggesting that these chromosomal regions may influence migraine pathology by regulating gene expression in specific neural tissues.

Notably, core genes such as CHD4 [35,36] and PPP1CC [37] play central roles in the PPI network and are implicated in critical biological processes, including epigenetic regulation, neural development, and synaptic plasticity. For instance, CHD4, as a core component of the chromatin remodeling complex, regulates epigenetic modifications that may influence neuronal differentiation and neuroinflammation [38,39], the latter being recognized as an important pathological hallmark of migraine [38]. Additionally, PPP1CC participates in neural precursor cell differentiation and synaptic transmission by modulating the phosphorylation status of cyclins [40]. Aberrant expression of PPP1CC may contribute to neuronal damage and synaptic remodeling in migraine patients [41]. Collectively, these detailed analyses of gene functions provide new insights into the neurobiological mechanisms underlying migraine.

In terms of drug targets and translational applications, we identified 13 druggable genes with significant causal relationships. Among these, protective genes such as CHD4 and HVCN1 are involved in anti-inflammatory and ion channel regulation pathways, while risk genes, such as GBA2, are implicated in lipid metabolism and cell adhesion processes. Based on the comprehensive SMR, colocalization, and PPI analysis results, the top five in the comprehensive score ranking are NR1D1, THRA, NCOR2, CHD4, and BACE2 (see Figure 11). It is particularly worth noting that the four genes, NR1D1, THRA, NCOR2, and CHD4, not only showed through the gene–drug interaction analysis that they had a highly reliable association with first-line migraine treatment drugs [42,43] (PPI score > 0.7), but there is also a potential synergistic effect with known migraine therapeutic targets, such as the CGRP pathway and the 5-HT receptor.

Furthermore, drug repurposing analysis identifies existing marketed drugs (e.g., GSK4112, Dextrothyroxine Sodium, benzbromarone, vorinostat) as potential migraine treatments, offering theoretical support. The corresponding candidate drugs converge three core pathophysiological axes in migraine (see): (1) vascular dysregulation: TGF-β superfamily modulation (e.g., TGFB3) to balance tone, suppress pro-inflammatory (TGF-β1/2), and enhance anti-inflammatory (TGF-β3) signaling, improving endothelial function. (2) Neuroimmune inflammation: NRF2/HO-1 or complement inhibition (e.g., C5a blockers) to curb microglial activation and IL-1β release, blocking CSD cascades. (3) Neuronal hyperexcitability: BACE2 (Aβ inhibition + Nav1.6 stabilization) and nAChR modulation (e.g., varenicline) reduce CSD susceptibility. Multi-target agents (varenicline, SR9011) synergize via peripheral/central regulation. Epigenetic (HDAC inhibitors) and circadian (REV-ERBα agonists) mechanisms innovate beyond traditional paradigms. Supplementary Table S10

However, BACE2's controversial role and pan-TGF-β inhibitor toxicity require further study. And the current candidate drugs are primarily prioritized based on gene–drug interaction scores. Additional functional experiments are required to verify their effects on core migraine-related pathways (e.g., CGRP signaling and ion channel regulation) and to evaluate their therapeutic potential.

Given that differences in outcome databases may affect the results, we employed the most recent and large-scale GWAS datasets to enhance the reliability of our findings. Compared with previous studies [44], we not only observed the repetition of some results, but also identified new potential targets, including REV3L and BACE2. Both this study and Zhang et al. (2024) [44] identified TGFB3 as a significant gene in both blood and brain tissue, and associated with lipid metabolism and insulin-like growth factor pathways [45,46]. Xiong et al., 2024 [47], using proteomic MR, discovered associations between proteins such as FCAR and UBE2L6 and migraines, and reported protein interactions between HBQ1 and topiramate (a preventive drug), but did not identify consistent results across different tissues. Similarly, we found that BACE2, CHRNB1, and TGFB3 had protein interactions with topiramate, but not such interactions with HBQ1. Although NCOR2 showed heterogeneity in different MR validation methods in this study, previous population studies have confirmed its association with migraine or epileptia-related phenotypes [48]. Our study shares a similar methodology to that of Sun et al. (2024) [49]. However, our study was based on large-scale data and systematically screened drug genes and candidate drugs, whereas Sun et al.'s (2024) [49] focused on identifying a single target, GSTM4, by integrating eQTL, pQTL, and GWAS data. In contrast to our study, which employed SMR and HEIDI as primary analytical methods, Sun et al. (2024) [49] used them only for validation purposes following initial target identification.

This study has several methodological strengths. We integrated multi-source GWASs and eQTL data, combined with SMR, multi-method MR validation, and strict horizontal pleiotropy control (p < 5 × 10⁻⁸), effectively enhancing the reliability of causal inference. However, due to the reliance on publicly aggregated data, differences in sample size, racial representativeness, and inconsistent data quality control standards may affect the validity of IVs [50]. For instance, stricter LD pruning or sample overlap control may lead to a reduction in instrumental variables, thereby affecting the significance of the results [51]. In addition, the sample size of this study is relatively large. Since Cochran's Q test may show significant results due to the "large sample power" [52,53], we do not take the test results as the sole criterion. This might lead us to include more significant results. Furthermore, although confounding factors were excluded through the PheWAS database, the biological functions of some genetic variations have not been fully elucidated, and the risk of horizontal pleiotropy may still exist. In the future, it is necessary to further verify the robustness of causal associations based on larger-scale multi-ethnic cohorts, combined with fine phenotypic stratification and functional experiments [54].

Complementing these mechanistic considerations is the need to address the limited generalizability of our findings, namely the overreliance on European-ancestry cohorts. Approximately 85–90% of GWASs for common traits, including migraine, focus on European populations, with cross-ancestry (African, East Asian) analyses hindered by small sample sizes and data heterogeneity. While key associated genes (e.g., "TGFB3", "CHRNB1") localize to conserved or low-variant regions, future work must integrate non-European cohorts via international collaborations and leverage cross-ethnic functional annotations (gnomAD, TOPMed) to evaluate locus-specific effects, strengthening the global relevance of our conclusions.

Moreover, our current druggability filter, focused on genes with established experimental or structural evidence, may overlook hidden opportunities. Of the 60 candidate genes, 47 could exhibit "cryptic" druggability via uncharacterized pockets or emerging modalities. To address this, we plan to systematically re-evaluate these targets using AlphaFold2-based conformational sampling to uncover novel binding sites. Concurrently, we will explore PROTACs (proteolysis-targeting chimeras) and RNA-targeting strategies to expand the druggability landscape beyond traditional pocket-dependent targets, thereby broadening therapeutic possibilities.

While expanding druggability horizons is critical, translating these targets into actionable therapies demands rigorous mechanistic validation. Although multi-omics (GWAS-eQTL-PheWAS) links genes like "NR1D1/CHD4" to migraine-related pathways (circadian disruption, neuroinflammation), it will be necessary in the future to conduct in vitro verification (e.g., neuronal excitability assays) to directly verify their causal role in migraine pathology mechanism. Moreover, further studies should be conducted to MR static nature inherently limits our ability to capture the "dynamic" gene expression shifts across migraine's distinct phases. Prospective, time-stratified sampling of attacks remains logistically arduous and ethically complex, and publicly available resources like UK Biobank or GTEx lack annotations for migraine phases. To bridge this gap, future studies should integrate time-resolved transcriptomic profiling from migraine animal models paired with targeted clinical sampling during premonitory/ictal phases to dissect these patterns.

Another layer of complexity arises from the pleiotropic nature of many candidate genes and their associated drugs (e.g., "NR1D1", "CHD4", "BACE2"), which risk off-target effects. To mitigate this, we propose a multi-level optimization framework: (1) tissue-specific targeting (e.g., CGRP/TRPV1 nanoparticles or cell-type promoter-driven AAVs) to confine action to migraine-relevant tissues (trigeminal nerves, hypothalamus); (2) combinatorial therapy (synergistic drug pairing or dose fractionation) to enhance efficacy while reducing systemic toxicity; and (3) precision dosing aligned with circadian biology (e.g., REV-ERBα agonists) or real-time biomarkers (e.g., CGRP monitoring) to match migraine's variability.

Taken together, these limitations underscore a clear roadmap for future research. Based on our findings, the following directions are recommended for future research:(1)Conduct cross-ancestry validation studies to explore the impact of gene-environment interaction on migraine;(2)Verify the functions of key genes (such as MICU1, UFL1, LY6G5C, PPP1CC) in the pathology of migraine through CRISPR/Cas9 or animal models;(3)Develop targeted regulatory strategies for genes such as NR1D1, THRA, NCOR2, CHD4, and evaluate their roles in neuroinflammation and synaptic plasticity;(4)Design multi-target combination therapies tailored to the different phases of migraine (e.g., acute vs. preventive), and integrate real-world adverse reaction data from monotherapies [55] to optimize treatment safety and efficacy.

In summary, this study provides a comprehensive analysis of the genetic-epigenetic molecular network of migraine. We identify multiple potential drug targets and repurposable drugs, offering new insights into the mechanisms of migraine and potential avenues for translational therapy [40]. Future research combining experimental verification and multi-omics analysis will be crucial to further elucidate migraine pathogenesis and promote precise diagnosis and treatment strategies.

This study is based on multiple publicly available databases, which are summarized in Supplementary Table S10. As shown in Figure 12, this study included relatively comprehensive cis-eQTL genetic variants as the IVs for gene expression. We performed SMR analysis for 13 regions in the brain and whole blood. We used the GTEx V8 (n = 114–670), PsychENCODE (n = 1387) project and eQTLGen (n = 31,684) summarized data as the main analysis. Due to the differences of various brain tissues, we also conducted subgroup analysis using BrainMeta v1 eQTL (n = 1194) and BrainMeta v2 (n = 2865) to supplement and verify the results. Multiple MR analysis was conducted on the key genes obtained from main analysis. To further examine the pleiotropy, linkage disequilibrium (LD), and heterogeneity at the target gene level, we utilized colocalization analysis, PPI, and GO enrichment analysis to evaluate the associations and mechanisms of action among genes. Finally, the clinical applicability of the candidate drug targets was evaluated through PheWAS.

To identify common SNPs (MAF ≥ 1%) associated with the expression of migraine drug target genes in whole blood and the brain, we extracted publicly available eQTL data. The genetic tools for brain target gene expression originated from cis-eQTL summary data in the prefrontal cortex of the PsychENCODE project (n = 1387) and from whole blood and 13 brain regions of the GTEx V8 project (n = 114–670). The genetic tool for whole blood target gene expression was derived from the eQTLGen Consortium (https://www.eqtlgen.org/↗), which integrated 37 datasets covering 31,684 blood and peripheral blood mononuclear cell samples and reported 16,989 genes expressed in whole blood. All SNPs used in the above eQTL analysis originated from cis-regulatory regions (within 1 Mb of each gene), with a default p-value threshold of 5 × 10⁻⁸ (data sources are summarized in Supplementary Table S11).

The aggregated statistical data of migraine outcomes were derived from the study by Hautakangas Heidi et al. [12]. This study integrated 873,341 individuals of European ancestry (102,084 cases and 771,257 controls) from five sets of European ancestry studies (see Supplementary Table S10).

We employed the SMR method to evaluate the association between target gene expression levels (whole blood and brain eQTL) and migraine risk (GWAS). SMR is a Mendelian randomization tool based on aggregated data, which is used to explore causal relationships between genetic variations and phenotypes [16]. We used SMR Version 1.3.1 [16] with default settings, including a cis-eQTL p-value threshold of <5 × 10⁻⁸, MAF > 0.01, and excluded SNPs in strong LD (r² > 0.9) with top eQTLs. Additionally, we removed SNPs with weak LD (r² < 0.05) with top eQTLs. SMR results were derived from the most significant SNP per gene. To account for multiple testing, we used the Benjamini–Hochberg method [56] to correct for multiple comparisons and identify statistically significant associations. The HEIDI test was additionally conducted to assess horizontal pleiotropy. For genes that reached the corrected significance threshold (FDR < 0.05) and passed the HEIDI test (p_HEIDI > 0.05), the SMR site map was generated using the method described on the SMR webpage (https://yanglab.westlake.edu.cn/software/smr/#Overview↗, accessed on 13 April 2025). The effect direction on migraine was determined based on b_SMR, and the odds ratio (OR) and 95% confidence interval (CI) for the expression level of each significantly associated gene and its association with migraine risk were calculated, respectively.

The robustness of causal relationships was tested by excluding reverse causality, horizontal pleiotropy, and uncertainty in colocalization. Reverse causality was excluded using the Steiger directionality test (p < 0.05). Horizontal pleiotropy was examined by testing the association of other genes within a 2 megabase (Mb) range of each genetic IV, and SMR analysis was used to determine if these adjacent genes were associated with migraine risk. Bayesian colocalization analysis was performed to assess whether two traits share a common causal variant, using the "coloc" package with default parameters (https://github.com/chr1swallace/coloc↗, accessed on 18 April 2025). This analysis calculates the posterior probabilities of five hypotheses (H0-H4) regarding whether a single causal variant is shared between the two traits. Hypothesis 4 (PPH4) proposes that both traits are associated with the genetic variant, and the associations are driven by the same causal variant. We used the coloc.abf algorithm to identify genes with PPH4 > 0.80 as determined by the algorithm.

The key genes obtained from SMR analysis were subjected to multiple MR analyses. Sensitivity analyses included horizontal pleiotropy, heterogeneity, and leave-one-out sensitivity tests. Horizontal pleiotropy was tested using MR–Egger regression [16,57]; significant intercept terms in MR–Egger analysis indicated horizontal pleiotropy [58]. Heterogeneity was quantified by calculating the Cochran's Q statistic; a p-value < 0.05 for the Cochran's Q test indicated significant heterogeneity [59]. Outliers detected by the MR pleiotropy residual sum and outlier test (MR-presso) were excluded, and the remaining IVs were reanalyzed [58]. The leave-one-out test examined the potential impact of a single SNP on causal effect estimation by sequentially excluding each SNP, calculating the meta-effect of the remaining SNPs, and assessing whether a single SNP drove the causal association [58].

To conduct a more comprehensive summary analysis of eQTL data in brain tissue, SMR analysis was continued using the BrainMeta v1 eQTL summary data (n = 1194) by Qi et al. (2018) [60]. BrainMeta v1 eQTL is a group from MeC on GTEx Brain (GTEx Consortium 2017 Nature) [61], CMC [62], and ROSMAP (eQTL data from the meta-analysis of Nat Neurosci et al., 2017) [63]. Only SNPs within a 1Mb distance from each probe were available. In addition, chromosomal SMR subgroup analysis was conducted using BrainMeta v2 eQTL, and these data were from Qi et al. (2022) [64]. Cis-eQTL analysis of 16,704 eGenes in 2865 cerebral cortex samples from 2443 unrelated individuals of European origin Summary statistics (https://yanglab.westlake.edu.cn/software/smr/#eQTLsummarydata↗, accessed on 1 April 2025) was used to supplement the validation master analysis results.

The data of druggable genes mainly came from the study of Su et al. (2023) [32], integrating a total of 5883 druggable genes reviewed by DGIdb v4.2.0 (https://www.dgidb.org/↗) and Finan et al. (2017) [65]. In addition, 2532 overlapping druggable genes integrated by Kang-Fu Yin et al. (2024) [33] were also utilized for refinement and verification.

To conduct a more comprehensive screening of candidate genes, we used an online database resource, Pharos (https://pharos.nih.gov/↗), to evaluate the potential druggability of undruggable genes. To identify the role of key druggable genes in migraine, PPI was conducted via the STRING database. First-line intervention drugs for the attack and remission periods were obtained from clinical guidelines related to adult migraines [42,43], and the target genes associated with these different drug types were identified from the DrugBank database. Enrichment analysis was applied with significance thresholds of p ≤ 0.05 and Q ≤ 0.05 [66]. Candidate drug prediction was performed using DGIdb 5.0 (https://dgidb.org↗).

PheWAS analysis of candidate genes was conducted using the AstraZeneca PheWAS portal (https://www.azphewas.com↗), a publicly available repository of gene–phenotypic associations. These data were generated using sequencing data and phenotypic data collected in the UK. Based on the results of the phenotypic analysis, potential side effects or horizontal pleiotropy of the key targets were identified, which strengthened the validity of our analysis. The threshold for identifying genes with horizontal pleiotropy (i.e., associations with multiple phenotypes) in the PheWAS Manhattan plot was set at −log₁₀(p-value) < 8.