Authors: Jaana Simola, Felix Siebenhühner, Vladislav Myrov, Katri Kantojärvi, Tiina Paunio, J. Matias Palva, Elvira Brattico, Satu Palva
Categories: Article, Biological sciences, Cognitive neuroscience, Neuroscience
Source: iScience
Neuronal oscillations, their inter-areal synchronization, and scale-free dynamics constitute fundamental mechanisms for cognition by regulating communication in neuronal networks. These oscillatory dynamics have large inter-individual variability that is partly heritable. We hypothesized that this variability could be partially explained by genetic polymorphisms in neuromodulatory genes. We recorded resting-state magnetoencephalography (MEG) from 82 healthy participants and investigated whether oscillation dynamics were influenced by genetic polymorphisms in catechol-O-methyltransferase (COMT) Val^158^Met and brain-derived neurotrophic factor (BDNF) Val^66^Met. Both COMT and BDNF polymorphisms influenced local oscillation amplitudes and their long-range temporal correlations (LRTCs), while only BDNF polymorphism affected the strength of large-scale synchronization. Our findings demonstrate that COMT and BDNF genetic polymorphisms contribute to inter-individual variability in neuronal oscillation dynamics. Comparison of these results to computational modeling of near-critical synchronization dynamics further suggested that COMT and BDNF polymorphisms influenced local oscillations by modulating the excitation-inhibition balance according to the brain criticality framework.
**Subject ** Biological sciences, Neuroscience, Cognitive neuroscience
Neuronal oscillations and their inter-areal synchronization play fundamental mechanistic roles in cognitive functions, regulating neuronal processing and communication across distributed brain areas (Fries, 2015; Singer, 1999; Womelsdorf et al., 2014). In humans, dynamics of neuronal oscillations and their long-range phase synchronization, observed non-invasively in magneto- (MEG) and electroencephalography (EEG) recordings, serve fundamental roles in a variety of sensory, motor, and cognitive functions (Fell and Axmacher, 2011; Jensen et al., 2014; Klimesch et al., 2007; Marzetti et al., 2019; Palva and Palva, 2012; Samaha et al., 2020; Siegel et al., 2012; Thut et al., 2012). During rest, brain activity involves frequency-specific inter-areal correlation structures throughout the 1–100 Hz range (Arnulfo et al., 2020; Betti et al., 2018; Brookes et al., 2011; Foster et al., 2016; Hipp and Siegel, 2015; Marzetti et al., 2013; Siebenhühner et al., 2020; Zhigalov et al., 2017). The intrinsic spatial organization of spontaneous resting-state functional connectivity (FC) observed with MEG (Betti et al., 2018; Oswald et al., 2017; Siebenhühner et al., 2020), intra-cranial EEG (Arnulfo et al., 2020), and functional magnetic resonance imaging (fMRI) (Hahn et al., 2020; Vidaurre et al., 2021) predicts individual variability in task performance. Resting-state FC networks reflect individual trait-like behavior (Seitzmana et al., 2019) that could arise from the underlying individual brain structural organization (Cabral et al., 2014; Mostame and Sadaghiani, 2020) and neuromodulation (McCormick et al., 2020; Pfeffer et al., 2021; van den Brink et al., 2018; van den Brink et al., 2019).
Amplitude envelopes of fast (>1 Hz) oscillations are also characterized by slow (0.1–1 Hz) and infra-slow (0.01–0.1 Hz) fluctuations (Hiltunen et al., 2014; Hipp and Siegel, 2015; Monto et al., 2008). These fluctuations are scale-free and exhibit power-law autocorrelations in their strength known as long-range temporal correlations (LRTCs) (Linkenkaer-Hansen et al., 2001; Monto et al., 2008; Palva et al., 2013; Poil et al., 2012; Thiery et al., 2018; Zhigalov et al., 2015). In addition to instantaneous oscillation dynamics, inter-individual variability of LRTCs also predicts variability in task-related behavioral LRCTs (Palva et al., 2013; Smit et al., 2013). Scale-free LRTCs are hallmarks of systems with critical or near-critical dynamics, suggesting that the brain operates near a critical phase transition between disorder and order, i.e., between inadequate and excessive synchronization, respectively (Cocchi et al., 2017), which enables high flexibility and variability.
A large fraction of the variability in fast and slow oscillation dynamics could be heritable. MEG spectral power shows similarities between siblings (Leppäaho et al., 2019; Salmela et al., 2016). Twin studies indicate that frontal midline theta activity is genetically related to symptoms and response time variability in attention-deficit/hyperactivity disorder (McLoughlin et al., 2014) and that peak frequencies (Van Pelt et al., 2012) as well as LRTCs (Linkenkaer-Hansen et al., 2007) of local cortical oscillations have a genetic basis, so that up to 60% of the variance in neuronal LRTCs (Linkenkaer-Hansen et al., 2007) and up to 85% of the variance in cognitive ability (Posthuma et al., 2001) can be attributed to genetic factors. However, what these genetic factors are and how they could contribute to the inter-individual variability of synchronization dynamics, including oscillation amplitudes, large-scale synchronization, and LRTCs, have remained central unresolved questions.
We hypothesized that polymorphism in genes regulating neuromodulation is one genetic factor underlying variability in oscillatory dynamics. Compelling evidence indicates that especially the neuromodulatory effects of the catecholamines dopamine (DA) and noradrenaline (NE) (Cools, 2019; Pfeffer et al., 2021; Robbins and Arnsten, 2009; van den Brink et al., 2018) and the indolamine serotonin (5-hydroxytryptamine; 5-HT) (Harris and Thiele, 2011; Lucki, 1988; Štrac et al., 2016) contribute to variability in neuronal dynamics (Ferguson and Cardin, 2020; Harris and Thiele, 2011; Lee and Dan, 2012; van den Brink et al., 2019).
Brain DA and NE levels are regulated by their catabolism via COMT enzyme (Männistö and Kaakkola, 1999; Parkin et al., 2018) whose activity is influenced by a common SNP in the COMT gene resulting in variation substituting valine (Val) to methionine (Met) at codon 158 Val^158^Met (captured as rs4680). COMT Val^158^ homozygotes have higher COMT activity and therefore lower brain DA and NE levels than Met^158^ homozygotes, with activity for heterozygotes falling in the middle (Chen et al., 2004). A considerable fraction of the variance in serotonin levels, on the other hand, is explained by Val^66^Met (rs6265) polymorphism, resulting in a valine (Val) to methionine (Met) substitution at codon 66, of the gene coding for brain-derived neurotrophic factor (BDNF), which influences 5-HT receptor binding (Fisher et al., 2017) and the activity-dependent secretion of BDNF (Egan et al., 2003). Convincing evidence emerging from numerous fMRI studies indicates that COMT Val^158^Met (rs4680) polymorphism is associated with modulation of signals especially in the prefrontal cortex (PFC) and with performance of cognitive functions dependent on PFC (Barnett et al., 2007; Meyer-Lindenberg et al., 2005; Mier et al., 2010; Papenberg et al., 2020; Robbins and Arnsten, 2009), while BDNF Val^66^Met polymorphism is implicated in many aspects of behavior and cognition (Beste et al., 2010; Matsuo et al., 2009; Notaras et al., 2015). In addition to modulation of 5-HT levels, activity-dependent release of BDNF also influences the maintenance and formation of neuronal networks (Park and Poo, 2013; Zagrebelsky and Korte, 2014), thereby influencing cognition and behavior over lifetime.
Although the effects of COMT Val^158^Met and BDNF Val^66^Met polymorphisms on slow neuronal dynamics have been well-studied with fMRI, their influence on fast neuronal oscillations has remained uncharted. We hypothesized that COMT Val^158^Met and BDNF Val^66^Met polymorphisms could influence inter-individual variability in oscillation amplitudes, phase synchronization, and their LRTCs. Importantly, critical dynamics in the brain are thought to be controlled primarily by the excitation-inhibition (E/I) ratio so that critical dynamics emerge at balanced E/I while excess of inhibition or excitation leads to sub- or super-critical dynamics, respectively (Plenz and Thiagarajan, 2007; Poil et al., 2012; Shew et al., 2009). We thus further hypothesized that these polymorphisms would specifically exert their influence on critical brain dynamics via changes in the E/I ratio.
Genome-wide association studies (GWAS) have demonstrated the clinical relevance of genes on functional and structural MRI data (Pereira et al., 2018) and for local oscillatory activity in the brain using sensor-level EEG and MEG (Salmela et al., 2016; Smit et al., 2018). Since our specific interest here was on the neuromodulation-related genetic influences of polymorphisms in COMT Val^158^Met and BDNF Val^66^Met candidate genes on neuronal oscillations, we did not use a GWAS approach (which would require a huge sample size N >1000) nor polygenic risk scores (which do not exist for oscillation dynamics). To test our hypotheses, we measured resting-state brain dynamics with MEG and estimated local oscillation amplitudes, large-scale phase synchronization, and LRTCs from individually source-reconstructed MEG data. We then investigated their correlation with the COMT Val^158^Met and BDNF Val^66^Met polymorphisms estimated from SNP analysis (Figure 1) and compared the results with predictions from computational modeling. In line with the predictions of the brain criticality framework and computational modeling, we found that the COMT Val^158^Met and BDNF Val^66^Met polymorphisms influenced local oscillation amplitudes and LRTCs, while BDNF Val^66^Met polymorphism also influenced the strength of large-scale synchronization.
Figure 1 A schematic of study outline and rationale(A) Genotyping of the participants (N = 82) was used to identify the COMT Val^158^Met (rs4680) and BNDF Val^66^Met (rs6265) genetic polymorphisms.(B) Broadband MEG data were recorded from each participant with 306-channel MEG during ∼8-min resting state.(C) MEG data were source-reconstructed using minimum-norm estimation with individual MRI-based head and cortical surface anatomy and collapsed to cortical parcels (upper panel). Parcel-wise broadband signals (shown for inferior frontal gyrus, iFG, in second panel) were filtered to narrowband signals (signal and amplitude envelope shown here for the iFG in third panel). From the narrowband signal, instantaneous phase was extracted (shown for iFG and posterior central sulcus, poCS, in lower panel). Mean amplitudes and LRTCs were estimated for each parcel from amplitude envelope time series, and inter-parcel phase synchronization for all parcel pairs from phase time series.(D) In group analyses, we assessed the anatomical distribution of mean amplitudes, LRTCs, and synchronization within the polymorphism groups and identified statistically significant differences between these groups (Figure 2, Figure 3, Figure 4).
We first evaluated mean oscillation amplitudes between 3 and 60 Hz averaged across all cortical parcels separately for participants grouped according to the COMT Val^158^Met (Val/Val, Val/Met, and Met/Met) and the BDNF Val^66^Met (Val/Val and the combined Val/Met and Met/Met) polymorphisms. We first confirmed that age or the gender ratios did not differ significantly between the polymorphism groups (Tables 1 and 2). We then confirmed the validity of using the canonical frequency bands theta (θ, 3–7 Hz), alpha (α, 8–14 Hz), beta (β, 14–30 Hz), and gamma (γ, 30–60 Hz) by an analysis of spatial similarity across frequencies using Louvain clustering (Figure S1).
As observed commonly in resting state, oscillation amplitudes peaked in the α band in all COMT and BDNF polymorphism groups (Figure 2A). This was similarly evident when mean oscillation amplitudes were averaged within the canonical frequency bands (Figure S2). We found no influence of COMT Val^158^Met and BDNF Val^66^Met polymorphisms on the mean oscillation amplitudes in the canonical bands using univariate ANOVAs that controlled for age effects (p > 0.0125, Bonferroni-corrected across frequency bands) (Figure S2). We repeated the same analysis at the level of individual wavelet frequencies and found no significant effects there either (Figure 2A).
Figure 2 COMT and BDNF polymorphisms influence oscillation amplitudes(A) Oscillation amplitudes averaged across parcels for the COMT (Val/Val, purple, n = 18; Val/Met, turquoise; n = 48; and Met/Met, orange, n = 16) and BDNF polymorphism groups (Val/Val, purple, n = 66; Met carriers, orange, n = 16). Shaded areas represent 95% bootstrapped confidence intervals of the group means.(B) The fraction P of parcels with a significant positive or negative pairwise difference (shown above and below zero line, respectively) between the COMT and BDNF polymorphisms groups in oscillation amplitudes (Mann–Whitney-U test, p < 0.05). The light gray area indicates the 5% alpha level for findings expected by chance under the null hypothesis. The dark gray area indicates the Q-level which defines the chance level for false positives at p
q< 0.05 in any of the 26 frequencies.(C) Cortical parcels where the amplitudes differed significantly between polymorphism groups (indicated above the brain surfaces) in the θ (3–7 Hz), β–γ (20–60 Hz), and α to low-β (7–18 Hz) frequency ranges in at least one of the individual frequencies. The color scale indicates the max-normalized amplitude differences for COMT (red) and BDNF (blue) polymorphisms independently, so that for each effect the value 1 represents the maximum observed effect, and so that parcels where the effects are co-localized have purple hues. Abbreviations are formed from a combination of i) Areas: a = anterior, ca= calcarine, d = dorsal, i = inferior, la = lateral, m = middle, mrg = marginal, opc = opercular, p = posterior, s = superior, tr = triangular and ii) Lobes: C = central/cingulate, F = frontal, P = parietal, O = occipital, T = temporal, and iii) G = Gyrus, S = Sulcus. Additional subC = subcentral gyrus/sulcus, preCN = precuneus.
As the polymorphism effects could also be anatomically confined, we tested significant pairwise differences between polymorphism groups for each wavelet frequency and for each individual parcel. For COMT, the amplitudes were stronger (Mann-Whitney-U test, p < 0.05, FDR and Q-level corrected) for Val/Met heterozygotes than for either Met or Val homozygotes in the θ band (3–7 Hz) in up to 13% and 6% of cortical parcels, respectively (Figure 2B), whereas in the high-β to γ (20–60 Hz) bands, amplitudes were larger for Met homozygotes in 12% of parcels compared to Val/Met, and 6% compared to Val/Val. For BDNF, oscillation amplitudes were stronger for the Val homozygotes than for Met carriers in θ and in the α to low-β bands (up to 13% of parcels). We next visualized significant differences in oscillation amplitudes between the polymorphism groups in the cortical anatomy (Figure 2C). Stronger θ and α-low-β band oscillation amplitudes for the COMT Val/Met heterozygotes were localized to the central sulcus (CS), anterior, mid-, and posterior cingulate cortex (aCC, mCC, and pCC), inferior frontal gyrus (iFG), right-hemispheric middle frontal gyrus (mFG), and inferior parietal gyrus (iPG), inferior part of the insula (iINS), and precuneus (preCN) (Figure 2C) these parcels mostly belonging to the frontoparietal network (FPN) and default mode networks (DMN) (Figure S3). Of these, iPG and left iINS and left CS showed also larger amplitudes for BDNF Val homozygotes than for Met carriers. In the β–γ bands, larger amplitudes for COMT Met homozygotes compared to Val/Met heterozygotes were observed in superior precentral sulcus (spreCS), corresponding to frontal eye fields and in temporal cortices while stronger amplitudes for the Val homozygotes were found in aCC and pCC. For BDNF Val homozygotes, oscillation amplitudes were larger than in Met carriers in iFS, parieto-occipital sulcus (POS), and postcentral sulcus (poCS) in θ to low-β bands. COMT and BDNF jointly affected θ and α to low-β band oscillation amplitudes in the CS and iPG.
We quantified long-range temporal correlations (LRTCs) in oscillation amplitude fluctuations with detrended fluctuation analysis (DFA) (Hardstone et al., 2012) and estimated the mean DFA scaling exponents in all wavelet frequencies (3–60 Hz). DFA exponents were first averaged across all parcels separately for each COMT Val^158^Met and BDNF Val^66^Met polymorphism group for each frequency (Figure 3A). Using univariate ANOVAs that controlled for age, we observed a difference of group mean whole-brain DFA exponents between COMT polymorphism in the canonical γ band [F(2, 75) = 3.80, p = 0.027, ηp^2^ = 0.092] and for BDNF, Val homozygotes showed significantly higher exponents than Met carriers in the canonical α band [F(1, 76) = 4.19, p = 0.044, ηp^2^ = 0.052] (Figure S2). However, neither of these effects were significant after Bonferroni correction over frequency bands (p > 0.0125).
Figure 3 Influence of COMT and BDNF polymorphisms on LRTCs(A) Detrended fluctuation analysis (DFA) exponents averaged across parcels for the COMT and BDNF polymorphism groups. Shaded areas represent 95% bootstrapped confidence intervals of group means. Black bar with asterisk denotes the range of frequencies with significant differences between the COMT and BDNF polymorphism groups (univariate ANOVA, with age as covariate, ∗: p < 0.05, controlled for multiple comparisons with cluster-based permutation statistics).(B) Fractions (P) of parcels with a significant positive or negative difference between COMT and BDNF polymorphism groups in the DFA exponents (Mann–Whitney-U test, p < 0.05), as in Figure 2B.(C) Cortical areas where significant differences between polymorphism groups in the DFA exponents were found in the θ, β–γ, and α– low-β bands (Mann–Whitney U test, p < 0.05, corrected). Colors and abbreviations as in Figure 2C.
At the level of individual wavelet frequencies, the age-controlled univariate ANOVAs revealed for COMT a cluster of 3 consecutive frequencies from 35 to 45 Hz that remained significant after controlling for multiple comparisons with cluster-based permutation statistics (Fcluster = 11.64 >Fperm = 6.25 at α = 0.05). For BDNF, a significant cluster was found between 7 and 9 Hz (Fcluster = 13.08 >Fperm = 8.53 at α = 0.05) (Figure 3A).
Similar to oscillation amplitudes, COMT Val^158^Met and BDNF Val^66^Met polymorphisms had a more robust effect on DFA exponents at the resolution of cortical parcels indicating large anatomical heterogeneity. DFA exponents were greater for COMT Val/Met heterozygotes than for Met homozygotes (Mann-Whitney-U test, p<0.05, FDR and Q-level corrected) in the θ band (up to 20% of parcels) and in the β–γ band (up to 40%) (Figure 3B). Furthermore, θ band DFA exponents were greater for Val/Val than for Met/Met (16%) and β–γ DFA exponents were greater in Val/Met than for Val/Val (26%) (Figure 3B). For BDNF, a very robust effect of polymorphism was found for DFA exponents (Mann-Whitney-U test, p<0.05, FDR and Q-level corrected). These were greater for Val homozygotes than Met carriers in the θ band (35%) and in α band (68%) (Figure 3B). The localization of significant differences in DFA exponents between polymorphism groups (Mann-Whitney-U test, p<0.05, FDR and Q-level corrected, Figure 3C) revealed wide-spread group differences. In the θ band, large differences between COMT Val/Met heterozygotes and Met homozygotes and between BDNF Val homozygotes and Met carriers were co-localized to preCS, mCC, and mFG (Figure 3C). DFA exponents in θ band were also greater for COMT Val/Met group in INS, pCC, and in β–γ band particularly in parcels belonging to somatomotor (SM) network and DMN (see Figure S3). BDNF Val homozygotes had greater LRTC exponents in α and low-β bands than Met carriers in nearly all parcels (Figure 3C). Generally, the effects of COMT and BDNF polymorphisms on LRCTs were widespread, in line with prior work and with the idea that LRTCs reflect global fluctuations in neuronal activity across the brain’s modular structure (Monto et al., 2008; Palva et al., 2013; Zhigalov et al., 2017).
Finally, for each polymorphism, we co-localized its effects on both oscillation amplitudes and DFA exponents. We found that joint effects on amplitudes and exponents were strongest in the parcels of PFC and PPC, and CS and anterior cingulate structures (Figure S4).
We next set out to assess the influence of COMT Val^158^Met and BDNF Val^66^Met polymorphisms on inter-areal phase synchronization of neuronal oscillations. Phase synchrony is a key characteristic of individual cortical dynamics and has fundamental influence on behavior (Brookes et al., 2011; Fell and Axmacher, 2011; Fries, 2015; Palva and Palva, 2012; Singer, 1999; Vidaurre et al., 2021). We estimated phase synchrony among all cortical parcels with the weighted phase-lag index (wPLI) that is insensitive to linear source mixing (Palva et al., 2018) and computed the mean synchronization over all parcel pairs for each of the COMT Val^158^Met and BDNF Val^66^Met polymorphism groups (Figure 4A). Both within the canonical frequency bands and in individual wavelet frequencies, we did not find any significant differences (Kruskal-Wallis test, p > 0.05) among the COMT Val^158^Met polymorphism groups (Figures 4A and S2). Within canonical bands, we found that mean synchronization for BDNF Val homozygotes was significantly stronger (Kruskal-Wallis test, H = 6.15, p = 0.013) than for Met carriers in the α band (Figure S2). In individual wavelet frequencies, Val homozygotes exhibited stronger synchronization (Kruskal-Wallis test, p<0.05) in two frequency clusters, from 5 to 6.6 Hz and from 8.9 to 10.1 Hz (Figure 4A), which both remained significant after cluster-based permutation correction (Hcl,1 = 14.81 and Hcl,2 = 10.57 vs Hperm = 7.14 at α = 0.05).
Figure 4 BDNF polymorphisms influence large-scale phase synchronization(A) The mean strength of phase synchronization averaged across the whole brain for the COMT and BDNF polymorphism groups with 95% bootstrapped confidence intervals. Black bars with asterisks denote the frequencies where a significant group difference was found (Kruskal-Wallis, ∗: p < 0.05, controlled for multiple comparisons with cluster-based permutation statistics).(B) The differences in the mean phase synchrony (Δ-wPLI) between BDNF Val homozygotes and Met carriers within and between functional systems averaged across hemispheres and wavelet frequencies in canonical frequency bands. Red color indicates stronger synchronization for BDNF Val homozygotes and blue for Met carriers. Stars denote the brain system pairs where the group difference was significant (Kruskal-Wallis, ∗: p < 0.05, uncorrected, ∗∗: p < 0.01, Benjamini-Hochberg corrected). Abbreviations: DAN = Dorsal attention network, DMN = Default Mode network, FPN = Frontoparietal network, Lim = Limbic system, SM= Somatomotor network, VAN = Ventral attention network, Vis = Visual network.
To acquire neuroanatomical insight into the patterns of inter-areal synchronization in functional brain systems, we collapsed the parcel-wise synchronization estimates within functional systems (Yeo et al., 2011) (Figure S3) and canonical frequency bands (Figure S1). Synchronization between or within the subsystems did not differ significantly among the COMT polymorphisms (Figure S5). However, BDNF Val homozygotes exhibited stronger synchronization than Met carriers in the α band between and within nearly all subsystems (Figures 4B and S5, Kruskal-Wallis test, corrected with Benjamini-Hochberg). We also computed synchronization separately for each hemisphere and between hemispheres. Interestingly, group differences in α band synchronization were widespread and more robust for left- and inter-hemispheric connections than for the right-hemispheric connections (Figure S6). Additionally, β-band synchronization was stronger for BDNF Val homozygotes than Met carriers within and between several left-hemispheric networks, including FPN, dorsal attention network (DAN), ventral attention (VAN), and Limbic networks (Lim) and DMN (Figure S6).
We next set out to address the effects of COMT Val^158^Met and BDNF Val^66^Met polymorphisms in the context of brain criticality. Brain criticality is assumed to be regulated by a control parameter (K) that reflects the gross excitation/inhibition (E/I) ratio (Deco et al., 2021; Shew et al., 2009; Stephani et al., 2020). The human brain is thought to normally operate in a slightly subcritical regime (Priesemann et al., 2014), where LRTCs exhibit large inter-individual variability (Linkenkaer-Hansen et al., 2001; Palva et al., 2013; Smit et al., 2013; Zhigalov et al., 2015). We hypothesized that COMT Val^158^Met and BDNF Val^66^Met polymorphisms would bias the E/I balance and hence the “control parameter” that tunes the distance to criticality. The carriers of different variants should therefore be at different positions away from the critical point.
To assess the relationships of the observable brain dynamics measures with the control parameter, we implemented a nested Kuramoto model with separable local network (cortical parcel) and large-scale (whole-brain network) levels (Siebenhühner et al., 2020) of synchronization dynamics exhibiting a transition from sub- to super-critical dynamics (Figure 5A). We used the model to predict oscillation synchronization dynamics for two hypothetical cohorts with different values of the control parameter (K) regulating their position in the critical regime. Cohort A had lower and cohort B had higher control parameter (K) value positioning cohort B closer to the critical regime (Figure 5A) associated with stronger local amplitudes, LRTCs, and large-scale synchronization (Figure 5B). The model further predicted larger mutual correlations between DFA and oscillation amplitudes as well as between DFA and phase synchronization for cohort B with higher K (Figure 5C).
Figure 5 Model and experimental evidence for modulation of brain criticality(A) Computational model of local and inter-areal synchronization dynamics in the human connectome shows that both local amplitude (blue) and inter-areal synchronization (green) increase monotonically with increasing control parameter K and exhibit a phase transition between low- (subcritical) and high-synchronization (super-critical) regimes. Conversely, the DFA scaling exponent (yellow), a measure of long-range temporal correlations (LRTCs), peaks at the critical point. The model predicts that two cohorts (A and B, blue and red) with simulated subjects operating in different parts of the subcritical regime should exhibit group differences in amplitude, LRTCs, and inter-areal synchronization.(B) Local amplitude A, DFA exponents, and inter-areal phase synchronization (PLV) in modeled data for two cohorts A and B. All metrics are greater for the cohort B that is designed to operate closer to the critical point (see panel A) than cohort A. The box ends indicate the lower quartile (Q1, 25^th^ percentile) and upper quartile (Q3, 75^th^ percentile) across cortical regions/parcels, notches indicate the median, and whiskers indicate the range of Q1 – 1.5 ∗ IQR and Q3 + 1.5 ∗ IQR, where IQR is the inter-quartile range (Q3-Q1).(C) Correlations between amplitude and DFA exponents, and between PLV and DFA, in modeled data are also larger for cohort B than the more subcritical cohort A.(D) Correlations of local amplitudes and DFA exponents for COMT and BDNF polymorphism groups in MEG data in θ band.(E) Correlations between large-scale synchrony (wPLI) and DFA exponents for BDNF in MEG data differ between polymorphism groups in θ, but not in α band.
To investigate whether COMT and BDNF polymorphisms could influence synchronization dynamics via modulation of K, in line with the model predictions, we computed these dynamics separately for each polymorphism group in the canonical frequency bands (Figure S2). In line with model predictions, the elevated θ band amplitudes and LRTCs in COMT Val/Met group (Figure S2) were associated with stronger correlation of local oscillation amplitudes and LRTCs (Figure 5D, Val/Met > Met/Met in 95% and Val/Met > Val/Val in 84% of permutations). For BDNF, elevated θ band amplitudes and LRTCs in Val homozygotes (Figure S2) were associated with stronger correlation of local oscillation amplitudes and LRTCs (Figure 5D, Val > Met in 92% of permutations, for γ bands see Figure S7). In contrast, elevated α band global synchronization in BDNF Val homozygotes (Figure S2) was not associated with elevated correlation of synchronization with DFA (Figure 5E) whereas in θ band the correlation between synchronization and LRTCs was higher for BNDF Val homozygotes (Figure 5E, Val > Met in 98% of permutations).
These findings thus provide converging evidence that COMT Val^158^Met and BDNF Val^66^Met polymorphisms influence local oscillation dynamics as predicted by the framework of brain criticality under the hypothesis that these polymorphisms influence the control parameter of cortical critical dynamics. In contrast, the most robust effects of BDNF on large-scale synchronization were not in line of the model predictions, suggesting that these were mediated mainly by other factors than the control parameter of brain criticality.
Inter-individual variability in spontaneous oscillatory dynamics influences individual behavioral performance and cognitive ability (Donhauser and Baillet, 2020; Palva and Palva, 2012, 2018; Siebenhühner et al., 2020; Vidaurre et al., 2021). While the causes of this variability are diverse, it has been shown in human MEG and EEG data that part of it is attributable to genetics (Leppäaho et al., 2019; Linkenkaer-Hansen et al., 2007; Van Pelt et al., 2012). While macroscopic oscillations can be influenced through regulation of dynamic circuit motifs on a cellular level (Leppäaho et al., 2019; Van Pelt et al., 2012; Womelsdorf et al., 2014), there has been little research on how the genetic biases of cellular mechanisms translate into alterations of macroscopic brain dynamics observable in electrophysiological neuroimaging. We tested here a hypothesis that the COMT Val^158^Met and BDNF Val^66^Met genetic polymorphisms contribute to the variability in emergent local (oscillation amplitudes and LRTCs) and global (inter-areal phase synchronization) dynamics of neuronal oscillations. We showed that both COMT and BDNF polymorphisms influenced local oscillation amplitudes and LRTCs. However, only BDNF polymorphism impacted the strength of global synchronization. Using computational modeling, we furthermore showed that the influence of COMT and BDNF on local oscillation dynamics is achieved in line with brain criticality framework. This suggests that COMT and BNDF polymorphisms influence synchronization dynamics via modulations of the brain's E/I balance in line with the established role of COMT and BDNF in the regulation of neuronal activity via DA, NE, and 5-HT.
Prior observations about COMT and BDNF polymorphism effects on electrophysiological activity are scarce and mainly restricted to EEG and MEG sensor level (Bodenmann et al., 2009; Bonetti et al., 2021; Garcia-Garcia et al., 2017; Rodríguez-Rojo et al., 2018; Roy et al., 2020). In these studies, COMT Val homozygotes were found to have lower α power than Met homozygotes (Bodenmann et al., 2009), and COMT polymorphism was found to influence phase-resetting of γ oscillations in attention (Garcia-Garcia et al., 2017). Furthermore, BDNF has been implicated in resting-state oscillation amplitude differences in delta (0.5–3.5 Hz), α, and β bands but not in θ band amplitudes (Roy et al., 2020) and in reduced γ band FC between anterior and posterior regions in Val/Met compared to Val/Val participants (Rodríguez-Rojo et al., 2018). We used here a data-driven approach for mapping the effects of COMT and BDNF polymorphism on both instantaneous oscillation amplitudes and their LRTCs from source-reconstructed MEG data. We found anatomically specific θ−α band oscillations amplitudes to be greater for COMT Val/Met heterozygotes than for COMT Val or Met homozygotes, and θ and α oscillation amplitudes to be greater for BDNF Val homozygotes than for Met carriers. Similar to local oscillation amplitudes, we also found that COMT and BDNF polymorphisms modulated LRTCs in oscillation amplitude dynamics and hence the variability of the slow fluctuations in oscillation strengths. Anatomically specific LRTC scaling exponents were larger for COMT Val/Met heterozygotes than for Val or Met homozygotes in θ and γ bands and larger for BDNF Val homozygotes than for Met carriers in α-low-β bands, demonstrating that near-critical dynamics of brain oscillations is influenced by COMT and BDNF polymorphisms.
Our results provide robust evidence for COMT polymorphism influencing local amplitude dynamics. As COMT Val/Met heterozygotes have intermediate COMT enzymatic activity, and therefore intermediate brain DA and NE levels compared to Val and Met homozygotes (Chen et al., 2004), and as brain dynamics operating near criticality have been proposed to be beneficial for individuals (see below), our results are consistent with previous findings showing that the impact of neuromodulation on cognition follows an “inverted-U-shape” function where both excessively low and excessively high neuromodulation impair performance (Cools and D’Esposito, 2011; Farrell et al., 2012; Harris and Thiele, 2011; Meyer-Lindenberg et al., 2005; Papenberg et al., 2020; Vijayraghavan et al., 2007).
The effects of COMT Val^158^Met polymorphism were localized particularly to FPN and DMN, in agreement with the role of DA in cognitive control (Cools and D’Esposito, 2011; Ott and Nieder, 2019), as well as with previous fMRI studies in which the association between the polymorphism and neuronal activity was predominantly observed in PFC (Meyer-Lindenberg et al., 2005; Mier et al., 2010; Papenberg et al., 2020; Robbins and Arnsten, 2009) which has low DA reuptake and higher COMT levels. BDNF Val homozygotes exhibited stronger α band amplitudes bilaterally in frontoparietal regions and parcels of the visual system. In line with previous observations (Palva et al., 2013), differences in LRTC exponents were widespread, suggesting that they reflect state-dependent changes affecting most of the brain. Overall, these results support the idea that COMT and BDNF influence local oscillation dynamics in a variety of cortical regions important for cognition via DA, NE, and 5-HT.
We found that BDNF Val^66^Met but not COMT Val^158^Met polymorphism influenced large-scale global inter-areal phase synchronization. Global α band synchronization was stronger for BDNF Val homozygotes than Met carriers in all functional subsystems except within the visual system. BDNF Val homozygotes also had stronger θ synchronization within DAN and stronger left-hemispheric β synchronization between various subsystems, particularly involving DAN and DMN. This robust modulation of synchronization dynamics by BDNF polymorphism could contribute to significant differences across individuals in behavior and cognition, which have been shown to be strongly influenced by spatial and spectral characteristics of oscillatory synchronization (Betti et al., 2018; Mostame and Sadaghiani, 2020; Siebenhühner et al., 2020). However, large-scale synchronization was not influenced by COMT Val^158^Met polymorphism, although it has been shown that neuropharmacologically elevated levels of catecholamines are related to fMRI-based changes in functional connectivity (van den Brink et al., 2018). These differences might be explained by different effects of catecholaminergic levels on BOLD signal and electrophysiological connectivity, or by COMT polymorphism having more subtle effects on neuronal dynamics compared to neuropharmacological interventions.
Near-critical synchronization dynamics have been proposed to enable transient and flexible computations that are beneficial for individuals by maximizing information storage and transmission (Cocchi et al., 2017; Daffertshofer et al., 2018; Deco and Jirsa, 2012; Kinouchi and Copelli, 2006; Palva and Palva, 2018; Priesemann et al., 2014; Shew and Plenz, 2013). This notion has been supported experimentally by behavioral data (Simola et al., 2017) and in vitro studies (Shew et al., 2009). Here, we show that genetic polymorphisms in COMT and BDNF influenced individual near-critical synchronization dynamics and proximity to the critical point. The systems-level E/I ratio (or gain) is thought to be the primary control parameter for the regulation of critical brain dynamics (Shew et al., 2009; Stephani et al., 2020). We postulated that the observed differences in neuronal synchronization dynamics are due to neuromodulation-related shifts in the system’s net E/I ratio. Comparison of experimental and computational model data showed that COMT Val^158^Met and BDNF Val^66^Met polymorphisms influenced local oscillatory dynamics and their correlations in the θ band in line with the brain criticality framework. More specifically, COMT Val^158^Met heterozygotes and BDNF Val^66^ homozygotes had larger amplitudes and LRTC exponents and stronger correlations between these than the other polymorphism groups in θ band, indicating that they are associated with a higher value of the underlying control parameter, i.e., with greater net excitatory drive or neuronal gain (Ferguson and Cardin, 2020). This is in accordance with the framework in which neuromodulatory drive influences whole-brain neural states (Braun et al., 2021; Munn et al., 2021).
In contrast, we found no evidence that the influence of BDNF Val^66^Met polymorphism on global α-band synchronization could be explained within the brain criticality framework. This finding implies that the strongest effects of BDNF polymorphisms on global synchronization might not be due to a shifts in global excitation, but to other mechanisms such as the influence of BDNF on maintenance, maturation, and formation of neuronal networks (Park and Poo, 2013; Zagrebelsky and Korte, 2014) putatively via receptor TrkB (neurotrophic receptor tyrosine kinase 2) signaling (Castrén and Monteggia, 2021; Winkel et al., 2021). These could impact large-scale network synchronization via changes in structural connectivity (Mostame and Sadaghiani, 2020) or via changes in brain gray matter whose density has been shown to be associated with neuroreceptor and neurotransporter availabilities (Manninen et al., 2021) and whose thickness is associated with the structure of oscillatory networks (Mahjoory et al., 2020).
Our data showed that COMT Val^158^Met and BDNF Val^66^Met polymorphisms contribute to the variability of oscillation dynamics, their LRTCs, and synchronization, and add on to the scarce data on the genetic basis of oscillatory neuronal dynamics. In addition to COMT and BDNF polymorphisms, many other genetic factors are likely to contribute to oscillation dynamics, as observed for the cortical structure (Grasby et al., 2020), but this remains to be charted in future studies. Because oscillation amplitudes, LRTCs, and global synchronization dynamics are predictive of behavioral variability over both long (Palva et al., 2013; Smit et al., 2013) and sub-second time scales (Palva and Palva, 2012), COMT Val^158^Met and BDNF Val^66^Met polymorphisms may underlie individual differences in trait-like behavior and cognitive performance via their effect on critical oscillatory dynamics and synchronization. Differences in these dynamics caused by genetic polymorphisms are likely to influence individuals’ cognitive and mental development and may predispose them to specific brain diseases. In line with the latter idea, alterations in the brain's E/I ratio have been suggested to predispose individuals to dementia and Alzheimer disease (Pusil et al., 2019). Our results thus point in new directions for the investigation of how genetic factors could translate to behavioral variability via modulations of the systems-level neuronal dynamics.
We did not screen participants for their polymorphisms in COMT and BDNF genes prior to the MEG recordings. Due to this, one limitation of our study is the unequal group sizes, particularly for BDNF polymorphism for which the Val/Val group was the largest while there were fewer participants for the Val/Met and only 2 participants for the Met/Met group.
Further studies are also required to verify that the influence of BDNF and COMT polymorphisms on oscillation dynamics is mediated via shifts in excitability. This could be done for example by pharmacological intervention.
Requests for further information and resources should be directed to and will be fulfilled by the lead contact, Satu Palva (satu.palva@helsinki.fi).
Ethical restrictions apply to data and original data cannot be shared on a public server. Data underlying figures, statistics, and main conclusions will be uploaded to Dryad server upon acceptance.
Data from 82 healthy human volunteers (age 18 to 55 years old; mean 29 years; 6 left-handed; 44 female) were collected for this study. This sample size corresponds to a meta-analysis of neuroimaging studies for COMT effects (Mier et al., 2010) that concluded based on pooled effect size estimates that 62 is the sample size necessary to achieve 80% power to detect association with PFC activation at an α-level of 0.05 (two-tailed). The study was performed in accordance with the Declaration of Helsinki and with permission of the Ethical Committee of the Helsinki University Central Hospital. All participants gave a written informed consent prior to the recordings.
Our sample for the COMT gene consisted of 18 Val/Val, 48 Val/Met, and 16 Met/Met individuals (Table 1). A one-way ANOVA confirmed that the COMT polymorphism groups did not differ in terms of age [F(2, 46) = 0.05, p = 0.984] and the duration of the resting state recording [F(2, 46) = 0.32, p = 0.726]. A chi-square test further confirmed that the percentages of female and male participants did not differ between the COMT polymorphism groups [X^2^(2) = 2.11, p = 0.349]. The BDNF sample consisted of 66 Val/Val, 14 Val/Met, and 2 Met/Met individuals (Table 2). As the Met/Met polymorphism had a low frequency, we combined Met carriers (i.e., Val/Met and Met/Met polymorphisms) into one group (see also Beste et al., 2010). An independent samples t-test (2-tailed) confirmed that the BDNF polymorphisms did not differ in terms of age [t(32) = 0.20, p = 0.847] or the length of the resting state recording [t(32) = 0.76, p = 0.444]. A chi-square test showed no differences in the percentages of female and male participants between the BDNF polymorphisms [X^2^(1) = 0.00, p = 0.995].
Eyes-open resting-state brain activity was recorded from all participants with 306-channel MEG (Vectorview, Elekta-Neuromag, LtD) at a sampling rate of 600 Hz for ∼8 min (duration 7.8 ± 2.9 min, mean ± SD Tables 1 and 2). T1-weighted anatomical magnetic resonance imaging (MRI) MP-RAGE scans were obtained for each participant at a resolution of 1× 1 × 1 mm using a 1.5 T MRI scanner (Siemens, Germany).
Blood samples for the collection of DNA were obtained from all subjects before each MEG session and DNA was extracted from blood samples according to standard procedures at the Finnish National Institute for Health and Welfare. DNA samples were genotyped at Estonian Genome Center using Infinium PsychArray-24 v1.1 (Illumina). Quality control (QC) was performed with PLINK 1.9 (www.cog-genomics.org/plink/1.9/). Markers were removed for missingness (> 5%) and Hardy-Weinberg equilibrium (p < 1 x 10^−6^). Individuals were checked for missing genotypes (> 5%), relatedness (identical by descent calculation, PI_HAT > 0.2) and population stratification (multidimensional scaling). The dataset collected for this study is part of a broader project “Tunteet” that includes also other datasets and paradigms (*e.g.,*Bonetti et al., 2021; Kliuchko et al., 2018).
Signal space separation method (tSSS) (Taulu et al., 2005) with the Maxfilter software (Elekta Neuromag, Helsinki, Finland) was used to suppress extracranial noise, interpolate bad channels and to co-localize the recordings in signal space individually for each participant. Next, independent component analysis (ICA, MATLAB toolbox Fieldtrip (Oostenveld et al., 2011)) was used to extract and exclude signal components that were correlated with eye movement, blink and cardiac artifacts. The preprocessed MEG data were filtered into 26 logarithmically spaced frequency bands, fmin = 3 Hz, fmax = 60 Hz, using complex Morlet Wavelets with parameter m = 5 (making a compromise between time and frequency resolution). After filtering, the time-series data were decimated to a frequency-dependent sampling rate of between 2 and 5 times the center frequency. FreeSurfer software (http://surfer.nmr.mgh.harvard.edu/) was used for volumetric segmentation of the MRI data, surface reconstruction, flattening, cortical parcellations and neuroanatomical labeling with the Destrieux atlas. We then iteratively split the largest parcel in the atlas in two, until we arrived at a resolution of 400 parcels. MNE software (https://mne.tools/stable/index.html) (Gramfort et al., 2014) was used for creating cortically constrained source models and for the preparation of the forward and inverse operators (Dale et al., 2000) using noise covariance matrices (NCMs) obtained from broad-band filtered data (125–195 Hz, with 150 Hz notch). MNE inverse operators were then computed for all wavelet frequencies and used to project the sensor-space data into source space. The source models had dipole orientations fixed to the pial surface normal, yielding 5000–8000 source vertices per hemisphere. As in previous studies (Korhonen et al., 2014; Siebenhühner et al., 2020), source-vertex time-series were then collapsed into parcel time series with individually source-reconstruction-accuracy- (fidelity-) optimized collapse operators. To enhance the possibility of detecting true connections amongst the spurious connections, we used an atlas of 400 parcels adapted from the Destrieux atlas by iteratively splitting parcels along their most elongated axis, using the same parcel splits for all subjects (Palva et al., 2010, 2011; Rouhinen et al., 2020). The 400-parcel data were collapsed to 148 parcels before statistical analyses to reduce the influence of the large inter-individual function-anatomical variability.
The parcel-wise narrow-band filtered time-series were used for cortex-wide mapping of induced oscillation amplitudes (A) and LRTCs in oscillation amplitude envelopes, as well as for mapping cortex-wide phase-synchronization networks. Data were averaged for each cortical parcel and frequency separately for subjects in each COMT and BDNF polymorphisms group. Oscillation amplitudes were extracted from Morlet-wavelet filtered complex-valued time-series while LRTC exponents of the amplitude envelopes were quantified with detrended fluctuation analysis (DFA) (Hardstone et al., 2012). The power-law scaling exponent β was defined as the slope of linear regression of the function F(Δt) plotted in log–log coordinates, estimated using a least-squares algorithm. To estimate cortex-wide synchronization connectomes, phase synchronization between all 400 parcels was computed using the weighted phase-lag index (wPLI) (Siebenhühner et al., 2020; Vinck et al., 2011). The wPLI is insensitive to artificial connections arising as direct effects of zero-phase lagged linear signal mixing that is a major issue in connectivity analysis using MEG/EEG data (Palva and Palva, 2012). To reduce the amount of spurious interactions, we removed from the analysis parcels with fidelity < 0.2, which led to the exclusion of 9.27% of possible parcel pairs (Korhonen et al., 2014; Siebenhühner et al., 2020). The mean connection strengths were obtained by averaging the strength of all remaining edges for each subject and then over subjects in each COMT and BDNF polymorphism group. For functional subsystem analysis, all pairwise interactions were averaged across the frequencies in one of the bands defined by clustering analyses (Figure S1) and then collapsed into the 7 functional systems of the Yeo parcellation (Yeo et al., 2011) (Figure S3).
We used a nested Kuramoto model (Siebenhühner et al., 2020) to investigate the covariance and mutual correlations of local synchronization (node/parcel amplitude), inter-areal synchronization, and local LRTCs in a heuristic model of synchronization dynamics. The model was adapted from a conventional Kuramoto model so that it comprised 100 nodes (corresponding to 100 cortical parcels) and each node was modeled by a Kuramoto model of 500 oscillators. The model was defined so that for each node k, the temporal evolution of the phase θ of each oscillator i was described
where ωi is the oscillation frequency of an oscillator, K is the local coupling (control) parameter, φl is the average phase of the node l, Wk,j is the connection strength between the node k and other nodes j, and N(0,σ) is Gaussian noise. We set all oscillators’ frequencies to a constant value (10 Hz). W was approximated from the log-scaled structural connectome of the Schaefer parcellation with 100 parcels (Schaefer et al., 2018) multiplied by a global coupling (control) parameter L at L = 0.5. Node time series mimicking the parcel time series in MEG data were obtained by averaging the complex time series of the oscillators in the node. Local synchronization was measured by the absolute value of the complex node time series (Kuramoto order) that is comparable to the amplitude of MEG parcel time series. We predicted the synchronization dynamics of neuronal oscillations for two hypothetical cohorts that differed in their control parameters regulating their position in the subcritical regime. (A, control parameter from 8.25–9.25; and B, control parameter from 8.75–9.75) We estimated pairwise synchronization here with the Phase Locking Value (PLV), since, unlike in MEG data, in modeled data there is no volume conduction. Further, in the absence of conduction delays, the coupling function used in the Kuramoto model leads exclusively to near-zero-lag phase synchrony, which is not detectable with wPLI, necessitating the use of PLV. To test the model predictions, we estimated the correlation between amplitudes and DFA exponents, and between synchronization and DFA exponents, across all nodes, with Spearman’s correlation coefficient.
Before statistical testing, time-series were collapsed into the 148 parcels of the original Destrieux atlas (Destrieux et al., 2011) to limit the number of comparisons and to decrease the influence of large function-anatomical variability. Group-level upper and lower confidence limits (2.5% and 97.5%) were computed for each polymorphism group with a bootstrapping approach, using N = 1000 resamplings with replacement of the subjects in the group.
To validate the use of canonical frequency bands, we first estimated spatial similarity across frequency-bands separately for oscillation amplitudes and DFA exponents using machine learning by Louvain community detection (Blondel et al., 2008) with gamma value of 1.09. This analysis confirmed the presence of canonical frequency bands of theta (θ, 3–7 Hz), alpha (α, 8–14 Hz), beta (β, 14–30 Hz), and gamma (γ, 30–60 Hz) based on spatial similarity (Figure S1). The effects of COMT and BDNF polymorphisms on oscillation amplitudes and DFAs over the whole brain were tested separately for each Morlet frequency (Figures 2A and 3A), as well as for the canonical frequency bands θ, α, β and γ (Figure S2), using univariate ANOVA with the polymorphisms as between-participants factors and age as a covariate. Partial eta squared (ηp^2^) is reported as a measure of effect size.
In order to control for multiple comparisons at the level of individual wavelet frequencies, we used cluster-based non-parametric permutation testing as described by Maris and Oostenveld (2007). Briefly, we clustered consecutive wavelet frequencies which were significant in the ANOVA, summed F-values within clusters, and computed the F-value permutation distribution (based on random assignment of participants to polymorphism groups, N = 1000) for the cluster with the largest summed F. Any cluster for which the summed F-value was larger than the 95^th^ percentile of this permutation distribution was deemed significant.
We then used the Mann-Whitney-U test (p < 0.05) to test pairwise differences between the polymorphism groups separately for each wavelet frequency (Figures 2B and 3B) and for each of the 148 parcels (Figures 2C and 3C). Multiple comparisons were accounted for in two steps. First, we used false discovery rate reduction (FDR) by discarding as many least-significant observations as were predicted to be false by the alpha-level (5% at α = 0.05) for each time-frequency bin. Second, we estimated a threshold Q to define the joint probability, pq< 0.05, for residual fraction of significants that could arise by chance in any of the frequencies (see dark gray area in Figures 2 and 3B). This threshold gives an additional estimate of expected false positives.
For inter-areal synchronization, confidence intervals were computed with a bootstrapping approach (N = 1000). As the synchronization data violated the normality assumption, we tested, for the canonical frequency bands and each Morlet wavelet frequency, whether there was a significant difference between polymorphism groups with a non-parametric Kruskal-Wallis test (p < 0.05) (Figures S2 and 4A). We used the same cluster-based non-parametric permutation testing approach to control for multiple comparisons at the wavelet frequency level as we had used for amplitude and DFA. Connection density (K) of the mean connection strengths was used to index the proportion of significant edges from all possible interactions. Differences between the polymorphism groups in connectivity within and between functional subsystems (Yeo et al., 2011) (Figure S3) were then tested for statistical significance with the Kruskal-Wallis test (p < 0.05) and corrected for multiple comparisons with Benjamini-Hochberg (Figure 4B).
Finally, we estimated the correlation coefficient between amplitudes and DFA exponents and between synchronization and DFA exponents, for each parcel and wavelet frequency, across subjects within the polymorphism groups, using Spearman’s rank correlation test. Values were then averaged across all parcels and within frequency bands. Confidence limits from 2.5% to 97.5% were obtained from bootstrapping (N = 1000) as above. Significance testing for difference between polymorphism groups was done with permutation statistics (N = 1000) where the subjects were randomly assigned, without replacement, to original-sized polymorphism groups.
This work was supported by the Jane and Aatos Erkko Foundation, by the Academy of Finland (SA 1266402, 1267030 to S.P., SA 1266745, 1296304 to J.M.P., and SA 1294761 to J.S.).
Conceptualization, S.P. and J.M.P.; Methodology, F.S., S.P., J.M.P., and V.M.; Software, F.S., J.M.P., and V.M.; Formal Analysis (MEG), J.S., F.S., and V.M.; Formal Analysis (Genotyping), K.K. and T.P.; Investigation, E.B.; Writing – Original Draft, J.S., S.P., F.S., J.M.P., and E.B.; Writing – Review & Editing, J.S., F.S., S.P., J.M.P., and V.M.; Visualization, J.S., F.S., and V.M.; Supervision, S.P., J.M.P., and E.B.; Funding Acquisition, E.B., J.M.P., S.P., T.P., and J.S.
The authors declare no competing interests.
Published: September 16, 2022
Jaana Simola, Email: jaana.simola@helsinki.fi.
Satu Palva, Email: satu.palva@helsinki.fi.