Introduction

Neural oscillations in the mammalian brain are thought to play an important role in coordinating neural activity across different brain regions, allowing for the integration of sensory information, the control of motor movements, and the maintenance of cognitive functions (Pfurtscheller and Da Silva, 1999; Caplan et al., 2003; Buzsaki and Draguhn, 2004; Jensen and Mazaheri, 2010; Giraud and Poeppel, 2012; Schalk, 2015; Fries, 2015). Detecting neural oscillations is important in neuroscience as it helps unravel the mysteries of brain function, understand brain disorders, investigate cognitive processes, track neurodevelopment, develop brain-computer interfaces, and explore new therapeutic approaches. Thus, detecting and analyzing the “when”, the “where”, and the “what” of neural oscillations is an essential step in understanding the processes that govern neural oscillations.

For example, detecting the onset and offset of a neural oscillation (i.e., the “when”) is necessary to understand the relationship between oscillatory power/phase and neural excitation, an essential step in explaining an oscillation’s excitatory or inhibitory function (Pfurtscheller and Da Silva, 1999; Canolty et al., 2006; Klimesch et al., 2007; Haegens et al., 2011; de Pesters et al., 2016). Localizing the brain area or layer that generates the oscillation (i.e., the “where”) provides neuroanatomical relevance to cognitive and behavioral functions (Buzsaki and Draguhn, 2004; Miller et al., 2010). Lastly, determining the oscillation’s fundamental frequency (i.e., the “what”) indicates underlying brain states (Penfield and Jasper, 1954; Buzsaki and Draguhn, 2004). Together, the “when”, the “where”, and the “what” can be seen as the fundamental pillars in investigating the role of oscillations in interregional communication throughout the brain (Fries, 2015). These fundamental pillars can also provide insight into the functional purpose (i.e., the “why”), underlying mechanisms (i.e., the “how”), and pathologies (i.e., the “whom”) of neural oscillations (Buzsaki and Draguhn, 2004; Buzsaki, 2006).

The detection of neural oscillations has historically been extensively studied in the frequency-(Wen and Liu, 2016; Donoghue et al., 2020; Ostlund et al., 2022), time-(Hughes et al., 2012; Gips et al., 2017), and time-frequency domains (Chen et al., 2011; Wilson et al., 2022; Neymotin et al., 2022). With the notable exception of Gips et al. 2017, these studies assume that neural oscillations are predominantly sinusoidal and stationary in their frequency. However, there is an increasing realization that neural oscillations are actually non-sinusoidal and exhibit spurious phase-amplitude coupling (Belluscio et al., 2012; Cole et al., 2017; Scheffer-Teixeira and Tort, 2016; Gips et al., 2017; Donoghue et al., 2022). A recent review paper on methodological issues in analyzing neural oscillations (Donoghue et al., 2022) identified determining the fundamental frequency of non-sinusoidal neural oscillations as the most challenging problem in building an understanding of how neural oscillations govern interactions throughout the brain.

Fast Fourier Transform (FFT) is the most commonly used method to detect neural oscillations. The FFT separates a neural signal into sinusoidal components within canonical bands of the frequency spectrum (e.g., theta, alpha, beta). The components of these canonical bands are typically considered to be functionally independent and involved in different brain functions. However, when applied to non-sinusoidal neural signals, the FFT produces harmonic phase-locked components at multiples of the fundamental frequency. While the asymmetric nature of the fundamental oscillation can be of great physiological relevance (Mazaheri and Jensen, 2008; Cole et al., 2017; Donoghue et al., 2022), its harmonics are considered to be an artifact produced by the FFT that can confound the detection and physiological interpretation of neural oscillation (Belluscio et al., 2012; Donoghue et al., 2022).

An example of an unfiltered electrocorticographic recording from auditory cortex (Figure 1A) illustrates the non-sinusoidal nature of neural oscillations. The associated FFT-based power spectrum (Figure 1B) exhibits multiple peaks over 1/f noise even though only one oscillatory signal is visibly present in the time domain signal. Whether the peaks over 1/f at 12 and 18 Hz, are harmonics of 6 Hz oscillations or independent oscillations remains unknown. This ambiguity affects the ability to accurately and efficiently identify neural oscillations and understand their role in cognition and behavior. For this illustrative example of non-sinusoidal neural oscillation, we used a phase-phase coupling analysis (Belluscio et al., 2012) to determine whether the exhibited 18 Hz beta oscillation is a harmonic of the 6 Hz theta oscillation. This analysis confirmed that the beta oscillation was indeed a harmonic of the theta oscillation (Figure 1E and F). In marked contrast, for a sinusoidal neural oscillation, a phase-phase coupling analysis could not fully ascertain whether the oscillations are phase-locked and thus are harmonics of each other (Figure 1G-L). This ambiguity, combined with the exorbitant computational complexity of the entailed permutation test and the requirement to perform the analysis across all cross-frequency bands over all channels and trials render phase-phase coupling impracticable for determining the fundamental frequency of neural oscillations in real-time and, thus, the use in closed-loop neuromodulation applications.

Examples of non-sinusoidal and sinusoidal neural oscillations recorded from the human auditory cortex.

Detecting the presence, onset/offset, and fundamental frequency of non-sinusoidal oscillations is challenging. This is because the power spectrum of the non-sinusoidal theta-band oscillation (A) exhibits multiple harmonic peaks in the alpha and beta bands (B). The peaks of these harmonics are also exhibited in the time-frequency domain (C). To determine whether these peaks are independent oscillations or harmonics of the fundamental frequency, we tested whether fundamental theta oscillation and potential beta-band harmonic oscillations exhibit a 1:3 phase-locking (D-F), i.e., whether the beta-band oscillation is a true 3rd harmonic of the fundamental theta-band oscillation. In our test, we found that the theta-band oscillation was significantly phase-locked to the beta-band oscillation with a 1:3 ratio in their frequencies (F). This means that the tested theta- and beta-band oscillations are part of one single non-sinusoidal neural oscillation. We applied the same statistical test to a more sinusoidal neural oscillation (G). Since this neural oscillation more closely resembles a sinusoidal shape, it does not exhibit any prominent harmonic peaks in the alpha and beta bands within the power spectrum (H) and time-frequency domain (I). Consequently, our test found that the phase of the theta-band and beta-band oscillations were not phase-locked (J-L). Thus, this statistical test suggests the absence of a harmonic structure.

In this study, we aim to define the principle criteria that characterize a neural oscillation and to synthesize these criteria into a method that accurately determines the duration (“when”), location (“where”), and fundamental frequency (“what”) of non-sinusoidal neural oscillations. For this purpose, we introduce the Cyclic Homogeneous Oscillation (CHO) detection method to identify neural oscillations using an auto-correlation analysis to identify whether a neural oscillation is an independent oscillation or a harmonic of another oscillation. Auto-correlation is a statistical measure that assesses the degree of similarity between a time series and a delayed version of itself.

Thus, auto-correlation can explain the periodicity of a signal without assuming that the signal is sinusoidal. Further, the peaks in the output of the auto-correlation function indicate the fundamental frequency of the neural oscillation. As shown in Figure 2, irrespective of the shape of neural oscillation (Figure 2A and C), the fundamental frequency can be determined from the positive peak-to-peak intervals (see Figure 2B and D). Despite auto-correlation being a well-known method to identify the fundamental frequency of a signal, its application to neural oscillations has been impeded by the requirement to accurately determine the onset and offset of the oscillation.

Using auto-correlation to determine the fundamental frequency of non-sinusoidal and sinusoidal neural oscillations recorded from the human auditory cortex.

(A) Temporal dynamics of non-sinusoidal (left) and sinusoidal (right) neural oscillation and (B) their auto-correlation. The periodicity of peaks in the auto-correlation reveals the fundamental frequency of the underlying oscillation. Asymmetry in peaks and troughs of the auto-correlation is indicative of a non-sinusoidal oscillation.

To overcome this limitation, we combine the auto-correlation method with the Oscillation Event (OEvent) method (Neymotin et al., 2022) to determine the onset/offset of oscillations. In this approach, OEvent determines bounding boxes in the time-frequency domain that mark the onset and offset of suspected oscillations. Each bounding box is generated by identifying a period of significantly increased power from averaged power spectrum. To further improve OEvent, we replaced the empirical threshold that identifies bounding boxes in the time-frequency domain with a parametric threshold driven by an estimation of the underlying 1/f noise (Donoghue et al., 2020), as shown in Figure 3A.

Procedural steps of CHO.

(A) First, to identify periodic oscillations, CHO removes the underlying 1/f aperiodic noise in the time-frequency space and generates initial bounding boxes of candidate oscillations. (B) In the second step, CHO rejects bounding boxes that exhibit less than two oscillatory cycles. (C) In the final step, CHO limits the analysis to only those bounding boxes that exhibit the same frequency in the time-frequency map and auto-correlation. Each remaining bounding box is characterized by onset/offset, frequency range, center frequency, and number of cycles.

Furthermore, we improved OEvent to reject any short-cycled oscillations that could represent evoked potentials (EP), event-related potentials (ERP), or spike activities, as shown in Figure 3B. In general, EPs or ERPs in neural signals generate less than two cycles of fluctuations. Large-amplitude EPs, ERPs, and spike activities can result in spurious oscillatory power in the frequency domain (de Cheveigné and Nelken, 2019; Donoghue et al., 2020, 2022).

In the final step, we determine the oscillation’s periodicity and fundamental frequency by identifying positive peaks in the auto-correlation of the signal. As shown for a representative oscillation in Figure 3C, the center frequency of the highlighted bounding box is 24 Hz, but the periodicity of the underlying raw signal does not match the calculated fundamental frequency of 7 Hz. Consequently, this bounding box at 24 Hz will be rejected. Finally, we merge those remaining bounding boxes that neighbor each other in the frequency domain and overlap more than 75% (Neymotin et al., 2022) in time.

In summary, the presented CHO method identifies neural oscillations that fulfill the following three criteria: 1) oscillations (peaks over 1/f noise) must be present in the time and frequency domains; 2) oscillations must exhibit at least two full cycles; and 3) oscillations must have auto-correlation. These criteria are supported by studies in the neuroscience literature (Buzsaki and Draguhn, 2004; Niedermeyer and da Silva, 2005; Buzsaki, 2006; Cohen, 2014; de Cheveigné and Nelken, 2019; Donoghue et al., 2020, 2022). The synthesis of these criteria into the presented method allows us to detect and identify non-sinusoidal oscillations and their fundamental frequency. This is because criteria #1 (i.e., the presence of an oscillation) and #2 (i.e., the length of the oscillation) identify potential oscillations, which are then tested to be fundamental oscillations using an auto-correlation analysis using criteria #3 (i.e., the periodicity of an oscillation).

To verify and validate CHO, we applied the above-presented principle criteria on simulated non-sinusoidal signals and human electrophysiological signals, including electrocorticographic (ECoG) signals recorded from the lateral brain surface, electroencephalographic signals (EEG) recorded from the scalp, and local field potentials recorded from the hippocampus using stereo EEG (SEEG). We further validated our approach by comparing CHO to other commonly used methods.

To determine the spectral accuracy in detecting the peak frequency of non-sinusoidal oscillations, we compared CHO to established methods, including the fitting of oscillations using 1/f (FOOOF, also known as specparam, Donoghue et al. 2020), the OEvent method (Neymotin et al., 2022), and the Spectral Parameterization Resolved in Time (SPRiNT, Wilson et al. 2022) methods. Moreover, to determine the spectro-temporal accuracy in detecting both the peak frequency and the onset/offset of non-sinusoidal oscillations, we compared CHO with the OEvent method.

The selection of FOOOF, SPRiNT, and OEvent is based on their fundamental approaches. To the best of our knowledge, FOOOF is the most representative method for detecting the peak frequency of neural oscillations. SPRiNT expands the FOOOF method into the time-frequency domain, and OEvent can determine the onset/offset of the detected oscillations.

Results

The following sections describe the results of our study: The first section presents simulation results by comparing the accuracy of CHO with that of existing methods in detecting non-sinusoidal oscillations. The second section reports physiological results by comparing the accuracy of CHO with that of established methods in detecting oscillations within in-vivo recordings.

Synthetic results

To determine the specificity and sensitivity of CHO in detecting neural oscillations, we applied CHO to synthetic non-sinusoidal oscillatory bursts (2.5 cycles, 1–3 seconds long) convolved with 1/f noise, also known as pink noise, which has a power spectral density that is inversely proportional to the frequency of the signal. As shown in Figure 4, we generated 5s-long 1/f signals composed of pink noise and added non-sinusoidal oscillations of different lengths (one cycle, two-and-a-half cycles, 1s-duration, and 3s-duration). The rightmost panel of Figure 4A shows two examples of non-sinusoidal oscillations (two-and-a-half cycles and 2s-duration) along with their power spectra. As can be seen in Figure 4A, longer non-sinusoidal oscillations exhibit stronger harmonic peaks.

Performance of CHO in detecting synthetic non-sinusoidal oscillations.

(A) We evaluated CHO by verifying its specificity, sensitivity, and accuracy in detecting the fundamental frequency of non-sinusoidal oscillatory bursts (2.5 cycles, 1–3 seconds long) convolved with 1/f noise. (B-D) CHO outperformed existing methods in detecting the fundamental frequency of non-sinusoidal oscillation (FOOOF: fitting oscillations one over f (Donoghue et al., 2020), OEvent (Neymotin et al., 2022): Oscillation event detection method, and SPRiNT (Wilson et al., 2022): Spectral Parameterization Resolved in Time) in specificity and accuracy, but not in sensitivity. CHO exhibited fewer false-positive and more true-negative detections than existing methods. (C) However, at SNR-levels of alpha oscillations found in EEG and ECoG recordings (i.e., -7 dB and -6 dB, respectively), the sensitivity of CHO in detecting the peak frequency of non-sinusoidal oscillation is comparable to that of SPRiNT. (D) This means that the overall accuracy of CHO was higher than that of existing methods. (E-G) CHO outperformed existing methods in detecting the fundamental frequency and onset/offset of non-sinusoidal oscillation. (F) Similar to the results shown in (C) CHO can effectively detect the fundamental frequency and onset/offset for more than half of all oscillations at SNR-levels of alpha oscillations found in EEG and ECoG recordings.

Figure 4—figure supplement 1. SNR Histograms of EEG and ECoG.

Figure 4—figure supplement 2. Synthetic sinusoidal oscillations.

Our results in Figure 4B-D demonstrate that CHO outperforms conventional techniques in specificity and accuracy for detecting the peak frequency of non-sinusoidal oscillations. High specificity depends on high true-negative and low false-positive rates. For conventional methods, we expected harmonic oscillations to increase the false-positive rate and one-cycled oscillations to decrease the true-negative rate. As expected, conventional methods detected harmonic and one-cycled oscillations as true oscillations. For example, the average specificity of SPRiNT was below 0.3, which was significantly lower than the robust specificity of CHO across the entire range of SNR. We also observed that CHO requires a higher SNR to detect the presence of oscillations. Sensitivity depends on the true-positive and the false-negative rates. We found existing methods to be overly sensitive in detecting the presence of oscillations. At the same time, this severely limits their specificity and, thus, their ability to accurately detect the presence and frequency of an oscillation. Based on our physiological datasets, we found the average SNR of oscillations in EEG and ECoG to be -7 dB and -6 dB, respectively (Figure 4—figure Supplement 1). When tested at these physiologically-motivated SNR levels, and found that the sensitivity of CHO is comparable to that of SPRiNT. Overall, when considering the accuracy combined with specificity and sensitivity, CHO outperformed all other methods in detecting the peak frequency of non-sinusoidal oscillations at the physiologically motivated SNR levels.

In addition to determining the accuracy in detecting the presence of oscillations and determining their peak frequency, we also determined the accuracy of all methods in detecting the onset and offset of oscillations. This comparison is limited to OEvent because FOOOF and SPRiNT methods cannot determine the onset and offset of short oscillations. In this analysis, CHO outperformed the OEvent method in specificity but not sensitivity, as shown in Figure 4E-G. Specifically, we found performance trends similar to those in our previous simulation result (Figure 4B-D). Thus, CHO outperforms conventional techniques in specificity for detecting both the peak frequency and on-set/offset of oscillations.

Empirical results

We further assessed CHO by testing it on electrophysiological signals recorded from human subjects. Specifically, we evaluated CHO on electrocorticographic (ECoG, x1–x8, 8 subjects) and electroencephalographic (EEG, y1–y7, 7 subjects) signals recorded during the pre-stimulus period of an auditory reaction time task. Furthermore, we also evaluated CHO on signals recorded during resting state from cortical areas and hippocampus using ECoG (ze1–ze8, N=6) and stereo EEG (zs1–zs6, 6 subjects).

Electrocorticographic (ECoG) results

In the auditory reaction time task, we expected to observe neural low-frequency oscillations during the pre-stimulus period within task-relevant areas, such as the auditory and motor cortex. As we expected, we found alpha and beta oscillations within these cortical areas. We compared the topographic distribution of the oscillations detected by FOOOF with those detected by CHO. As shown in Figure 5A for one representative subject, FOOOF detected the presence of alpha, and beta oscillations within temporal and motor cortex. In contrast, while CHO also detected alpha oscillations in temporal and motor cortex, it only detected beta oscillations in motor cortex. We found this pattern to be consistent across subjects, as shown in Figure 5B and Figure 5—figure Supplement 1.

Validation of CHO in detecting oscillations in ECoG signals.

(A) We applied CHO and FOOOF to determine the fundamental frequency of oscillations from ECoG signals recorded during the pre-stimulus period of an auditory reaction time task. FOOOF detected oscillations primarily in the alpha- and beta-band over STG and pre-motor area. In contrast, CHO also detected alpha-band oscillations primarily within STG, and more focal beta-band oscillations over the pre-motor area, but not STG. (B) We investigated the occurrence of each oscillation within defined cerebral regions across eight ECoG subjects. The horizontal bars and horizontal lines represent the median and median absolute deviation (MAD) of oscillations occurring across the eight subjects. An asterisk (*) indicates statistically significant differences in oscillation detection between CHO and FOOOF (Wilcoxon rank-sum test, p<0.05 after Bonferroni correction).

Figure 5—figure supplement 1. ECoG results using FOOOF and CHO for all subjects.

We compared neural oscillation detection rates between CHO and FOOOF across eight ECoG subjects. We used FreeSurfer (Fischl, 2012) to determine the associated cerebral region for each ECoG location. Each subject performed approximately 400 trials of a simple auditory reaction-time task. We analyzed the neural oscillations during the 1.5-second-long pre-stimulus period within each trial. CHO and FOOOF demonstrated statistically comparable results in the theta and alpha bands despite CHO exhibiting smaller median occurrence rates than FOOOF across eight subjects. Notably, within the beta band, excluding specific regions such as precentral, pars opercularis, and caudal middle frontal areas, CHO’s beta oscillation detection rate was significantly lower than that of FOOOF (Wilcoxon rank-sum test, p < 0.05 after Bonferroni correction). This suggests comparable detection rates between CHO and FOOOF in premotor and Broca’s areas, while the detection of beta oscillations by FOOOF in other regions, such as the temporal area, may represent harmonics of theta or alpha, as illustrated in Figure 5A and B. Furthermore, FOOOF exhibited a higher sensitivity in detecting delta, theta, and low gamma oscillations overall, although both CHO and FOOOF detected only a limited number of oscillations in these frequency bands.

Electroencephalographic (EEG) results

We expected that the EEG would exhibit similar results as seen in the ECoG results. Indeed, the EEG results mainly exhibit alpha and beta oscillations during the pre-stimulus periods of the auditory reaction time task, as shown in Figure 6. Specifically, FOOOF found alpha oscillations in mid-frontal and visual areas and beta oscillations throughout all areas of the scalp. In contrast, CHO found more focal visual alpha and pre-motor beta. Furthermore, the low gamma oscillations detected by CHO were also more focal than those detected by FOOOF. We found these results to be consistent across subjects (see Figure 6B, C and Figure 6—figure Supplement 1).

Validation of CHO in detecting oscillations in EEG signals.

(A) We applied CHO and FOOOF to determine the fundamental frequency of oscillations from EEG signals recorded during the pre-stimulus period of an auditory reaction time task. FOOOF primarily detected alpha-band oscillations over frontal/visual areas and beta-band oscillations across all areas (with a focus on central areas). In contrast, CHO detected alpha-band oscillations primarily within visual areas and detected more focal beta-band oscillations over the pre-motor area, similar to the ECoG results shown in Figure 5. (B) We investigated the occurrence of each oscillation within the EEG signals across seven subjects. An asterisk (*) indicates statistically significant differences in oscillation detection between CHO and FOOOF (Wilcoxon rank-sum test, p<0.05 after Bonferroni correction). CHO exhibited lower entropy values of alpha and beta occurrence than FOOOF across 64 channels. (C) We compared the performance of FOOO and CHO in detecting oscillation across visual and pre-motor-related EEG channels. CHO detected more alpha and beta oscillations in visual cortex than in pre-motor cortex. FOOOF detected alpha and beta oscillations in visual cortex than in pre-motor cortex.

Figure 6—figure supplement 1. All EEG results using FOOOF and CHO.

We assessed the difference in neural oscillation detection performance between CHO and FOOOF across seven EEG subjects. We used EEG electrode locations according to the 10-10 electrode system (Nuwer, 2018) and assigned each electrode to the appropriate underlying cortex (e.g., O1 and O2 for the visual cortex). Each subject performed 200 trials of a simple auditory reaction-time task. We analyzed the neural oscillations during the 1.5-second-long pre-stimulus period. In the alpha band, CHO and FOOOF presented statistically comparable outcomes. However, CHO exhibited a greater alpha detection rate for the visual cortex than for the pre-motor cortex, as shown in Figure 6B and C. The entropy of CHO’s alpha oscillation occurrences (3.82) was lower than that of FOOOF (4.15), with a maximal entropy across 64 electrodes of 4.16. Furthermore, in the beta band, CHO’s entropy (4.05) was smaller than that of FOOOF (4.15). These findings suggest that CHO may offer a more region-specific oscillation detection than FOOOF. As illustrated in Figure 6C, CHO found fewer alpha oscillations in pre-motor cortex (FC2 and FC4) than in occipital cortex (O1 and O2), while FOOOF found more beta oscillations occurrences in pre-motor cortex (FC2 and FC4) than in occipital cortex. However, FOOOF found more alpha and beta oscillations in visual cortex than in pre-motor cortex. Consistent with ECoG results, FOOOF demonstrated heightened sensitivity in detecting delta, theta, and low gamma oscillations. Nonetheless, both CHO and FOOOF identified only a limited number of oscillations in delta and theta frequency bands. Contrary to the ECoG results, FOOOF found more low gamma oscillations in EEG subjects than in ECoG subjects.

Onset and offset of neural oscillations

So far, we have established that CHO can localize beta rhythms within pre-motor cortex in EEG and ECoG. Here, we are interested in determining the accuracy of the onset/offset detection of neural oscillations. For this purpose, we tested whether CHO, applied to signals recorded from auditory cortex during an auditory reaction-time task, can accurately detect the transition between resting and task periods. Specifically, we expected CHO to detect the offset times of neural oscillations after the stimulus onset (i.e., a beep tone that remained until a button was pressed). Based on the principle of event-related de-/synchronization (ERD/ERS, Pfurtscheller and Da Silva 1999), cortical neurons may be de-synchronized to process an auditory stimulus. As shown in Figure 7, CHO successfully detected offset times of 7 Hz neural oscillations. During the pre-stimulus period, the distribution of the onset time remains uniform, reflecting the subject waiting for the stimulus. In contrast, after the stimulus onset, the distribution of onset times becomes Gaussian, reflecting the variable reaction time to the auditory stimulus. Of note, the detection of onset times peaks 950 ms post-stimulus, which occurs significantly later than the button press that happens 200 ms post-stimulus (Figure 7B).

Application of CHO in determining the spatiotemporal characteristics of neural oscillations in ECoG signals during a reaction-time task.

(A) We selected those cortical locations (red) from all locations (black) that exhibited a significant broadband gamma response to an auditory stimulus in a reaction-time task. (B) In this task, the subjects were asked to react as fast as possible with a button press to a salient auditory stimulus. (C-D) Onset and offset times of detected neural oscillations. Fundamental oscillations were centered around 7 Hz (left histogram). Onset and offset times during pre-stimulus period exhibited a uniform distribution, indicating that 7 Hz oscillations randomly started and stopped during this period. A trough in the onset and a peak in the offset of 7 Hz oscillations is visible from the histograms, indicating a general decrease of the presence of neural oscillations immediately following the auditory stimulus. The subjects responded with a button press within 200 ms of the auditory stimulus, on average. The prominent peak in the offset and onset of oscillations at 300 ms and 950 ms post-stimulus, respectively, indicates a suspension of oscillations in response to the auditory stimulus, and their reemergence after the execution of the button press behavior.

Similar to the distribution of onset times, the distribution of offset times remained uniform throughout the pre-stimulus period. After stimulus onset, the distribution becomes Gaussian, with a peak of offset detections at 300 ms post-stimulus, or 200 ms post-response (i.e., the button press) (Figure 7C).

In summary, this means that, on average, the detected 7 Hz oscillations de-synchronized 250 ms and synchronized 900 ms, post-stimulus, respectively.

Stereoelectroencephalographic (SEEG) results

We also investigated neural oscillations within the hippocampus. Specifically, we were interested in the frequency and duration of hippocampal oscillations, which are known to be non-sinusoidal and a hallmark of memory processing (Buzsaki, 2006; Lundqvist et al., 2016). Using the CHO method, we plotted a representative example of detected hippocampal fast theta bursts (Lega et al., 2012; Goyal et al., 2020), as shown in Figure 8. As expected, the non-sinusoidal alpha-band oscillations also resulted in harmonic oscillations in the beta band, which, while not clearly visible in the power spectrum (Figure 8B), can be clearly seen in the time-frequency analysis (Figure 8D and Figure 8E). In contrast to the ECoG and EEG results, the frequency of beta-band oscillations in the hippocampus exhibited a frequency close to the alpha-band (8–14 Hz). CHO found primarily alpha-band oscillations in the hippocampus (see Figure 8—figure Supplement 1). When comparing the consistency between CHO and FOOOF across hippocampal locations, CHO exhibits more specific results with less overlap between alpha and beta locations and almost no detection in the low-gamma band (30–40 Hz). For example, subject zs4 in Figure 8—figure Supplement 1 shows alpha and beta locations mutually supplement each other when using CHO but not when using the FOOOF method. However, we did not find a statistically significant difference between CHO and FOOOF due to the small number of subjects and variable electrode locations within hippocampus across the six SEEG subjects.

Application of CHO in determining the fundamental frequency and duration of hippocampal oscillations in SEEG signals during resting state.

(A) We recorded hippocampal oscillations from one representative human subject implanted with SEEG electrodes within the left anterior hippocampus. (B) Power spectrum (blue) and 1/f trend (red) for one electrode within the anterior-medial left hippocampus (red dot in A). The power spectrum of a 10-second-long hippocampal signal indicates the presence of neural activity over a 1/f trend across a wide frequency band up to 30 Hz. (C) In marked contrast to the relatively unspecific results indicated by the power spectrum, CHO detected several distinct hippocampal fast theta bursts. (D) This detection is based on first denoising the power spectrum using 1/f fitting (principle criterion #1 of CHO), which yields initial bounding boxes that include short-cycled oscillations and harmonics. (E) The auto-correlation step then successfully removes all short-cycled oscillations and harmonics, with only those bounding boxes remaining that exhibit a fundamental frequency.

Figure 8—figure supplement 1. All results from six SEEG subjects using the FOOOF and CHO methods.

Frequency and duration of neural oscillations

Here, we are interested in identifying the predominant frequency and duration of neural oscillations for specific brain areas during the resting state. For this purpose, we first determined the specific Brodmann area of each recording electrode using an intracranial electrode localization tool, Versatile Electrode Localization Framework (VERA, Adamek et al. 2022). Next, we investigated electrodes belonging to the primary auditory cortex (i.e., Brodmann areas 41 and 42), as shown in Figure 9A. We found that 7 and 11 Hz oscillations were the predominant neural oscillations for electrodes near the primary auditory cortex. The average duration of an 11 Hz oscillation was 450 ms. Next, our results for primary motor cortex (i.e., Brodmann area 4) showed that 7 Hz was the pre-dominant oscillation frequency in the motor cortex with 450 ms duration on average, as shown in Figure 9B. We found that motor cortex exhibits more beta-band oscillations (around 500 ms duration) than the auditory cortex. Next, Broca’s area exhibited characteristics similar to those of the motor cortex, however, with a predominant beta-band frequency of 17 Hz, which is lower than the 22 or 24 Hz oscillations found in the motor cortex (Figure 9C). Lastly, using SEEG electrodes, we investigated neural oscillations within the human hippocampus (Figure 9D). This analysis showed that 8 Hz was the predominant oscillatory frequency in the hippocampus with a 450 ms duration on average. During the resting state, neural alpha- and beta-band oscillations within the hippocampus were shorter than in the motor cortex (p<0.05, Wilcoxon rank sum test, N=6).

Application of CHO in determining the fundamental frequency and duration of neural oscillations in auditory cortex, motor cortex, Broca’s area, and hippocampus during resting state.

This figure presents the distribution of detected oscillations in a 2-dimensional frequency/duration histogram and projected onto frequency and duration axes. The red line indicates the rejection line (less than two cycles). (A) In primary auditory cortex (Brodmann area 41/42), the most dominant frequency and duration in the auditory cortex was 11 Hz with 450 ms duration. (B) The primary motor cortex’s most dominant frequency was 7 Hz with 450 ms duration, but more beta rhythms were detected with >500 ms duration than in auditory cortex. (C) Broca’s area exhibits similar characteristics to that of motor cortex, but dominant beta-band oscillations were found to be less present than in motor cortex. (D) Hippocampus primarily exhibits 8 Hz oscillations with 450 ms duration.

Discussion

Our novel CHO method demonstrates high precision and specificity in detecting neural oscillations in time and frequency domains. The method’s specificity enables the detailed study of spatiotemporal dynamics of oscillations throughout the brain and the investigation of oscillatory biomarkers that index functional brain areas.

High specificity for detecting neural oscillations

In our simulation study, CHO demonstrated high specificity in detecting both the peak and on-set/offset of neural oscillations in time and frequency domains. This high specificity directly results from the three criteria we established in this study. The first criterion was that neural oscillations (peaks over 1/f noise) must be present in the time and frequency domain. The 1/f trend estimation served as a threshold to reject aperiodic oscillatory power in the neural signals (Donoghue et al., 2020).

Next, the second condition was that oscillations must exhibit at least two complete cycles. This condition distinguishes periodic oscillations from evoked/event-related potentials (EP/ERP) and spike artifacts. EP/ERP have spectral characteristics that are similar to those of theta or alpha frequency oscillations. To discriminate EP/ERPs from genuine oscillations, we reject them if they don’t exhibit peaks over 1/f or if they have fewer than two cycles.

The third and final condition is that oscillations should share the same periodicity as their auto-correlation. This is because positive peaks in the auto-correlation can identify the oscillation’s fundamental frequency even if it is non-sinusoidal. The bounding boxes help us to identify possible onsets/offsets of neural oscillations. Moreover, calculating the auto-correlation of the raw signals within a bounding box provides the true periodic frequency of the raw signal. We then reject any bounding boxes for which the periodicity of the raw signal is not in alignment with the true periodic frequency revealed by the auto-correlation. This third condition is important in rejecting harmonic peaks over 1/f noise in the frequency domain. Furthermore, it is also effective in rejecting spurious oscillations, which are broadly generated by spike activities in the frequency domain (de Cheveigné and Nelken, 2019).

To calculate the auto-correlation, we first needed to determine the onset/offset of the potential oscillations. The first and second criteria serve as a triage in finding the onset/offset of genuine oscillations. Thus, these three principle criteria were essential to reject aperiodic harmonic oscillations and increase CHO’s specificity in detecting both the peak frequency and the onset/offset of non-sinusoidal oscillations. We also evaluated CHO on purely sinusoidal oscillations (see Figure 4—figure Supplement 2). The results of this analysis show that even in the absence of any asymmetry in the oscillations, CHO still outperforms existing methods in specificity. It further shows that the sensitivity increases with increasing SNR. Even though this analysis is based on synthetic sinusoidal oscillations, our results demonstrated that existing methods are susceptible to noise which results in the detection of spurious oscillations. However, as expected, both FOOOF and SPRiNT methods exhibited reasonable specificity when applied to sinusoidal signals.

Focal localization of beta oscillations

Beta oscillations occur within the 13–30 Hz band throughout various brain regions, including the motor cortex. In the motor cortex, beta oscillations are thought to be involved in motor planning and execution. Studies have shown that beta oscillations increase and decrease in power during movement preparation and movement execution, respectively (Pfurtscheller and Da Silva, 1999; Jenkinson and Brown, 2011; Doyle et al., 2005; Senkowski et al., 2006). In our empirical results based on the presented ECoG dataset, CHO found focal beta oscillations to occur within pre-motor and frontal cortex prior to the button response, as shown in Figure 5. These findings were consistent across subjects. Conventional methods found alpha and beta oscillations in the auditory cortex, while CHO found only select beta oscillations. This suggests that most of the beta oscillations detected by conventional methods within auditory cortex may be simply harmonics of the predominant asymmetric alpha oscillation. Along the same line, conventional methods found beta and low gamma oscillations in pre-motor and frontal areas, while CHO found predominantly beta oscillations. This suggests that low gamma oscillations detected by conventional methods are harmonics of beta oscillations.

In the EEG results, CHO found focal visual alpha and motor beta oscillations, while the FOOOF found frontal and visual alpha and beta oscillations across broad scalp areas, as shown in Figure 6. In contrast to the ECoG results, neither CHO nor FOOOF auditory found alpha oscillations within the temporal areas. This is interesting as FOOOF exhibits a better sensitivity than CHO and suggests that auditory alpha rhythms may be difficult to observe in EEG. Similar to the ECoG results, our analysis confirmed that non-sinusoidal alpha and beta oscillations generate harmonic oscillations in both beta and low gamma in EEG. This shows that our CHO method, which has a high specificity, can detect focal motor beta oscillations.

Harmonic oscillations in human hippocampus

Recent studies suggest that the frequency range of hippocampal oscillations is wider than previously assumed (<40 Hz in Cole and Voytek 2019, or 3–12 Hz in Li et al. 2022) and that it does not match the conventional frequency range of theta/alpha rhythms (Buzsaki, 2006). This realization stems from the recognition that neural oscillations are non-sinusoidal, and thus require a wide frequency band to be fully captured (Cole and Voytek, 2019; Donoghue et al., 2022). Adopting a wider frequency band provides more frequency options in fitting the non-sinusoidal shape of brain waves. The recognition of the need to expand the frequency band within oscillation analysis is not limited to the hippocampus. Our ECoG and EEG results show that harmonics can occur in any brain area and frequency band because neural oscillations are inherently non-sinusoidal. A recent study showed that the phase of wide-band oscillations could better predict neural firing (Davis et al., 2020).

CHO can determine the fundamental frequency of non-sinusoidal oscillations when applied within a wide-band analysis, as shown in Figure 8E. Moreover, CHO provides onset/offset and the frequency range of an oscillation, allowing us to investigate non-sinusoidal features, such as the degree of asymmetry and amplitudes of trough/peak (Cole and Voytek, 2019).

Identifying onset/offset of neural oscillations and its application

Although the frequency of neural oscillation has been extensively investigated, the onset/offset and duration of neural oscillations have remained elusive. Using CHO, the onset/offset, and duration of neural oscillations can be revealed, as shown in Figure 7 and Figure 9. Knowing the onset/offset and duration of a neural oscillation is essential for realizing closed-loop neuromodulation. This is because neuromodulation may be most efficient when electrical stimulation is delivered phase-locked to the underlying ongoing oscillation (Chen et al., 2011; Cagnan et al., 2017, 2019; Zanos et al., 2018; Shirinpour et al., 2020). For example, deep-brain stimulation in phase with ongoing oscillation can reduce the stimulation necessary to achieve the desired therapeutic effect (Cagnan et al., 2017, 2019). This improved efficiency in delivering the stimulation therapy reduces power consumption and thus enhances the battery life of the implanted system (Chen et al., 2011). Longer battery life means fewer battery changes (which require surgical procedures), or for rechargeable systems, fewer recharging sessions (which require the user’s attention). Realizing phase-locked neuromodulation requires detecting the duration of an ongoing oscillation with high specificity and delivering the electrical stimulation at a predicted oscillation phase. The detection and identification with high specificity thus enable neuromodulation applications that depend on phase-locked electrical stimulation.

Moreover, the temporal precision of CHO in detecting neural oscillations can improve the effectiveness of neurofeedback-based systems. For example, a neurofeedback system may provide targeted feedback on the magnitude of the user’s alpha oscillation to improve attention and in turn improve task performance. For this purpose, the system must detect the frequency, onset/offset, and duration of the user’s alpha oscillation with high specificity. High specificity requires distinguishing other oscillations and artifacts from true physiological alpha-band oscillations. The identification of true neural oscillations with the high specificity of CHO thus enables targeted neurofeedback applications to enhance or restore task performance.

Illuminating the when, where, what, why, how, and whom of neural oscillations

In our study, we focused on the temporal dynamics (“when”), spatial distribution (“where”), and fundamental frequency (“what”) of neural oscillations. However, fully understanding the role of neural oscillations in cognition and behavior also requires investigating their underlying mechanisms (“how”), functional purpose (“why”), and pathologies (“whom”).

Temporal Dynamics – the “when”

CHO demonstrated high specificity in detecting the onset and offset of fundamental non-sinusoidal oscillations (see Figure 4E). Using CHO, our study revealed the temporal dynamics of oscillations within the temporal lobe in an auditory reaction-time task. We identified the onsets and offsets of 7 Hz oscillations and, thus, the boundaries in oscillatory activity between resting and task engagement. Our results show a rapid decrease in oscillatory activity for the duration of the auditory stimulus, followed by a rapid reemergence of the oscillatory activity following the cessation of the auditory stimulus (see Figure 7C and D). These results shed light on the temporal dynamics of neural oscillatory activity in cognitive processes and how the brain adapts to environmental stimuli.

Spatial Distribution – the “where”

CHO revealed the spatial distribution of neural oscillations in EEG, SEEG, and ECoG recordings. The spatial distribution of fundamental neural oscillations, and their absence during task engagement, can reveal underlying shared functional organization. CHO can be applied to a wide range of neu-roimaging techniques such as EEG, MEG, ECoG, and SEEG to elucidate the involvement of different brain regions in various cognitive functions. For example, using CHO, our study found focal specific alpha oscillations over occipital (visual) cortex in EEG and focal beta oscillations over parietal (motor) cortex in ECoG. These results demonstrate the utility of CHO in precisely mapping the spatial distribution of neural oscillations across the brain, and in revealing shared functional organization of brain networks.

Fundamental Frequency – the “what”

CHO revealed the fundamental frequencies of asymmetric neural oscillations recorded from the scalp, auditory cortex, motor cortex, Broca’s area, and hippocampus. Distinct brain states can be identified based on the fundamental frequency of their underlying neural oscillation. CHO showed high specificity in determining the fundamental frequency of synthetic non-sinusoidal oscillations (see Figure 4B). When applied to ECoG and SEEG signals, CHO revealed distinct fundamental frequencies of oscillations found within auditory cortex, motor cortex, Broca’s area, and hippocampus (see Figure 9). CHO can be applied in real time to detect the fundamental frequency and the onset/offset of neural oscillations. Characterizing neural oscillations in real time can make transitions in brain states observable to the investigator. For example, investigators can characterize brain dynamics during wakefulness, sleep, or specific cognitive tasks by tracking changes in oscillatory activity during different behavioral states. This information provides insights into the brain’s adaptability and flexibility in response to internal and external cues and could inform closed-loop neuromodulation.

Underlying Mechanisms – the “how”

Accurate detection of neural oscillations aids in deciphering the underlying mechanisms governing their generation and synchronization. In our study, we focused on determining the temporal dynamics, spatial distribution, and fundamental frequency of neural oscillations. The results of our study, and more specifically the CHO method itself, provide a methodological foundation to systematically study oscillatory connectivity and traveling oscillations throughout cortical layers and brain regions to create insights into unraveling the generating mechanism of neural oscillations. The information gained from such studies could create a better understanding of neural circuitry at the network level and could inform computational models that help refine our knowledge of the complex mechanisms underlying brain function.

Functional Purpose – the “why”

Neural oscillation detection plays a crucial role in uncovering the functional significance of oscillatory activity. In our study, CHO detected focal alpha oscillations over occipital (visual) cortex in EEG and focal beta oscillations over parietal (motor) cortex in ECoG during the pre-stimulus period of an auditory reaction-time task (see Figure 5 and Figure 6). The presence of these oscillations during the pre-stimulus period implicates visual-alpha and motor-beta oscillations in inhibition. We found the same inhibitory oscillatory phenomenon over the auditory cortex, however, with a fundamental frequency of 7 Hz, indicating functional independence between inhibitory oscillations found in visual, motor, and auditory cortex (see Figure 7C and D). The approach presented in this study could be expanded to studying attention, memory, decision-making, and more by correlating neural oscillations with specific cognitive processes. Further, applying cross-frequency and phase-amplitude coupling analysis to oscillations detected by CHO could illuminate the role of neural oscillations in facilitating information processing and communication between brain regions.

Pathologies – the “whom”

Detecting and characterizing neural oscillations has significant implications for the study of neurological and psychiatric disorders. For example, recent studies reported that patients affected by severe Parkinson’s disease exhibited more asymmetry between peak and trough amplitudes in beta oscillations (Cole et al., 2017; Jackson et al., 2019). The high specificity demonstrated by CHO in detecting asymmetric neural oscillations could benefit the investigation of neural pathologies. Specifically, CHO could improve the quality of asymmetry measurements by providing onset/offset detection of the beta oscillations with high specificity. Abnormalities in neural oscillations are often associated with various pathologies. Detecting and characterizing aberrant oscillatory patterns could lead to identifying biomarkers for specific disorders and insights into their underlying mechanisms. These advancements could aid the development of targeted therapies and treatments for these conditions.

Illuminating neural oscillations

Overall, developing a reliable neural oscillation detection method is crucial for advancing our understanding of brain function and cognition. The presented CHO method opens up new avenues of research by contributing to the investigation of temporal dynamics, spatial distribution, brain states, underlying mechanisms, functional purpose, and pathologies of neural oscillations. Ultimately, a comprehensive understanding of neural oscillations will deepen our knowledge of the brain’s complexity and pave the way for innovative approaches to treating neurological and psychiatric disorders.

Limitations

The results of this study show that our CHO method favors specificity over sensitivity when SNR is low. More specifically, CHO exhibited a low sensitivity due to the high false-negative rate in a low-SNR environment. This means that even though there are oscillations present in the recorded signals, CHO cannot detect them when they are drowned in noise. To investigate whether this is an issue in real-world applications, we determined the averaged SNR of alpha oscillations in EEG (-7 dB) and ECoG (-6 dB). Based on our evaluation of synthetic data, we found that at these physiologically-motivated SNR levels, CHO can detect 50–60% of all true oscillations. This sensitivity could be further improved by averaging across spatially correlated locations, e.g., within the hippocampus.

One potential approach to reducing the dependency of sensitivity on SNR is to apply a wavelet transform in the estimation of the time-frequency map of the signal. Wavelet transform can better capture short cycles of oscillations. Currently, CHO uses a Hilbert transform method rather than Wavelet or short-time fast Fourier transform (STFFT) because it is easy to implement in MATLAB and provides better control over the spectral shape (i.e., better accuracy in detecting peak frequency of oscillations, Cohen 2014). Despite the theoretical advantages of wavelet over Hilbert transform, in developing our CHO method, we found no significant differences when we used different approaches to estimate the time-frequency map. This finding is further supported by a comparative study shown by Bruns in 2004. However, because our CHO method is modular, the FFT-based time-frequency analysis can be replaced with more sophisticated time-frequency estimation methods to improve the sensitivity of neural oscillation detection. Specifically, a state-space model (Matsuda and Komaki, 2017; Beck et al., 2022; Brady and Bardouille, 2022; He et al., 2023) or empirical mode decomposition (EMD, Fabus et al. 2022; Quinn et al. 2021) may improve the estimation of the auto-correlation of the harmonic structure underlying non-sinusoidal oscillations. Furthermore, a Gabor transform or matching pursuit-based approach may improve the onset/offset detection of short burst-like neural oscillations (Kuś et al., 2013; Morales and Bowers, 2022).

Another avenue to improve the sensitivity of CHO is to modify the third criterion to better distinguish neural oscillations from background noise. When we performed each detection step within CHO, as shown in Figure 3, we captured oscillations in a low-SNR situation. However, applying the third criterion rejected many possible bounding boxes. Thus, developing a better conceptual framework to reject harmonic peaks in the spectral domain may decrease the false-negative rate and, in turn, increase the sensitivity in low-SNR situations.

Another limitation of this study is that it does not assess the harmonic structure of neural oscillations. Thus, CHO cannot distinguish between oscillations that have the same fundamental frequency but differ in their non-sinusoidal properties. This limitation stems from the objective of this study, which is to identify the fundamental frequency of non-sinusoidal neural oscillations. Over-coming this limitation requires further studies to improve CHO to distinguish between different non-sinusoidal properties of pathological neural oscillations. The data that is necessary for these further studies could be obtained from the wide range of studies that have linked the harmonic structures in the neural oscillations to various cognitive functions (van Dijk et al., 2010; Schalk, 2015; Mazaheri and Jensen, 2008) and neural disorders (Cole et al., 2017; Jackson et al., 2019; Hu et al., 2023). For example, Cole and Voytek 2019 showed that a harmonic structure of beta oscillations can explain the degree of Parkinson’s disease, and Hu et al. 2023 showed the number of harmonic peaks can localize the seizure onset zone.

Conclusions

Neural oscillations are thought to play an important role in coordinating neural activity across different brain regions, allowing for the integration of sensory information, the control of motor movements, and the maintenance of cognitive functions. Thus, better methods to detect and characterize neural oscillations, especially those that are asymmetric, can greatly impact neuroscience. In this study, we present Cyclic Homogeneous Oscillation (CHO) as a method to reveal the “when”, the “where”, and the “what” of neural oscillations. With this method, we overcome the confounding effect of detecting spurious oscillations that result from harmonics of the non-sinusoidal neural oscillations (Donoghue et al., 2022). In our study, we demonstrate that solving this problem yields scientific insights into local beta oscillations in pre-motor areas, the onset/offset of oscillations in the time domain, and the fundamental frequency of hippocampal oscillations. These results demonstrate the potential for CHO to support closed-loop neuromodulation (brain-computer interfaces and neurofeedback) and neural oscillation detection systems to implement various neurological diagnostic and therapeutic systems and methods.

Methods and Materials

Electrophysiological data

Eight human subjects implanted with ECoG electrodes (x1–x8, 4 females, average age = 41±14) participated in an auditory reaction time task at the Albany Medical Center in Albany, New York. The subjects were mentally and physically capable of participating in our study (average IQ = 96±18, range 75–120, Wechsler 1997). All subjects were patients with intractable epilepsy who underwent temporary placement of subdural electrode arrays to localize seizure foci before surgical resection. The implanted electrode grids were approved for human use (Ad-Tech Medical Corp., Racine, WI; and PMT Corp., Chanhassen, MN). The platinum-iridium electrodes were 4 mm in diameter (2.3 mm exposed), spaced 10 mm center-to-center, and embedded in silicone. The electrode grids were implanted in the left hemisphere for seven subjects (x1, x3, x6, and x7) and the right hemisphere for five subjects (x2, x4, x5, and x8). Following the placement of the subdural grids, each subject had postoperative anterior-posterior and lateral radiographs and computer tomography (CT) scans to verify grid location. These CT images, in conjunction with magnetic resonance imaging (MRI), were used to construct three-dimensional subject-specific cortical models and derive the electrode locations (Coon et al., 2016).

A further seven healthy human subjects (y1–y7, all males, average age = 27±3.6) served as a control group for which we recorded EEG while performing the same auditory reaction time task. These subjects were fitted with an elastic cap (Electro-cap International, Blom and Anneveldt 1982) with tin (Polich and Lawson, 1985) scalp electrodes in 64 positions according to the modified 10-20 system (Acharya et al., 2016).

In addition, six human subjects implanted with ECoG electrodes (ze1–ze6, 1 female, mean age 46, range between 31 and 69) participated in resting state recording at the Albany Medical Center in Albany, New York. All six subjects had extensive electrode coverage over the lateral STG. Patients provided informed consent to participate in the study, and additional verbal consent was given prior to each testing session. The Institutional Review Board at Albany Medical Center approved the experimental protocol. Electrodes were comprised of platinum-iridium and spaced 3-10mm (PMT Corp., Chanhassen, MN).

Lastly, six human subjects implanted with SEEG electrodes (zs1–zs6, three females, average age = 46±16.6) participated in resting state recordings at the Barnes Jewish Hospital in St. Louis, Missouri. All subjects were patients with intractable epilepsy who underwent temporary placement of subdural electrodes to localize seizure foci prior to surgical resection. All subjects provided informed consent for participating in the study, which was approved by the Institutional Review Board of Washington University School of Medicine in St. Louis.

The implanted SEEG electrodes were approved for human use (Ad-Tech Medical Corp., Racine, WI; and PMT Corp., Chanhassen, MN). The platinum-iridium electrodes were 2 mm in length (0.8 mm diameter) and spaced 3.5–5 mm center-to-center. Following the placement of the stereo EEG electrodes, each subject had postoperative anterior-posterior and lateral radiographs and computer tomography (CT) scans to verify electrode location. These postoperative CT images, in conjunction with preoperative magnetic resonance imaging (MRI), were used to construct three-dimensional subject-specific cortical models and derive the electrode locations (Coon et al., 2016).

Data collection

We recorded EEG, ECoG, and SEEG signals from the subjects at their bedside using the general purpose Brain-Computer Interface (BCI2000) software (Schalk et al., 2004), interfaced with eight 16-channel g.USBamp biosignal acquisition devices (for EEG), one 256-channel g.HIamp biosignal acquisition device (g.tec., Graz, Austria, for ECoG), or one Nihon Kohden JE-120A long-term recording system (Nihon Kohden, Tokyo, Japan, for SEEG) to amplify, digitize (sampling rate 1,200 Hz for EEG and ECoG and 2,000 Hz for SEEG) and store the signals. To ensure safe clinical monitoring of ECoG signals during the experimental tasks, a connector split the cables connected to the patients into a subset connected to the clinical monitoring system and a subset connected to the amplifiers.

Task

The subjects performed an auditory reaction task, responding with a button press to a salient 1 kHz tone. For their response, the subjects used their thumb contralateral to their ECoG implant. In total, the subjects performed between 134 and 580 trials. Throughout each trial, the subjects were first required to fixate and gaze at the screen in front of them. Next, a visual cue indicated the trial’s start, followed by a random 1–3 s pre-stimulus interval and, subsequently, the auditory stimulus. The stimulus was terminated by the subject’s button press or after a 2-s time out, after which the subject received feedback about his/her reaction time. This feedback motivated the subjects to respond as quickly as possible to the stimulus. We penalized subjects with a warning tone to prevent false starts if they responded too fast (i.e., less than 100 ms after stimulus onset). We excluded false-start trials from our analysis. We were interested in this task’s auditory and motor responses in this study. This required defining the onset of these two responses. We time-locked our analysis of the auditory response to the onset of the auditory stimulus (as measured by the voltage between the sound port on the PC and the loudspeaker). For the motor response, we time-locked our analysis to the time when the push button was pressed. To ensure the temporal accuracy of these two onset markers, we sampled them simultaneously with the EEG/ECoG signals using dedicated inputs in our biosignal acquisition systems. We defined baseline and task periods for the auditory and motor response. Specifically, we used the 0.5-s period prior to the stimulus onset as the baseline for the auditory response and the 1-s to 0.5-s period prior to the button press as the baseline for the motor response. Similarly, we used the 1-s period after stimulus onset as the task period for the auditory response and the period from 0.5-s before to 0.5-s after the button press as the task period for the motor task.

Data pre-processing

As our amplifiers acquired raw, unfiltered EEG/ECoG/SEEG signals, we removed any offset from our signals using a 2nd-order Butterworth highpass filter at 0.05 Hz. Next, we removed any common noise using a common median reference filter (Liu et al., 2015). To create the common-mode reference, we excluded signals that exhibited an excessive 60 Hz line noise level (i.e., ten times the median absolute deviation). To improve the signal-to-noise ratio of our recordings and to reduce the computational complexity of our subsequent analysis, we downsampled our signals from 1200 Hz or 2000 Hz to 400 Hz or 500 Hz, respectively, using MATLABs “resample” function, which uses a polyphase antialiasing filter to resample the signal at the uniform sample rate.

Phase-phase coupling

To demonstrate phase-locking, as illustrated between theta and beta oscillations in Figure 1E and Figure 1K, we utilized the n:m phase-phase coupling method described in Belluscio et al. 2012. Specifically, we calculated the “mean radial distance”: , where j indexes the samples in time, and N represents the number of samples (epoch length in seconds × sampling frequency in Hz). Rn:mequals 1 when Δϕnm(tj) is constant for all time samples tj, and 0 when Δϕnm is uniformly distributed. Of note, Δϕnm(tj) equals nϕf1 (tj) - mϕf2 (tj), with f1 and f2 being two different frequency bands.

A novel oscillation detection method

We propose a novel method based on principle criteria to identify neural oscillations’ when, where, and what. The principle criteria are as follows: 1. Oscillations (peaks over 1/f noise) must be present in the time and frequency domain. 2. Oscillations must exhibit at least two full cycles. 3. The periodicity of an oscillation is the fundamental frequency of the oscillation. The procedural steps of CHO adhere to these principle criteria, as shown in Figure 3 and Algorithm 1. First, we apply a time-frequency analysis to determine power changes for each frequency component over time. To measure the significant spectral power increase over the time domain, we use the 1/f fitting technique as the principal threshold. In other words, the proposed method only considers those oscillations that emerge above the underlying 1/f noise. Thus, any oscillation with smaller power than 1/f noise is not considered to be an oscillation. To accomplish this, we subtract the underlying 1/f noise within the time-frequency domain. Specifically, we divide the time domain into four periods and estimate the minimum 1/f aperiodic fit across these periods (see Line 5 in Algorithm 1). After the subtraction of the underlying 1/f noise, we calculate the averaged power difference between the signal and the 1/f noise (named sigma). If the spectral power exceeds two times sigma, we consider the oscillation to exhibit significant power above the 1/f noise (see Line 12 in Algorithm 1). Next, we cluster time points with significant power over 1/f noise to generate initial bounding boxes as shown in Figure 3A; this idea is adopted from a previous study (Neymotin et al., 2022) (see Lines 10–20 in Algorithm 1).

Next, as the second principle criterion, we only consider those oscillations that exhibit at least two full cycles. This restriction allows CHO to distinguish oscillations from confounding event-related potentials (ERPs) or evoked potentials (EPs). In general, the frequency characteristics of those potentials often overlap with neural oscillations (e.g., theta power of ERPs and theta power of theta rhythm). However, ERPs or EPs never exhibit more than two cycles (see Line 23 in Algorithm 1). Therefore, we reject those bounding boxes that exhibit less than two cycles. An example is shown in Figure 3B.

Algorithm 1

CHO detection method

Lastly, we calculate the periodicity of an oscillation using an autocorrelation analysis to determine the fundamental frequency of the oscillation. Non-sinusoidal signals are known to exhibit harmonics in the frequency domain, significantly increasing the false-positive detection rate —the confounding factor addressed by CHO’s third criterion. The power spectrum of the non-sinusoidal oscillations has additional harmonic peaks over 1/f noise, even though the periodicity of the signal does not match the harmonic peak frequency. Therefore, the positive peaks of the oscillation’s autocorrelation represent the oscillation’s periodicity and fundamental frequency (see Figure 2). As shown in Figure 3C, the center frequency of the bounding box is 24 Hz, but the periodicity of the raw signal within the bounding box does not match 24 Hz. Consequently, this bounding box will be rejected (see Line 30 in Algorithm 1). Finally, the method merges those remaining bounding boxes that neighbor each other in the frequency domain and that overlap more than 75% in time (Neymotin et al., 2022).

The MATLAB code that implements CHO and sample data is available on GitHub (https://github.com/neurotechcenter/CHO).

Tradeoffs in adjusting the hyper-parameters that govern the detection in CHO

The ability of CHO to detect neural oscillations and determine their fundamental frequency is governed by four principal hyper-parameters. Adjusting these parameters requires understanding their effect on the sensitivity and specificity in the detection of neural oscillations.

The first hyper-parameter is the number of time windows (N in Line 5 in Algorithm 1), that is used to estimate the 1/f noise. In our performance assessment of CHO, we used four windows, resulting in estimation periods of 250 ms in duration for each 1/f spectrum. A higher number of time windows results in smaller estimation periods and thus minimizes the likelihood of observing multiple neural oscillations within this time window, which otherwise could confound the 1/f estimation. However, a higher number of time windows and, thus, smaller time estimation periods may lead to unstable 1/f estimates.

The second hyper-parameter defines the minimum number of cycles of a neural oscillation to be detected by CHO (see Line 23 in Algorithm 1). In our study, we specified this parameter to be two cycles. Increasing the number of cycles increases specificity, as it will reject spurious oscillations. However, increasing the number also reduces sensitivity as it will reject short oscillations.

The third hyper-parameter is the significance threshold that selects positive peaks within the auto-correlation of the signal. The magnitude of the peaks in the auto-correlation indicates the periodicity of the oscillations (see Line 26 in Algorithm 1). Referred to as “NumSTD,” this parameter denotes the number of standard errors that a positive peak has to exceed to be selected to be a true oscillation. For this study, we set the “NumSTD” value to 1. Increasing the “NumSTD” value increases specificity in the detection as it reduces the detection of spurious peaks in the auto-correlation. However, increasing the “NumSTD” value also decreases the sensitivity in the detection of neural oscillations with varying instantaneous oscillatory frequencies.

The fourth hyper-parameter is the percentage of overlap between two bounding boxes that trig-ger their merger (see Line 31 in Algorithm 1). In our study, we set this parameter to 75% overlap. In-creasing this threshold yields more fragmentation in the detection of oscillations, while decreasing this threshold may reduce the accuracy in determining the onset and offset of neural oscillations.

Validation on synthetic non-sinusoidal oscillations

While empirical physiological signals are most appropriate for validating our method, they generally lack the necessary ground truth to characterize neural oscillation with sinusoidal or non-sinusoidal properties. To overcome this limitation, we first validated CHO on synthetic non-sinusoidal oscillatory bursts (2.5 cycles, 1–3 seconds long) convolved with 1/f noise to test the performance of the proposed method.

As shown in Figure 4, we generated five second-long periods comprised of 1/f noise (i.e., pink noise). We added non-sinusoidal oscillations with different amplitudes and lengths. The amplitudes of non-sinusoidal oscillations vary between 5 and 20 microvolts, while the pink noise remains at 10 microvolts in amplitude. The signal-to-noise (SNR) was calculated by the snr() function in the Signal Processing Toolbox of MATLAB, which determines the signal-to-noise ratio in decibels of the non-sinusoidal burst by computing the ratio between summed squared magnitudes of the oscillation and the pink noise, respectively. We simulated ten iterations for each amplitude. For each iteration, we tested four different lengths of non-sinusoidal oscillations (one cycle, two-and-a-half cycles, one second, and three seconds long).

We generated non-sinusoidal oscillations by introducing asymmetry between the trough and peak periods of sinusoidal waves. To generate this asymmetric nature of an oscillation, we applied a 9:1 ratio between trough and peak amplitudes, as shown in an example of Figure 4A. To smooth the onset and offset of the non-sinusoidal oscillations, we used Tukey (tapered cosine) window function with a 0.40 ratio for the taper section (Bloomfield, 2004). Of note, the smaller the Turkey ratio within the taper section, the higher the occurrence of high-frequency artifacts.

To evaluate the performance of CHO, we calculated the specificity and sensitivity of CHO in detecting non-sinusoidal oscillations. High specificity depends on high true-negative and low false-positive detection rates. In contrast, high sensitivity depends on high true-positive and low false-negative detection rates. In this simulation, we expected harmonic oscillations to increase the false-positive detection rate, and one-cycled oscillations to decrease the true-negative detection rate within conventional methods. Thus, harmonic oscillations and one-cycled oscillations decrease the specificity, not sensitivity.

For evaluating the performance of each method in determining the fundamental frequency of the oscillations, we defined an accurate detection as one that exhibited a difference between the ground truth peak frequency and detected frequency of less than 1.5 Hz. Furthermore, to evaluate the performance of each method in detecting the onset/offset of the oscillations, we calculated the correlation between the envelope of the ground truth oscillation and the detected oscillation. We defined those onset/offset detections as accurate if the correlation was positive and the p-value was smaller than 0.05.

Acknowledgements

This work was supported by the National Institutes of Health (NIH) grants R01-MH120194, R01-EB026439, U24-NS109103, U01-NS108916, U01-NS128612, P41-EB018783, the McDonnell Center for Systems Neuroscience and Fondazione Neurone.

Authors’ Contributions

Conceptualization: HC, PB. Methodology: HC, MA, JTW, PB. Data Curation: PB. Formal analysis: HC. Visualization: HC, MA. Funding acquisition: JTW, PB. Writing – original draft: HC, PB. Writing – review & editing: HC, MA, JTW, PB. All authors read and approved the final version of the manuscript.

Competing Interests

One U.S. patent (Provisional Application Serial No.63/326,257) related to systems and methods for detection of neurophysiological signal oscillations described in this manuscript was filed on March 31, 2022. The inventors/contributors of this patent involve some of the manuscript authors, including HC, MA, JTW, PB.

Data Availability

Datasets may be provided to interested researchers upon reasonable request to the corresponding author.

Code Availability

The Matlab code and sample data used for CHO are available at https://github.com/neurotechcenter/ CHO.

SNR Histograms of EEG (A) and ECoG (B).

Synthetic sinusoidal oscillations.

ECoG results using FOOOF and CHO for all subjects.

Results from seven EEG subjects using the FOOOF and CHO methods.

All results from six SEEG subjects using the FOOOF and CHO methods.