A total of 25 DBS implantation patients were examined. In these patients we simultaneously recorded from ensembles of up to 23 well-isolated neurons from either Vim/Vop or STN, depending on the site of electrode location. Recording sessions varied substantially in terms of duration and target acquisition rate, as limited by individual patient pathology and motivation. Neurons from these subcortical areas were classified by oscillatory firing patterns and tuning to target, movement, direction, and tremor. Moreover, neuronal ensemble data served as input for an offline linear prediction model to reconstruct cursor position. Finally, neuronal pairs were analyzed for evidence of functional synchrony.
STN cells (N=168) exhibited a higher (p<0.01, Mann-Whitney U test) mean firing rate than Vim/Vop cells (N=83): 15.8 ± 1.95 Hz and 11.7 ± 1.02 Hz, respectively (mean ± 1 SE in both cases). In both subcortical areas we found substantial populations of oscillatory neurons, as well as neurons strongly tuned to target, movement, direction, and tremor. Furthermore, neurons in both subcortical areas tended to show tuning to multiple parameters () rather than belonging to disjoint sets. For example, the number of Vim/Vop cells tuned to both target and tremor was higher than would be expected under statistical independence (two-tailed Fisher’s exact test, p<0.05). At the ensemble level, a substantial number of analyzed cell pairs were found to exhibit synchrony. Curiously, all tremor-associated neurons exhibited synchrony within the recorded neuronal ensemble.
Pairwise Classifications for Single Units
Neuronal tuning to target and movement
Both Vim/Vop and STN neurons represented target appearance and movement onset (). Of all tested single units tested, 29.2% of 168 Vim/Vop cells and 22.9% of 83 STN cells were found to be tuned to target appearance. Both of these percentages represent significant populations (Binomial test, p<<0.001 in both cases). shows example PETHs for three highly responsive neurons.
Behavioral tuning of subcortical neurons
Figure 3 Example PETHs. (a) Strongly tuned units to target appearance; (i) Vim/Vop cell, Patient M, 465 trials; (ii) Vim/Vop cell, Patient M, 375 trials; (iii) STN cell, Patient H, 310 trials. (b) Strongly tuned units to movement time; (i) Vim/Vop cell, Patient (more ...)
As explained above, trials with reaction times outside the 200–1000 ms range were discarded for movement tuning, and only sessions with at least 50 valid trials were subjected to further statistical analysis. Because of the additional reaction-time criterion, fewer single units were analyzed for movement tuning than for target tuning. Of these, 34.7% of 75 Vim/Vop cells and 42.3% of 26 STN cells were found to be tuned to movement. Both of these percentages represent statistically significant populations (Binomial test, p<<0.001 in both cases). shows example PETHs for three highly responsive neurons.
We also found a strong positive correlation between the strength of target appearance tuning and that of movement tuning, for both Vim/Vop and STN cells. When controlling for the number of session trials, target tuning strength significantly predicted movement tuning strength (β=0.70, p<<0.001 for Vim/Vop; β=0.71, p<<0.001 for STN). This result is consistent with , which indicates that a larger-than-expected number of neurons in both subcortical areas were tuned to both target and movement.
A portion of the correlation between target tuning and movement tuning may be explained by a tight temporal offset between target appearance and movement time. However, visual inspection of some tuned units indicated a clear decoupling of the neural encoding of target appearance and movement. Sorted raster plots from example neurons are shown in . From these, we derived the color maps in , each showing two clear bands of increased spike density. In both panels, the vertical bands are independent of movement time and are clearly related to target appearance (about 450 ms post-appearance). The diagonal bands have a near-unity slope, indicating a clear time-locked relationship between neuronal activity and movement time. For both units, the second peak in firing rate occurred about 300 ms after the defined movement time. From these data it can be concluded that these neurons were tuned to both target appearance and movement; they modulated their firing rates in relation to both events.
Figure 4 Separation between single unit response to target appearance and movement time. Time along the X axis is relative to target appearance. Panels (a) and (c) show spike raster plots relative to target appearance, with individual trials sorted by movement (more ...)
Modest differences were seen in the aggregate response patterns of Vim/Vop and STN cells classified as responsive to either target or movement (). Both cell types exhibited a mean response that peaked following target appearance (); STN cells peaked later on average. The mean response of the Vim/Vop cells peaked immediately before movement while that of the STN cells peaked concurrently with movement (). Differences in the relative lags for Vim/Vop and STN neuronal activation likely reflect the position of thalamic and STN neurons in the network hierarchy of motor control (Marsden et al., 2001
; Guillery and Sherman, 2002
; Gradinaru et al., 2009
). The motor regions of the thalamus are more involved with intention, with signals arriving prior to motor cortex activation, whereas the collaterals from motor cortex to STN deliver signals at the time of motor activation.
Figure 5 Mean normalized PETHs for all responsive units (p<0.05) from both subcortical areas using two event triggers: (a) target appearance, and (b) movement time. Reported time is relative to the event trigger. Prior to aggregation, individual PETHs (more ...)
For both aggregates, the differences between the mean Vim/Vop and STN responses were statistically significant (χ2 test, p<<0.001). However, similar proportions of Vim/Vop and STN cells were tuned to target appearance; the same was also true for movement tuning (two-tailed Fisher’s exact test, p>0.05 in both cases).
Another metric of interest for the behavioral responsiveness of subcortical neurons was directional tuning; gives the directional tuning results for both single units and multiunits. Of all tested single units, 25.3% of 75 Vim/Vop cells and 19.2% of 26 STN cells were found to exhibit directional tuning. Both of these percentages represent statistically significant populations (Binomial test, p<0.001 for Vim/Vop, p<0.01 for STN). shows example PETHs for three strongly tuned neurons. Similar proportions of Vim/Vop and STN cells were tuned to direction (two-tailed Fisher’s exact test, p>0.05).
Despite the clear separation in the neuronal response to leftward and rightward movements, note the transient regions of convergence that occurred in . In , for example, the neuronal responses to each direction converged just prior to movement. For many tuned neurons in both Vim/Vop and STN, the degree of directional modulation varied throughout the temporal window.
For both Vim/Vop and STN cells, we found strong positive correlations between directional tuning strength and the strength of both target tuning and movement tuning. When controlling for the number of session trials, target tuning strength significantly predicted directional tuning strength (β=0.16, p<0.05 for Vim/Vop; β=0.38, p<0.01 for STN). Similarly, movement tuning strength significantly predicted directional tuning strength (β=0.24, p<0.01 for Vim/Vop; β=0.43, p<0.05 for STN). The latter finding is consistent with the Vim/Vop pairwise classification result in .
Properties of multiunits
Whereas single units are identifiable as distinct neurons, a multiunit is likely comprised of distant neurons with lower SNR spike profiles. From , it can be seen that a significant population of analyzed multiunits were tuned to target, movement, and direction (Binomial test, p<0.01 in all cases). A substantial number of these tuned multiunits were found on the same recorded channel as tuned sorted units. Furthermore, when controlling for the number of session trials, the target tuning strength of single units significantly predicted the target tuning strength of same-channel multiunits (β=0.18, p<0.01). This confirms the presence of correlated tuning in nearby neurons. Thus, a substantial amount of encoded information was present in subcortical multiunits, arguing for the potential inclusion of these signals in future analyses of ensemble activity.
In order to identify potentially pathological neurons within the recorded subcortical populations, we analyzed the tremor sensitivity of single units using the discussed peri-event phase histogram (PEPH) approach; the results are given in . Of all single units tested, 12.4% of 169 Vim/Vop cells and 15.9% of 82 STN cells were found to be correlated to observable hand tremor. Both of these percentages represent statistically significant populations (Binomial test, p<0.001 in both cases). shows example PEPHs for three strongly tremor-sensitive neurons. These results demonstrate that for highly tuned units, the dependence of spike rate on tremor phase remained stable throughout the recording session (), even if the mean firing rate varied substantially . Similar proportions of Vim/Vop and STN cells were tuned to tremor (two-tailed Fisher’s exact test, p>0.05).
Tremor tuning of subcortical neurons
Figure 6 Example peri-event phase histograms (PEPHs) triggered on hand tremor phase, for strongly tremor tuned units. (a) Vim/Vop cell, Patient M, 4478 tremor periods. (b) Vim/Vop cell, Patient M, 1259 tremor periods. (c) Vim/Vop cell, Patient M, 4642 tremor periods. (more ...)
For Vim/Vop cells (but not STN cells), we found a positive correlation between the strength of tremor tuning and that of directional tuning. When controlling for the number of session trials, directional tuning strength significantly predicted tremor tuning (β=0.34, p<0.05). However, we found no relationship between tremor tuning and either undirected target or movement tuning (p>0.05 for all cases).
However, these results do not distinguish whether these tremor tuned neurons are involved in a pathological mechanism that causes tremor or merely reflect somatosensory signals indicative of tremor.
To explore neuronal oscillations in Vim/Vop and STN and their relationship to patient pathology, we inspected the autopower spectra of single unit spike trains for strong frequency peaks, yielding peak frequency and SNR (). Of all tested single units with SNR > 2, the distribution of peak frequencies showed a clear bimodal distribution with a border between low- and high-frequency oscillations at about 2.5 Hz ().
Figure 7 (a) Distribution of peak frequencies for the spike train autopower spectra of 397 analyzed single units. Units with peak frequency above 10 Hz (17.4% of all units) are not shown in this figure; they are distributed with near uniformity in the 10–25 (more ...)
Spike train spectra from Vim/Vop and STN neurons also tended to possess large amounts of energy at low frequencies, suggestive of 1/f (pink) noise. This power-law distribution has been described for cortical neurons as a stochastic process (Davidsen and Schuster, 2002
), but may also be related to slow modulations of patient attention and arousal. The observed distribution may explain the 2.5 Hz trough () that separates neurons dominated by 1/f noise from those exhibiting strong oscillations in the tremor-relevant frequency range (2.5–7.5 Hz). Only neurons with sufficient power within this frequency range were eligible to be classified as oscillatory.
The oscillatory classification results are shown in ; 21.5% of 274 Vim/Vop cells and 17.9% of 123 STN cells were classified as oscillatory. No difference was seen in the proportions of oscillatory Vim/Vop and STN cells (two-tailed Fisher’s exact test, p>0.05). shows example interspike interval (ISI) plots for three highly oscillatory cells. Note that all three ISI histograms exhibit some degree of bimodality, indicative of periodic bursting behavior.
Example ISI histograms for highly oscillatory neurons. (a) Vim/Vop cell, Patient W. (b) STN cell, Patient V. (c) Vim/Vop cell, Patient J.
shows the smoothed autopower spectra of spike trains for all analyzed single units, with each individually normalized horizontal trace corresponding to a distinct unit. From this figure, one can visually identify some of the highly oscillatory units as well as observe the congruity between multiple units from the same patient. The difference between the mean normalized spectra for Vim/Vop and STN cells () is statistically significant (χ2 test, p<<0.001). From , it is clear that STN cells tended to concentrate power at a lower frequency (3 Hz rather than 4 Hz).
Figure 9 (a) Smoothed autopower of spike trains for all sorted units with sufficient spike count. The autopower spectra (determined by Welch's method) of each horizontal trace, corresponding to a distinct unit, has been individually smoothed (using a 0.5 Hz sliding (more ...)
The pairwise classification results in reject the notion that oscillatory neurons and behaviorally tuned neurons form disjoint sets. Furthermore, we found no relationship between spike autopower peakedness and the strength of any of the three (target, movement, direction) behavioral tuning metrics (p>0.05 for all cases, for both Vim/Vop and STN). The lack of a clear anticorrelation suggests that the sets of behavioral neurons and oscillatory neurons are far from disjoint. Instead, they appear to exist as overlapping populations.
Our next analysis intended to uncover a relationship between strong oscillatory neuronal patterns and observable hand tremor. However, we did not find any clear relationship. Linear regression analysis revealed no relationship between peak frequency (2.5–7.5 Hz range) of spike train autopower spectra and corresponding hand acceleration autopower spectra (p>0.05 for both Vim/Vop and STN). No relationship was found between the sharpness of the two spectra for Vim/Vop neurons (p>0.05), but we did observe a marginally significant positive relationship for STN neurons (β=0.72, p=0.038). shows overlaid spectra for the spike train autopower and hand acceleration autopower of three highly oscillatory units. For all three cells (representative of the population as a whole), the peak frequencies do not coincide. On the other hand, we did find a marginally significant (β=0.24, p=0.055) correlation between spike autopower peakedness and tremor tuning strength for Vim/Vop cells (p>0.1 for STN cells). These findings call into question the presumed causal linear relationship between the two, suggesting the possibility of an elusive nonlinear relationship.
Figure 10 Comparison of spike train autopower and hand acceleration autopower for three highly oscillatory units. (a) Vim/Vop cell, Patient W. (b) STN cell, Patient V. (c) Vim/Vop cell, Patient J. Note that for all three cells, peak frequency does not coincide. (more ...)
Our heterodyne decoding analysis further explored this relationship by applying a nonlinear frequency shifting approach to the autopower spectra. Of 239 analyzed single units, 13.3% of Vim/Vop cells and 17.3% of STN cells were found to be tremor associated via heterodyne decoding (). Both of these percentages represent statistically significant populations (Binomial test, p<<0.001), but the difference between them is not significant (two-tailed Fisher’s exact test, p>0.05).
Heterodyne decoding is a more sensitive method for detecting tremor correlations than spectral peak analysis if the tremor signal is wide-bandwidth or prone to phase changes. Indeed, this schema may better serve to explain the relationship between the oscillatory activity of neurons and observed tremor. In fact, this can explain the similarity in the proportions of tuned neurons in and . Furthermore, Vim/Vop firing indicated a strong positive correlation between tremor tuning strength, identified using PEPHs, and heterodyne tremor tuning strength (β=0.31, p<0.001). This relationship was marginally significant in STN cells (β=0.28, p=0.086). These findings are consistent with the pairwise classification results in , which indicated a higher than expected joint classification for the two tremor tuning analyses for Vim/Vop cells.
Efficacy of neuronal recordings for kinematic predictions
We also performed offline predictions of cursor motion using the recorded ensembles. The correlation coefficient (mean ± 1.98 SE) for each of the sessions is shown in . Although the predictions varied greatly across sessions and patients, the results compared favorably with our previous study (Patil et al., 2004
). The best session for each subcortical area (Vim/Vop, STN) was chosen for further analysis, and neuron dropping curves were generated for these two sessions and fitted to a hyperbolic function (Wessberg et al., 2000
). Extrapolation of the hyperbolic fit produced estimates of the approximate ensemble sizes required to achieve R2
=0.9: 106 Vim/Vop neurons or 397 STN neurons.
Figure 11 Dependence of offline BMI predictions on neuron ensemble size. Each data point corresponds to a recording session. Correlation coefficient (R) is indicated as mean ± 1.98 SE. The best session for each subcortical area (Vim/Vop, STN) was chosen (more ...)
We analyzed neuronal synchrony in pairs of sorted units and investigated how its prevalence varied across subcortical areas; the results are given in . Using the cross-correlation approach, 43.0% of Vim/Vop pairs and 25.8% of STN pairs were found to be significantly synchronous. Both of these percentages represent significant populations (Binomial test, p<<0.001 for both cases). shows example cross-correlation plots for three highly synchronous pairs, while shows the normalized JPSTH for the same three pairs. Whereas clearly shows temporal synchronization along the diagonal (and off-diagonals), the same result is not visually discernible in .
Figure 12 (a) Example plots of cross-correlation coefficient for three highly synchronously neuron pairs. Bootstrapped simulations of coefficient shown in red. (b) Example normalized JPSTHs for the same three neuron pairs. (i) Pair of Vim/Vop cells, Patient I; (more ...)
We found highly significant differences between the Vim/Vop and STN in terms of the proportions of synchronous pairs. A significantly higher proportion of Vim/Vop pairs were synchronous than STN pairs (two-tailed Fisher’s exact test, p<<0.001).
It has been reported that the level of tremor in parkinsonian patients is positively correlated to the degree of pairwise synchrony among STN cells (Levy et al., 2000
). To test the relationship between tremor tuning and local synchrony, we compared the subpopulation of both Vim/Vop and STN neurons which were synchronous with at least one other neuron in their respective ensembles to the subpopulation of neurons tuned to hand tremor (PEPH method). Only neurons fulfilling the criteria of both individual analyses were considered. The results are shown in . The observed proportions are significantly different (two-tailed Fisher’s exact test, p<0.01), indicating a clear interaction between tremor tuning and local synchrony. Only units synchronized to at least one other unit were tuned to tremor, whereas no unsynchronized units were tuned to tremor.
Comparison of synchrony (cross-correlation method) and tremor tuning (PEPH method)