Genetic variability within molecular core clock genes in cluster headache

Abstract

Background / Hypothesis

Headache attacks are reported to occur with circadian rhythmicity by 2/3 of individuals with cluster headache, hinting to potential dysfunctions of the molecular clock. The aim of this study was to investigate the contribution of genetic markers in core clock genes, alone or in combinations, to the genetic risk for cluster headache.

Methods

Seven markers in core clock genes BMAL1 and NPAS2 were genotyped using TaqMan qPCR in 707 individuals with cluster headache and 682 controls. Genetic data from eleven additional markers located in six other core clock genes (CLOCK, CRY1-2, and PER1-3) was obtained from the database of the Centre for Cluster Headache at Karolinska Institutet. Genotype analysis was applied for risk analysis for combinations of up to three markers. For validation we used a cluster headache cohort from the National Hospital for Neurology and Neurosurgery, London, UK.

Results and interpretation

Single marker analysis of the newly genotyped markers in BMAL1 and NPAS2 found rs3789327 and rs3768984 more frequently among individuals with cluster headache (p < 0.01 and p < 0.05 respectively). Multiallelic analysis revealed that the median number of risk alleles was eight for controls and nine for individuals with cluster headache, which justifies the analysis of combinations of markers. The analysis of combinations of up to three markers identified 258 out of 897 combinations to be associated with significant risk. Further investigation starting from the function of genes in the molecular clock transcription-translation feedback loop (TTFL) found that 80% of the combinations had >50% markers located in the positive arm of the TTFL. The comparison between Swedish and UK cohorts identified 39 concordant combinations, which confirmed the risk associated with rs3768984 (NPAS2), as well as the enrichment in markers in BMAL1, CLOCK, and NPAS2 in combinations associated with significant risk.

Conclusion

Our data points to molecular clock dysfunction playing a central role in the manifestation of cluster headache.

Graphical abstract

This is a visual representation of the abstract.

Keywords

circadian sleep

Introduction

Cluster headache (CH) is the most prevalent form of trigeminal autonomic cephalalgia with an occurrence of 1/1000.^1,2 Cluster headache attacks often occur during specific time points of the day, with a higher prevalence around 2–4 a.m.² Additionally, individuals with CH report a rhythmicity to the active bouts that respect a certain seasonality.^2,3 Finally, the hypothalamus, responsible for many complex behaviors such as sleep and arousal, is known to be activated during attacks.⁴ The hypothalamus contains the suprachiasmatic nucleus, known as the master circadian pacemaker, which helps entrain the clock genes to external stimuli.⁵ Altogether, this suggests that the underlying mechanisms of CH are somehow connected to the circadian clock.⁶

In mammals, the circadian clock allows the synchronization of behavior and body physiology to light-dark phases following a 24h-cycle.⁵ The molecular clock is an intracellular mechanism dictating circadian rhythm of cells. It is very complex and involves several transcription-translation feedback loops (TTFL). The core mechanism relies on the following genes: Brain and muscle ARNT-like (BMAL1 also known as ARNTL1), Circadian locomotor output cycles kaput (CLOCK), Neuronal PAS domain protein 2 (NPAS2), three Period genes (PER1, PER2, PER3) and two Cryptochrome genes (CRY1, CRY2).⁷ CLOCK and NPAS2 are paralogs and can both dimerize with BMAL1.⁸ This creates an activator complex that can bind to E-box elements, DNA response elements that can act as promoters or enhancers of target genes. Among genes that regulate cell metabolism are the family genes of Period and Cryptochrome.⁹ These genes encode repressive transcription factors that dimerize and bind to BMAL1/CLOCK, interrupting its transcriptional activity. The PER/CRY complex is gradually degraded by ubiquitylation, lifting the inhibition of BMAL1/CLOCK which allows the complex's transcriptional activity to start again.⁷ This alternating transcription creates a cycle that is tightly regulated by internal and environmental factors and involves different feedback loops that help fine tune the molecular clock.

Despite these diverse regulatory mechanisms, the molecular clock can be altered and has been suggested to contribute to several diseases ranging from sleep disturbances to psychiatric disorders.^10–13 In addition to sharing a circadian phenotype, both sleep and psychiatric disorders have been linked to CH as co-morbidities, and may therefore share certain genetic risk-factors with CH.^6,14 Indeed, it has been demonstrated that sleep is altered in individuals with CH,^6,15 and that CH is genetically correlated with mood disorders.¹⁶ Although a formal association has only been made with depression and attention-deficit/hyperactive-disorder, bipolar disorder (BD) is also often discussed in relation to CH and individuals with CH have a higher prevalence of BD than the general population (6.6% vs. 2.4%).¹⁷ Moreover, a study comparing BD and CH with healthy controls transcriptomes pointed out similar gene dysregulation profiles between CH and BD.¹⁸

Previously, genetic markers in clock genes have been studied and observed to be significantly more present in CH population. This has been the case for CLOCK and CRY,^19,20 though no association could be found for PER.^21,22 In this study, we complement the investigation of the core feedback loop of the molecular clock by genetic analysis of single nucleotide polymorphisms (SNPs) in BMAL1 and NPAS2.

On that premise, BMAL1 and NPAS2 SNPs were chosen from published articles on the criteria of having previously been associated with sleep variations (rs10766071, rs1562438, and rs3768984) or psychiatric disorders (rs11123857, rs11541353, rs3789327 and rs6738097). These SNPs were investigated for possible association with CH in a Swedish cohort using association analysis, whereafter a combinatorial multiallelic analysis was conducted including SNPs from all core clock genes studied in our cohort to date. The purpose of the multiallelic analysis was to investigate the joint influence of core clock genes on CH and further replicate the analysis in an independent cohort.

Material and methods

Material

Blood samples for the genotyping in the present study were collected from 682 healthy controls and 707 individuals with CH residing in Sweden, after written informed consent was obtained. Study participants with CH were recruited at the Karolinska University Hospital between 2014 and 2024. Controls were recruited at the Karolinska University Hospital, as anonymous healthy blood donors aged 18–60 years. Blood donors are suitable controls for rare disorders such as CH because they constitute a random cross section of the healthy population and there is no selection bias in the recruitment. All CH diagnoses were made by a neurologist following the International Classification of Headache Disorder (ICHD-III beta). The demographic information is available in Table 1. This study is reported in accordance with the STROBE guidelines for observational study.

Table 1.

Demographic characterization of cluster headache (CH) participants and controls.

	CH participants	Controls
Number of individuals (n)	707	682
Age (years)	52.0 ± 15.5	n/a
Age at onset (years)^a	32.0	-
Sex (female) % (n)	32.4 (229)	44.3 (302)
Chronic CH % (n)^b	10.9 (77)	-
Diurnal rhythmicity % (n)^c	68.5 (440)	-

CH = Cluster headache, ^a,b,c Information was available for a subset of study participants: ^a Age at onset (n = 591), ^bCH type (n = 703), ^c Diurnal rhythmicity of attacks (n = 642).

Ethical permit was approved by the Regional Ethical Review Board in Stockholm, (diary number 2014/656–31/4). All experiments were conducted in accordance with the declaration of Helsinki adopted by the World Medical Association regarding human samples.

Genotyping

Genotyping of seven SNPs in BMAL1 and NPAS2 was performed on DNA extracted from blood samples using standard protocols. Briefly, 25 µL of assay and 475 µL of Master mix (either TaqMan Genotyping Master Mix, or TaqPath ProAmp Master Mix (Thermo Fischer Scientific, Waltham, USA)) was added to 5 ng of dry DNA for the BMAL1 and NPAS2 assays (Table 2). Samples were genotyped using a 7500 FAST Real-Time instrument or a QuantStudio 5 Real-Time PCR System from Applied Biosystems by Thermo Fisher Scientific. Data acquisition and analysis was made with the Design & Analysis Software 2.7.0 supplied by Applied Biosystems (Thermo Fisher Scientific). The total call rate for all SNPs was >95%.

Table 2.

Genotyping results of SNPs in BMAL1 and NPAS2.

Gene	SNP	Assay ID	MA	MAF		Power of 80%	Position	n	OR (95% CI)	P_corr-value
Gene	SNP	Assay ID	MA	C	CH	Power of 80%	Position	n	OR (95% CI)	P_corr-value
BMAL1	rs10766071	C___31248742_20	G	0.010	0.022	0.217 < OR > 2.282	chr11:13291666	1312	2.054 (1.10–3.85)	0.175
	rs1562438	C___2160492_10	T	0.324	0.302	0.701 < OR > 1.390	chr11:13342653	1346	0.899 (0.76–1.06)	1
	rs3789327	C___2160503_20	G	0.437	0.503	0.734 < OR > 1.359	chr11:13363769	1265	1.312 (1.12–1.54)	0.006**
NPAS2	rs11123857	C__30726372_10	G	0.292	0.271	0.697 < OR > 1.394	chr2:100987350	1369	0.909 (0.77–1.07)	1
	rs11541353	C___2153849_10	T	0.181	0.155	0.697 < OR > 1.394	chr2:100977729	1276	0.835 (0.68–1.03)	0.651
	rs3768984	C___2153804_10	C	0.200	0.242	0.687 < OR > 1.407	chr2:100934525	1359	1.298 (1.08–1.56)	0.042*
	rs6738097	C__28972034_10	C	0.134	0.153	0.606 < OR > 1.522	chr2:100962080	1368	1.164 (0.94–1.44)	1

Analysis was made using logistic regression with sex as a covariate. C = Control, CH = Cluster headache, CI = Confidence interval, MA = Minor allele, MAF = Minor allele frequency, OR = Odds ratio, P_corr = P-value after Bonferroni correction, SNP = Single nucleotide polymorphism, *P-value < 0.05, **P-value < 0.01.

Combinatorial multiallelic analysis data

Genotype information for SNPs in core clock genes was compiled from genotyping data from the present study (BMAL1 and NPAS2) as well as from previously published data (CLOCK, CRY1, CRY2, PER1, PER2 and PER3) from our Swedish cohort (Online Supplementary Table S1).^19,20,22 For validation of the combinatorial multiallelic analysis we extracted genetic information from genome wide association study (GWAS) data from a CH cohort from the UK, described previously by O’Connor et al.²³ Out of the 18 SNPs genotyped in the Swedish discovery cohort, 17 were present in the UK replication cohort, while rs10766071 in BMAL1 was replaced by another SNP in high linkage disequilibrium (LD) present on the UK array; rs74921564 (GRCh37 11:13329090, D’ = 1, R² = 1 in CEU). Linkage disequilibrium was evaluated using LDLink/Proxy.²⁴ The inclusion criteria for the two cohorts were (1) each subject must have information on at least 50% of the SNPs included in analysis; and (2) for each SNP there must be information for at least 2/3 subjects (see Online Supplementary Table S1 and S2). The application of selection criteria yielded 774 controls and 705 individuals with CH in the discovery (Swedish) cohort; and 5558 controls and 908 individuals with CH in the replication (UK) cohort. To facilitate the presentation and interpretation of the results, the SNPs were sorted by gene, then by location on chromosome, while the genes were grouped by contribution to the positive or negative arms of the transcription-translation feedback loop (TTFL).

The risk allele for each SNP was defined as the marker with higher frequency in CH cases (regardless if it was a wildtype or mutant allele; see Online Supplementary Table S2). Combinations of SNPs were analyzed using the carrier status for subjects, i.e., “carrier” corresponded to subjects which were either heterozygous or homozygous for the risk allele. We used a brute force approach, with independent testing for all combinations of up to three SNPs, which yielded 987 combinations as follows: 18 single SNPs; 153 pairs; and 816 triplets. The contribution of individual SNPs in combinations surviving false discovery rate (FDR) correction was analyzed using posterior inclusion probability (PIP), with a prior inclusion probability of 1/18 (as expected for uniform sampling). SNPs with PIP > 1/18 were labelled as “enriched” and were therefore considered as relevant for the set of combinations yielding significant risk.

For validation of combinations, the labelling of risk alleles in the discovery cohort was applied to the replicationcohort, and the parallel testing of combinations of up to three SNPs was performed as described above. The subset of combinations surviving FDR correction in both discovery and replication cohorts was further analyzed for identification of enrichment for individual SNPs.

Statistical analyses

All SNPs were tested for Hardy-Weinberg equilibrium and allelic association analysis was made using logistic regression with sex as a covariate using PLINK1.9.²⁵ None of the SNPs were in linkage disequilibrium. P-values were corrected for multiple testing for seven tests. Haplotype analysis was performed using Haploview 4.1. P-values were rectified using a permutation test using 10,000 random permutations. The statistical analysis was made using Rstudio 4.3.2 (RStudio Team 2020). Visual representation of the data was made using the ggplot2 package and Matlab version R2025a (The Mathworks, Natick, MD, USA). Power calculations were made using PS Power and Sample Size Calculation program Version 3.0.²⁶ For all genotyped SNPs, we used MAFs from the ALFA Allele Frequency project. The risk for each combination was evaluated using chi-square statistics, followed by Benjamini-Hochberg correction for FDR (threshold set to 0.25; range for sensitivity analysis: 0.05–0.25, see Online Supplementary Table S3).

Results

Association analysis of BMAL1 and NPAS2

All SNPs analyzed were in Hardy-Weinberg equilibrium, except for rs11123857 (P-value = 0.013) in NPAS2. Three SNPs in BMAL1 were analyzed; rs10766071, rs1562438 and rs3789327. The minor allele of rs10766071 has a very low minor allele frequency (MAF) of 3% in the European population. In our cohort, the minor allele was more often present in the CH group compared to controls (2.18% vs 1.05%). However, this difference was not significant after correction for multiple testing (P_corr-value = 0.17), (Table 2). The analysis of rs1562438 showed no difference between the two populations (Table 2). The minor allele of rs3789327 was more common among study participants with CH than in the control group (P_corr-value = 0.006, OR = 1.312, 95%CI = 1.12–1.54) (Table 2). Four NPAS2 SNPs were studied: rs11123857, rs11541353, rs3768984 and rs6738097. No significant difference was observed in the distribution of rs11123857, rs11541353 and rs6738097 (Table 2). The minor allele of rs3768984, was significantly more prevalent in the CH group (MAF = 24%) than in the control group (MAF = 20%), P_corr-value = 0.042, OR = 1.298, 95%CI = 1.08–1.56.

Haplotype analysis of BMAL1 and NPAS2

Haplotype analysis revealed one association with BMAL1 (Table 3). The haplotype, “A-T-A” was associated with a decreased risk for CH, and the complementary haplotype “C-G-C” showed a trend for association with increased risk for CH. These data confirm the results from the allelic association analysis, as the “C-G-C” haplotype involves the minor allele of the two SNPs (rs10766071 and rs3789327) that were more common in individuals with CH. No specific haplotype in NPAS2 was associated with CH (Table 3).

Table 3.

Results from the haplotype analysis.

Gene	Haplotype	Control n	Control %	Case n	Case %	χ ² (df = 2 (BMAL1), 3 (NPAS2))	P-value	P_corr -value
BMAL1	A-C-A^#	549	0.411	529	0.386	1.729	0.189	0.646
	A-C-G	341	0.255	402	0.293	4.977	0.026*	0.107
	A-T-G	234	0.175	258	0.188	0.763	0.382	0.907
	A-T-A	198	0.148	153	0.112	7.962	0.005**	0.020*
	G-C-G	9	0.007	24	0.017	6.243	0.012*	0.052
NPAS2	A-T-C-A^#	603	0.443	601	0.426	0.785	0.376	0.965
	C-T-C-A	193	0.141	236	0.167	3.478	0.062	0.326
	A-T-T-G	199	0.146	175	0.124	2.746	0.097	0.472
	A-T-C-G	129	0.095	126	0.089	0.25	0.617	0.999
	A-C-C-A	106	0.078	121	0.086	0.532	0.466	0.989
	C-C-C-A	38	0.028	54	0.038	2.353	0.125	0.572
	A-C-C-G	32	0.023	32	0.022	0.015	0.904	1.000
	C-T-T-G	23	0.017	27	0.019	0.122	0.726	1.000
	A-T-T-A	21	0.015	15	0.011	1.216	0.270	0.881

Order of SNPs in the haplotype for BMAL1: rs10766071 – rs1562438 – rs3789327, for NPAS2: rs3768984 – rs6738097 – rs11541353 – rs11123857. ^#Wildtype haplotypes. χ ² = Chi-square, P_corr-value = P-value after 10,000 permutations, *P-value < 0.05, **P-value < 0.01, ***P-value < 0.001.

Combinatorial multiallelic risk analysis

Genotype data from this study was combined with raw data from three previous studies on CLOCK, CRY and PER genes, genotyped in the same Swedish cohort and reanalyzed for the purpose of performing a combinatorial multiallelic risk analysis.^19,20,22 The total number of included participants from our prior circadian gene studies and additional genotyped samples from the discovery cohort are specified in Online Supplementary Table S1.^19,20,22 All UK participants from a prior GWAS were included as a replication cohort.²³ Risk alleles for each SNP in the extended population were labelled based on differential allele frequency in the control and CH groups (Figure 1; see also Online Supplementary Table S2) in the discovery cohort. Virtually all subjects carried at least three risk alleles (99.9% of controls and 100% of CH group), and 90% of subjects carried at least five risk alleles (Figure 2A). The number of risk alleles per subject was approximately normally distributed (Figure 2A), with a higher number of risk alleles per subject in individuals with CH (median: 9) as compared to controls (median: 8) (p < 0.05, two-sample t-test with unequal variance). The analysis of risk for CH associated with combinations of up to three risk alleles identified 258 combinations which survived FDR correction (Figure 2B). Five SNPs appear significant independently (Figure 2C), all located in core clock genes included in the positive feedback arm (i.e., BMAL1, CLOCK, NPAS2) of the molecular clock TTFL. Combinations of two or three risk alleles yielded generally higher OR than any single SNP included (Figure 2C; see also Online Supplementary Figure S1). Notably, rs11541353 (NPAS2) was also included in combinations with OR below the OR as single risk allele (Figure 2C). The analysis of rate of inclusion for individual risk alleles identified three enriched SNPs (rs1562438, rs11123857, and rs8192440) in addition to the SNPs with significant risk independently. The preliminary analysis of change in risk from single risk alleles to pairs (Online Supplementary Figure S2) indicated a consistent increase in OR when information from an additional SNP is included in the analysis. Therefore, we further analyzed all combinations surviving FDR correction to investigate patterns of associations in relation to protein function. The purpose of the analysis was to investigate potential functional implications of combination clustering, namely whether the combinations tend to include predominantly risk alleles located in core clock genes in the positive or negative arms of the molecular clock TTFL. The parallel visualization of number of combinations including specific pairs of risk alleles and the change in OR for significant combinations of three risk alleles (Figure 3) highlights “driver pairs”, namely combinations of two SNPs which carry significant risk in combinations with another SNP. Driver pairs were distributed in quadrants defined by the contribution of genes to either arm of the TTFL: top left – driver pairs in the positive TTFL arm; top right and bottom left (symmetrical across the main diagonal) – drive pairs made of SNPs in both positive and negative arms; bottom right – driver pairs in the negative arm of the TTFL. The driver pairs were distributed evenly across the quadrants (>80% of all driver pairs appeared in significant combinations. However, the analysis of impact of adding the third SNP (Figure 3) showed that additional SNPs had rather modest effect on driver pairs in the positive arm of the TTFL (top left quadrant: numerous significant combinations, but small change in OR), but large impact when the driver pairs included SNPs in the negative arm of the TTFL (bottom quadrants: relatively lower number of combination with large impact on OR).

Figure 1.

Distribution of risk alleles in the discovery cohort. (A) Raster plot illustration of carrier status for risk alleles (white indicates non carrier). (B) Risk allele frequency in controls (left, in purple) and cluster headache (right, in orange).

Figure 2.

Multiple risk allele distribution. (A) Distribution of number of risk alleles per subject. Note that it is rather uncommon for a subject to carry less than three alleles (∼1.4%). Dashed lines depict median number of risk alleles per subject. CH = cluster headache (B) Volcano plot depicting the risk for developing cluster headache for combinations of up to three risk alleles. Note the high numbers of combinations of risk alleles which yield significant risk (red dots) as compared to individual risk alleles (blue dots). (C) Illustration of difference between risk increase due to individual single nucleotide polymorphisms (SNPs) (blue dots) vs. combinations including each specific SNP (red dots). Note that it is very rare to have a SNP included in combinations with lower OR than the OR of the SNP alone (e.g., rs11541353). The rate of inclusion of individual SNPs in combinations yielding significant risk increase was quantified using the posterior inclusion probability (PIP, bottom panel). PIP higher than expected by random sampling (1/18) largely coincides with SNPs which increase the risk independently (indicated by asterisks). However, SNPs not significant individually emerge as frequently included in combinations yielding significant risk increase (rs1562438, rs8192440, and rs11123857).

Figure 3.

Analysis of combinations of up to three risk alleles. The single nucleotide polymorphisms (SNPs) are displayed in the same order as in Figure 2, functional labelling of core clock genes displayed for simplicity: green – positive feedback arm; red – negative feedback arm. Left panel: occurrence of specific driver pairs of SNPs in significant combinations. Right panel: median change in odds ratio (OR) when a third risk allele is added to the driver pair panel) associated with large change in OR (intense color in right panel) indicate weak driver pairs, where adding a third risk allele yields a major change in OR (e.g., #4 in NPAS2 and #3 in BMAL1).

To validate the findings in the Swedish cohort we used data extracted from a GWAS study on CH from the National Hospital for Neurology and Neurosurgery, London, UK and headache-free controls from the Wellcome Trust Case Control Consortium (the 1958 birth cohort), and the National Blood Survey cohort.²³ The risk allele labelling from the discovery cohort applied on the replication cohort yielded good concordance in allele frequency across the cohorts (Figure 4A). The number of risk alleles per subject was higher than in the discovery cohort, with a median number of 10 per subject in both groups and 90% of subjects carrying at least eight risk alleles (Online Supplementary Figure S3A). The analysis of risk for single SNPs highlighted concordant increased risk for the mutated allele (C variant) for rs3768984 (NPAS2) (Figure 4B). The analysis of combinations of up to three risk alleles identified 205 combinations surviving FDR correction (Online Supplementary Figure S3B). The comparison between combinations surviving FDR correction in both cohorts identified 39 concordant combinations (Figure 4C). The analysis of individual risk allele contribution to concordant combinations confirmed the enrichment for six SNPs in the discovery cohort, located mostly in clock genes in the positive arm of the TTFL (Figure 4D). Notably, rs3768984 (NPAS2) was included in 31 concordant combinations (Figure 4D; see also Online Supplementary Figure S4). The distribution of PIP apparently favored risk alleles located in BMAL1, CLOCK and NPAS2 in both discovery and replication cohorts. Lastly, the combinations were classified in four categories based on the proportion of SNPs located in the positive arm of the molecular clock TTFL (Figure 4E). The analysis of composition of combinations found that SNPs in the positive arm of the molecular clock TTFL account for at least 50% of the risk alleles (i.e., 1 out 2, or 2 out of 3) in >80% of concordant combinations, which is consistent with the proportion in the discovery cohort (Figure 4E; see also Online Supplementary Figure S4). A sensitivity analysis with increasingly stringent FDR threshold (0.25–0.05) yielded increasingly higher enrichment in combinations dominated by SNPs in the positive arm of the TTFL in the discovery cohort, while the enrichment remained stable in the concordant combinations (Online Supplementary Table S3).

Figure 4.

Validation of single nucleotide polymorphism (SNP) combinations using replication cohort. (A) Correlation between frequency of risk alleles (as defined in the discovery cohort) in the two cohorts. Alleles with frequency over 0.5 depict SNPs where the minor/mutated allele is protective. (B) Comparison between risks for single SNPs. (C) Concordant SNP combinations. Note that 36 out of 39 combinations with significant risk in both cohorts increase the risk to develop cluster headache. (D) Probability of inclusion in concordant combinations. Individual SNPs which appear more often than expected by random selection (enriched) are highlighted by darker colored columns. (E) Functional classification of SNP combinations based on transcription-translation feedback loop (TTFL) arms. Note that SNP combinations located in the positive arm of the TTFL dominate in both discovery and replication cohorts, as well as in the concordant combinations.

Considering these findings, we asked whether perceived diurnal rhythmicity of CH attacks is related to the genetic background of the subjects. The comparison between individuals with CH reporting diurnal rhythmicity and controls found largely similar patterns as when the entire cohort was analyzed (see Online Supplementary Figure S5A and B), and the composition of combinations of risk alleles surviving FDR correction displayed a similar bias towards SNPs located in clock genes in the positive arm of the molecular clock TTFL. Remarkably, the direct comparison between individuals with CH with and without reported diurnal rhythmicity yielded no significant differences.

Discussion

Our results complement previous studies investigating clock genes and their potential implication in CH. Among the new SNPs in BMAL1 and NPAS2 investigated in this study, two (rs3789327 and rs3768984) were found to be associated with CH after correction for multiple testing and were further highlighted in a haplotype analysis. The combinatorial multiallelic risk analysis including all additional core clock gene SNPs previously investigated in our Swedish CH cohort brings new perspectives on clock genes and their potential contribution to CH.^19,20,22 Briefly, SNPs located in core clock genes involved in the positive feedback arm of the TTFL are enriched in individuals with CH. In addition, they seem to act as drivers overriding combinations with other SNPs and show consistent association with CH, as verified in a replication cohort.

Sleep disturbance can be a manifestation of circadian disruption, and the phenotype has been consistently documented in experimental models where gene expression was silenced, such as BMAL1-knockout macaques or NPAS2-deficient mice.^27,28 In human populations, sleep disturbances have also been found to be associated with specific SNPs. For instance, rs3768984, in NPAS2, was significantly associated with changes in sleep and wake onset times. Interestingly, similar phenotypes have been reported recently in individuals with CH, namely longer sleep latency and lower sleep efficiency.¹⁵ Sleep variations were also observed in remission phases, which strongly suggest that the circadian/sleep phenotype is driven by underlying biological mechanisms, and not by the headache attacks, therefore hinting at a genetic contribution.¹⁵ Circadian disruption has also been associated with psychiatric disorders. Genetic manipulations of BMAL1 in animal models yield behaviors classically linked to psychiatric disorders. Depending on deletion type and/or location of the variation, these consequences resemble symptoms of BD, schizophrenia, and depression.^10,12,29 In the present study rs3789327, previously associated with BD,³⁰ was identified as a CH risk-allele and was confirmed by haplotype analysis.

The analysis of risk associated with single SNPs has some limitations. Firstly, no genetic marker in core clock genes is known to be pathogenic for CH, and the functional impact of the analyzed markers are not known. When the present results were analyzed together with genotype data from previous studies on the same cohort in the combinatorial multiallelic risk analysis,^19,20,22 we found that all study participants carried at least three risk alleles (in both the control and CH groups). Consequently, the lack of association of clock-related genes in prior GWASs may be explained, at least in part, by the focus on the initial individual SNP analysis instead of multiallelic or combinatorial analysis.¹⁶ This raised the question whether combinations of risk alleles may yield patterns meaningful from a functional perspective. The conspicuous enrichment of risk alleles located in clock genes in the positive arm of the TTFL points towards functional alterations in timekeeping in individuals with CH. Specifically, the timekeeping function following circadian entrainment is supposedly altered, which can explain the increase in headache attacks around the equinoxes, i.e., when daylight length changes fastest.^7,31 In addition, the circannual pattern of CLOCK expression was found to be blunted in CH subjects as compared to controls.³² At a molecular level, effects downstream of BMAL1 activation include the regulation of melatonin synthesis. The circadian rhythm of melatonin secretion is controlled by the rate-limiting enzyme aralkylamine N-acetyltransferase (AA-NAT), with peak expression at night.³³ AA-NAT gene expression is regulated by BMAL1/NPAS2 heterodimers binding to an E-box motif.³⁴ In individuals with CH, the concentration of melatonin metabolites in urine is lower than in controls regardless of disease activity, while serum melatonin has been reported to be lower during cluster periods than in remission.^35,36

We hypothesize a mechanism linked to the circadian rhythm which modulates the internal environment in a way that allows the attacks to occur and may constitute an interesting pharmacological target. In the case of CH, commonly used preventive treatments interact with many different molecular targets or receptors. The molecular clock is therefore an alternative activity with treatments interacting on different levels. One example is the most used preventive treatment: verapamil, a calcium channel blocker.³⁷ Calcium signaling dictates a variety of cellular signaling pathways, which includes the molecular clock.^9,38 Different publications show the effect of verapamil on circadian modulation^39,40 although no clear molecular mechanism has been uncovered yet.⁴¹ Lithium is prescribed as a last resort in chronic CH patients unresponsive to other treatments due to high toxicity. Interestingly, the pharmacodynamic profile of lithium is consistent with reinforcing the positive feedback arm of the TTFL. The phosphorylation of BMAL1 by glycogen synthase kinase 3 beta (GSK3b) initiates ubiquitylation and subsequent proteasomal degradation.^42,43 Lithium inactivates GSK3b,⁴⁴ which leads to enhanced coupling between the feedback loops within the TTFL and yields more robust oscillations in clock gene expression.^42,45 One can further speculate that preventing proteasomal degradation of cytosolic BMAL1 by melatonin may stabilize the molecular clock machinery.⁴⁶ However, evidence for positive effects of melatonin in CH patients is limited to date.^37,47,48

The main strength of this study is its large sample size and the well characterized cohort consisting of Swedish study participants with verified CH diagnoses. By contrast, controls consist of healthy blood donors, which is a limitation as some of them could potentially have CH or develop CH in the future. However, considering the prevalence of CH, the number of undetected cases in the control group should be minimal. Since the genotyping was executed on a Swedish cohort, a replication study was performed in an independent CH cohort from the UK which further strengthened our results. Nevertheless, it would be valuable to investigate if these results can be replicated in cohorts from other geographical regions, and to further expand our scope to include other genetic and functional markers in these core-clock genes.

Conclusion

More studies are needed to understand the circadian component of CH using both molecular and physiological approaches. Here we have used genetic screening as a tool and show a link between two new candidate genes and CH, both have a major contribution to the molecular mechanisms driving the circadian rhythm. It remains to be seen if the circadian component is a driver in the pathophysiological events leading to CH, in which case the molecular clock might constitute a therapeutic target in CH.

Article highlights

Genetic markers in BMAL1 and NPAS2 were found to be associated with an increased risk for cluster headache

The combinations of risk alleles carried by individuals with cluster headache are enriched in the positive feedback arm of the TTFL

The analysis of risk associated with combination of non-pathogenic genetic markers can be used for identifying new treatment targets

Supplemental Material

sj-docx-1-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-docx-1-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Supplemental Material

sj-png-2-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-png-2-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Supplemental Material

sj-png-3-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-png-3-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Supplemental Material

sj-png-4-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-png-4-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Supplemental Material

sj-png-5-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-png-5-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Supplemental Material

sj-png-6-cep-10.1177_03331024261440136 - Supplemental material for Genetic variability within molecular core clock genes in cluster headache

Supplemental material, sj-png-6-cep-10.1177_03331024261440136 for Genetic variability within molecular core clock genes in cluster headache by Clémence Deborgies Sanches, Stefan Spulber, Felicia Jennysdotter Olofsgård, Carmen Fourier, Anna Sundholm, Maria Lantz, Christina Sjöstrand, Elisabet Waldenlind, Anna Steinberg, Henry Houlden, Manjit Matharu, Caroline Ran and Andrea Carmine Belin in Cephalalgia

Footnotes

Acknowledgements

We thank Ann-Christin Karlsson for help with study participant recruitment.

ORCID iDs

Clémence Deborgies Sanches

Stefan Spulber

Felicia Jennysdotter Olofsgård

Carmen Fourier

Anna Sundholm

Maria Lantz

Christina Sjöstrand

Elisabet Waldenlind

Anna Steinberg

Henry Houlden

Manjit Matharu

Caroline Ran

Andrea Carmine Belin

Ethical considerations

Ethical permit was approved by the Regional Ethical Review Board in Stockholm (registration number 2014/656 − 31/4).

Consent to participate

Blood samples for the genotyping in the present study were collected from 682 healthy controls and 707 individuals with CH residing in Sweden, after written informed consent was obtained.

Author contributions

Conceptualization: C.D.S, C.R., A.C.B. Data curation: C.D.S, S.S. Formal analysis: C.D.S., S.S. Funding acquisition; A.C.B. Investigation: C.D.S, S.S., F.J.O, C.F, A.Su., M.L., C.S., E.W., A.St, H.H., M.M. Methodology: C.D.S, S.S., C.R. Project administration: A.C.B. Resources: C.D.S., S.S., F.J.O., C.F., A.Su., M.L., C.S., E.W. A.St, H.H., M.M. Software: C.D.S., S.S. Supervision: C.R., A.C.B. Validation: C.R., A.C.B. Visualization: C.D.S., S.S., C.R., A.C.B. Writing – original draft: C.D.S., S.S., C.R., A.C.B. Writing – review & editing: C.D.S., S.S., F.J.O., C.F., A.Su., M.L., C.S., E.W., A.St, H.H., M.M., C.R., A.C.B.

Funding

This research work was supported by the Mellby Gård Foundation, the Swedish Brain Foundation (FO2025-0312), Swedish Research Council (2025-02368), Karolinska Institutet Research Funds (2024-02448) Region Stockholm (ALF project) (FoUI-1019690) and The Centre for Neuromusculoskeletal Restorative Medicine (CNRM) (SS).

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Data availability statement

Raw data files will be made available after reasonable request to the corresponding author.

Supplemental material

Supplemental material for this article is available online.

References

Headache Classification Committee of the International Headache Society (IHS). The international classification of headache disorders, 3rd edition. Cephalalgia 2018; 38: 1–211.

Fourier

Ran

Steinberg

, et al. Sex differences in clinical features, treatment, and lifestyle factors in patients with cluster headache. Neurology 2023; 100: e1207–e1220.

Benkli

Kim

Koike

, et al. Circadian features of cluster headache and migraine: a systematic review, meta-analysis, and genetic analysis. Neurology 2023; 100: e2224–e2236.

May

Bahra

Büchel

, et al. Hypothalamic activation in cluster headache attacks. Lancet 1998; 352: 275–278.

Welsh

Takahashi

Kay

. Suprachiasmatic nucleus: cell autonomy and network properties. Annu Rev Physiol 2010; 72: 551–577.

Pilati

Torrente

Alonge

, et al. Sleep and chronobiology as a key to understand cluster headache. Neurol Int 2023; 15: 497–507.

Kim

Lazar

. Transcriptional control of circadian rhythms and metabolism: a matter of time and space. Endocr Rev 2021; 41: 707–732.

DeBruyne

Weaver

Reppert

. CLOCK and NPAS2 have overlapping roles in the suprachiasmatic circadian clock. Nat Neurosci 2007; 10: 543–545.

Hastings

Maywood

Brancaccio

. Generation of circadian rhythms in the suprachiasmatic nucleus. Nat Rev Neurosci 2018; 19: 453–469.

10.

Bengesser

Reininghaus

Lackner

, et al. Is the molecular clock ticking differently in bipolar disorder? Methylation analysis of the clock gene ARNTL. World J Biol Psychiatry 2018; 19: S21–S29.

11.

Johansson

Willeit

Smedh

, et al. Circadian clock-related polymorphisms in seasonal affective disorder and their relevance to diurnal preference. Neuropsychopharmacology 2003; 28: 734–739.

12.

Christiansen

Bouzinova

Fahrenkrug

, et al. Altered expression pattern of clock genes in a rat model of depression. Int J Neuropsychopharmacol 2016; 19: pyv135.

13.

Walker

Walton

DeVries

, et al. Circadian rhythm disruption and mental health. Transl Psychiatry 2020; 10: 28.

14.

Robbins

. The psychiatric comorbidities of cluster headache. Curr Pain Headache Rep 2013; 17: 13.

15.

Ran

Jennysdotter Olofsgård

Steinberg

, et al. Patients with cluster headache show signs of insomnia and sleep related stress: results from an actigraphy and self-assessed sleep study. J Headache Pain 2023; 24: 14.

16.

Winsvold

Harder

AVE

Ran

, et al. Cluster headache genomewide association study and meta-analysis identifies eight loci and implicates smoking as causal risk factor. Ann Neurol 2023; 94: 713–726.

17.

Robbins

. The bipolar spectrum in migraine, cluster and chronic tension headache. Eur Neurol Rev 2008; 3: 123.

18.

Costa

Squassina

Piras

, et al. Preliminary transcriptome analysis in lymphoblasts from cluster headache and bipolar disorder patients implicates dysregulation of circadian and serotonergic genes. J Mol Neurosci 2015; 56: 688–695.

19.

Fourier

Ran

Zinnegger

, et al. A genetic CLOCK variant associated with cluster headache causing increased mRNA levels. Cephalalgia 2018; 38: 496–502.

20.

Fourier

Ran

Sjöstrand

, et al. The molecular clock gene cryptochrome 1 (CRY1) and its role in cluster headache. Cephalalgia 2021; 41: 1374–1381.

21.

Ofte

Tronvik

Alstadhaug

. Lack of association between cluster headache and PER3 clock gene polymorphism. J Headache Pain 2016; 17: 18.

22.

Jennysdotter Olofsgård

Ran

Fourier

, et al. PER Gene family polymorphisms in relation to cluster headache and circadian rhythm in Sweden. Brain Sci 2021; 11: 1108–1117.

23.

O’Connor

Fourier

Ran

, et al. Genome-wide association study identifies risk loci for cluster headache. Ann Neurol 2021; 90: 193–202.

24.

Machiela

Chanock

. LDlink: a web-based application for exploring population-specific haplotype structure and linking correlated alleles of possible functional variants. Bioinformatics 2015; 31: 3555–3557.

25.

Purcell

Neale

Todd-Brown

, et al. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet 2007; 81: 559–575.

26.

Dupont

Plummer

. Power and sample size calculations A review and computer program. Control Clin Trials 1990; 11: 116–128.

27.

Qiu

Jiang

Liu

, et al. BMAL1 Knockout macaque monkeys display reduced sleep and psychiatric disorders. Natl Sci Rev 2019; 6: 87–100.

28.

Dudley

Erbel-Sieler

Estill

, et al. Altered patterns of sleep and behavioral adaptability in NPAS2-deficient mice. Science 2003; 301: 379–383.

29.

Mansour

Wood

Logue

, et al. Association study of eight circadian genes with bipolar I disorder, schizoaffective disorder and schizophrenia. Genes Brain Behav 2006; 5: 150–157.

30.

Maciukiewicz

Dmitrzak-Weglarz

Pawlak

, et al. Analysis of genetic association and epistasis interactions between circadian clock genes and symptom dimensions of bipolar affective disorder. Chronobiol Int 2014; 31: 770–778.

31.

Riederer

Schankin

. Circannual worsening in cluster headache correlates with daylight change. Headache 2024; 64: 1336–1338.

32.

Oliveira

Neves-Costa

Pedroso

, et al. CLOCK gene circannual expression in cluster headache. Cephalalgia 2024; 44: 3331024241303592.

33.

Simonneaux

Ribelayga

. Generation of the melatonin endocrine message in mammals: a review of the complex regulation of melatonin synthesis by norepinephrine, peptides, and other pineal transmitters. Pharmacol Rev 2003; 55: 325–395.

34.

Chong

Bernard

Klein

. Characterization of the chicken serotonin N-acetyltransferase gene: activation via clock gene heterodimer/e box interaction. J Biol Chem 2000; 275: 32991–8.

35.

Waldenlind

Gustafsson

Ekbom

, et al. Circadian secretion of cortisol and melatonin in cluster headache during active cluster periods and remission. J Neurol Neurosurg Psychiatry 1987; 50: 207–213.

36.

Liampas

Siokas

Brotis

, et al. Meta-analysis of melatonin levels in cluster headache—review of clinical implications. Acta Neurol Scand 2020; 142: 356–367.

37.

Wei

Goadsby

. Cluster headache pathophysiology — insights from current and emerging treatments. Nature Reviews Neurology 2021; 17: 308–324.

38.

Colwell

. Linking neural activity and molecular oscillations in the SCN. Nat Rev Neurosci 2011; 12: 553–569.

39.

Burish

Han

Mawatari

, et al. The first-line cluster headache medication verapamil alters the circadian period and elicits sex-specific sleep changes in mice. Chronobiol Int 2021; 38: 839–850.

40.

Barloese

Haddock

Lund

, et al.

Chronorisk in cluster headache: a tool for individualised therapy?

Cephalalgia 2018; 38: 2058–2067.

41.

Petersen

Barloese

MCJ

Snoer

, et al. Verapamil and cluster headache: still a mystery. A narrative review of efficacy, mechanisms and perspectives. Headache 2019; 59: 1198–1211.

42.

Sahar

Zocchi

Kinoshita

, et al. Regulation of BMAL1 protein stability and circadian function by GSK3β-mediated phosphorylation. PLoS One 2010; 5: e8561.

43.

Paul

Johnson

Jope

, et al. Disruption of circadian rhythmicity and suprachiasmatic action potential frequency in a mouse model with constitutive activation of glycogen synthase kinase 3. Neuroscience 2012; 226: 1–9.

44.

Klein

Meltont

. A molecular mechanism for the effect of lithium on development. Proc Natl Acad Sci U S A 1996; 93: 8455–8459.

45.

Sawai

Okamoto

Muranaka

, et al. In vivo evaluation of the effect of lithium on peripheral circadian clocks by real-time monitoring of clock gene expression in near-freely moving mice. Sci Rep 2019; 9: 10909.

46.

Vriend

Reiter

. Melatonin feedback on clock genes: a theory involving the proteasome. J Pineal Res 2015; 58: 1–11.

47.

Smedfors

Jennysdotter Olofsgård

Steinberg

, et al. Use of prescribed and non-prescribed treatments for cluster headache in a Swedish cohort. Brain Sci 2024; 14: 48.

48.

Brandt

Doesborg

PGG

Haan

, et al. Pharmacotherapy for cluster headache. CNS Drugs 2020; 34: 171–184.