Abstract
Several lines of inquiry have separately identified beta oscillations, synchrony, waveform shape, and phase-amplitude coupling as important but sometimes inconsistent factors in the pathophysiology of Parkinson's disease. What has so far been lacking is a means by which these neurophysiological parameters are interrelated and how they relate to clinical symptomatology. To clarify the relationship among oscillatory power, bursting, synchrony, and phase-amplitude coupling, we recorded local field potentials/electrocorticography from hand motor and premotor cortical area in human subjects with c (N = 10) and Parkinson's disease (N = 22) during deep brain stimulator implantation surgery (14 females, 18 males). We show that motor cortical high beta oscillations in Parkinson's disease demonstrate increased burst durations relative to essential tremor patients. Notably, increased corticocortical synchrony between primary motor and premotor cortices precedes motor high beta bursts, suggesting a possible causal relationship between corticocortical synchrony and localized increases in beta power. We further show that high beta bursts are associated with significant changes in waveform shape and that beta-encoded phase-amplitude coupling is more evident during periods of high beta bursting. These findings reveal a deeper structure to the pathologic changes identified in the neurophysiology of Parkinson's disease, suggesting mechanisms by which the treatment may be enhanced using targeted network synchrony disruption approaches.
SIGNIFICANCE STATEMENT Understanding Parkinson's disease pathophysiology is crucial for optimizing symptom management. Present inconsistencies in the literature may be explained by temporal transients in neural signals driven by transient fluctuations in network synchrony. Synchrony may also act as a unifying phenomenon for the pathophysiological observations reported in Parkinson's disease. Here, simultaneous recordings from motor cortices show that increases in network beta synchrony anticipate episodes of beta bursting. We furthermore identify beta bursting as being associated with changes in waveform shape and increases in phase-amplitude coupling. Our results identify network synchrony as a driver of various pathophysiological observations reported in the literature and account for inconsistencies in the literature by virtue of the temporally variable nature of the phenomenon.
- deep brain stimulation
- electrocorticography
- motor cortex
- phase amplitude coupling
- synchrony
- waveform analysis
Introduction
Parkinson's disease (PD) has been strongly linked to pathologic beta oscillations (12–35 Hz) in cortico-subcortical circuits (Kühn et al., 2006; Bronte-Stewart et al., 2009; Eusebio and Brown, 2009; Eusebio et al., 2011; Whitmer et al., 2012; Brittain et al., 2014; Steiner et al., 2017). In particular, Parkinson's disease has been associated with increases in beta power in cortical (He et al., 2017) and subcortical (Steiner et al., 2017) structures, beta synchrony between (Cagnan et al., 2015) and within (Cagnan et al., 2015) subcortical structures, subthalamic nucleus (STN) beta bursting (Tinkhauser et al., 2017b), alterations in beta waveform shape (Cole et al., 2016) and increased beta-encoded phase-amplitude coupling (PAC) in the motor cortex (de Hemptinne et al., 2013).
Despite these findings, the literature contains numerous inconsistencies. For example, there are reports both of increased (He et al., 2017) and unaltered (Crowell et al., 2012) motor cortical beta power in Parkinson's disease when compared with normal controls. Similarly, while motor cortical PAC has been reported as increased in Parkinson's disease by several groups (de Hemptinne et al., 2013; Kondylis et al., 2016), it is not consistently increased on an individual basis (Miocinovic et al., 2015). Likewise, STN beta band power has been reported as elevated in as few as 50% of Parkinson's disease subjects in some cohorts (Steigerwald et al., 2015). Moreover, studies have reported both prokinetic (Feurra et al., 2011) and antikinetic (Eusebio et al., 2008; Pogosyan et al., 2009) effects of beta oscillations. Therapeutic administration of levadopa has likewise been observed to both decrease (Brown and Marsden, 1999) and increase (Sörös et al., 2017) cortical beta power.
These inconsistent findings may be accounted for, at least in part, by variability in the temporal domain of these signals, which has only recently been appreciated. Absolute beta power, measured over seconds to minutes, obscures details of the moment-to-moment variation in beta amplitude (beta bursting), which could be a critical aspect in the functional role of neural oscillations. Similarly, it has been shown that PAC is temporally dynamic (Malekmohammadi et al., 2015), suggesting that the observed absence of motor cortical PAC in single subjects may be related to temporal fluctuations in this signal.
A model to unify these disparate findings has yet to be developed and tested. Synchrony has been shown to play an important role in the cortical pathophysiology of Parkinson's disease (Kühn et al., 2009). Corticocortical beta coherence correlates with Parkinson's disease duration and severity and is reduced by levodopa in proportion to clinical improvements (Pollok et al., 2013). The role of synchrony in beta amplification has been described in more detail by Cagnan et al. (2015). Taking advantage of dual electrode placement in the STN and internal globus pallidus (GPi), the authors report that certain phase alignments between beta oscillations of the two nuclei favor transient increases in beta power, analogous to beta bursting. Degree of amplitude enhancement is proportional to the duration of time the two nuclei remained locked in a beta-enhancing phase. Furthermore, levodopa treatment reduced how often amplitude-enhancing phase relationships occurred (Cagnan et al., 2015).
Dynamic changes in beta power are a key regulator of movement in the basal ganglia, and evidence indicates that the temporal variability of the STN beta amplitude envelope correlates inversely with Parkinson's disease symptomatology (Little et al., 2012). Tinkhauser et al. (2017b) have shown that the duration of beta bursts in the STN is linked to clinical scores of Parkinson's disease severity. Their results also demonstrated a reduction in beta burst duration with levodopa therapy. Finally, and perhaps most intriguingly, the authors showed that interhemispheric measures of beta phase synchrony are increased during temporally overlapping beta bursts compared with bursts that did not overlap (Tinkhauser et al., 2017b). Subsequent studies have shown that targeting beta bursts in the STN using adaptive (“closed loop”) deep brain stimulation (DBS) can specifically curtail beta burst duration (Tinkhauser et al., 2017a) and has greater symptomatic benefit compared with continuous DBS (Little et al., 2013).
The foregoing studies by Tinkhauser et al. (2017b), and Cagnan et al. (2015) have shown that phase synchrony is associated with increased beta amplitude. Synchrony has also been linked to changes in waveform shape (Sherman et al., 2016), and waveform shape, in turn, has been closely associated with PAC (Cole et al., 2016). This article investigates the relationship between beta synchrony, amplitude, waveform shape, and PAC in an attempt to relate these variables under a unified explanatory framework. We acknowledge the argument that in signal analysis, amplitude and phase are intimately linked, providing a challenge to comparisons of phase measures across different amplitude ranges. We justify our analysis approach on the basis that theoretical studies suggest that phase estimation is minimally affected by amplitude decreases unless these are close to zero for prolonged periods (Cohen, 2014).
Materials and Methods
Surgical procedure.
The study was approved by the Institutional Review Board at the University of California, Los Angeles, and all participants provided written informed consent. All DBS implantations were carried by a single surgeon (N.P.), and all participants were implanted solely on the grounds of clinical necessity. Implantation targets were determined according to the standard operating protocol. We studied a cohort of 32 individuals undergoing DBS for essential tremor (ET; 10 patients; 3 male, 7 female; mean age, 58.9 years, SD, 15.2 years) or Parkinson's disease (22 patients; 15 male, 7 female; mean age, 64.4 years; SD, 7.72 years). We use subjects with essential tremor for controls as deep brain stimulation surgery provides access to similar cortical areas but these patients do not exhibit the same clinical symptoms of rigidity and bradykinesia as seen in patients with Parkinson's disease. Leads were implanted in bilateral GPi for the Parkinson's disease cohort and ventral intermediate nucleus of thalamus for the essential tremor cohort (Table 1, patient demographic details). Written informed consent was obtained for insertion of an eight contact subdural electrocorticography (ECoG) strip with 4 mm platinum-iridium electrodes and 1 cm spacing (AdTech Medical) through a burr hole fashioned for insertion of the DBS lead. Parkinson's medications were withdrawn at least 12 h before recording as per the standard clinical protocol.
Patient demographics
Experimental protocol.
Recordings were made at a sampling rate of 2400 Hz, bandpass filtered at 0.1 and 1000 Hz using two designated recording g.USBamp amplifiers (GTech), and stored using BCI2000 V3 software. A single scalp reference was used for all electrodes. All analysis was performed offline using MATLAB version 2016a (MathWorks) in conjunction with the Chronux toolbox (Bokil et al., 2010). Electrode locations were determined postoperatively using preoperative MRI and computed tomography (CT), intraoperative fluoroscopy, and postoperative CT using a technique described in depth previously (Randazzo et al., 2016; Malekmohammadi et al., 2018). Briefly, the cortical surface was reconstructed using Freesurfer (Dale et al., 1999), while OsiriX software was used to build a reconstruction of skull and stereotactic frame reference points (Ratib and Rosset, 2006). DBS leads and ECoG contacts were then mapped onto the model using the intraoperative fluoroscopy and postoperative CT in a purpose-built MATLAB graphic user interface (GUI; Fig. 1A). The first contact to be sited unambiguously anterior to the central sulcus was considered to lie over the motor cortex. Premotor cortex (PM) was determined to be the contact immediately anterior to primary motor cortex (M1; Fig. 1B). A map of all contacts used for analysis has been projected onto a 3D reconstruction of a standardized cortical surface (Fig. 1E,F).
A, B, Using postoperative CT and intraoperative fluoroscopy (A), the locations of the ECoG contacts (M1, green dot; PM, magenta dot; unused contacts, blue dots) are mapped onto a reconstructed cortical surface (B). The first contact immediately anterior to the central sulcus (CS; red line) was taken to be motor cortex (green dot), and the contact anterior to that was taken to be PM (magenta dot). Blue dots indicate unused contacts on the ECoG strip. C, Raw trace from one subject (blue top trace), with the high beta filtered trace below (blue bottom trace). Episodes of beta bursting when the beta analytic envelope exceeds the 75th percentile (black dotted line) are highlighted in magenta on both traces. Burst duration shown in red. D, Waveform analysis. Top black trace shows the filtered signal with green and magenta dots indicating peaks and troughs, respectively. The peaks and troughs are assigned to corresponding locations on the raw signal in the blue trace below. The inset shows a magnified segment of the trace with brackets at the peaks and troughs showing the points used for calculation of peak and trough sharpness measures. E, F, Lateral (E) and superior (F) views of electrode positions across the cohort (M1, blue dots; PM, red dots) projected onto a standard Montreal Neurologic Institute (MNI) 3D cortical reconstruction. Electrodes not used in the analysis are faded gray dots. PM electrodes were 1 cm anterior to the electrode labeled M1.
Patients were woken from anesthesia after the formation of burrholes in locations determined according to preoperative clinical stereotactic planning. Bilateral DBS leads were sited, and patients were allowed to recover from the propofol anesthesia for a minimum of 30 min before clinical testing of stimulation parameters was conducted with a clinical neurologist and a technical team. Once satisfactory lead placement was confirmed both physiologically and clinically, recordings were conducted in a quiet setting. The subjects were instructed to relax for epochs of 30 s while keeping their eyes open and refraining from moving their body. A minimum of five and a maximum of six periods of 30 s rest recording (150 s total) were collected and analyzed for each subject.
Data analysis.
Data were analyzed in a monopolar configuration (without rereferencing) to preserve the phase of oscillations spatially proximal to each source. Subsequent analyses using the debiased weighted phase lag index (dwPLI) accounted for zero phase difference oscillations that may represent global signals rather than spatially separated yet synchronized oscillations, the latter of which is the phenomenon of interest.
Preprocessing consisted of 60 Hz line noise removal using sine wave fitting over a moving window (rmlinesmovingwinc function; window, 3 s; step, 2 s; Chronux Toolbox). Removal of artifacts used previously described techniques (AuYong et al., 2018), including using a custom GUI inspection of time–frequency spectra and manual removal of artifact. Automatic artifact removal used full-wave rectification followed by calculation of the first derivative and filtering using a five-point median filter. An artifact was defined as a region in which the first derivative exceeded 5 SDs for the subject calculated across the entirety of their recording. Removed segments of artifact were replaced by linear interpolation derived from data 2 ms before and after the artifact. Removed artifact constituted episodes of high-amplitude activity, excessive harmonics, or time series with high rates of voltage change.
Power was calculated using Welch's estimate (Welch, 1967) using a 1 s window with zero overlap. Power spectra between 4 and 400 Hz were normalized relative to the total power between 4 and 400 Hz. This accounted for variation in the total power across individuals and provided for a more standardized comparison across individuals.
Periods of beta bursting were determined by filtering each subject's data using a zero-phase FIR (finite impulse response) bandpass filter [cutoff frequencies: low beta, 12–20 Hz; high beta, 21–35 Hz; eeglab function eegfilt; order = 3*fix(Fs/locutoff)]. An analytic envelope was extracted using a Hilbert transform, and bursting was defined as any period of oscillation amplitude that crossed the 75th percentile power for that subject within a particular band and persisted for at least 50 ms (Fig. 1C; Tinkhauser et al., 2017a,b). After filtering, the phase and amplitude segments of the signals were sorted into burst and nonburst segments. To construct burst duration across a spectrum of frequencies (see Fig. 2B), we filtered each signal as described above in bands of 5 Hz width centered on frequencies 4–40 Hz in 1 Hz increments. For each band, bursts were extracted and the mean burst duration was calculated. The mean burst duration was calculated for high beta between 21 and 35 Hz by taking the mean burst duration across this frequency range for each subject. For analysis of burst duration characteristics only, all episodes of amplitude crossing the 75th percentile threshold were included. We identified periods of beta bursting in the raw signal by indexing back to the unfiltered signal in subsequent analyses.
Statistical comparisons of burst and preburst synchrony assessed with PSI and dwPLI
Synchrony was assessed using two related techniques, phase synchrony index (PSI) and dwPLI. Custom scripts were developed for the calculation of each synchrony measure and calculated using methods detailed previously (Vinck et al., 2011; Cohen, 2014). Both are purely phase-based measures of connectivity and thus relatively resistant to any confounding effects of amplitude that may be peculiar to episodes of bursting. These two scores are complimentary as clustering at zero or π phase angles is not recognized as increased synchrony using dwPLI alone, and a consistent phase lag or phase lead (positive or negative phase difference relationship) is not revealed by using PSI alone. To ascertain the temporal relationship between increases in synchrony and increases in power, synchrony (PSI and dwPLI) at burst onset was calculated for each subject and plotted alongside normalized power at burst onset. Maximum values of the first derivative for each subject's curve and the time of peak first derivative relative to burst onset were calculated for each subject. To assure that our measures of synchrony increase were not solely due to the threshold adopted (point of maximum derivative), we used an alternative measure as verification. For each subject, the mean across the preburst and burst period (−50 and +50 ms with respect to burst onset at 0 ms) was taken. The time point at which this mean was met or exceeded constituted our alternative threshold.
To test whether our increased synchrony was merely a result of increasing amplitude during beta bursts, we generated synthetic signals with no inherent synchrony between signals. This approach allowed us to explore how low amplitudes and transient increases in amplitude (as in beta bursting) might compromise synchrony measurements. To perform this in silico analysis, we generated 400 synthetic signals of 10 s duration using sine waves with frequencies in the 1–20, 21–35, 36–50, and 51–150 Hz ranges. Each signal was composed of the 20 sine waves (5 from each of the four frequency bands) randomly selected from the frequency bands. Each sine wave was modified with noise [a vector of randomly generated values on the open interval (0, 0.05)], a randomized phase shift, and a direct current offset of −0.5 a.u. All other signal parameters and analysis (sampling rate, filter parameters, and PSI calculations) were identical to those used our in vivo analysis. Testing was performed on 200 pairs of signals generated in the above manner. We elected to modulate only one of our signal amplitudes as our in vivo analysis was concerned with episodes of bursting only at the motor cortex (regardless of amplitudes at the “premotor” cortex). Signal 2 was thus modulated by taking the dot product of the raw synthetically generated signal with a Gaussian kernel of equal length having 1 at the zenith, a mean of 5 s, and an SD of 0.5 s. This produced a signal that varied in amplitude over the course of 10 s from near zero amplitude to equal amplitude with signal 1, before decreasing again to near zero amplitude at the end of the 10 s (see Fig. 4D, blue trace). Maximal PSI levels for each iteration were defined as the mean PSI within ±0.5 s (±1 SD) of the maximal amplitude at the 5 s point in the signal. The time point of synchrony increase was defined as the point at which PSI levels between signal 1 and signal 2 reached or exceeded the maximal amplitude threshold, and the time point of synchrony decrease was defined as the time at which synchrony dropped below this level.
Waveform analysis is an evolving field and attempts to characterize the nonsinusoidal features of neural oscillations that have been made previously (Cole et al., 2017). Here we adopted the procedure used by Cole et al. (2017) and focused on waveform changes in the high beta band because it was in this band that we have identified significant changes in burst duration specific to Parkinson's disease. First, we used a zero phase FIR to bandpass filter our signal (high beta: 20–35 Hz; order = 360; eeglab function eegfilt) and identified the zero crossings of the first derivative corresponding to peaks and troughs. In the raw signal, peaks corresponded to the extrema of the signal between two troughs identified in the filtered signal. Likewise, troughs were the extrema of the raw signal between two peaks in the filtered signal. To constitute a valid peak or trough, the extrema had to be accompanied by a net positive gradient 2 ms before the peak and a net negative gradient 2 ms after the peak. Likewise, troughs were defined only in regions of the raw trace where a negative gradient preceded the trough for 2 ms and a positive gradient followed the trough for 2 ms. Figure 1D shows a representative sample of the waveform analysis illustrating where peaks and troughs fall on the raw signal with the high beta and high beta analytic envelope superimposed.
PAC was calculated according to methods described in detail previously (Malekmohammadi et al., 2015). In brief, we used the modulation index (MI) of Tort et al. (2010) and parameterized with reference to previous PAC analyses (Malekmohammadi et al., 2015). Data were filtered using FIR filtering (high pass then low pass) at frequencies of 8–38 Hz in 2 Hz steps with a 2 Hz bandwidth. For optimal detection of PAC relationships, we used variable bandwidths for amplitude encoding set to be twice the center frequency of the phase-encoding band (Aru et al., 2015). Amplitude encoding bands were centered on frequencies between 40 and 200 Hz in 5 Hz steps. We extracted phases and amplitudes using Hilbert transforms on our phase and amplitude encoding signals respectively. Twelve bins of 30° width were used to categorize phase and derive phase-amplitude histograms.
Statistics.
All statistical testing was conducted in MATLAB using the Statistics and Machine Learning Toolbox. Power differences and mean burst durations were compared between the essential tremor and Parkinson's disease cohorts using a subject-based average over the high beta band followed by Student's t test. To compare the distributions of beta bursting across cohorts, we compiled a probability distribution for each subject according to the proportion of the whole signal that fell into duration bins of 50 ms width. Duration bins were from 50 to 750 ms in 50 ms steps, the final bin was defined as >750-∞ ms. A log-logistic distribution was fitted to the resulting distribution for each subject, and a μ value (mean duration of the fitted curve) was calculated. This distribution was selected as our burst duration data are continuous and non-negative, and exhibit a strong positive skew. Group duration distributions were compared using a two-sample Kolmogorov–Smirnov (K–S) test, and μ values were then compared across groups using a two-sample t test. Correlations with clinical measurements were calculated using a Spearman's ρ test of significance. Correlational analyses were conducted across the whole cohort, and where correlations against more than one clinical score aggregation method were conducted, the α value for significance was adjusted for multiple comparisons using Bonferroni correction (Fig. 2D–F). Unified Parkinson Disease Rating Scale (UPDRS) scores for the ET cohort were taken to be zero. To compare synchrony across cohorts and conditions, we performed a two-way ANOVA on mean synchrony outcomes (dwPLI/PSI) for each subject. Violin plots were delineated using 100 query points between the maximum and minimum of the dataset based on an optimized Gaussian kernel (MATLAB function boxplot; https://uk.mathworks.com/matlabcentral/fileexchange/46545-alternative-box-plot). Circles on the plots indicate points lying outside Q1 −1.5*IQR or Q3 + 1.5*IQR (where Q1, Q3, and IQR are the first quartiles, third quartiles and the interquartile range, respectively).
For clinical scores, we used the UPDRS, part III (UPDRS-III) subscore item 3.3 to obtain a rigidity score for the upper limb contralateral to the ECoG strip. Measures of Parkinson's disease symptomatology were correlated with mean burst duration and mean burst PSI for each subject using a nonparametric Spearman's ρ test.
To compare results from the waveform analysis, sharpness and steepness ratios were derived for each subject in each condition by dividing the mean peak sharpness by the mean trough sharpness for the former and the mean rising steepness by the mean falling steepness for the latter. Differences between conditions and groups were assessed using a multivariate ANOVA (MANOVA) to determine whether group means differed. To compare changes across cohorts and conditions with respect to waveform sharpness and steepness ratios, we generated a vector in the 2D plane of the sharpness/steepness ratio space that joined each individual's ratio from nonbursting epochs to the corresponding ratio during bursts. The lengths and angles of vectors could be compared using a t test for vector lengths and a parametric Watson–Williams multisample test for equal means from the CircStat Toolbox (Berens, 2009). This is in effect a one-way ANOVA adapted for circular data. To assess the statistical significance of PAC, the MI was calculated for the original data and then, to create a null distribution, was also calculated for 500 surrogate datasets composed of phase-shuffled original data. MI scores for original data and surrogate datasets were transformed into z scores, and an MI for the original dataset was considered significant if it lay 1.96 SDs above the mean of the surrogate dataset, corresponding to a p ≤ 0.025 significance level. Only MI scores exceeding this significance level were included in further analyses.
PAC results were compared in the low and high beta bands by calculating a mean MI for each subject over frequency ranges of 12–20 Hz (low beta) and 20–35 Hz (high beta) phase-encoding frequencies with 50–200 Hz amplitude encoding frequencies for both bands, as used in previous studies (de Hemptinne et al., 2013). Mean values of MI were compared using a Mann–Whitney U test for differences between cohorts and to compare bursting and nonbursting conditions within each cohort.
Results
Cohort demographic comparison
The ages of the ET cohort were not normally distributed by Anderson–Darling (AD) testing (AD test statistic = 0.81, p = 0.022). Accordingly, a nonparametric Mann–Whitney U test was used to compare the ages of the PD and ET cohorts. There were no significant differences in the ages of the two cohorts (Mann–Whitney U test z value = 0.77, p = 0.44). Table 1 details patient demographics. Our PD cohort included four tremor-dominant individuals (Table 1, subjects 13, 15, 20, and 30) and 18 akinetic-rigid individuals as defined by scoring systems published previously (Lian et al., 2019).
High beta power and bursts in Parkinson's disease
Analysis of power spectra demonstrated significantly greater high beta power at rest in the Parkinson's disease cohort compared with the essential tremor cohort (mean essential tremor = −20.99; 95% CI, ±1.53; mean Parkinson's disease = −17.45; 95% CI, ±0.67; two-tailed t test: t = −2.33, p = 7.6 × 10−6; Fig. 2A), and this increased high beta power was mirrored by an increased burst duration for Parkinson's disease compared with essential tremor in the high beta band (Fig. 2B). Anderson–Darling tests of normality found that the burst duration data for the Parkinson's disease cohort was not normally distributed (test statistic = 0.87, p = 0.021), and so nonparametric testing was used to compare the cohorts (mean high beta burst duration essential tremor = 143.24 ± 8.20 ms; mean high beta burst duration for Parkinson's disease = 156.34 ±7.40 ms; Mann–Whitney U test z value = −2.38, p = 0.017). To further characterize differences in high beta burst duration between the cohorts, we calculated the distributions of burst durations. A preliminary two-sample K–S test showed that the burst duration distributions differed significantly between the two groups (two-sample K–S test, K–S value = 0.18, p = 0.001; Fig. 2C). We then modeled high beta burst durations using a log-logistic distribution for each subject to characterize further the differences between the cohorts. This distribution was selected as it is suitable for continuous, non-negative data with positive skew. Burst analysis closely parallels common applications of this distribution in survival and time-to-failure analyses. In this distribution, fitted μ values (normally the logarithm of the distribution mean) have been transformed by taking the exponential so that they represent the mean of the fitted distribution in milliseconds. Sigma values (σ is the scale parameter of the log-logistic distribution) represent the concentration of burst durations either closer to zero (low σ values) or further from zero (high σ values). Fitting of the log-logistic distribution for the population yielded a log likelihood of −75,376.8 (estimated μ = 5.21, SE = 0.007; estimated σ = 0.457, SE = 0.0035). Both the μ and σ values differed for essential tremor and Parkinson's disease (mean μ essential tremor = 88.5 ms, 95% CI, ±2.03 ms; mean μ value Parkinson's disease = 108.4, 95%, CI ± 2.1; Mann–Whitney U z value = −2.70, p = 0.007; Fig. 2C, inset; the mean σ for essential tremor = 0.42, 95% CI, ±0.017; mean σ for Parkinson's disease = 0.47, 95% CI, ±0.018; AD test statistic = 1.26, p = 0.002, Mann–Whitney U test z value = −3.27, p = 0.001), suggesting that the Parkinson's disease cohort had a tendency to stay in high beta bursts for greater durations compared with the essential tremor cohort (increased μ) and that increased mean burst duration was due to a general shift toward increased frequencies of higher duration bursts rather than occasionally more extreme values (increased σ).
A, Normalized power spectrum comparisons show that Parkinson's disease has a higher resting power in the high beta band compared with essential tremor. Inset, High beta band demarked by black broken lines and statistical comparisons by cohort of average high beta powers. ET cohort shown in black and PD cohort shown in red. B, The same region of the power spectrum (high beta) also demonstrates higher average beta burst duration. Inset, Statistical comparison of cohort mean high beta burst durations. Color conventions as in A. C, For each subject, the distribution of high beta burst durations can be represented as a log-logistic distribution, where mean cohort distributions are plotted in black (essential tremor) and red (Parkinson's disease). Fitting a μ value (mean duration) to the distribution allows us to compare the duration of high beta bursting directly across cohorts. Asterisks indicate the results of the K-S test implying that the distributions differ significantly. Parkinson's disease shows a more rightward skew in their distribution when compared with essential tremor (higher σ values, main plot). The inset demonstrates that the parameterized distributions for each cohort (μ values) also differ between the two groups with the Parkinson's disease cohort having a higher mean μ value compared with essential tremor. Error bars represent 95% confidence intervals. D–F, Correlations between μ and three UPDRS-III scores; a total of three different μ/UPDRS-III comparisons were corrected for using the Bonferroni method. D, Mean high beta burst duration (μ) correlates with severity of Parkinson's disease as quantified by the aggregate UPDRS-III contralateral upper limb score (the sum of rigidity and bradykinesia scores for the contralateral upper limb only) and survived Bonferroni correction. Red line shows robust least squares fit. E, μ showed a positive correlation with aggregate contralateral hemibody (arm and leg combined) bradykinesia/rigidity scores (p = 0.016, r = 0.42), which survived Bonferroni correction. F, Although there was a positive correlation between μ and rigidity score for the contralateral upper limb, it did not survive Bonferroni correction (p = 0.018, r = 0.42).
Finally, we compared the μ values for each subject in the cohort with the contralateral upper extremity aggregate bradykinesia and rigidity score on the UPDRS-III. This showed a positive relationship between ratings of Parkinson's disease symptomatology and burst duration (μ) at rest (Bonferroni-corrected p = 0.029, Spearman's r = 0.45; Fig. 2D). Here we highlight the relationship between μ and aggregated contralateral upper limb rigidity/bradykinesia because our ECoG placement targeted the more medially placed arm area of M1. In total, we compared μ with three different measures of hemibody UPDRS-III scores (Fig. 2D–F), all of which showed a trend toward a positive correlation and two of which survived Bonferroni correction (Fig. 2D,E). Comparison of high beta burst duration (μ) with aggregated hemibody (upper and lower limb) rigidity/bradykinesia scores revealed a significant positive correlation (Bonferroni corrected, p = 0.048; Spearman's r = 0.42; Fig. 2E). Our μ scores showed a trend toward positive correlation with rigidity alone (UPDRS-III 3.3); however, this third correlation did not survive Bonferroni correction (Bonferroni corrected, p = 0.054; Spearman's r = 0.42; Fig. 2F).
Beta bursts and corticocortical synchrony
It has been previously established that synchrony between M1 and the STN is increased during beta burst activity (Tinkhauser et al., 2018). Here we show corticocortical coupling is likewise increased during episodes of beta bursting. We also demonstrate that synchrony during bursting does not differ between essential tremor and Parkinson's disease despite there being a higher level of synchrony before bursting in the Parkinson's disease group. Finally, we demonstrate that M1 -PM synchrony increases precede M1 beta power increases and that the average synchrony during beta bursting correlates significantly with clinical features of Parkinson's disease.
All PSI values met the Anderson–Darling test criteria for normality (all p values > 0.09), and so comparisons were performed using t tests. Before beta bursting, the Parkinson's disease cohort demonstrates greater M1-premotor PSI compared with the essential tremor cohort (mean preburst PSI for essential tremor = 0.7; 95% CI, ±0.027; mean preburst PSI for Parkinson's disease = 0.76; 95% CI, ±0.027; two-tailed t test, t = −3.21; Bonferroni corrected, p = 0.0125; Fig. 3A). In both groups, PSI increased during bursts (mean preburst PSI essential tremor = 0.69; 95% CI, ±0.03; mean burst PSI essential tremor = 0.86; 95% CI, ±0.038; two-tailed t test, t = −18.6; Bonferroni corrected p = 6.87 × 10−8; mean preburst PSI Parkinson's disease = 0.76, 95% CI ± 0.03, mean burst PSI, Parkinson's disease = 0.91; 95% CI, ±0.02; two-tailed t test, t = −22.3; Bonferroni corrected p = 1.8 × 10−15), but, most importantly, PSI did not differ between the groups during bursts (mean essential tremor burst PSI = 0.86; 95% CI, ±0.04; mean Parkinson's disease burst PSI = 0.91; 95% CI, ±0.02; two-tailed t test, t = −2.38; Bonferroni corrected p = 0.097). The pattern of results for dwPLI closely paralleled those for PSI; however, dwPLI did not meet the test of normality and so was analyzed using nonparametric Mann–Whitney tests (Table 2, Fig. 3B).
A, PSI is higher during high beta bursting than immediately before episodes of high beta bursting. Immediately before high beta bursts the Parkinson's disease cohort (red) have a significantly higher PSI when compared with the essential tremor cohort (gray). Both cohorts have significantly higher PSI during high beta bursts compared with immediately before high beta bursts; however, there is no difference between cohorts in the average PSI during bursting. B, An identical pattern of results was identified for dwPLI analysis (see also Table 2). C, Average high beta envelope (blue) and PSI (red) aligned to burst onset (black broken line at 0 s). Ticks show the point of maximum derivative calculated for each subject for each measure. Shading indicates 95% confidence intervals. Broken lines in color indicate cohort mean times of amplitude (blue) and PSI (red) increases. D, The same plot as C but for dwPLI (magenta line) showing the temporal relationship between the maximum derivatives of dwPLI and amplitude (blue line). Conventions are as in C. E, A violin plot comparing average times of maximum derivatives for amplitude (blue violin), PSI (red violin), and dwPLI (magenta violin). Both increases in synchrony measures significantly precede increases in amplitude. F, Scatter plot showing that the mean bursting PSI correlates with symptomatic severity as assessed by UPDRS-III contralateral upper limb bradykinesia/rigidity scores. *p ≤ 0.05, ***p ≤ 0.001. N.S., Not Significant.
It has been shown that phase synchrony can precede power increases by several milliseconds (Womelsdorf et al., 2007). To attempt to discern whether synchrony increases were preemptory of power increases, we evaluated the maximum gradient of power and synchrony increases relative to burst onset, see Figure 3, C and D. This analysis was performed on the whole cohort together to demonstrate a general point regarding synchrony and amplitude. Time distributions of maximal amplitude did not meet the Anderson–Darling test of normality, and so comparisons were made using nonparametric tests. We found that the time point of the maximal rate of increase in preburst synchrony (both PSI and dwPLI) preceded the time point of the maximum rate of power increase (both measured relative to bursting onset, denoted as time 0). Of the three measures compared (PSI, dwPLI, and power), PSI showed the earliest increase and occurred significantly earlier than power increases (mean time of greatest power increase = −2.68 ms; 95% CI, ±1.61 ms; mean time of greatest PSI increase = −20.52 ms; 95% CI, ±1.61 ms; Mann–Whitney U test z value = 6.57; p = 4.96 × 10−11; Fig. 3E), with dwPLI increasing maximally shortly after PSI and also showing significant anticipation of power increases (mean time of greatest power increase = −2.68 ms; 95% CI, ±1.61 ms; mean time of greatest dwPLI increase = −13.26 ms; 95% CI, ±6.06 ms; Mann–Whitney U test z value = 3.70; p = 2.20 × 10−4). Thus, for PSI which shows the earliest rise in synchrony, there is a temporal lead with reference to the time of maximal power increase. This temporal lead can be calculated as −20.52 ms − (−2.68 ms) = −17.84 ms, or a PSI increase 17.84 ms in advance of the power increase. In light of this temporal precedence, we suggest that these findings could reflect a causal relationship with inter-regional (PM–M1) synchrony increases causing increases in locally measured oscillatory power (M1).
Although synchrony measures during bursts did not significantly differ between the two groups, the difference between essential tremor and Parkinson's disease groups approached significance, and there was still intersubject variability. We therefore explored the relationship between the magnitude of synchrony and clinical symptoms. Average PSI during bursts was positively correlated with contralateral upper limb bradykinesia/rigidity scores across all subjects (Spearman's r = 0.40, p = 0.023; Fig. 3F).
Two approaches were used to assess the robustness of our findings. First, we developed a different method for assessing the timing of increase in synchrony to demonstrate that the results are independent of the criterion used for determining synchrony increases. Second, we sought to show that observed increases in synchrony were not merely necessitated by an increase in amplitude leading to improved accuracy of phase position and, hence, to increased synchrony.
The point of the maximal synchrony increase could be a suboptimal means for determining whether synchrony increases precede amplitude bursts. We therefore repeated our analysis using the time point at which synchrony increases beyond the mean taken across the preburst and burst periods. This approach yielded results identical to our initial analysis, and full details are depicted in Figure 4A–C. PSI increased before amplitude (mean time of amplitude increase to threshold = −0.77 ms; 95% CI, ±2.86 ms; mean time of PSI increase to threshold = −8.57 ms; 95% CI, ±1.50 ms; Mann–Whitney U test z value = 6.44, p = 1.13 × 10−10). Similarly, for dwPLI the time point at which the threshold was exceeded preceded the time point at which the threshold was crossed for the normalized amplitude (mean time of amplitude increase to threshold = −0.77 ms; 95% CI, ±2.86 ms; mean time of dwPLI increase to threshold = −9.28 ms; 95% CI, ±4.02 ms; Mann–Whitney U test z value = 5.13, p = 2.88 × 10−7).
A, Using a distinct criterion for PSI increase produces results identical to those in prior analyses. To check the validity of our burst synchrony analysis (Fig. 3C–E), we repeated the analysis replacing the criterion for a synchrony increase with an alternative measure. The time point of an amplitude or synchrony increase was determined to have occurred when the synchrony or amplitude measurement met or exceeded the mean taken across the entire preburst/burst episode. All graph conventions are identical to Fig. 3C–E. This analysis replicated previous findings to show that the measured increases in synchrony (Fig. 3C, PSI, red line; Fig. 3D, dwPLI, magenta line) preceded the maximal rate of increase in the amplitude (Fig. 3C,D, blue line; A, B). Ticks on the graphs indicate the time points of amplitude threshold crossings (Fig. 3C,D, blue ticks on amplitude traces; A, B) and time points of synchrony threshold crossings (respectively, Fig. 3C,D, red ticks for PSI, magenta ticks for dwPLI; A, B, red ticks for PSI, magenta ticks for dwPLI). C, A comparison of threshold-crossing time points across amplitude and synchrony measures. Amplitudes (blue trace) superseded the mean close to the burst onset (0 ms), whereas synchrony measures increased to supersede thresholds several milliseconds before burst onset (red and magenta violins). D, Varying the amplitude of one signal has little effect on the measurement of PSI. We generated two synthetic signals, signal 1 (black trace) and signal 2 (blue trace), with no inherent synchrony and modulated signal 2 by taking the dot product of the raw signal with a Gaussian kernel having 1 at the zenith, a mean of 5 s and an SD of 0.5 s. Our results show that PSI increased and decreased steeply as the amplitude-modulated signal 2 reached 0.4% and 0.39% of signal maximum for PSI (D, red trace) at onset and offset, respectively (D, black dotted lines show the mean time point at which PSI reached maximal levels). Blue ticks on the signal 2 trace indicate the PSI onset and offset time points at which individual iterations of the experiment reached threshold synchrony levels and dropped below threshold levels of synchrony. Threshold PSI levels for each iteration were defined as the mean PSI within ±0.5 s (±1 SD) of the maximal amplitude at the 5 s point in the signal. ***p < = 0.001.
To show that our synchrony findings were not necessitated by amplitude increases during bursting, we generated 400 synthetic signals composed of random frequencies plus noise and gradually modulated the amplitude of one signal while measuring PSI between pairs of signals (Fig. 4D). This analysis shows that maximal PSI onset is achieved when the modulated signal reaches 0.4% of maximum amplitude (PSI “onset” point; 95% CI, ±0.029%) and synchrony measurement again becomes compromised when the signal reaches 0.39% of maximum amplitude (PSI “offset” point; 95% CI, ±0.027%). These points of the signal corresponded to signal envelope amplitudes of 0.0017 a.u. (mean envelope amplitude for PSI onset point) and 0.0016 a.u. (mean envelope amplitude for PSI offset point). Across the entirety of our in vivo recordings, the minimum amplitude of the high beta envelope was 0.0058 μV, which is more than three times the level at which we found PSI was compromised in our in silico experiments. This result suggests that low amplitudes alone cannot account for the burst-related changes in synchrony we have observed. Furthermore, our burst analyses are conducted using data surrounding the 75th percentile of signal amplitude meaning that our analyses were unlikely to have included epochs of extremely low amplitude. This offers some assurance that our results are not necessitated by the a priori relationships between phase estimation and signal amplitude.
High beta synchrony and waveform shape
Plotting peak/trough sharpness and rising/falling steepness ratios across conditions is one means by which to characterize changes in overall waveform shape (Fig. 5A). This ratio measure has the advantage of negating the potentially confounding effects of power differences between bursting and nonbursting periods as both peak and trough sharpness would be expected to be similarly increased by increases in amplitude excursion during bursting episodes. The same argument holds for the ratio measurement of rising and falling steepness. There was a significant shift from a relatively steeper rising phase and relatively sharper peak during nonburst periods (Fig. 5A, black “X”) toward a relatively steeper falling phase and relatively sharper trough during bursting (Fig. 5A, red “X; MANOVA, d = 1, p = 1.2 × 10−6). Nonburst to burst waveform changes can be decomposed into vectors with length and direction to quantify the nature of waveform changes (Fig. 5B,C). Figure 5C illustrates that there was a significant difference between the essential tremor and Parkinson's disease cohorts in the extent (as measured by Euclidean distance on the two-dimensional waveform shape plane; Fig. 5C, thick black and red vectors), but not the direction (angle) of this burst-dependent waveform change (ET mean Euclidean distance = 0.30; 95% CI, ±0.099; mean PD Euclidean distance = 0.15; 95% CI, ±0.033; two-tailed Student's t test, t value = 3.88, p = 5.24 × 10−4). This is reflected in the larger vector size for ET (Fig. 5C, solid black vector) compared with PD (Fig. 5C, solid red vector). Our postulation that changing waveform shape reflects changes in synchrony among the underlying neural elements would predict that a difference in waveform changes could be linked to the symptomatology of Parkinson's disease. In support of this prediction, we find that the Euclidean distance of waveform shape is inversely correlated with clinical measures of Parkinson's disease severity (Spearman's r = −0.40, p = 0.022; Fig. 5D). Across the cohort a smaller Euclidean distance (change in waveform shape) between nonbursting and bursting waveforms was associated with higher scores of Parkinson's disease clinical symptomatology. Despite differences in the extent of waveform change and a correlation between UPDRS-III scores and the magnitude of waveform change, we found no significant difference between the cohorts when comparing waveform shape loci during both bursting and nonbursting conditions (Fig. 5E). An attempt to replicate previous findings correlating sharpness ratios with contralateral rigidity scores (Cole et al., 2016) did not reveal a significant relationship in our cohort (Spearman's r = −0.12, p = 0.52; Fig. 5F).
A, waveform ratio plots show that during nonbursting periods (black markers) the waveform is predominantly steeper on the rising phase and has a sharper peak, as evidenced by the concentration at the origin. During bursting (red markers), the waveforms predominantly occupy the lower left quadrant showing a steeper falling phase and a sharper trough. B, Each subject had their waveform shift quantified as a directional vector to facilitate the comparison of waveform changes in terms of magnitude and direction. C, There was no significant difference between essential tremor and Parkinson's disease in the direction (vector angle) of waveform shift, but there was a significant difference in the vector-based magnitude (vector length) of their waveform change as seen by comparison of solid black (ET) and solid black (PD) vector lengths. D, Magnitude of waveform change showed a significant negative correlation with Parkinson's disease symptomatology. E, There was no significant difference in waveform shape between cohorts in either nonburst (black markers) or burst conditions (red markers). In this diagram, the waveform centroid of the essential tremor cohort during the nonburst condition is marked with a black X, and the waveform centroid of the Parkinson's disease cohort is marked with an O. The same markers in red represent the waveform centroids of the cohorts during bursting conditions. F, We found no significant correlation between peak sharpness across the whole recording and contralateral upper limb rigidity scores contrary to findings in previous articles (Cole et al., 2017).
PAC is specific to beta bursting
PAC was calculated for bursting and nonburst periods separately. Anderson–Darling tests of normality showed that the data for PD cohort PAC was not normally distributed (low beta PD nonburst MI AD test statistic = 2.32, p = 5.00 × 10−4; low beta PD burst MI AD test statistic = 2.15, p = 5.00 × 10−4), and comparison of the data was therefore made using nonparametric testing. Figure 6A shows the group average PAC for nonbursting and bursting segments across all subjects. There were significant burst versus nonburst differences in both the low and high beta band-encoded PAC (Fig. 6A,B). Low beta phase-encoded PAC increased significantly during bursts compared with nonbursting (mean low beta nonburst MI = 1.20 × 10−4; 95% CI, ±6.15 × 10−3; mean low beta burst MI = 0.0013; 95% CI, ±5.43 × 10−5; Mann–Whitney U test, z value = −4.78, p = 1.70 × 10−3; Fig. 6A,B). Likewise, high beta phase-encoded PAC showed similar results with a pronounced elevation in PAC during bursts compared with nonburst periods (mean high beta nonburst MI = 6.96 × 10−5; 95% CI, ±4.14 × 10−5; mean high beta bursting MI = 6.96 × 10−4; 95% CI, ±2.88 × 10−4; Mann–Whitney U z value = −4.76, p = 3.93 × 10−6; Fig. 6A,B). There was no significant difference for MI in the essential tremor and Parkinson's disease cohorts during nonburst or burst periods in the high or low beta phase-encoded bands.
A, Across the whole cohort, PAC during nonburst periods is minimal compared with PAC observed during bursts. White and magenta broken lines indicate the regions of PAC used to calculate the statistical comparisons shown in B. White broken lines delineate low beta PAC area and magenta broken lines delineate high beta PAC area. B, For both low beta and high beta-encoded PAC, there was a significant increase in MI during bursting compared with nonbursting epochs. C, Average bursting high beta phase-encoded PAC correlated with the magnitude of waveform change during high beta bursting. ***p < = 0.001.
To investigate the relationship between cortical PAC and waveform shape, we compared cortical high beta PAC with magnitude of changes in high beta waveform shape. High beta PAC data did not meet the test of normality (Anderson–Darling test statistic = 2.93, p = 5 × 10−4), and so a Spearman's nonparametric test was used to assess correlation. Our data show a significant positive correlation between the extent of waveform shape change (Euclidean distance on the two-dimensional waveform shape plane) and average high beta phase-encoded PAC (Spearman's r value = 0.52, p = 0.0026; Fig. 6C).
Discussion
Our results demonstrate a temporally dynamic pattern of beta bursting that links dynamic network synchrony with locally observed phenomena of increased high beta power and high beta phase-encoded PAC. Here we have provided evidence linking high beta corticocortical synchrony, motor cortex high beta amplitude increases, high beta waveform changes, and high beta-encoded PAC. Previous studies have shown that in Parkinson's disease, the STN beta burst durations and amplitude are modulated by treatment with levodopa (Tinkhauser et al., 2017b). Decreases in burst duration were also correlated with symptomatic improvement (Tinkhauser et al., 2017b). The same study showed that longer beta bursts are associated with increases in local and interhemispheric phase synchronization. Our results demonstrate that cortical synchrony plays a similar role in cortical beta bursting; our measures of synchrony were elevated across the M1 and PM cortices during beta bursting. Importantly, we demonstrate that bursts are preceded by increases in synchrony. This is consistent with the prior report by Cagnan et al. (2015) of amplitude-enhancing phase relationships between the globus pallidus and STN. However, we now demonstrate that the phenomenon of phase synchronization resulting in amplitude enhancement (or bursting) occurs at other nodes in the motor network and is not specific to the subcortical basal ganglia. Notably, we also show that higher corticocortical synchronization during bursts and longer duration bursts both correlate with aggregate measures of contralateral Parkinson's disease symptomatology. Our measures of synchrony, burst duration, and waveform shape change all correlated with contralateral upper limb rigidity alone, but the correlation did not survive Bonferroni correction. The increased strength of relationship among synchrony, burst duration, waveform change, and aggregate measures of rigidity/bradykinesia could indicate that these processes play a generalized role in underpinning the clinical features of Parkinson's disease. These results together suggest that hypersynchrony across the corticobasal ganglia circuit could be driving increases in beta amplitude, waveform shape changes, and motor cortical PAC. Present forms of adaptive DBS stimulate at high frequency only when the amplitude of beta rises above a preset threshold, thus effectively shortening burst duration and decreasing burst amplitude (Tinkhauser et al., 2017a). That we observe >17 ms of increasing PSI preceding cortical high beta burst spotlights an anticipatory biomarker that could be used for adaptive DBS.
Increased beta power, previously reported in the motor cortex, may be attributable to the longer duration of beta bursts, which we report here. Indeed, the mean duration of bursts in this spectral region correlates with the severity of contralateral limb rigidity/bradykinesia, suggesting that the prolongation of these transient increases in high beta power is playing a pathophysiological role in Parkinson's disease. Our results could help explain why Parkinson's disease-related elevation of beta power in the motor cortex (He et al., 2017) is an inconsistent finding (Crowell et al., 2012). Total high beta power could remain unchanged despite changes in the underlying dynamics of high beta bursting. Given the dynamicity of high beta bursts, the underlying phenomenon may not always be reflected in grosser time-aggregated measurement of overall power. Intuitively, this would necessitate either less frequent beta bursting, reduced amplitude of beta bursts, or both of these phenomena. Of note, previous studies have found that beta bursting in M1 is not significantly different in terms of duration or frequency when Parkinson's disease and dystonia cohorts are compared (Wang et al., 2018). Further investigations are necessary to evaluate the relationship between bursts and disease symptomatology across movement disorder diagnoses. It may be the case that increased burst duration is not a feature unique to Parkinson's disease pathophysiology, just as increased PAC is not unique to Parkinson's disease (Kondylis et al., 2016).
We have linked transient episodes of high beta power both to increases in corticocortical synchrony and to changes in waveform shape. Waveform shape change has been shown here to relate to Parkinson's disease symptomatology and to high beta-encoded PAC in the cortex. Waveform analysis is an emerging field of study but has already made contributions to the study of Parkinson's disease pathophysiology (Cole et al., 2016). Cole et al. (2016) analyzed the motor cortical signal during DBS and found that STN DBS had the effect of “smoothing” the nonsinusoidal features of beta oscillations in M1. The authors postulated that this reflected a decrease in the synchrony of synaptic inputs reflecting reduced coherence in the corticobasal ganglia circuits. Here we show that waveform characteristics indeed change relative to the bursting state of the underlying signal and corticocortical synchrony. We furthermore demonstrate that this change in waveform shape is impaired in the PD cohort and the degree of impairment correlates with clinical symptoms of Parkinson's disease. Nonburst epochs appear to favor a higher peak-to-trough ratio and rising-to-falling steepness ratio, whereas bursting periods of the signal show the opposite. We propose that waveforms are dynamically changing between these two bursting and nonbursting states and that temporal collapsing of waveform features neglects important features specific to bursting and nonbursting epochs. Cole et al. (2016) found that the M1 sharpness ratio was positively correlated with patient rigidity scores. We did not replicate this finding with our methodology; however, we instead find that synchrony, waveform shape change, and mean burst duration are correlated with composite contralateral upper limb rigidity/bradykinesia (Figs. 2D,E, 3F, 5D). We further note that on close inspection of previous findings it appears that a large part of the previous correlations of Cole et al. (2016) could have been carried by a small number of data points, possibly representing outliers.
PAC has been strongly implicated in Parkinson's disease pathophysiology and is often increased in Parkinson's disease despite not being found to be elevated in all subjects with Parkinson's disease (de Hemptinne et al., 2013). PAC has been shown by our analysis to be largely a phenomenon confined to beta bursts. It has also been related to waveform shape with higher peak/trough ratios being linked with increased beta–gamma PAC (Cole et al., 2016). Our findings may explain why not all Parkinson's disease subjects are found to have elevated PAC. Once we acknowledge that beta-to-broadband gamma PAC is not uniform throughout the recording, we are presented with the possibility that certain recording episodes may fail to capture these important epochs of coupling, thus explaining why some patients fail to demonstrate increased PAC. Our findings establish a correlation between the magnitude of waveform change and levels of high beta phase-encoded PAC. PAC has been proposed to act as an inhibitor of cortical movement plans (Yanagisawa et al., 2012), and so changes in the dynamics of this phenomenon have great potential to act as direct or proximal causes of rigidity observed in Parkinson's disease. According to our hypothesis outlined above, and supported by the evidence presented, synchrony across cortical circuits could be driving changes in oscillation amplitude, shape, and cross-frequency oscillatory couplings.
Many questions remain unanswered in light of these findings. Attempts must be made to determine whether bursts change during stimulation. We predict on the basis of our current findings that during DBS burst duration will decrease and the Euclidean distances of waveform shape change will increase. This would account for the modulatory effect of DBS on PAC observed previously (de Hemptinne et al., 2015) and would be predicted by the correlation we observe between waveform shape changes and PAC (Fig. 6C). We have proposed that burst dynamics could also account for the individual differences seen in both power changes and PAC increases observed in Parkinson's disease. Further work is required to establish whether burst duration is indeed the primary explanatory factor for these intersubject variations. Finally, it is important to establish how burst duration, synchrony, waveform shape, and PAC change during movement. We propose that stimulation will have similar effects on these parameters to those observed during movement, explaining why patients obtain symptomatic benefit from the changes induced by DBS.
Footnotes
The authors declare no competing financial interests.
This work was supported by the National Institute of Biomedical Imaging and Bioengineering (Grant K23-EB-014326), by the National Institute of Neurological Disorders and Stroke (Grant R01-NS-097782), and by philanthropic support from the Casa Colina Center for Rehabilitation. A.B.O. is a recipient of the US State Department Fulbright Science and Technology Scholarship and received additional funding from the UCLA Neuroscience Interdepartmental Program. M.M. also was supported by a postdoctoral fellowship from American Parkinson Disease Association. We thank Joni Ricks-Oddie at the University of California, Riverside, for specialist statistical advice. We also thank William Speier of the Pouratian Laboratory, University of California, Los Angeles, for help and advice.
- Correspondence should be addressed to Andrew B. O'Keeffe at andrewbokeeffe{at}yahoo.com
