|Home | About | Journals | Submit | Contact Us | Français|
In Parkinson’s disease (PD), striatal dopamine denervation results in a cascade of abnormalities in the single unit activity of downstream basal ganglia nuclei that include increased firing rate, altered firing patterns, and increased oscillatory activity. However, the effects of these abnormalities on cortical function are poorly understood. Here, in humans undergoing deep brain stimulator implantation surgery, we utilize the novel technique of subdural electrocorticography in combination with subthalamic nucleus (STN) single unit recording to study basal ganglia-cortex interactions at the millisecond time scale. We show that in patients with PD, STN spiking is synchronized with primary motor cortex (M1) local field potentials in two distinct patterns: First, STN spikes are phase-synchronized with M1 rhythms in the theta, alpha, or beta (4-30 Hz) bands. Second, STN spikes are synchronized with M1 gamma activity over a broad spectral range (50-200 Hz). The amplitude of STN spike-synchronized gamma activity in M1 is itself rhythmically modulated by the phase of a lower frequency rhythm (phase-amplitude coupling), such that “waves” of phase-synchronized gamma activity precede the occurrence of STN spikes. We show the disease specificity of these phenomena in PD, by comparison with STN-M1 paired recordings performed in a group of patients with a different disorder, primary cranio-cervical dystonia. Our findings support a model of the basal ganglia-thalamocortical loop in PD in which gamma activity in primary motor cortex, modulated by the phase of low frequency rhythms, drives STN unit discharge.
Much work has been done to define the electrophysiology of the basal ganglia nuclei in humans with movement disorders. In Parkinson’s disease (PD) patients off dopaminergic medications subthalamic nucleus (STN) neurons have increased neuronal firing rates (Hutchison et al., 1998; Steigerwald et al., 2008; Schrock et al., 2009), oscillatory firing patterns in the theta (4-8 Hz), alpha (8-12 Hz) and beta (12-30 Hz) ranges (Rodriguez-Oroz et al., 2001; Levy et al., 2002b), and exaggerated synchronization to neighboring STN units and to STN local field potentials (LFPs) in the beta range (Levy et al., 2002b; Kuhn et al., 2005; Weinberger et al., 2006; Moran et al., 2008). Excessive basal ganglia neuronal spike synchronization is also well established in animal models of parkinsonism (Nini et al., 1995; Goldberg et al., 2004; Costa et al., 2006).
The role of the motor cortex in parkinsonian pathophysiology is less clear. Scalp electroencephalography and invasive LFP studies in humans without movement disorders show that beta frequency oscillations are normally present in motor cortex (Pfurtscheller et al., 1997; Crone et al., 1998; Miller et al., 2007). In rodent models, cortical beta rhythms are excessively synchronized to basal ganglia spike discharge (Mallet et al., 2008; Brazhnik et al., 2012), but it is not known if this is true in humans with PD. Further, interactions between basal ganglia spiking and higher frequency cortical activity have not been explored. The broad range of gamma frequencies between 50 and 200 Hz, often referred to as “broadband gamma” (Ossandon et al., 2011), is probably a surrogate measure of underlying asynchronous spiking (Manning et al., 2009), and task-related changes in broadband gamma are thought to reflect localized cortical function (Lachaux et al., 2012). Further, broadband gamma activity in primary motor cortex is elevated in PD (Crowell et al., 2012), as well as excessively coupled to the phase of low frequency cortical rhythms (phase-amplitude coupling) (de Hemptinne et al., 2012). Thus, analysis of the interaction of basal ganglia spiking with cortical broadband could yield great insight into basal ganglia-cortex dynamics in PD.
To address these questions, we recorded primary motor cortex (M1) arm area LFPs simultaneously with STN unit discharge, in awake patients undergoing deep brain stimulator (DBS) implantation surgery. M1 LFPs were recorded using subdural electrocorticography (ECoG), a technique that combines excellent spatial and temporal resolution with sufficient signal amplitude to resolve high frequency broadband activity, in addition to low frequency rhythms, relatively free of artifact. To provide a comparison group, we performed similar recordings in patients with primary cranio-cervical dystonia who had minimal symptoms in the contralateral arm. We tested the hypothesis that in PD, STN single unit discharge at rest is excessively synchronized to cortical oscillations, to cortical broadband gamma activity, and to epochs of cortical phase-amplitude coupling. Our findings support a model of the basal ganglia-thalamocortical loop in PD in which cortical 4-30 Hz rhythms are phase locked to gamma activity, and in which phase modulated “waves” of gamma activity in motor cortex drive STN hyperactivity.
Study subjects were recruited from a population of movement disorders patients scheduled to undergo DBS implantation at one of two campuses: the University of California at San Francisco (UCSF), or the San Francisco Veterans Affairs Medical Center (SFVAMC). Subjects had a diagnosis of idiopathic PD, or primary cervical or cranio-cervical dystonia, confirmed by a movement disorders neurologist (JLO or NBG). The dystonia patients included in this study were simultaneously participating in a clinical trial to evaluate the effectiveness of STN DBS in dystonia (Ostrem et al., 2011). Informed consent was obtained prior to surgery under a protocol approved by the UCSF/SFVAMC Institutional Review Board, according to the declaration of Helsinki. During the consent process, it was explained to each prospective subject that the temporary intra-operative placement of the cortical recording strip was performed solely for research purposes. Potential study subjects were characterized within 60 days prior to surgery using the following standardized rating scales: For PD patients, the Unified Parkinson’s Disease Rating Scale part III (UPDRS-III) following withdrawal of anti-parkinsonian medications for 12 hours; for dystonia patients, the Burke Fahn Marsden Dystonia Rating Scale (BFMDRS) movement score and the Toronto Western Spasmodic Torticollis Rating Scale (TWSTRS) severity score. For PD patients, preoperative UPDRS-III subscores for the hemibody contralateral to brain recordings were used in some analyses for correlations with physiological data. However, tremor was also scored for each individual recording as “present” or “absent” intra-operatively based on off-line visual inspection of limb electromyography (EMG) and accelerometry (detailed further below). Preoperative medications in PD were summarized as levodopa-equivalent doses (LED) using the following conversion factors: ropinirole × 20; pramipexole × 100; levodopa with decarboxylase inhibitor × 1; controlled release levodopa with decarboxylase inhibitor × 0.7; levodopa with decarboxylase and COMT inhibitor × 1.3 (Wenzelburger et al., 2002).
Subject inclusion criteria were: age 21-75, normal brain magnetic resonance imaging (MRI) examination, sufficient disease severity in the setting of optimal medical management to justify treatment by DBS, and ability to cooperate during awake neurosurgery. Further criteria for the primary dystonia group: dystonia affecting predominantly cervical or cranio-cervical muscles with minimal arm involvement, and no treatment with injected botulinum toxin for three months prior to surgery. 29 PD patients and 6 dystonia patients met study criteria and were included.
All subjects underwent typical procedures for planning and stereotactic surgical placement of deep brain stimulator electrodes into the STN (Starr et al., 2002). On the surgical planning MRI, the STN target was identified as a signal hypointensity on T2-weighted fast spin echo imaging lateral to the anterior border of the red nucleus, typically 11-13 mm from the midline. For targeting M1 for placement of the ECoG array, we identified a point on M1, 3 cm from the midline, based on anatomic identification of the central sulcus (Figure 1A). The intention was to target the arm-related area of M1, slightly medial to the “hand knob” (Yousry et al., 1997), in the same parasagittal plane as the typical surgical entry site for placement of STN DBS electrodes. A radio-opaque marker was stereotactically placed on the scalp over this location. When bilateral DBS implantations were planned, cortical recording was performed on one side only, typically the side with the clearest anatomic demarcation of the central sulcus. After drilling of frontal burr holes, and opening of the dura mater, we placed a 6-contact subdural ECoG strip (3mm contacts, 1cm spacing) (Ad-Tech, Racine, WI, or Integra, Plainsboro NJ) on the brain surface and directed it posteriorly in a parasagittal plane under fluoroscopic control to provide contacts covering both M1 and primary sensory cortex (S1). The locations of the burr holes were determined solely by the selection of the safest entry point for the intended DBS trajectory, and no additional skull or scalp exposure was needed for ECoG strip placement. After placement of the ECoG strip and guide tube for microelectrode recording, the burr hole was sealed with a fibrin sealant.
ECoG contact locations were confirmed anatomically using intraoperative CT (iCT) (O-arm, Medtronic, Inc.) (Shahlaie et al., 2011) or lateral fluoroscopy. iCT images were computationally fused to the preoperative planning MRI using standard surgical planning software (Framelink 5.1, Medtronic Inc., Minneapolis), and windowed so as to visualize each ECoG contact with respect to underlying MRI anatomy (Figure 1B). For cases where lateral fluoroscopy only was used, a lateral X-ray documented the location of the ECoG contact in the anterior-posterior direction with respect to the radio-opaque scalp marker.
For physiological confirmation of electrode location, we recorded somatosensory evoked potentials (SSEP) generated by median nerve stimulation. The stimulation parameters were: pulse frequency = 2Hz, pulse width = 200 μsec, pulse train length = 160, current 25-35 mAmp (typically 50% higher than the threshold for visible thumb twitch). Signals recorded from the ECoG strip during the SSEP were sampled at 5000 Hz, band-pass filtered (1-500Hz), and amplified × 7000. Signal averaged SSEPs for each adjacent contact pair (computationally remontaged) were then visually inspected to determine the M1 location by reversal of the N20 waveform. For each bipolar pair, the more posterior contact was used as the “active” electrode while the more anterior was the “reference”. The posterior contact of the most posterior bipolar contact pair that showed a negative going waveform was considered to be the contact best localized to M1 (Figure 1C). There was high concordance between anatomic and physiologic determinations of contact position with respect to central sulcus.
Cortical LFPs and STN units were recorded in a state of alert rest, using the Guideline 4000 system (FHC Inc, Bowdoin, ME) or the Alpha Omega Microguide Pro (Alpha Omega Inc, Nazareth, Israel). These are customized 14-channel data collection systems approved by the United States Food and Drug Administration for human use. Five-channel cortical LFP recordings were obtained from each of the five most posterior ECoG contacts (1-5), each differentially referenced to the most anterior contact (6). A needle electrode in the scalp served as the ground. LFP signals were band-pass filtered 1-500 Hz, amplified × 7000 and sampled at 1000 Hz (Guideline 4000 system) or 1500-3000 Hz (Alpha Omega system). A 60 Hz notch filter was employed to reduce line artifact. Data collected with the Alpha Omega system were downsampled to 1000 Hz (matlab function resample). All subjects had been sedated with propofol for the initial surgical exposure, but data collection was performed at least 30 minutes after stopping propofol. This is sufficient time for all neuronal effects of this agent to be eliminated (Fechner et al., 2004; Raz et al., 2008). Simultaneous EMG (biceps brachii, extensor carpi radialis, and gastrocnemius) and arm accelerometry were performed to confirm the absence of voluntary movement, and presence or absence of tremor, during recording. All raw, unfiltered EMG and accelerometry recordings were visually inspected off line, and if 4-6 Hz oscillations were present in at least one EMG or accelerometry trace for a minimum continuous duration of 10 seconds, the corresponding single unit/M1 LFP recording was scored as “tremor present”; otherwise it was scored as “tremor absent”. When present, 4-6 Hz oscillations typically were seen during most or all of the recording, and scoring was based on concurrence between two examiners. Because tremor can occur transiently, some recordings from a single subject could be scored as “tremor present” while others were “tremor absent”.
STN neuronal recordings were obtained using glass-coated platinum/iridium microelectrodes with impedance 0.4 –1.0 MOhm at 1000 Hz (FHC Inc., or Alpha Omega Inc.) as previously described (Starr et al., 2002). Microelectrodes were advanced into the brain using a motorized or manual microdrive (Alpha Omega Inc., or Elekta Inc., Stockholm). In a typical surgical case, one to two microelectrode penetrations were made serially through the STN on each side, separated by 2–3 mm. Signals were band-pass filtered (300 Hz to 3 kHz), amplified, played on an audio monitor, displayed on an oscilloscope, and digitized (20-kHz sampling rate). Cells were recorded at approximately every 300–800 μm along each trajectory, and each recording lasted 30-120 seconds. Neurons were screened for movement-related activity based on audible changes in action potential discharge evoked by passive (investigator-initiated) contralateral limb movements (i.e., shoulder, elbow, wrist, hip, knee, and ankle joints). Proprioceptive responsiveness of a neuron in relation to movements of one or more joints was determined by concurrence between the examiner and one or more operating room staff based on audiovisual assessments of the response. Only neurons recorded in the dorsal 4 mm of the nucleus, in a region where movement related activity was detected, and in which a minimum of 1000 spikes were digitized, were used for subsequent analysis.
Most STN units (89%) were recorded on the same trajectory as the final DBS electrode placement. Postoperative MRI, or CT computationally fused to the preoperative MRI (Framelink 5.1 software, Medtronic, Inc, Minneapolis MN) (Figure 1D), was performed to verify appropriate localization of the DBS electrode in the STN.
Digitized spike trains were imported into off-line spike-sorting software (Plexon, Dallas, TX) for discrimination of single populations of action potentials by principal components analysis (Adamos et al., 2008). The software generated a record of spike times (subsequently reduced to millisecond accuracy) for each action potential waveform detected. The first 1000 spike times were used to calculate discharge rate and oscillatory activity in the 0- to 200-Hz range. Analyses were performed using Matlab software (The MathWorks, Natick, MA). Neuronal data were accepted as single unit data in this study only if action potentials could be discriminated with a high degree of certainty, as indicated by the presence of a clear refractory period of at least 3 msec in the interspike interval (ISI) histogram. Neurons whose action potential morphology varied greatly in synchrony with the cardiac cycle were excluded. Oscillations in the spike train at 0-50 Hz were evaluated using the “global spike shuffling” method (Rivlin-Etzion et al., 2006) to eliminate the artifactual autocorrelations that arise from the neuronal refractory period. Spike time stamps were converted to a data stream consisting of 1-ms bins in which the occurrence or absence of a spike was represented by 1 or 0, respectively, in that bin. A 2,048-point fast Fourier transform with Hanning window was used, resulting in a spectral resolution of 0.5 Hz. Similar analysis was performed on “control” data in which ISIs had been randomly shuffled 100 times (Rivlin-Etzion et al. 2006). Statistically significant peaks in the spike train data (after normalization with the spectrum of the shuffled data) were determined by using the 300- to 500-Hz part of the spectrum as the control segment and its standard deviation was used as a measure of random fluctuations in the spectrum. Each frequency point between 0.25 and 200 Hz was then checked for deviation from the expected power, at a significance level of p= 0.0025, after correction for multiple (100) comparisons. Only 4-50 Hz oscillations were plotted as few occurred outside of this range.
For analysis of cortical LFPs, the raw strip recordings were computationally re-referenced to adjacent bipolar pairs by subtraction of the common reference. The M1 contact was referenced to either the contact immediately anterior or the contact immediately posterior to it, depending on which pair showed a larger spike-timed average waveform (described further below). The M1 LFP power spectrum was calculated using the Welch method (Matlab function pwelch) with a 512 point FFT, for 1.96 Hz resolution, using the entire LFP recording (duration of recordings 34.3 to 140.0 sec).
The first 1000 consecutive spike times of each recording were used for all analyses of spike-cortex synchronization. To provide a sensitive measure of spike-LFP interactions and to calculate latencies between the timing of spikes and peaks in cortical oscillatory activity, we analyzed the spike-timed average (STA) of the M1 LFP. The STA was computed in a 1 second window centered on the time of occurrence of STN spikes. We used the STA to calculate an STN spike-M1 LFP oscillation modulation index MIO as follows (also illustrated in Figure 3 B-D):
MIO values >3 were considered to be statistically significant. This threshold resulted in the occurrence of false positives (recordings in which the reverse time STA analytic amplitude Mrev also crossed the 3.0 SD level at some time within the 400 msec STA test window) at a frequency of <5%. For all STAs that exceeded the threshold of significance, we tabulated: 1) the frequency of peak spectral power, calculated as described in step 1 above (Figure 3B). 2.) The instantaneous phase value of the band-pass filtered STA at spike time (time zero) (Figure 3C). 3.) The latency from time zero to the time of the peak of the amplitude envelope of the band-pass filtered STA (Figure 3C). 4.) The duration of the interval over which the amplitude envelope of the STA crossed and stayed above the threshold of 3 standard deviations (Figure 3D). Of note, the exact time at which an STA amplitude envelope crosses its significance threshold is in part related to the frequency of the underlying oscillation, not just the strength of the STN-M1 interaction. Therefore, analyses of latency and duration of spike-LFP synchronization were only used as general descriptors, not for detailed comparison between subject groups.
To confirm the overall frequency characteristics of spike-M1 synchronization obtained by analyses of STN spike-timed averages of the M1 LFP, we also calculated STN spike-cortical LFP coherence using traditional methods for determining coherence between a point process and a field potential, described in Halliday et al. (Halliday et al., 1995).
Since low frequencies (<50 Hz) predominate in the LFP power spectrum and its STA, we utilized additional methods to study synchronization between STN spikes and cortical broadband gamma activity. First, we visualized spike-broadband gamma interactions by plotting the spike-timed scalogram of the M1 LFP, in 1 second windows centered on the time of occurrence of STN spikes, from 0 to 300 Hz, using a wavelet decomposition (Figure 4A). A Morlet wavelet of 5 cycles in width was convolved with the M1 LFP time series to estimate an amplitude and phase of the signal at each frequency for every point in time.
Next, we calculated an STN spike-M1 gamma modulation index MIgamma using similar methods as above for MIO, but first isolating the 50-200 Hz portion of the cortical LFP spectrum (also explained in Figure 4 B-E). Prior to signal averaging, we band-pass filtered the raw cortical LFP at 50-200 Hz to isolate broadband gamma (100 tap FIR filter, implemented using matlab function filtfilt), and calculated the amplitude envelope of broadband gamma time series, from the magnitude of its analytic amplitude. This time series is referred to as gamma_AE. We computed the spike-timed average (STA) of gamma_AE, in 1 second windows centered on the time of occurrence of STN spikes (N=1000 spikes). This STA waveform is indicative of the degree to which cortical broadband gamma activity is synchronized to the timing of STN spikes. We noted that, when spike-synchronized peaks in broadband gamma occurred, the amplitude of cortical broadband gamma was invariably modulated by the phase of a low frequency rhythm between 4 and 30 Hz (Figure 4B). We therefore defined the spike-gamma modulation index MIgamma in a manner that accounted for this behavior, using the amplitude envelope of the STA of gamma_AE. The method for determining MIgamma was identical to steps 1-5 described above for calculating the spike-M1 oscillation modulution index MIO, except that the analysis utilized the STA of gamma_AE rather than the STA of the unfiltered M1 LFP.
MIgamma values > 3 were considered to be statistically significant. For all spike-LFP pairs with a significant MIgamma, we calculated: 1) the frequency of peak spectral power of the STA of gamma_AE, (Figure 4B). 2.) The instantaneous phase value of the band-pass filtered STA of gamma_AE at spike time (time zero) (Figure 4D). 3.) The latency from time zero to the time of the peak of the amplitude envelope of the band-pass filtered STA of gamma_AE (Figure D). 4.) The duration of the interval over which the amplitude envelope of the STA of gamma_AE crossed and stayed above the threshold of 3 standard deviations (Figure 4E).
For the subset of recordings that had modulation indices >3 and that had >2000 consecutive spikes available for analysis, we checked the effect of spike number (500, 1000, or 2000) on the modulation indices MIO and MIgamma. For most recordings, spike-M1 synchronization reached statistical significance using a spike number of 1000, justifying the use of trains of 1000 spikes for most analyses in the study.
Parameters calculated from the STN spike-cortical LFP synchronization analyses were summarized by their means and standard deviations, or medians and ranges depending on the normality of the data (tested using Matlab function lillietest). Of note, for most study subjects, multiple STN-M1 paired recordings were sampled, from different locations within the motor territory of the STN (Tables 1 and and2).2). We did not wish to assume a priori that different recordings from the same subject were statistically independent. Therefore, for comparisons between disease groups (Figure 5A), results from standard statistical tests (Wilcoxon rank sum test) were also checked using a general estimation equations method (GEE, SPSS software), with subject number as a within-subjects variable to account for potential clustering of within-subjects data (Hanley et al., 2003). Similarly, correlations of physiological measures with symptom severity (Figures 3E and and4F)4F) were performed both with linear regression (assuming independence of all samples), and using a GEE model that did not assume sample independence within-subjects. Phase distributions were analyzed for their difference from a uniform distribution using the Rayleigh test (Moran et al., 2008). Categorical data were analyzed using chi square or Fisher exact tests.
Characteristics of the study subjects, including preoperative medications, are provided in Tables 1 and and2.2. There were 29 subjects in the PD group, and six in the dystonia group. Primary dystonia patients had mainly cranio-cervical involvement, with minimal symptoms in the arm contralateral to the side of recording in most subjects. Two dystonia subjects had action induced arm tremor or writer’s cramp, resulting in nonzero BFMDRS scores for the contralateral arm (Table 2), but did not have tremor or dystonic posturing of limbs at rest.
In PD patients, a total of 116 simultaneous recordings of STN single units with M1 LFPs were analyzed. 46 of these recordings occurred during contralateral arm tremor, while 70 occurred during no detectable contralateral arm tremor. The spike firing rate (mean +/- STD) was 35.6 +/- 14.1 Hz. Twenty-eight units had significant oscillations in the spike train in the 4-50 Hz range (Figure 2A and C). For subsequent analyses, oscillations were further subcategorized as beta frequency (12-30 Hz), theta-alpha frequency (4-12 Hz, which includes tremor frequency and its first harmonic), or neither (Table 3). The distribution of single unit oscillation frequencies was bimodal, with one beta cluster and one theta-alpha cluster, consistent with prior reports (Levy et al., 2000; Rodriguez-Oroz et al., 2001; Moran et al., 2008). Two units oscillated in both frequency bands (Figure 2A, middle example). Oscillations near tremor frequency were found in some recordings in the absence of detectable tremor during the recording, but were more likely to occur when tremor was present (details in Table 3).
Power spectra of the M1 LFP (Figure 2B) were inspected for peaks (local maxima) in the beta and tremor bands. All but one M1 recording had a spectral peak in the beta band, while 68 also had a separate peak in the theta or alpha frequency bands (Figure 2B, D). The presence or absence of contralateral limb tremor during the recording was not predictive of the presence or absence of a tremor frequency peak in the M1 LFP power spectrum (Figure 2D and Table 3).
Synchronization of STN single units to cortical rhythms was examined primarily by analysis of STN spike-timed averages of the cortical LFP. The STA waveform is indicative of the degree to which cortical oscillations are phase synchronized to the timing of STN spikes. An illustrative example is shown in Figure 3A-B. When compared across all five cortical LFP recording pairs (Figure 3B), spike-cortex synchronization was always maximal in a contact pair that covered M1. The methods of quantifying the frequency, phase, and magnitude (modulation index MIO), of the STN spike synchronization to M1 LFP oscillations, are shown in Figure 3B-D. The magnitude of spike-M1 LFP synchronization correlated with the severity of contralateral bradykinesia (using preoperative UPDRS-III bradykinesia scores, Figure 3E), and this was true whether or not within-subjects recordings were treated as statistically independent (p-value for linear correlation of independent samples of .015 with Person’s r of .224, versus p value from the GEE model of .001 with a beta coefficient of 0.37).
For those spike-M1 pairs with significant modulation (defined as MIO>3, 56 recordings), the latency, frequency and phase relationships between spike and LFP are shown in Figures 3F-H. Three recordings had significant spike-LFP modulation at two different frequencies. The median onset of spike-LFP synchrony was 103 msec prior to the spike time (p<.001 for difference from zero, nonparametric sign test) and the median time of peak modulation (defined as shown in Figure 3C) was 18.5 msec prior to the spike time (not significantly different from zero by the sign test) (Figure 3F). The distribution of STN-M1 LFP synchronization frequencies, derived from the power spectra of the spike-timed averages, showed clusters in the theta and beta bands. (Figure 3G). STN-M1 LFP synchronization at or near tremor frequency (4-12 Hz) was much more prominent when limb tremor was present during the recording (Table 3). Similar results were found when STN spike–M1 LFP coherence was analyzed by the Halliday method rather than by utilizing spike time averages of M1 (data not shown). The distribution of phases for STN-M1 synchronization was significantly different from a uniform distribution (p=.001, Rayleigh test), with a mean phase angle of 1.6 radians, corresponding to the falling slope of the bandpass filtered M1 LFP sinusoid (Figure 3C, H). STN spike synchronization to the M1 LFP was not restricted to neurons that had prominent arm related activity during intraoperative somatosensory examination, but was also found for neurons whose discharge was primarily modulated by leg movement or had no detectable response to passive movement (p=0.83, chi square) (Figure 3I), notwithstanding the fact that cortical recordings were made exclusively over the arm territory of M1.
The use of the spike-timed LFP average for analysis of spike-M1 synchrony (Figure 3) mainly explores the low frequency range of the LFP power spectrum (<50 Hz), as these frequencies predominate in the LFP and in its STA. However, we were also interested in synchronization between STN spikes and cortical broadband gamma activity (50-200 Hz), as the latter is thought to reflect underlying cortical spiking activity (Manning et al., 2009), and therefore offers a means of exploring the relationship between STN spiking and cortical spiking. To visualize STN spike-M1 gamma interactions, we plotted spike-timed scalograms, illustrated for one representative example in Figure 4A. In the scalogram, spectral power is normalized separately at each frequency (to the mean of its wavelet convolution), allowing visual examination of gamma frequencies that are too small to observe in the spike-timed LFP average illustrated in Figure 3.
Inspection of spike timed-M1 scalograms revealed synchronization of broadband gamma with STN spikes, in a highly distinctive pattern: “waves” of increased M1 broadband gamma activity, alternating with troughs of broadband gamma. The occurrence of this pattern in the scalogram implies that there are epochs of cortical phase-amplitude coupling whose phase frequency is itself synchronized to the timing of STN spikes. The highest amplitude gamma waves typically preceded the time of occurrence of STN spikes. This oscillating pattern of broadband gamma activity had an amplitude envelope typically in the beta, alpha, or theta range. The frequency range of broadband gamma that was most prominently synchronized to STN spiking was 50-200 Hz (Figure 4A), so this region of the spectrum was utilized for further quantitative analysis.
Figure 4B-E shows the method of quantifying the parameters of STN spike–M1 gamma synchronization, including the frequency and phase of the amplitude envelope for gamma activity, the latency from peak gamma modulation to STN spike time, and the overall strength of the spike-gamma interaction (modulation index MIgamma). The magnitude of spike-M1 gamma synchronization did not correlate with lateralized total preoperative UPDRS scores or with contralateral bradykinesia scores, but did have a modest negative correlation with contralateral tremor scores (Figure 4F). These results were invariant to the statistical method used for correlation (treating all recording as independent samples, p=.038 with r of -.193, or treating within-subjects recordings as non-independent samples in a GEE model, p=.043 and beta coefficient of -.226, for the inverse correlation with tremor). For those spike-M1 pairs with significant modulation of M1 gamma (defined as MIgamma >3, 46 recordings), the latency, frequency and phase relationships for synchronization of STN spikes to the M1 gamma amplitude envelope are shown in Figures 4 G-I. The median onset of spike-gamma synchrony was 107.1 ms prior to the spike time (p<.001 for difference from zero, nonparametric sign test) and the median time of peak modulation (defined as shown in Figure 4D) was 48 msec prior to the spike time (p=.05 for difference from zero, sign test, Figure 4G). The distribution of frequencies for the spike-synchronized M1 gamma envelope, derived from the power spectra of the spike-timed average of gamma activity, shows that the gamma amplitude envelope oscillated most prominently in the theta and beta bands (Figure 4H). Spike-gamma envelope synchronization at any frequency, and especially at beta frequencies, was more likely to occur in the absence of contralateral tremor (Table 3). This contrasts with the relationship between spike-LFP synchronization and contralateral tremor (Table 3), and suggests a relationship between symptom profile and the type of synchronization (MIO versus MIgamma), not just the frequency, of basal-ganglia cortex interactions. The distribution of phases for STN-M1 gamma envelope synchronization showed variable phase relationships that did not differ from a uniform distribution (p>.05, Rayleigh test) (Figure 4I). STN spike synchronization to M1 gamma was not restricted to neurons that had prominent arm related activity during intraoperative somatosensory examination (p=0.72, chi square, Figure 3J).
Although normal control data for invasive recordings cannot be obtained from humans, we were able to assess the disease specificity of our findings in PD subjects by comparison with patients with primary cranio-cervical dystonia (Figure 5), who also underwent placement of stimulators into the STN (Ostrem et al., 2011). Twenty-six spike-LFP pairs from six dystonia patients were available for analysis. For STN synchronization to M1 oscillations in the 4-50 Hz range, the magnitude of synchronization (MIO) was less than in PD (p=.005 for Wilcoxon rank sum test and .009 for GEE model, Figure 5A, left). Six spike-LFP pairs in dystonia patients had significant synchronization (MIO>3), in the theta, alpha, or beta ranges (Figure 5B, left), but the proportion of all recordings with significant beta range spike-LFP synchronization was less than for PD (4 out of 26 in dystonia versus 41 out of 116 in PD, p=.05, chi square). For STN synchronization to M1 gamma activity, six paired recordings in dystonia also showed “waves” of amplitude modulated broadband activity, slightly preceding the time of occurrence of STN spikes, one of which occurred at two different frequencies of the sinusoidal amplitude envelope (Figure 5B, right). However, the pattern of spike-gamma synchronization was strikingly different from that seen in PD, in that the frequency of the sinusoidal amplitude envelope for broadband gamma was never in the beta band in dystonia patients but was commonly in the beta band in PD patients (0 out of 26 in dystonia versus 23 out of 116 in PD, p=.008 by Fisher exact test).
Two forms of spike-M1 synchronization were studied here: spike synchronization to M1 rhythms (4-50 Hz), and spike synchronization to M1 gamma (50-200 Hz). Several lines of evidence suggest that these represent distinct phenomena. First, across all recordings with a significant MIO or significant MIgamma, the value of these modulation indices were not correlated (Figure 6A and B). Second, although the spike-gamma modulation had an oscillating amplitude envelope (Figure 4 A-C), the frequency of the amplitude envelope for the STN spike-M1 gamma modulation did not correlate with the frequency for spike-M1 rhythm modulation (Figure 6B). Thus, the mechanisms for generating these two modes of synchronization may be at least partly independent.
To examine cortex-basal ganglia interaction in humans with movement disorders, we recorded M1 local field potentials simultaneously with STN single unit activity in patients undergoing STN DBS implantation for Parkinson’s disease, as well as in a comparison group with primary cranio-cervical dystonia. LFPs represent summed, synchronized subthreshold activity in pre and postsynaptic elements near the recording electrode, and are increasingly used in studies of neuronal synchronization in the motor system in normal as well as disease states (Murthy and Fetz, 1996; Goldberg et al., 2004; Hammond et al., 2007). In general, the phase of LFP oscillations biases the probability of spike firing. In the normal state, this modulation is dynamic and task specific, representing a flexible mechanism whereby LFP oscillations can link functionally related brain areas (Fries, 2005; Leventhal et al., 2012).
Excessive neuronal oscillatory activity within specific basal ganglia nuclei is a well-documented feature in parkinsonian animal models as well as in human PD (Dostrovsky and Bergman, 2004; Hammond et al., 2007). Similarly, simultaneously recorded neurons within the subthalamic nucleus or internal globus pallidus have oscillations that are strongly synchronized (Nini et al., 1995; Levy et al., 2002a; Hanson et al., 2012), and basal ganglia neurons tend to be synchronized to LFPs recorded in the same structure (Moran et al., 2008) (Kuhn et al., 2005). However, the relation of basal ganglia oscillations to cortical activity is less well studied. In humans, interactions between cortex and basal ganglia at fast time scales have primarily been studied by recording basal ganglia LFPs in conjunction with scalp electroencephalography (Williams et al., 2002; Fogelson et al., 2006; Lalo et al., 2008; Eusebio et al., 2009) or magnetoencephalography (Litvak et al., 2011). These studies have elucidated specific patterns of cortex-basal ganglia coherence, but in general have not utilized nonparkinsonian comparison groups.
Here, we provide the first demonstration of the synchronization of basal ganglia spiking with cortical oscillations in humans with PD. This behavior is relatively specific to the parkinsonian state, as it is less prominent in a different movement disorder of basal ganglia origin, primary craniocervical dystonia. We show that synchronization occurs across a broad range of low frequency rhythms (theta, alpha and beta). Spike-LFP synchronization at tremor (theta) frequency is much more likely when tremor is present during the recordings (Figure 3G). However, tremor frequency oscillations in spike discharge and M1 LFPs do occur in the absence of tremor (Figure 2), suggesting that neuronal oscillatory phenomena at tremor frequency are not exclusively a consequence of tremor-frequency peripheral input. It is possible that tremor frequency discharges at different points in the circuit represent a “parkinsonian endophenotype”, but that the clinical manifestation of tremor specifically requires synchronization of these oscillations between basal ganglia and cortex.
Our findings support the evolving view that excessive synchronization between basal ganglia and cortex is a prominent feature of PD. We provide evidence that the parkinsonian state alters the statistical relationship between LFP phase and spike timing in global brain networks, not just in local basal ganglia or cortical circuits (Goldberg et al., 2004; Gatev and Wichmann, 2009; Li et al., 2012). Since many brain functions appear to depend on this mechanism for task performance (Fries, 2005), its disruption in a critical motor circuit could results in widespread impairment of motor functions. Our findings in humans are broadly consistent with several rodent studies of parkinsonism showing that dopamine denervation produces excessive 20-40 Hz synchronization of basal ganglia spiking with motor cortical LFP oscillations (Walters et al., 2007; Mallet et al., 2008; Galati et al., 2009; Brazhnik et al., 2012). The strong involvement of primary motor cortex in abnormal synchronization in humans has therapeutic implications, supporting the possibility of interrupting pathological synchronization at the level of the cortex (Drouot et al., 2004; Wu et al., 2007; Li et al., 2012), or of using a cortically based signal detector as a driver of basal ganglia deep brain stimulation in a closed loop paradigm (Rosin et al., 2011).
Our most novel finding is that STN spiking is also synchronized to M1 broadband gamma activity. Spike time averaged cortical gamma activity has a highly distinctive pattern, with “waves” of alternating high and low gamma amplitude, whose peak amplitude typically precedes the timing of STN spikes (median latency of -48 msec). These gamma “waves” have an oscillatory envelope in the theta, alpha or beta ranges, demonstrating an interaction between the phase of low frequency rhythms and the amplitude of broadband gamma.
This finding has important implications for the role of cross-frequency interactions in movement disorders. Previous studies of STN LFPs have called attention to the possibility of phase-amplitude interactions as a biomarker of PD (Lopez-Azcarate et al., 2010; Ozkurt et al., 2011). These studies showed coupling between the beta rhythm and a narrow bandwidth, very high frequency activity (250-300 Hz). Subsequently, we used M1 LFP recordings to demonstrate that PD patients have excessive coupling of the phase of the beta rhythm to broadband gamma amplitude (de Hemptinne et al., 2012). The present study directly links this abnormal cortical phase amplitude coupling to basal ganglia discharge, and indicates that epochs of M1 phase amplitude coupling statistically predict the timing of STN spiking. Further, we show a symptom and disease specific pattern of spike-synchronized gamma activity: Tremor severity is inversely related to the magnitude of spike-synchronized gamma activity in PD (at any frequency, Table 3 and Figure 4F), and beta phase modulation of spike-synchronized gamma activity occurred uniquely in the parkinsonian state, not in cranio-cervical dystonia patients. In the cortex, we did not observe any narrow bandwidth, very high frequency peaks of gamma activity in the 250-300 Hz range that were distinct from the broadband gamma range, as has been observed in the STN LFP (Foffani et al., 2003; Lopez-Azcarate et al., 2010; Ozkurt et al., 2011)
In classical models of PD, STN hyperactivity is thought to result from abnormal synaptic transmission through intrinsic basal ganglia pathways involving the striatum, primarily the substantia nigra compacta (SNc)-striatum-globus pallidus externa (GPe)-STN pathway (Bergman et al., 1990). Our findings suggest an additional mechanism for STN hyperactivity in PD. Since cortical broadband gamma activity is thought to reflect underlying action potential firing (Manning et al., 2009), we propose that M1 drives STN hyperactivity via the corticosubthalamic pathway (Nambu et al., 2002). The organization of M1 spiking activity into phase-modulated waves of high and low amplitude is ideally suited to drive STN activity, since strong phasic inputs to STN at beta envelope frequencies would drive STN spiking more efficiently than an asynchronous increase in corticosubthalamic activity (Beurrier et al., 1999; Plenz and Kital, 1999; Bevan et al., 2002), and at a faster rate than phasic inputs at lower envelope frequencies. The latency between M1 population spiking and STN spiking is longer than can be explained by fast monosynaptic neurotransmission. This may relate to the membrane dynamics of STN cells, in which repetitive slow depolarization is especially effective for triggering spikes in “burst mode”, thought to be characteristic of parkinsonism (Beurrier et al., 1999). A primary role for pathological activity in the corticosubthalamic pathway is also attractive in light of recent investigations in rodent models of PD, showing that modulation of the cortiocosubthalamic pathway is important for amelioration of motor deficits (Gradinaru et al., 2009; Li et al., 2012). Further, modeling work indicates that network oscillations in the basal ganglia-thalamocortical circuit are enhanced in the setting of excessive corticosubthalamic driving of intrinsic basal ganglia pathways (Holgado et al., 2010).
All recordings studied here were performed in a state of alert rest. We did not study patients during a movement task due to the extreme difficulty of holding stable unit activity for a long period of time during patient movement. We did not study patients in the on medication state, so it is not clear how dopamine replacement therapy would affect the observed STN-cortex interactions (Williams et al., 2002). Due to the smaller number of dystonia subjects implanted, the number of units in dystonia patients with sufficient spike numbers for comparison with PD was relatively small. Because our PD subjects had longstanding disease and had received chronic levodopa, our findings may represent compensatory changes rather than primary abnormalities, or may reflect long term effects of chronic levodopa treatment that do not “wash out” in 12 hours off of levodopa (Picconi et al., 2003).
Using the novel approach of simultaneous recording of M1 LFPs and STN single units in awake humans undergoing DBS implantation surgery, we studied basal ganglia – cortex interactions on a fast time scale. In humans, we confirm the evidence from animal models that basal ganglia unit discharge in the parkinsonian state is strongly synchronized to cortical oscillatory activity. We further show that STN spiking is synchronized with cortical broadband gamma, and that the latter occurs in a phase modulated pattern that begins prior to the occurrence of STN spikes. Abnormal cortical phase amplitude coupling may be an important mediator of network oscillations in the basal ganglia-thalamocortical circuit.
This work was supported by the National Institutes of Health [R01NS069779 to P.A.S.]; and the Dystonia Medical Research Foundation. We thank Leslie Markun and Marta San Luciano for assistance with clinical data and statistical analyses, and Coralie de Hemptinne for critical review of the manuscript.
Conflict of Interest: The authors declare no competing financial interests