## Abstract

The synchronization between neural oscillations at different frequencies has been proposed as a core mechanism for the coordination and integration of neural systems at different spatiotemporal scales. Because neural oscillations of different frequencies can only fully synchronize when their “peak” frequencies form harmonic relationships (e.g., *f*_{2} = *f*_{1}/2), the present study explored whether the transient occurrence of harmonic cross-frequency relationship between task-relevant rhythms underlies efficient cognitive processing. Continuous EEG recordings (51 human participants; 14 males) were obtained during an arithmetic task, rest and breath focus. In two separate experiments, we consistently show that the proportion of epochs displaying a 2:1 harmonic relationship between alpha (8–14 Hz) and theta (4–8 Hz) peak frequencies (i.e., alpha_{peak} ≈ 10.6 Hz; theta_{peak} ≈ 5.3 Hz), was significantly higher when cognitive demands increased. In addition, a higher incidence of 2:1 harmonic cross-frequency relationships was significantly associated with increased alpha–theta phase synchrony and improved arithmetic task performance, thereby underlining the functional relevance of this cross-frequency configuration. Notably, opposite dynamics were identified for a specific range of “nonharmonic” alpha–theta cross-frequency relationships (i.e., alpha_{peak}/theta_{peak} = 1.1–1.6), which showed a higher incidence during rest compared with the arithmetic task. The observation that alpha and theta rhythms shifted into harmonic versus nonharmonic cross-frequency relationships depending on (cognitive) task demands is in line with the notion that the neural frequency architecture entails optimal frequency arrangements to facilitate cross-frequency “coupling” and “decoupling”.

**SIGNIFICANCE STATEMENT** Neural activity is known to oscillate within discrete frequency bands and the interplay between these brain rhythms is hypothesized to underlie cognitive functions. A recent theory posits that shifts in the peak frequencies of oscillatory rhythms form the principal mechanism by which cross-frequency coupling and decoupling is implemented in the brain. In line with this notion, we show that the occurrence of a cross-frequency arrangement that mathematically enables coupling between alpha and theta rhythms is more prominent during active cognitive processing (compared with rest and non-cognitively demanding tasks) and is associated with improved cognitive performance. Together, our results open new vistas for future research on cross-frequency dynamics in the brain and their functional role in cognitive processing.

## Introduction

Decades of research suggest that the rhythmic synchronization of neural populations underlie neural communication and information transmission in the brain (Varela et al., 2001; Fries, 2005, 2015). Because neural oscillations at different frequencies have been shown to represent distinct functions (Kopell et al., 2010), anatomical sources (Debener et al., 2005) and scales of cortical integration (von Stein and Sarnthein, 2000; Canolty and Knight, 2010), the synchronization between these rhythms (i.e., “cross-frequency coupling”) has been proposed as a core mechanism for the coordination and integration of neural systems at different spatiotemporal scales (Canolty and Knight, 2010; Palva and Palva, 2018). In support for this notion, significant relationships have been demonstrated between behavioral performance in cognitive tasks and the degree of cross-frequency interactions quantified through phase–phase synchrony (Siebenhühner et al., 2016), amplitude–amplitude correlations (Popov et al., 2018), and phase–amplitude coupling (Axmacher et al., 2010; Lisman and Jensen, 2013).

Mathematically, two oscillators of different frequencies can only fully synchronize when forming harmonics (e.g., *f*_{2} = *f*_{1}/2). Unlike nonharmonic relationships (e.g., *f*_{2} = *f*_{1}/1.6), a harmonic frequency arrangement leads to a pattern of frequent and regular excitatory phase meetings, thereby enabling (cross-frequency) synchronization (Pletzer et al., 2010). Based on this mathematical reality, a theoretical framework was posited by Klimesch (2012, 2013) and Pletzer et al. (2010), hypothesizing that the specific frequency architecture of neural oscillations may approximate a geometrical series of exponent “2” (e.g., δ = 2.5, θ = 5, α = 10, β = 20, γ = 40, etc.) to facilitate “harmonic” relationships (and therefore, synchronization) between different rhythms. Accordingly, transient changes in the “peak” frequency of different rhythms are proposed to form a flexible mechanism to maximize cross-frequency coupling or decoupling (Klimesch, 2012). For example, the specific arrangement of the peak frequencies of two adjacent frequency bands (e.g., alpha and theta) will determine whether they will form a harmonic (e.g.10/5 Hz = 2) or nonharmonic (e.g., 8/5 Hz = 1.6) ratio and, consequently, whether their synchronization will be facilitated or prevented. In line with this hypothesis, previous research has shown that peak frequencies in different bands can transiently shift under different experimental conditions (Foffani et al., 2005; Haegens et al., 2014; Lowet et al., 2017; Mierau et al., 2017). To date, however, these shifts in peak frequencies have not been investigated in the context of facilitating or altering cross-frequency dynamics. Indeed, despite the intuitive appeal of interpreting the neural frequency architecture in terms of a geometrically organized system of harmonic oscillators (Pletzer et al., 2010; Klimesch, 2012, 2013), its experimental application and functional relevance remains unexplored.

In this study, continuous EEG recordings (19 electrodes; 51 participants) were obtained during an arithmetic task with an important working memory component and two non-cognitively demanding conditions (i.e., rest and breath focus) to assess whether transient harmonic relationships between the peak frequencies of two adjacent rhythms are of functional relevance for cognitive processing. Given the well established role of alpha–theta cross-frequency interactions in tasks involving working memory and executive control (Schack et al., 2005; Kawasaki et al., 2010; Dimitriadis et al., 2016; Akiyama et al., 2017; Li et al., 2017; Popov et al., 2018), we specifically focused on the occurrence of transient harmonic relationships between alpha (8–14 Hz) and theta (4–8 Hz) peak frequencies. Because, from a mathematical perspective, a 2:1 harmonic relationship between alpha and theta peaks is anticipated to facilitate cross-frequency coupling (allowing frequent and regular cross-frequency excitatory phase meetings), we here specifically tested (1) whether the transient occurrence of this relationship will be more prominent during the arithmetic task condition (compared with the non-cognitively demanding conditions), and (2) whether an increased occurrence of harmonic cross-frequency relationships will be predictive of improved arithmetic task performance.

## Materials and Methods

#### Participants

A total of 21 and 33 adults participated in Experiment 1 (5 males, age range 20–28 years) and Experiment 2 (9 males, age range 20–29 years), respectively. Participants were recruited through flyers and via social media from various campuses of the University of Leuven. Informed written consent was obtained from all participants before the study. Consent forms and study design were approved by the local Ethics Committee for Biomedical Research at the University of Leuven. Participants were compensated for their participation (10 € per hour). One participant of Experiment 1 was excluded from the final analyses because of technical problems during data acquisition. Two participants of Experiment 2 were excluded because of excessive EEG artifacts.

#### Experimental design

In Experiment 1, continuous EEG recordings were performed during the following three conditions: “arithmetic”, “rest”, and “breath focus” (with the order randomized across subjects). During the experimental arithmetic condition, participants sat on a chair with their eyes closed and were instructed to perform (in silence) a simple arithmetic task iteratively. Participants were required to sum two numbers (1 digit integers: *x*_{1}, *x*_{2}) to then sum the last number (*x*_{2}) and the obtained summation number (*x*_{1} + *x*_{2}) iteratively (i.e., *x*_{1} + *x*_{2} = *x*_{3}; *x*_{2} + *x*_{3} = *x*_{4,} etc.). When the summation exceeded 200, participants were required to press a key button with their right middle finger and report the obtained summation number. After the key press, two new numbers were given and the procedure was repeated until the participant completed ∼6 min of EEG recordings. The two numbers provided at the beginning of each trial were the same for all participants and they were always presented aurally by the experimenter. Note that, by design, the number of trials that were completed within these 6 min varied across participants (average number of trials = 4, range 3–7). For each participant, accuracy (100 − absolute difference from correct answer) and response time (time from aural presentation of stimuli until key press) were averaged across trials to compute a performance index per subject. During the rest condition participants sat with their eyes closed for a total duration of 6 min and were asked “to relax, and not to fall asleep”. Finally, during the breath focus condition, participants were explicitly asked to focus on their breathing while closing their eyes for a duration of 6 min, after the experimenter read out loud a short introductory text on meditative breathing techniques.

In Experiment 2, continuous EEG recordings were performed during similar conditions as those adopted in Experiment 1, but with the following adjustments. For the experimental arithmetic condition, the number of trials was now fixed to *N* = 5, to obtain an equal number of trials across participants (i.e., regardless of individual response time). As such, the total EEG recoding time was now variable across participants (mean = 11.52 min, SD = 4.63). Furthermore, in Experiment 2, the breath focus condition was slightly altered to match its motor requirements to the motor requirements of the arithmetic condition. Specifically, in the breath focus condition, participants were now instructed to explicitly focus on their breathing by silently counting their breaths and, importantly, to press a key button with their right middle finger every tenth exhalation. The breath focus condition was performed for a duration of ∼15 min.

#### EEG recordings and preprocessing

The Nexus-32 system and BioTrace software (Mind Media) were used to perform electroencephalography (EEG) recordings. Continuous EEG was recorded with a 21-electrode cap positioned according to the 10-20 system (MediFactory). Skin abrasion and electrode paste (Nuprep) were used to reduce the electrode impedances <10 kΩ during the recordings. The EEG signal was amplified using a unipolar amplifier with a sampling rate of 256 Hz.

Raw EEG data from both experiments are available at the open science framework webpage: https://osf.io/gh6q3/?view_only=2f3e881eec294739b1f2519d8c522bf9.

Preprocessing was performed with the MATLAB-based toolbox “Letswave 6” (MATLAB version r2017b). Raw EEG data were filtered using a 0.5–40 Hz bandpass filter (Butterworth filter order = 4) and were segmented into 0.5 s epochs to remove those epochs with amplitudes >100 mV. The remaining epochs were then concatenated and the continuous signals were mathematically re-referenced off-line to linked ears (average of A1 and A2). Next, independent component analysis (ICA) was used to manually select and remove components with spatial topographies characteristic of eye blinks and horizontal eye movements. For this purpose, ICA decomposition was performed with the RUNICA algorithm (as implemented in EEGLAB) using the SQUARE method. On average, two components were rejected per subject.

#### Transient peak detection and determination of cross-frequency relationships

The time-frequency representation of the data was obtained using short-term fast Fourier transform (STFFT) computed through the MATLAB toolbox Letswave 6 (Hanning window length of 1 s and sliding step of 25 bins; i.e., 90.23% overlap). A frequency precision of 10 points per frequency was chosen (i.e., from 0.1 Hz until 40 Hz in 400 frequency lines). Then, transient peak frequencies of the theta (4–8 Hz) and alpha (8–14 Hz) bands were detected for each (overlapping) 1 s epochs of transformed data using the find local maxima function implemented in MATLAB r2017b (i.e., *findpeaks*). With this algorithm, data samples that are larger than their two neighboring samples were identified as “local peaks” within the specified alpha and theta frequency ranges. For a limited number of epochs (5.23% of total number of epochs, across experiments, conditions, subjects, and electrodes), no clear peaks were detected in the theta (4–8 Hz) or alpha (8–14 Hz) bands and these epochs were therefore excluded from further analyses. When two or more peaks were detected in one frequency band, only the peak with highest amplitude was selected. Based on the identified peak frequencies, the numerical ratio of the alpha and theta peaks (peak-frequency_{alpha}/peak-frequency_{theta}) was calculated for each epoch and rounded to the first decimal place (e.g., 10/5 Hz = 2.0). Hence, the obtained ratio-values ranged between 1.1 and 3.4 in steps of 0.1. Finally, the proportion of epochs in which the alpha–theta peak ratio equaled 2.0 (henceforward termed “harmonic locking”) was determined for each electrode, subject, and condition. Figure 1*A* visualizes the transient (epoch-wise) variability of alpha and theta peak frequencies over time (i.e., 10 s) for an exemplary subject and electrode, as well as the transient numerical ratio over time. Figure 1*B* visualizes the frequency spectra of two exemplary epochs in which the identified alpha and theta peak frequencies formed a harmonic (2:1) versus a nonharmonic (1.6:1) cross-frequency relationship. Figure 1*C* visualizes the distribution of alpha and theta peak frequencies that yielded harmonic 2:1 cross-frequency relationships (across experiments, conditions, subjects, and electrodes). While distinct pairs of alpha and theta peak frequencies led to harmonic 2:1 relationships, it can be seen from the distribution that harmonic relationships were predominantly evident for alpha peak frequencies centered ∼10.6 Hz (1.33 SD) and theta peak frequencies centered ∼5.3 Hz (0.67 SD).

#### Statistical analysis

For each experiment (1 and 2), condition-related differences in the proportion of epochs that displayed 2:1 harmonic locking were examined using paired-sample *t* tests, controlling for the type I error rate arising from multiple comparisons across electrodes through nonparametric cluster-based permutation statistics (Maris and Oostenveld, 2007) as implemented in the MATLAB toolbox FieldTrip (Oostenveld et al., 2011). Data were randomly partitioned (Monte Carlo approximation; 1000 permutations) and the maximal cluster-level statistics (the sum of *t* values across spatially adjacent electrodes) were extracted for each random partition to estimate a “null” distribution of effect size. Then, the cluster-corrected *p* value was defined as the proportion of random partitions in the null distribution whose test statistics exceeded the one obtained for each significant cluster (*p* < 0.05) in the original (non-shuffled) data.

The same statistical procedure was used to assess whether interindividual differences in harmonic locking during arithmetic task were related to task performance (i.e., accuracy/response time). In this case, analyses were performed across experiments (combined dataset; *n* = 51) using Pearson's partial correlation coefficient (controlling for “experiment”) as the test statistic.

Finally, nonparametric permutation statistics were also used to explore the specificity of the 2:1 harmonic cross-frequency ratio for inducing condition-related differences (see Results, Exploring the ratio specificity of the condition effect). For this analysis, data from Experiments 1 and 2 were also combined (i.e., the common conditions arithmetic and rest). In short, the proportion of epochs displaying any of the other possible cross-frequency relationships were computed (i.e., ratios within a range of 1.1–3.4, with a step-size of 0.1; separately for each subject, condition, and electrode) and paired-sample *t* tests were performed comparing arithmetic versus rest conditions, separately for each electrode and cross-frequency ratio. Here, cluster-based permutation statistics were adopted to assess statistical significance while controlling for multiple comparisons across electrodes and cross-frequency ratios, i.e., identifying significant clusters based on spatial (electrodes) and cross-frequency ratio adjacency.

#### Secondary analyses

##### Mean alpha–theta peak frequencies and their relation to harmonic locking.

In the primary analyses, we assessed condition-related differences in 2:1 harmonic locking (% of epochs in which the alpha and theta peak frequencies display a harmonic ratio; i.e., alpha–theta = 2.0), and relationships between harmonic locking and arithmetic task performance. Secondary analyses were performed to explore whether similar condition-related differences and relationships with task performance are evident when mean alpha and theta peak frequencies are analyzed separately. Correlation analyses were also performed to assess whether interindividual differences in mean alpha and theta peak frequencies were related to harmonic locking. Mean peak frequencies of the alpha and theta band were estimated per subject, electrode, and condition by averaging the peak frequency of transiently detected peaks over time. Similar to the primary analyses, the significance of condition-related differences (arithmetic vs rest conditions; paired-sample *t* test) and the correlation analyses (Pearson's partial correlation controlling for experiment) was assessed using nonparametric cluster-based permutation statistics (see Statistical analysis).

##### Cross-frequency phase synchronization and its relationship to harmonic locking.

Because a harmonic relationship between alpha–theta rhythms mathematically enables stable and frequent phase meetings, it is expected that a greater incidence of harmonic locking is related to increased alpha–theta cross-frequency phase synchronization. Hence, secondary analyses were performed to assess whether interindividual differences in cross-frequency phase synchronization are related to interindividual differences in harmonic locking. In addition, we explored whether cross-frequency phase synchronization was significantly modulated during arithmetic task compared with rest and whether it was related to arithmetic task performance. Alpha–theta 2:1 phase synchrony was computed within each electrode as follows. First, the EEG signal was filtered for alpha (8–14 Hz) and theta (4–8 Hz) bands with a plateau-shaped zero phase digital filter with transition zones of 15% (Cohen, 2014). Second, Hilbert transform was applied to the filtered data to obtain phase angle time series. Third, phase locking value (PLV) time series were computed with a sliding window of 500 ms (for MATLAB code, see Scheffer-Teixeira and Tort, 2016). Finally, PLVs were averaged over time to obtain a single estimation of alpha–theta 2:1 phase synchrony per subject, electrode and condition. Similar to previous analyses, the significance of condition-related differences (between arithmetic and rest conditions; paired-sample *t* test) and the correlation analyses (Pearson's partial correlation controlling for experiment) was assessed using nonparametric cluster-based permutation statistics (see Statistical analysis).

##### Transient peak detection accounting for the 1/*f* trend in the EEG spectrum.

Neural oscillations are embedded in a background signal that displays a decrease in power with increasing frequency (i.e., 1/f background spectrum). To rule out the possibility that this spectral pattern biased transient peak detections toward lower frequencies, we repeated the analysis taking only into account those alpha and theta peak frequencies that exceed the 1/*f* background spectrum. For this purpose, we estimated the 1/*f* trend (separately for each subject, electrode) by fitting a line in log–log space to the EEG frequency spectrum (averaged over time) using the *robustfit* function in MATLAB (Whitten et al., 2011; for the implementation of this or similar approaches, see Caplan et al., 2015; Watrous et al., 2018). Next, transient peaks in the alpha and theta bands exceeding the estimated 1/*f* trend line were identified and harmonic relationships between alpha/theta peaks were computed (see Transient peak detection and determination of cross-frequency relationships). The statistical analysis was performed as described earlier (see Statistical analysis).

##### Assessing the possibility of artifactual harmonics due to non-sinusoidal signals.

Non-sinusoidal (e.g., spiky) signals can produce harmonic peaks in power spectra and have been highlighted as a source of “artifactual” cross-frequency coupling (Kramer et al., 2008; Scheffer-Teixeira and Tort, 2016; Cole and Voytek, 2017). In this way, if the identified harmonic alpha–theta peak frequencies were induced by non-sinusoidal properties of the slower (theta) oscillations, the power of the faster (alpha) oscillations should be tightly linked to the power of the slower (theta; Palva et al., 2005; Palva and Palva, 2018). To assess this possibility, we computed Spearman's rank-order correlations between the power time series of the most frequently observed harmonic frequencies in the alpha and theta bands. In particular, the power time series of three different pairs of harmonic frequencies (i.e., 5 and 10 Hz; 5.5 and 11 Hz; 6 and 12 Hz) were extracted through STFFT (see Transient peak detection and determination of cross-frequency relationships) and Spearman's rank-order correlations (as implemented in MATLAB r2017b) were performed for each pair of harmonic frequencies per subject and electrode. Finally, Spearman's ρ and *p* values were averaged across subjects and electrodes to obtain a single estimation of power correlation for each of the three alpha–theta frequency pairs.

## Results

### Condition-specific differences in alpha–theta harmonic locking

In both experiments, cluster-based permutation analyses revealed that the proportion of epochs displaying 2:1 harmonic locking was significantly higher during the arithmetic condition, compared with the rest and breath focus conditions (Fig. 2).

In Experiment 1, the comparison of harmonic locking between the arithmetic condition and the rest/breath focus conditions led to the identification of respectively two (*t*_{(19)} = 12.7; *p* = 0.002; *t*_{(19)} = 5.8; *p* = 0.019) and one (*t*_{(19)} = 7.56; *p* = 0.006) significant positive clusters (Fig. 2, top left and middle), indicating a significant increase in the proportion of harmonic locking during arithmetic task compared with rest (in electrodes T5, P3, Pz, and O_{2}) and breath focus conditions (in electrodes and T3, T5, and P3). No significant clusters were identified for the comparison between the rest and breath focus conditions.

In Experiment 2, the comparison of harmonic locking between arithmetic and rest/breath focus conditions led to the identification of one significant positive cluster covering all electrodes (Fig. 2, bottom left and middle; *t*_{(30)} = 60.9; *p* < 0.001; *t*_{(30)} = 49.62; *p* < 0.001). Similar to Experiment 1, no significant differences were evident between the rest and breath focus conditions.

### Relationship between alpha–theta harmonic locking and arithmetic task performance

Across Experiments 1 and 2 (*n* = 51), partial correlation analyses revealed a significant positive correlation between interindividual differences in arithmetic task performance (i.e., accuracy/response time) and harmonic locking during arithmetic task, thereby indicating that better performers tended to display a greater proportion of transient 2:1 alpha–theta cross-frequency relationships. As visualized in Figure 3*A*, the positive relationship between harmonic locking and task performance was evident in a posterior cluster encompassing electrodes Pz, P4, and O_{2} (mean *r* value = 0.41; *t*_{(50)} = 9.62; *p* = 0.003) as well as in a right frontotemporal cluster including electrodes T4 and F8 (mean *r* value = 0.42; *t*_{(50)} = 6.59; *p* = 0.004). Figure 3*B* visualizes the relationship within each electrode.

### Exploring the ratio specificity of the condition effect

Here we explore whether the reported condition-effect in harmonic locking (i.e., increased occurrence during arithmetic task) was specific to this frequency arrangement. To do so, condition-related differences (between arithmetic and rest conditions) were assessed for each cross-frequency ratio (i.e., ratios within a range of 1.1–3.4, with a step size of 0.1). Significant clusters based on spatial (electrodes) and cross-frequency ratio adjacency were identified using cluster-based permutation statistics.

Figure 4*A* visualizes the occurrence of each cross-frequency relationship (proportion of epochs averaged across electrodes) separately for arithmetic and rest conditions (error bars indicate SD across subjects). In Figure 4*B*, condition-related differences in proportions (between arithmetic, and rest) are visualized by plotting, for each ratio, the respective paired-sample *t* test values (higher *t* values indicate a higher incidence of a ratio in the arithmetic versus the rest condition when averaging across electrodes). Cluster-based permutation statistics identified a significant positive cluster within a range of cross-frequency ratios between 1.8 and 3.0 (*t*_{(50)} = 611.19; *p* < 0.01), indicating a significantly higher occurrence of these ratios during the arithmetic, compared with the rest condition (see gray area in Fig. 4*B*, gray area; note that the effect was maximal ∼ratio 2.0). Also a significant negative cluster was identified within a range of cross-frequency ratios between 1.1 and 1.6, indicating a significantly lower occurrence of these ratios during the arithmetic, compared with the rest condition (*t*_{(50)} = −350.24; *p* = 0.002). The spatial distribution of these clusters is visualized in the topographical heat maps in Figure 4*C*.

In sum, these exploratory analyses revealed that an increase in cognitive demands (i.e., during the arithmetic task) induced a significant reduction in the occurrence of a set of lower-end ratios (i.e., <1.7) in addition to a relative increase in the occurrence of cross-frequency ratios clustered ∼2.0 (i.e., 1.8–3.0; Fig. 4*A*,*B*).

### Secondary analyses

#### Mean alpha–theta peak frequencies and their relation to harmonic locking

Across experiments (*n* = 51), cluster-based permutation analysis revealed that the mean peak frequency in alpha band was significantly faster during the arithmetic task compared with rest (*t*_{(50)} = 76.68; *p* < 0.001; cluster across all electrodes), whereas the mean peak frequency in theta band was significantly slower during the arithmetic task compared with rest (*t*_{(50)} = −62.48; *p* < 0.001; cluster across all electrodes). Furthermore, interindividual differences in harmonic locking during arithmetic task were significantly correlated with a faster mean alpha peak frequency (mean *r* value = 0.54; *t*_{(50)} = 83.45; *p* < 0.001; cluster across all electrodes) and a slower mean theta peak frequency (mean *r* value = −0.59; *t*_{(50)} = −97.14; *p* < 0.001; cluster across all electrodes). However, unlike harmonic locking, interindividual differences in mean alpha and theta peak frequencies were not significantly correlated to arithmetic task performance (no significant clusters).

Together, these results suggest that although the incidence of harmonic locking was facilitated by an acceleration of alpha peak frequencies, and a deceleration of theta peak frequencies; it appears that enhanced task performance is predominantly related to the transient formation of harmonic relationship between the respective alpha–theta peak frequencies rather than to a deceleration or acceleration of the peak frequencies per se.

#### Cross-frequency phase synchronization and its relation to harmonic locking

Subjects with a higher proportion of harmonic locking (percentage of epochs for which the alpha and theta peak frequencies formed a 2:1 harmonic ratio) during arithmetic task displayed a higher phase locking value (mean *r* value = 0.67; *t*_{(50)} = 123.87; *p* < 0.001; cluster across all electrodes). Cross-frequency phase synchronization between alpha and theta oscillations was also shown to be significantly higher during the arithmetic compared with the rest condition (significant positive cluster encompassing all electrodes; *t*_{(50)} = 85.17, *p* < 0.001). However, interindividual differences in cross-frequency phase synchrony were not significantly correlated to arithmetic task performance (no significant clusters). Together, these analyses confirm that a greater incidence of harmonic relationships between alpha–theta peak frequencies is related to increased alpha–theta phase synchrony.

#### Transient peak detection accounting for the 1/*f* trend in the EEG spectrum

Secondary analyses showed that the main pattern of results remained unchanged when accounting for the 1/*f* trend in the EEG spectrum (see Materials and Methods). In particular, the comparison of harmonic locking between the arithmetic and rest condition led to the identification of a significant positive cluster (*t*_{(50)} = 72.52; *p* < 0.001; cluster across all electrodes; analysis performed in combined datasets). Also similar to the primary analysis, a significant positive correlation between harmonic locking during arithmetic task and performance was identified in a posterior (mean *r* value = 0.41, *t*_{(50)} = 9.62, *p* = 0.009; electrodes Pz, P4, and O_{2}) and a right temporal cluster (mean *r* value = 0.42, *t*_{(50)} = 6.59, *p* = 0.018; electrodes T4 and F8).

#### Assessing the possibility of artifactual harmonics due to non-sinusoidal signals

Spearman's rank-order correlations between the time courses of alpha and theta power at distinct pairs of harmonic alpha–theta frequencies were extremely low and not significant (5 and 10 Hz: *r* = 0.038; *p* = 0.16; 5.5 and 11 Hz: *r* = 0.039; *p* = 0.15; 6 and 12 Hz: *r* = 0.039; *p* = 0.16; average correlation across subjects, electrodes, and conditions). As such, it is anticipated that the identified harmonic relationships between alpha–theta peak frequencies are unlikely to reflect harmonic peaks induced by non-sinusoidal properties of a single oscillator (as in this case the power of the alpha and theta oscillations should have been tightly linked).

## Discussion

The present study sought to assess whether a transient 2:1 relationship between peak frequencies in the alpha (8–14 Hz) and theta (4–8 Hz) bands is of functional relevance during an arithmetic task with an important working memory component. In line with our hypothesis, we showed in two separate experiments that the proportion of epochs presenting harmonic locking (i.e., 2:1 ratio between alpha and theta peak frequencies) increased significantly during arithmetic task compared with rest and breath focus (with or without a motor response). Notably, subsequent correlational analyses demonstrated that a higher incidence of harmonic locking during the arithmetic condition was related to improved arithmetic task performance.

A harmonic relationship between two rhythms of different frequencies is a necessary condition for their synchronization to occur. Because brain rhythms are known to be non-stationary in terms of transiently changing their peak frequency over time (Foffani et al., 2005; Lowet et al., 2017; Mierau et al., 2017), it has been proposed that cross-frequency synchronization could be continuously enabled or prevented by the transient occurrence of harmonic versus nonharmonic relationships between their peak frequencies (Pletzer et al., 2010; Klimesch, 2012). Here we specifically quantified opportunities for cross-frequency synchronization by transiently assessing the occurrence of harmonic relationships between the peak frequencies of two adjacent rhythms. By doing so, we consistently demonstrated that an increase in cognitive demands (i.e., retention and manipulation of selective information during the arithmetic task) was accompanied by a significant increase in the incidence of harmonic 2:1 cross-frequency relationships between alpha and theta peak frequencies. Furthermore, a relationship with task performance was evident, indicating that an increased occurrence of harmonic locking (i.e., transient 2:1 ratios) was associated with improved arithmetic performance. These results are in line with previous evidence pointing to the importance of alpha–theta interactions in tasks requiring working memory and executive control (Kawasaki et al., 2010; Dimitriadis et al., 2016; Akiyama et al., 2017; Popov et al., 2018). In this respect, it has been suggested that alpha oscillations subserve short-term maintenance of information (Kawasaki et al., 2010; Chik, 2013; Akiyama et al., 2017) as well as reactivation of semantic knowledge (Klimesch et al., 1999), whereas theta oscillations underlie executive control (Schack et al., 2005; Lisman and Jensen, 2013). Based on our results, we propose that the integration between cognitive functions encompassed in alpha and theta rhythms is facilitated through 2:1 harmonic relationships, as this frequency arrangement would enable cross-frequency (phase) synchronization (and therefore information transmission; Fries, 2005, 2015) to occur between their underlying neural networks. In line with this notion, secondary analyses confirmed that harmonic relationships between alpha–theta rhythms mathematically enable stable and frequent phase meetings, by showing that a greater incidence of harmonic alpha–theta relationships was related to increased alpha–theta 2:1 phase synchrony. Together, and in line with recent theoretical accounts (Pletzer et al., 2010; Klimesch, 2012, 2013), our results suggest that changes in the peak frequency of brain oscillations during different tasks (Foffani et al., 2005; Haegens et al., 2014; Mierau et al., 2017) could represent a neural mechanism to enable phase coupling to occur between task-relevant rhythms through the formation of harmonic cross-frequency relationships.

Pletzer et al. (2010) previously demonstrated that the golden mean (i.e., 1.618…) constitutes the best possible ratio to avoid spurious phase synchronizations between adjacent rhythms and related to this property, it has been hypothesized that its occurrence in EEG cross-frequency interactions may reflect an optimal configuration for establishing an efficient decoupling mechanism (Pletzer et al., 2010; Klimesch, 2012, 2013). Indeed, in the theoretical accounts by Klimesch (2012, 2013) and Pletzer et al. (2010), the EEG frequency architecture was proposed to entail “optimal” frequency domains for facilitating cross-frequency “coupling” (i.e., based on harmonic numerical ratios), but also for facilitating controlled cross-frequency decoupling (i.e., based on nonharmonic numerical ratios approximating the golden mean 1.618…). In this view, the peak frequency within a frequency band is proposed to transiently shift to guarantee either maximal coupling or decoupling with neighboring frequency domains. For example, as underlined by Klimesch (2012), “with theta oscillating at a dominant frequency of 5 Hz, the peak frequency of the alpha band may transiently shift from 10 to 8 Hz, to obtain separation from theta (8/1.618 = 5 Hz), or may stay at 10 Hz to enable optimal coupling with theta”. Interestingly, our exploratory analysis revealed a proportionally high incidence (both during rest and arithmetic task performance) of alpha–theta cross-frequency relationships approximating the golden mean (i.e., 1.6; Fig. 4*A*), thereby providing initial experimental support to the notion that this configuration forms a prevalent physiological state within the EEG frequency architecture.

Furthermore, contrary to harmonic locking, it was shown that the occurrence of cross-frequency ratios around the golden mean (1.6) was higher during the rest compared with the arithmetic condition. This finding may provide initial support to the notion that nonharmonic cross-frequency relationships based on the golden mean reflect an optimal frequency arrangement for facilitating cross-frequency decoupling. During resting-state, the brain's neural circuitry may adjust its frequency architecture to a state in which cross-frequency decoupling is actively facilitated [i.e., to avoid spurious (unwanted) phase synchronization; Pletzer et al., 2010]. However, when cognitive task demands increase (as during the arithmetic task), the brain's neural circuitry may inversely need to adjust its frequency architecture to a state in which, on the one hand, cross-frequency decoupling (1.618 ratio) is prevented; and on the other hand, cross-frequency coupling (2.0 harmonic ratio) is facilitated (Pletzer et al., 2010; Klimesch, 2012, 2013). Note that prior studies investigating cross-frequency dynamics in the brain have largely overlooked the presence and relevance of nonharmonic cross-frequency relationships because cross-frequency phase synchrony can be only estimated between harmonic frequencies (Palva et al., 2005; Palva and Palva, 2018; Lobier et al., 2018). In this view, future studies are needed to further explore the occurrence of nonharmonic cross-frequency relationships and their functional relevance with respect to cross-frequency interactions and cognitive processing.

The here adopted method for assessing cross-frequency dynamics (i.e., based on the calculation of transient numerical ratios between peak frequencies of 2 brain rhythms) follows a holistic and probabilistic approach as it quantifies short-lived changes in the neural frequency architecture that would (mathematically) enable cross-frequency coupling or decoupling. Despite its simplicity, we anticipate that this novel approach to electrophysiological data can open new vistas for future basic and translational neuroscience research. Nevertheless, some important limitations of the adopted approach and obtained results pattern should be underlined. First, the detection of transient alpha and theta peak frequencies might be significantly improved by an a priori spatial localization of brain activity occurring at each frequency band (Haegens et al., 2014). This type of source localization analysis could clarify whether the occurrence of alpha–theta harmonic relationships differs significantly between different brain regions. Unfortunately, source localization could not be reliably performed in our data given the small number electrodes (Song et al., 2015), so these relevant questions relating to the spatial distribution of the occurrence of alpha–theta harmonic relationships remain open for future research. In the same way, future research may be warranted to further explore the influence of transient alpha–theta cross-frequency relationships on cognitive performance. In this context, neurofeedback and/or brain stimulation protocols can be envisaged to assess whether the up or down-training of specific cross-frequency relationships (i.e., harmonic vs nonharmonic) have a causal impact on cognitive performance.

In summary, the present study suggests that transient numerical ratios between the peak frequencies of two brain rhythms are of functional relevance for cognitive processing. Specifically, we showed that the transient occurrence of harmonic 2:1 cross-frequency relationships between alpha (8–14 Hz) and theta (4–8 Hz) oscillations (mathematically facilitating phase synchrony) increased during an arithmetic task and correlated to improved cognitive performance. In addition, exploratory analysis revealed that the increase in harmonic 2:1 cross-frequency relationships during arithmetic task occurred at the expense of a reduced occurrence of set of nonharmonic cross-frequency relationships that would mathematically hinder phase synchrony. Our findings are interpreted in the light of recent theoretical accounts (Pletzer et al., 2010; Klimesch, 2012, 2013) positing that transient peak frequency shifts could encompass a flexible mechanism to facilitate cross-frequency coupling (harmonic relationships) and decoupling (nonharmonic relationships) between different brain rhythms.

## Footnotes

This work was supported by the Branco Weiss fellowship of the Society in Science–ETH Zurich and by Grants from the Flanders Fund for Scientific Research (FWO projects KAN 1506716N and G079017N). We thank Dante Mantini for his advice, and Laura Timmermans and Phaedra Lebegge for their help during data collection.

The authors declare no competing financial interests.

- Correspondence should be addressed to Kaat Alaerts at Kaat.Alaerts{at}kuleuven.be