β Band Rhythms Influence Reaction Times

Abstract Despite their involvement in many cognitive functions, β oscillations are among the least understood brain rhythms. Reports on whether the functional role of β is primarily inhibitory or excitatory have been contradictory. Our framework attempts to reconcile these findings and proposes that several β rhythms co-exist at different frequencies. β Frequency shifts and their potential influence on behavior have thus far received little attention. In this human magnetoencephalography (MEG) experiment, we asked whether changes in β power or frequency in auditory cortex and motor cortex influence behavior (reaction times) during an auditory sweep discrimination task. We found that in motor cortex, increased β power slowed down responses, while in auditory cortex, increased β frequency slowed down responses. We further characterized β as transient burst events with distinct spectro-temporal profiles influencing reaction times. Finally, we found that increased motor-to-auditory β connectivity also slowed down responses. In sum, β power, frequency, bursting properties, cortical focus, and connectivity profile all influenced behavioral outcomes. Our results imply that the study of β oscillations requires caution as β dynamics are multifaceted phenomena, and that several dynamics must be taken into account to reconcile mixed findings in the literature.


Introduction
b Rhythms (;13-30 Hz) are traditionally associated with the sensorimotor system where they are prominent (Pfurtscheller and Lopes da Silva, 1999). Beyond this sensorimotor role, b has been implicated in a wide range of cognitive phenomena including visual perception (Piantoni et al., 2010;Kloosterman et al., 2015), language processing (Weiss and Mueller, 2012), working memory (Axmacher et al., 2008;Siegel et al., 2009), long-term memory (Hanslmayr et al., 2016), decision-making (Wimmer et al., 2016;Wong et al., 2016), and reward processing (Marco-Pallarés et al., 2015). In nonhuman primates, b was shown to reflect top-down attention (Buschman and Miller, 2007), and in rodents, b was linked to working memory (Parnaudeau et al., 2013;Bolkan et al., 2017). However, the functional role of b is still unclear (Engel and Fries, 2010;Kilavik et al., 2012), as some studies report decreased b with task engagement (Hanslmayr et al., 2009), suggesting an inhibitory function, while others report the opposite (Kornblith et al., 2016), suggesting an excitatory function. Similarly, on the neural level, there have been mixed and contradictory findings on the relationship between b and other neural measures such as firing rate (Rule et al., 2017) and blood oxygenation level-dependent (BOLD) activity (Hanslmayr et al., 2011).
Current accounts of b mechanism and function have tried to reconcile these findings (Engel and Fries, 2010;Spitzer and Haegens, 2017). One account states that b band activity is related to the maintenance of the current sensorimotor or cognitive state via a top-down mechanism (Engel and Fries, 2010). Our account suggests that b band activity is involved in (re)activating latent sensorimotor and cognitive states (Spitzer and Haegens, 2017). We further propose that several b rhythms co-exist, including functionally inhibitory b as predominantly observed in sensorimotor regions, and functionally excitatory b as observed throughout cortex. These different b rhythms possibly operate at different frequencies (Spitzer and Haegens, 2017). At the neurophysiological level, we posit that while b events are likely excitatory in nature, there are several biologically plausible ways they could lead to functional inhibition, for example by activating inhibitory neurons or saturating excitatory neurons (Shin et al., 2017;Spitzer and Haegens, 2017).
b Activity has been characterized and modeled as transient, high-amplitude events or "bursts," which can be detected at the single-trial level (Lundqvist et al., 2016;Sherman et al., 2016). b Bursts have been observed both focally (Bonaiuto et al., 2021) and as part of long-range communication between brain regions, where b band synchrony is assumed to facilitate interareal connectivity (Seedat et al., 2020). One property of b that has received little attention is instantaneous variability in its peak frequency (Cohen, 2014). Here, we asked how frequency shifts within the b band influence behavior.
The influence of b on behavioral outcomes might depend on several factors such as b power, frequency, bursting properties, cortical focus, and connectivity profile. In the current experiment we investigated the relationship between single-trial b activity and behavior, specifically reaction times. Since analyzing neural activity in a prestimulus or pretarget interval is a convenient method to uncover the influence of ongoing neural activity on subsequent behavior (Rassi et al., 2019b), we made use of magnetoencephalography (MEG) data recorded during an auditory sweep discrimination task. To test how the various characteristics of b relate to behavior, we analyzed reaction times as a function of pretarget b differences in power, shifts in frequency, and bursting profiles, within and between motor and auditory cortices.

Participants
We recorded MEG in 35 adult participants, 28 of which we included in our analyses (22 female; mean age = 22.86 years, SD = 2.84; three participants excluded because of excessively noisy MEG data and four because of near-chance performance on the task). The study was approved by the local ethics committee. All participants gave informed consent before the experiment and were given monetary compensation for their participation.

Auditory target discrimination task
The auditory target discrimination task consisted of five rhythmic blocks and five nonrhythmic blocks (60 trials per block). The order of the blocks was randomized. In the rhythmic blocks, four cue tones were presented, separated by 0.5 s. Following the rhythmic cue, a target tone was presented at 0.5, 1, 1.5, or 2 s (80% of trials) or at 0.75, 1.25, 1.75 s (20% of trials) after the onset of the last cue tone. In the nonrhythmic blocks, the cue tone was presented continuously for a period of 1.5 s, followed by a target that was presented with a flat probability distribution within a window of 0.5-2 s after cue offset. The task for the participant was to determine whether the target tone (a 40-ms chirp) went up or down in pitch. As we had previously shown the experimental cueing manipulation not to produce behaviorally different effects (Wilsch et al., 2020;Lin et al., 2021), here we pooled all trial types (Fig.  1a).

Stimuli
The cue tones had a pitch frequency of 400 Hz, a sample rate of 44,100 Hz, and a duration of 40 ms (rhythmic blocks) or 1.5 s (nonrhythmic blocks). We used a Hanning taper to remove sharp edges. The target tone consisted of 30 different frequencies randomly drawn from within 500-1500 Hz. The target was a frequency-modulated sweep created with the MATLAB function chirp and was either increasing or decreasing in pitch. The sound had a 10ms cosine ramp fading in and fading out to avoid onset and offset click perception. The resulting target tone had a sample rate of 44,100 Hz and a duration of 40 ms.
We normalized all sounds (using peak normalization) to the same sound pressure level. We individually adjusted target stimuli to participants' discrimination threshold, using a custom adaptive-tracking procedure aiming for a discrimination performance between 65% and 85% correct responses. The threshold was the slope of the pitch increase and decrease, measured as the range from lowest to highest frequency (starting point to end). We presented each participant with a pair of sounds ("up" and "down") consisting of the a priori randomly generated frequencies, modulated depending on their individual threshold, i.e., the 30 base-frequencies were the same for each participant, with the participant's individual threshold changing the start and end frequencies of the sounds.

Data acquisition
Whole-head MEG data were recorded at a 1200-Hz sampling rate with a 275-channel CTF MEG system with axial gradiometers (CTF MEG Systems, VSM MedTech Ltd.) in a magnetically shielded room. To monitor the participants' head movements online and for offline co-registration of anatomic landmarks, three fiducial coils were placed at the nasion and both ear canals. Anatomical, T1weighted MRI scans for source localization purposes were obtained in a separate session, using either a 1.5 or 3 T Siemens MRI system (Siemens). To co-register the MEG and MRI data, we additionally mapped the scalp with Polhemus 3D (Polhemus).

MEG and MRI preprocessing
We processed the MEG data offline with the Fieldtrip toolbox (Oostenveld et al., 2011). First, we down-sampled the data to a sampling frequency of 300 Hz. We then applied a notch filter at 50 Hz to remove line noise. Next, we segmented trials into 6-s segments starting 1 s before cue onset. We rejected bad channels (;5%) and bad trials (;10%) via visual inspection before independent component analysis (runICA as implemented in Fieldtrip), which was used to visually detect and remove components representing eye blinks and heartbeats.
We co-registered the MRI to the CTF coordinate system using the fiducial points and the mapped scalp surface, and segmented the MRI image with SPM8 (as implemented in Fieldtrip).

MEG source reconstruction
We used the obob_ownft toolbox for source reconstruction (https://gitlab.com/obob/obob_ownft). In order to model virtual sensors at the locations of maximum evoked activity in both the auditory and motor sources in the right and left hemispheres respectively, we used a linearly constrained minimum variance (LCMV) beamformer approach (Van Veen et al., 1997). We first constructed volume conduction models of the participants' brains using a single-shell model of their individual anatomic scans (Nolte, 2003), which we then used to compute leadfields for each of 3000 gridpoints. Using these leadfields, we computed common spatial filters for each participant using time windows that included a baseline period and the evoked responses.
For the auditory source, we used a time window of 100 ms centered at the peak of the individual auditory evoked response, time-locked to the onset of the auditory cue, and a 100-ms baseline window before cue onset. For the motor source, we used an activation time window of 100 ms centered at the peak of the individual motor response, timelocked to the button press, and a 100-ms baseline window before the activation window). We then normalized the difference of the sources of the pre and post windows and projected onto their co-registered anatomic scans (Fig. 1b). For visualization, we normalized each participant's brain to Montreal Neurologic Institute space. We then identified the location of maximum pre versus post differences in auditory and motor sources in the right and left hemispheres, respectively. Using the spatial filters for these positions, we then extracted the time series for these two virtual channels.

Data analysis
We performed all further data analysis using the Fieldtrip toolbox (Oostenveld et al., 2011) and custom MATLAB code. We time-locked the source-reconstructed a b Figure 1. Trial sequence and region of interest definition. a, After a variable baseline delay (1-2 s), an auditory cue lasting 1.5 s played, followed by a variable pretarget delay (0.5-2 s). This pretarget delay was our analysis window. After target onset, participants responded as fast as possible, indicating via button press whether the target tone shifted upward or downward in pitch. b, Regions of interest (ROIs) were defined as the source location with maximum evoked activation versus baseline, based on the evoked response to the auditory cue for the auditory cortex ROI (left panel) and based on the evoked response to the button press for the motor cortex ROI (right). Showing source reconstruction for one representative subject (with a 95%-maximum activity threshold applied for illustrative purposes).
Research Article: New Research signals from auditory and motor cortices to the onset of the target tones and analyzed a 700-ms pretarget interval. To counteract the 1/f effect in the data, we took the derivative of the time-series data (Thomson and Emery, 2014). Note that whether or not we removed the 1/f component had no influence on any of our results. However, b activity in the resulting flattened spectra was more visually salient, so we used those for visualization.

Spectral power
To compare pretarget b power with baseline (i.e., the precue period) activity, and to test the relationship between reaction times and b power, we extracted 700 ms of precue and pretarget data and computed single-trial Fourier spectra (0-30 Hz) with a fast Fourier approach and a Hanning taper, padded to 2 s for a frequency resolution of 0.5 Hz. We log-transformed the single-trial power data and extracted power in the b frequency range (13-30 Hz). For the reaction time contrast, we split the pretarget data along the median reaction time, and averaged the power spectra for faster and slower reaction times. We tested for group-level differences in both contrasts with a clusterbased permutation approach (Maris and Oostenveld, 2007), clustering across b frequencies. We estimated effect sizes by calculating Cohen's d based on the average difference in the data within a cluster (Meyer et al., 2021).
For time-resolved analyses, we performed time-frequency transformation based on multiplication in the frequency domain, using a sliding time window of 250 ms in steps of 20 ms from À750 to 1250 ms relative to target onset, in steps of 0.5 Hz between 13 and 30 Hz. We then averaged power within this window to obtain a per-region normalization factor and divided each time-frequency point by that factor. Finally, we averaged across the frequency dimension within the b band and extracted single-trial b time-courses for the 700-ms pretarget window.
To test the relationship between the b time course and reaction time, we z-scored the power values and reaction times, removed those with z-values above three and below À3, and used linear regression (reaction time = b power * slope 1 intercept), relating each single-trial, pretarget time point of b power with the subsequent reaction time on that trial. This provided a time course of regression slopes per participant. We then generated time-courses of slopes obtained by randomly shuffling the correspondence between power values and reaction times, and tested for group-level differences between the real and shuffled data with a cluster-based permutation approach (Maris and Oostenveld, 2007), clustering across the time dimension (À700-0 ms). We estimated effect sizes by calculating Cohen's d based on the average difference in the data within a cluster (Meyer et al., 2021). We also obtained the R 2 values associated with each of those regressions, took the maximum across time points per participant, and reported the average across participants.

Burst properties
To examine b burst properties in the source-reconstructed signals, we used the time-frequency representations as described above. We computed the mean and standard deviation of power within a trial for each frequency, and marked the time-frequency points that exceeded two standard deviations above the mean and that lasted at least the duration of one cycle (defined as 1/frequency). We zoomed in on the precue and pretarget delays, and based on temporal and spectral adjacency, we clustered the marked time-frequency points into burst events. We then extracted six parameters of interest from these burst events: for each trial, we counted the number of burst events. Focusing on the event that contained the time-frequency point with the highest power, we extracted the maximum power, the time point with maximum power, the frequency with maximum power, the frequency range, and the time range. This gave us single-trial estimates of burst properties during the precue and pretarget intervals.
To contrast precue and pretarget burst properties at the group level, we used paired t tests, and estimated effect sizes by calculating Cohen's d. To examine the relationship between reaction times and burst properties, we zscored the burst properties in the pretarget interval and reaction times, and used linear regression analysis to relate them. To test for group-level relationships, we used paired t tests contrasting the regression slopes against a shuffled distribution.

Instantaneous frequency
To investigate the time course of the peak b frequency in the source-reconstructed signals, we analyzed instantaneous frequency as detailed by Cohen (Cohen, 2014). Briefly, we band-passed the single-trial data within the b frequency range, applied the Hilbert transform, extracted the phase angle time series, took the temporal derivative, and applied ten median filters. This resulted in single-trial time-series of instantaneous frequency during the precue and pretarget interval.
To compare pretarget instantaneous frequency with baseline, and to test the relationship between pretarget b frequency and reaction times (based on median split), we used a cluster-based permutation approach (Maris and Oostenveld, 2007), and estimated effect size by calculating Cohen's d based on the average difference in the data within a cluster (Meyer et al., 2021). To relate single-trial time-series of instantaneous frequency with reaction times, we used the same regression approach detailed above (see Spectral power) to obtain a time course of regression slopes, and tested them at the group level with a cluster-based permutation approach (Maris and Oostenveld, 2007).
To test the relationship between reaction time and peak frequencies in the power spectra, we averaged power spectra separately for slower and faster trials (based on median split), and detected the peaks of maximum power within the b range. We then contrasted the peaks at the group level with a paired t test, and estimated effect sizes by calculating Cohen's d.

Connectivity
To estimate the connectivity between auditory and motor cortices, we used the Fourier coefficients that we obtained in the spectral power analysis. As connectivity measures are Research Article: New Research not resolved on single-trials, we estimated them after splitting the data along the median reaction time. We computed the pairwise phase consistency (Vinck et al., 2010), a biasfree method of rhythmic synchronization. We also computed bi-variate, nonparametric Granger causality (Dhamala et al., 2008a,b), which gave us separate estimates of the connection strengths from motor to auditory cortex and vice versa. We finally contrasted slower versus faster trials on the group level with a cluster-based permutation approach (Maris and Oostenveld, 2007), clustering across the b band. We estimated effect sizes by calculating Cohen's d based on the average difference in the data within a cluster (Meyer et al., 2021).

Data and code availability
The data and code supporting the findings of this study can be found here: https://doi.org/10.17605/OSF.IO/ 4AYS5.

Results
Pretarget versus precue power, frequency, and burst properties First, we examined pretarget b properties in relation to a baseline (i.e., precue) interval (Fig. 2). We found that in motor cortex, b power decreased, and in both motor and auditory cortex, b frequency increased from baseline to pretarget interval.
In auditory cortex, we observed the upward shift in b frequency primarily from 210 to 75 ms before target onset ( Fig. 2e; p = 0 (Fig. 2f). There were no differences in spectral power ( Fig. 2d; Extended Data Fig. 2-1, bottom row), maximum power of the bursts, their peak frequency, or their timing relative to stimulus onset. Next, we tested whether pretarget b properties related to reaction times using two complementary approaches: a median-split approach to relate b measures to slow versus fast reaction times, and a regression approach to relate single-trial b measures to reaction times. The two approaches yielded the same results: in motor cortex, slower reaction times were related with higher b power, while in auditory cortex, slower reaction times were related with higher b frequency.

Spectral power
In motor cortex (Fig. 3a-c), slower reaction times were preceded by higher b power. Splitting the power spectra across the median reaction time revealed the effect was driven by differences in the 20-to 26-Hz frequency range ( Fig. 3a; cluster-based permutation test across frequencies: p = 1e-5, d = 0.65). This effect was present throughout the whole pretarget interval when looking at the time-resolved power envelopes ( Fig. 3b; cluster-based permutation test across time: p = 6e-4, d = 0.58). A timeresolved single-trial regression approach confirmed this effect as well ( Fig. 3c; slope = 0.030; cluster-corrected p = 4e-4; R 2 = 0.016). To further characterize this difference, we zoomed in on the b bursting profile and found that slower reaction times were preceded by more bursts (slope = 0.040; p = 0.0017) with wider time spans (slope = 0.037; p = 0.0277) and wider frequency spans (slope = 0.039; p = 0.0143).
In auditory cortex, when splitting the data along the median reaction time and contrasting the power spectra, we found an effect opposite to that observed in motor cortex, such that faster reaction times were preceded by higher b power (Fig. 3d; cluster-corrected p = 0.016, d = À0.63), an effect driven by differences in the 17-to 19-Hz frequency range. However, pretarget b power was not related to reaction times when using the time-resolved regression approach (slope = À0.008, no clusters). Given the discrepancy in results between the two approaches, we further investigated the observed difference in auditory cortex as a possible shift in peak frequency.

b Frequency
In auditory cortex (Fig. 3d-f), slower reaction times were preceded by a higher peak b frequency when splitting the power spectra across the median reaction time and detecting participants' individual peak b frequencies (t (24) = 2.4, p = 0.026, d = 0.48). This effect was most pronounced around 660-450 ms before target onset when looking at instantaneous frequency ( Fig. 3e; cluster-based permutation test across time: p = 0.006, d = 0.64). A time-resolved single-trial regression approach confirmed the same result ( Fig. 3f; slope = 0.011; cluster-corrected p = 0.032; R 2 = 0.019). When zooming in on the peak burst frequencies, we found the same relationship again (slope = 0.029; p = 0.030). In motor cortex, b frequency was not related to reaction times (slope = À0.005; no significant clusters).

Connectivity
We then asked whether auditory-motor b connectivity was related to reaction times (using the median-split approach). Slower reaction times were preceded by increased b connectivity between auditory and motor cortices, as quantified with pairwise phase consistency (p = 0.027, d = 0.57). The difference was most prominent at frequencies from 19 to 20 Hz (Fig. 4a). We then used Granger causality to check the directionality of this effect. There were no differences in auditory-to-motor b connectivity ( Fig. 4b; no significant clusters), but slower reaction times were preceded by increased motor-to-auditory b connectivity (p = 0.047, d = 0.31), most prominently at frequencies from 20 to 21 Hz (Fig. 4c).

Discussion
In an auditory target discrimination task, we sought to uncover the relationship between reaction times and various characteristics of the b rhythm. We found that slower (as compared with faster) reaction times were preceded by increased b power in motor cortex, increased b frequency in auditory cortex, and increased motor-to-auditory connectivity in the b range. The results were robust across our analysis approaches. We used a regression approach to relate single-trial reaction times to b measures, as well as a median-split approach to relate slower versus faster trials to changes in b measures, with both approaches yielding the same pattern of results. We further analyzed b activity separately in a time-resolved manner, a frequency-resolved manner, and by characterizing its burst profile, with all approaches yielding the same pattern of results.
b Activity over somato-motor cortex is traditionally viewed as a component of the m rhythm (the other component being motor alpha activity). This b rhythm has been associated with slower movement (or lack of movement) and therefore thought to reflect cortical inhibition (Jasper and Penfield, 1949;Gastaut, 1952;Hari, 2006). More recently, b activity has been found to occur in "bursts" of high-amplitude activity, and its bursting properties have been linked to impaired performance in somatosensory perception and attention tasks (Shin et al., 2017). Our results are consistent with the literature on somato-motor b activity, as we observed increased motor b before slower responses, and we found this increased activity likely reflected an increased number of pretarget bursts with wider time and frequency ranges.
In auditory cortex, b activity has been observed in tandem with alpha activity, and assumed to serve a similarly inhibitory function (Weisz et al., 2011Leske et al., 2014). However, in the current dataset, pretarget auditory b power was not robustly related to reaction times. This null finding precludes us from drawing decisive conclusions about auditory b power modulations. Instead, we found that an upward shift in b frequency and increased connectivity with motor cortex were related to slower reaction times. This result can be interpreted in light of the frequency-matching notion (Lowet et al., 2017). That is, on trials with slow responses, the difference in peak frequencies between auditory and motor cortex is reduced, resulting in stronger inhibitory synchronization.
Based on this set of results, we speculate that the b rhythm potentially serves different functions and operates via different mechanisms in different cortical locations. On the one hand, a local change in b power could relate to the local excitability of a brain region in an "a-like" way (Spitzer and Haegens, 2017), meaning it could reflect inhibitory processes (Klimesch et al., 2007; Jensen and a b c Figure 4. Auditory-motor cortex connectivity. a, Pairwise phase consistency between auditory and motor cortices for slow versus fast reaction times. b, Granger causality from auditory to motor cortex for slow versus fast reaction times. c, Granger causality from motor to auditory cortex for slow versus fast reaction times. Shaded regions around the line graphs represent the SEM. Dotted lines represent significant clusters (p , 0.05).
Mazaheri, 2010). On the other hand, a change in b frequency could relate to b 's possible role in interareal synchrony (Bressler and Richter, 2015). For example, a shift in the frequency of b connectivity was recently shown to reflect categorical decisions (Rassi et al., 2022). It is also possible that b rhythms at different frequencies serve different functions. Although a crude distinction between "higher" and "lower" b across cortical locations has been previously made (Kopell et al., 2011), attempts to assign them different functional roles have had mixed success (Spitzer and Haegens, 2017). The nonstationarity of b rhythm (and cortical rhythms in general) frequencies across time has so far received little attention. Models of b function account for the potential of different b rhythms occurring at different frequencies, assuming different cortical locations or different generators within a location. But it is so far under-appreciated that a single rhythm can shift in frequency over time (Cohen, 2014;Rassi et al., 2019a). Frequency shifts according to task demands have been observed in human EEG/MEG data for the alpha rhythm (Haegens et al., 2014;Samaha and Postle, 2015;Mierau et al., 2017;Wutz et al., 2018), and in nonhuman primate LFP data for the b rhythm (Kilavik et al., 2012). It has also been reported that slower a rhythms correlate with slower responses across subjects (Surwillo, 1961), but to our knowledge the relationship between b frequency and reaction times has not yet been investigated. We here report the opposite relationship for the b rhythm, such that faster (auditory) b correlated with slower reaction times within subjects.
Beyond local b dynamics, b has also been shown to be involved in long-range communication between cortical sites (Seedat et al., 2020). Here, we found increased b connectivity between motor and auditory cortex, specifically in the direction of motor to auditory cortex, before slower (vs faster) responses. It is unlikely that this effect was confounded by the power difference in motor cortex as we used a phase-based connectivity measure, and in addition, there were no robust power differences in auditory cortex. This finding is in line with the notion of covert active sensing, where the motor system actively coordinates sensory systems (Schroeder et al., 2010). In addition, oscillatory bursts have been described as channels for selective communication between brain regions, via a mechanism called frequency-division multiplexing (Akam and Kullmann, 2014). In this view, bursts at different frequencies act as distinct channels to selectively transmit neural codes to networks (regions) with appropriate filter settings that can selectively read out the codes. This interpretation links our connectivity results with the perspective that b is occurring in bursts.
Finally, our results imply that the analysis of b oscillations requires caution as b dynamics are multifaceted phenomena. For example, it is possible that observed power modulations are better explained as frequency shifts (as is the case for our results). It is also possible that the b rhythm serves different functions (i.e., inhibitory or excitatory) depending on the cortical region where it is found or depending on whether it is local or interareal. Future investigations could focus on interareal variability in b peak frequency, for example in intracranial human electrophysiological recordings.