Model-Based Evaluation of Closed-Loop Deep Brain Stimulation Controller to Adapt to Dynamic Changes in Reference Signal

High-frequency deep brain stimulation (DBS) of the subthalamic nucleus (STN) is effective in suppressing the motor symptoms of Parkinson's disease (PD). Current clinically-deployed DBS technology operates in an open-loop fashion, i.e., fixed parameter high-frequency stimulation is delivered continuously, invariant to the needs or status of the patient. This poses two major challenges: (1) depletion of the stimulator battery due to the energy demands of continuous high-frequency stimulation, (2) high-frequency stimulation-induced side-effects. Closed-loop deep brain stimulation (CL DBS) may be effective in suppressing parkinsonian symptoms with stimulation parameters that require less energy and evoke fewer side effects than open loop DBS. However, the design of CL DBS comes with several challenges including the selection of an appropriate biomarker reflecting the symptoms of PD, setting a suitable reference signal, and implementing a controller to adapt to dynamic changes in the reference signal. Dynamic changes in beta oscillatory activity occur during the course of voluntary movement, and thus there may be a performance advantage to tracking such dynamic activity. We addressed these challenges by studying the performance of a closed-loop controller using a biophysically-based network model of the basal ganglia. The model-based evaluation consisted of two parts: (1) we implemented a Proportional-Integral (PI) controller to compute optimal DBS frequencies based on the magnitude of a dynamic reference signal, the oscillatory power in the beta band (13–35 Hz) recorded from model globus pallidus internus (GPi) neurons. (2) We coupled a linear auto-regressive model based mapping function with the Routh-Hurwitz stability analysis method to compute the parameters of the PI controller to track dynamic changes in the reference signal. The simulation results demonstrated successful tracking of both constant and dynamic beta oscillatory activity by the PI controller, and the PI controller followed dynamic changes in the reference signal, something that cannot be accomplished by constant open-loop DBS.

Introduction

Parkinson's disease (PD) is characterized by degeneration of dopaminergic neurons in the substania nigra pars compacta (SNc) resulting in motor symptoms including bradykinesia, rest tremor, postural instability, and rigidity (Davie, 2008; Jankovic, 2008). High-frequency deep brain stimulation (DBS) of the subthalamic nucleus (STN) or globus pallidus internus (GPi) is a well-established surgical therapy to treat the motor symptoms of PD (Krack et al., 2003; Rodriguez-Oroz et al., 2005; Odekerken et al., 2016). Current clinical DBS technology is open loop—stimulation is always on and the stimulation parameters are tuned periodically through manual adjustments by health care professionals. The process of selection of DBS parameters is challenging due to the large number of parameters (Kuncel and Grill, 2004). Therefore, the efficacy of current open-loop DBS may be suboptimal and patients can experience side effects, including speech deficits and cognitive dysfunction (Deuschl et al., 2006; Okun and Foote, 2010; Massano and Garrett, 2012; Cyron, 2016).

Recent clinical studies suggest that closed-loop DBS (CL DBS) may be more efficient at suppressing PD motor symptoms with reduced side effects as compared to continuous high-frequency STN DBS (Rosin et al., 2011; Carron et al., 2013; Hebb et al., 2014; Rossi et al., 2016). However, the design of CL DBS controllers comes with several challenges including selection of a feedback signal reflecting PD symptoms and the capacity of the controller to adapt to dynamic changes in the reference signal (Hebb et al., 2014; Arlotti et al., 2016a; Parastarfeizabadi and Kouzani, 2017). Concurrent neuronal recordings and behavioral assessments from PD patients and animal models of PD showed a strong correlation between beta band oscillations (13–35 Hz) and PD motor symptoms, especially bradykinesia (Zaidel et al., 2010; Jenkinson and Brown, 2011; Little and Brown, 2012; Hoang et al., 2017), and beta band activity may be an appropriate feedback signal for CL DBS. However, beta oscillations in the basal ganglia desynchronize in preparation and during voluntary movement (Levy et al., 2002; Brittain and Brown, 2014). Therefore, a fixed beta power reference may not be appropriate for control of DBS, and it may be beneficial to include in the controller design the ability to adapt to dynamic changes in the reference signal.

The objective of this study was to design a controller for CL DBS that can adapt to dynamic changes in the reference signal. We evaluated the performance of a proportional integral (PI) controller using a network model of the basal ganglia (BG) (Kumaravelu et al., 2016). The parameters of the PI controller were tuned by coupling a linear controlled auto-regressive (CAR) model with Routh-Hurwitz stability analysis. The PI controller was successful in adapting to dynamic changes in the reference signal, and such a control scheme may be suitable for implementation in CL DBS systems.

Methods

A block diagram of the proposed CL DBS framework is shown in Figure 1A. The signal power of model neuron activity in the beta band was used as the feedback signal y(k), and the error e(k) between the actual beta power and the desired beta power y_sp(k) was sent to the PI controller to calculate the stimulation frequency u(k). Thus, the PI controller calculated the DBS frequency according to the variation of beta oscillatory power. The calculated DBS frequency determined the time of the next stimulation pulse I_dbs(t) delivered to a biophysical network model of the parkinsonian cortex-basal ganglia-thalamus (CTx-BG-Th) network. The selection of appropriate PI controller parameters was required for the actual beta power to track dynamic variations in the desired power. Below we propose a stability analysis method to calculate automatically the PI parameters.

FIGURE 1

Figure 1. (A) The CL DBS framework. The spike times of model neurons in the GPi were calculated, and the beta band power of these spike times was used as the feedback signal y (k). The error term e (k) between the desired beta power y_sp (k) and actual value y (k) was input to the PI controller to calculate the stimulation frequency u (k). The stimulation signal I_dbs (t) delivered to the cortex-basal ganglia-thalamus network model was subsequently determined. (B) The transformed linear system of the CL DBS system. This transformed linear system was used to determine the appropriate parameters for the PI controller, and the PI parameters were constant once calculated.

Computational Model of the Cortex-Basal Ganglia-Thalamus Network

We used a model of the CTx-BG-Th network as a test bed to evaluate the performance of the closed-loop control scheme (Kumaravelu et al., 2016), and a implementation of this model in MATLAB can be downloaded from ModelDB (https://senselab.med.yale.edu/modeldb/). The CTx-BG-Th model included the cortex, striatum, STN, globus pallidus externus (GPe), GPi and a thalamic nucleus, and each region was comprised of 10 single-compartment Hodgkin-Huxley type neurons. In the original publication, the model was validated extensively, including matching the responses evoked in the basal ganglia by cortical stimulation in rats (Kita and Kita, 2011), model neuron firing rates and patterns that were consistent with parkinsonian rats (Mallet et al., 2008), and responses to STN DBS at different frequencies that matched those measured experimentally (McConnell et al., 2012; So et al., 2012). Model BG neurons exhibited exaggerated low-frequency oscillatory activity in the parkinsonian state compared to the healthy condition, similar to that seen in vivo. Since, beta oscillatory activity is well-correlated with PD symptoms (Leventhal et al., 2012; Stein and Bar-Gad, 2013), we chose the beta band (13–35 Hz) power present in the activity of the GPi neurons as the model-based proxy for symptoms (Brocker et al., 2013) to evaluate the effectiveness of the CL DBS controller. There is a strong correlation between single unit firing and LFPs in the beta band in the STN (Levy et al., 2002; Kühn et al., 2005) and in the GP (Goldberg et al., 2004), and the power spectrum calculated from the single unit spike times of GP neurons was correlated with motor symptoms of parkinsonism (McConnell et al., 2012). Simulations were implemented in MATLAB R2016a and equations were solved using the forward Euler method with a time step of 0.01 ms; spectral analyses were performed using the “mtspecgrampt” function of the Chronux neural signal analysis package (chronux.org) (sliding 1 s window, 0.1 s step size and [3 5] tapers (3 is the time-bandwidth product and 5 is the number of tapers)). The spectrum of all 10 GPi neurons spike time series was calculated using the multi-taper spectral estimation method.

Identification of Relationship Between Stimulation Frequency and Beta Band Power of GPi Model Neurons Spike Times

The oscillations within the CTx-BG-Th network were similar across the different parts of the loop (Kumaravelu et al., 2016), for STN, GPi, and GPe both single neuron and local field potentials (LFPs) exhibited excessive beta band oscillation in the PD state, while for thalamus and cortex single neuron oscillation were not dominant (Stein and Bar-Gad, 2013). The beta band power of GPi model neurons spike times was chosen to characterize the model state. The dynamics of the CTx-BG-Th network were highly non-linear and therefore it was inappropriate to use the linear PI controller to control directly the network model of PD. A linear model of the plant between the stimulation frequency and the beta band power of GPi model neuron spike times was first identified using a CAR model. The structure of a CAR model was

$\begin{array}{l}(1+a_1{z}^{-1}+a_2{z}^{-2}+\cdots +a_{n_a}{z}^{-n_a})y(k)\\ =(b_0+b_1{z}^{-1}+b_2{z}^{-2}+\cdots +b_{n_b}{z}^{-n_b})u(k)+\epsilon (k)& (1)\end{array}$

where z was the lag operator, u(k) was the input signal (stimulation frequency) and y(k) was the output signal (beta power of GPi model neuron spike times), n_b and n_a were the order of input and output sequences, respectively, and ε(k) was assumed to be white noise. The identification process included the following steps:

1. Collect input and output data from the CTx-BG-Th network model.

2. Estimate model parameters a₁ ⋯ a_na and b₀ ⋯ b_nb.

3. Choose appropriate order parameters n_a and n_b.

4. Quantify the prediction accuracy of the CAR model.

The identification accuracy of the CAR model was highly dependent on the input output data that were selected, because not all data provided an equal amount of information (Ljung, 1999). The designed stimulation sequence was delivered to the CTx-BG-Th model (in the open loop), and the corresponding output data (beta band power) was calculated. To obtain more informative input/output data to identify the CAR model, the frequencies (input data) of the stimulation waveform were chosen randomly between 5 and 200 Hz. Figure 2A illustrates the stimulation sequence from 12 to 16 s, illustrating that each frequency continued for 0.4 s to ensure that at least two pulses were delivered for each random frequency. The simulation duration was 400 s, resulting in responses to 1,000 frequency samples. The stimulation sequence was delivered to the computational model of the CTx-BG-Th network, and spiking activity was recorded from GPi model neurons. The time window used to bin the beta power of GPi spike times was sensitive to the temporal dynamics of beta power when the stimulation frequencies were randomly changed (Figure 2C). Differences in beta power across time window bins were compared using one-way ANOVA with post-hoc Tukey's honestly significant difference (HSD) test, and statistical significance was defined as α = 0.05. The beta power varied across different time window bins (F = 252.54, p < 0.0001). When the time window bin was larger than 0.1 s, the calculated beta power was no <1.6 times the value with time window bin equaled to 0.1 s. The choice of the short 0.1 s bin enabled capture of small dynamic changes in beta power, as our objective was to implement a controller that responded to such changes. Bin sizes of 0.2 s or longer did not reflect the dynamic variation of the beta power, as indicated by the invariance to bin size. Since each frequency was delivered for 0.4 s and the bin used to calculate beta power was 0.1 s, the beta power obtained in Figure 2B was the average of four values within 0.4 s.

FIGURE 2

Figure 2. The stimulation sequence (A) and beta band power (B) obtained from the CTx-BG-Th network model to train the CAR model. Only data from 12 to 16 s are presented to improve visualization. (A) The stimulation frequency was randomly selected from 5 to 200 Hz, and for each frequency the corresponding stimulation sequence lasted for 0.4 s. (B) Circles represent the beta power value in each 0.1 s, the collected beta power within 0.4 s was the average of four values. (C) Mean ± standard deviation of beta power of GPi model neuron spike times plotted as a function of time window bin (50 trials). The mean value of beta power varies across different time window bin values, and values not sharing the same letter were significantly different (p < 0.05, Tukey's HSD).

We used the recursive least squares (RLS) method (Ljung, 1999) to estimate the CAR model parameters. The CAR model was transformed into a standard LS form (Ljung, 1999),

$\begin{array}{l}y(k)=-a_{1}y(k-1)-a_{2}y(k-2)-\cdots -a_{n_a}y(k-n_a)\\ +b_{0}u(k)+b_{1}u(k-1)+b_{2}u(k-2)+\cdots \\ +b_{n_b}u(k-n_b)+\epsilon (k)={\phi }^T(k)\theta +\epsilon (k)& (2)\end{array}$

where $\phi \left(k\right)=[-y\left(k-1\right),\cdots ,-y\left(k-n_a\right),u\left(k\right),\cdots ,u(k-{n_b)]}^T$ was the known sequence of input and output data, and $\theta {=\left[a_1,a_2,\cdots ,a_{n_a},b_0,b_1,\cdots ,b_{n_b}\right]}^T$ was the vector of unknown model parameters. From Equation (2), the current value of the output signal was correlated with the past input and output signals as well as the current input signal. Then, unknown model parameters were estimated by the RLS method,

$\begin{array}{ll}{y}_e(k)={\phi }^T(k)\widehat{\theta }& (3)\end{array}$

where $\widehat{\theta }$ was the estimated parameter vector calculated using the following equations:

$\begin{array}{ll}\{\begin{array}{c}\widehat{\theta }\left(k\right)=\widehat{\theta }\left(k-1\right)+K\left(k\right)\left[y\left(k\right)-{\phi }^T\left(k\right)\widehat{\theta }\left(k-1\right)\right]\\ K\left(k\right)=\frac{P\left(k-1\right)\phi \left(k\right)}{1+{\phi }^T\left(k\right)P\left(k-1\right)\phi \left(k\right)}\\ P\left(k\right)=\left[I-K\left(k\right){\phi }^T\left(k\right)\right]P\left(k-1\right)\hfill \end{array}& (4)\end{array}$

The root mean square error (RMSE) between the actual output signal and the CAR model predicted output signal was used to quantify the prediction accuracy of the CAR model,

$\begin{array}{ll}{e}_{RMSE}=\sqrt{\frac{1}{N}{\sum }_{k=1}^N{(y(k)-y_e(k))}^2}& (5)\end{array}$

The e_RMSE declined as the CAR model order (n_a and n_b) was increased (Figure 3A). Since the purpose of the identified CAR model was to design the PI controller but not to substitute for the original CTx-BG-Th network model, we were not interested in higher-order dynamics. Akaike's information criterion (AIC) was used to select the model order (McQuarrie and Chih-Ling, 1998),

$\begin{array}{ll}\text{AIC}=\frac{2K-2L}{N}+\frac{2K(K+1)}{N-K-1}& (6)\end{array}$

where K = n_a + n_b + 1 was the number of parameters to be estimated, N was the length of predicted data, and $L=-\frac{N}{2}ln(2\pi )-\frac{N}{2}ln\left(\frac{{e_{RMSE}}^2}{N}\right)-\frac{N}{2}$. When n_a = 3 and n_b = 3 the valued of AIC was minimized, thus, the structure of the CAR model was

$\begin{array}{l}{y}_e(k)=-a_1{y}_e(k-1)-a_2{y}_e(k-2)-a_3{y}_e(k-3)\\ +b_{0}u(k)+b_{1}u(k-1)+b_{2}u(k-2)\\ +b_{3}u(k-3)& (7)\end{array}$

and the corresponding estimated CAR model parameters in each iteration are shown in Figure 3B.

FIGURE 3

Figure 3. (A) The relationship between the CAR model order parameters (n_a and n_b) and the RMS error (e_RMSE) between the actual output signal and the CAR model predicted output signal. (B) The estimated CAR model parameters across iterations to minimize e_RMSE.

Selection of PI Controller Parameters

Although a common Proportional-Integral-Differential (PID) controller has three control terms (P, I, and D), we only chose the P and I terms, because the D action is sensitive to the model prediction accuracy (Aström and Hägglund, 1995). With the selected CAR model, e_RMSE = 27.9, there were still prediction error, and the D term was not used due to these inaccuracies of the CAR model. The transformed system with the CAR model substituted for the network model was used to choose the parameters of the PI controller (Figure 1B).

The structure of a discrete PI controller was (Aström and Hägglund, 1995),

$\begin{array}{ll}u(k)=u(k-1)+k_p\left[e(k)-e(k-1)\right]+k_{i}e(k)& (8)\end{array}$

and the aim was to select the P term and I term coefficients, k_p and k_i. The Routh-Hurwitz stability criterion (Gopal, 2002) was used to calculate automatically the PI parameters, where the selected PI controller must ensure the stability of the system. The forward transfer function of this system (Figure 1B) is given by Equation (9).

$\begin{array}{l}G(z)=\frac{Y(z)}{E(z)}=\frac{Y(z)}{U(z)}·\frac{U(z)}{E(z)}=\frac{b_0{z}^3+b_1{z}^2+b_{2}z+b_3}{z^3+a_1{z}^2+a_{2}z+a_3}·\frac{(k_p+k_i)z-k_p}{z-1}\\ =\frac{b_0(k_p+k_i)z^4+\left[b_1(k_p+k_i)-b_0{k}_p\right]z^3+\left[b_2(k_p+k_i)-b_1{k}_p\right]z^2+\left[b_3(k_p+k_i)-b_2{k}_p\right]z-b_3{k}_p}{z^4+(a_1-1)z^3+(a_2-a_1)z^2+(a_3-a_2)z-a_3}& (9)\end{array}$

The closed-loop transfer function was

$\begin{array}{ll}\Phi (z)=\frac{G(z)}{1+G(z)}& (10)\end{array}$

The characteristic equation of this system was

$\begin{array}{l}D\left(z\right)=1+G\left(z\right)=\left[1+b_0\left(k_p+k_i\right)\right]z^4+[\left(a_1-1\right)\\ +b_1\left(k_p+k_i\right)-b_0{k}_p]z^3+[\left(a_2-a_1\right)+b_2\left(k_p+k_i\right)\\ -b_1{k}_p]z^2+[\left(a_3-a_2\right)+b_3\left(k_p+k_i\right)\\ -b_2{k}_p]z-a_3-b_3{k}_p=0& (11)\end{array}$

According to the Routh-Hurwitz stability criterion, we substituted z with w, where $z=\frac{w+1}{w-1}$, and the variable of the characteristic equation became w.

$\begin{array}{l}\text{D}\left(w\right)=m_4{\left(\frac{w+1}{w-1}\right)}^4+m_3{\left(\frac{w+1}{w-1}\right)}^3+m_2{\left(\frac{w+1}{w-1}\right)}^2\\ +m_1\left(\frac{w+1}{w-1}\right)+m_0=0& (12)\end{array}$

Combining Equations (11) and (12), m₄ = 1+b₀(k_p+k_i), m₃ = (_{a1 − 1)+b1}(k_p+k_i) − b₀k_p, m₂ = (a₂ − a₁)+b₂(k_p+k_i) − b₁k_p, m₁ = (_{a3 − a2)+b3}(k_p+k_i) − b₂k_p, m₀ = −a₃ − b₃k_p. Then multiplying both sides of Equation (12) by (w − 1)⁴, such that, ${(w-1)}^{4}D(w)=n_4{w}^4+n_3{w}^3+n_2{w}^2+n_{1}w+n_0=0$, that is,

$\begin{array}{ll}{D}_1(w)=n_4{w}^4+n_3{w}^3+n_2{w}^2+n_{1}w+n_0=0& (13)\end{array}$

where n₄ = m₀+m₁+m₂+m₃+m₄, n₃ = −4m₀−2m₁+2m₃+4m₄, n₂ = 6m₀−2m₂+6m₄, n₁ = −4m₀ + 2m₁ − 2m₃ + 4m₄, n₀ = m₀ − m₁ + m₂ − m₃ + m₄.

The stability of this system was equivalent to the following conditions:

$\begin{array}{ll}{n}_i>0(i=0,1,2,3,4),n_3{n}_2>n_4{n}_1,n_3{n}_2{n}_1>n_4{n}_1^2+n_3^2{n}_0& (14)\end{array}$

Combining Equations (11)–(13), n_icould also be described as a function of k_p and k_i, and to ensure that all conditions in Equation (14) were satisfied, we chose k_p = 0.80, k_i = 0.05.

Closed-Loop Frequency Modulation

Considering the established physiological responses to different pulse repetition frequencies of DBS (Birdno and Grill, 2008), we constrained the calculated stimulation frequency to between 5 and 200 Hz. When the calculated frequency was larger than 200 Hz, it was set to 200 Hz; when the calculated frequency was <5 Hz, it was set to 5 Hz.

$\begin{array}{ll}u\left(k\right)=\{\begin{array}{cc}5& u\left(k\right)<5\hfill \\ u\left(k-1\right)+k_p\left[e\left(k\right)-e\left(k-1\right)\right]+k_{i}e\left(k\right)& 5\le u\left(k\right)\le 200\hfill \\ 200& u\left(k\right)>200\hfill \end{array}& (15)\end{array}$

The stimulation frequency was calculated using the PI controller, which required knowledge of the beta power at the kth and (k-1)th time points. The beta power of the kth time point was calculated from (t−0.1) s to t s, the beta power of the (k-1)th time point was calculated from (t₁-0.1) s to t₁ s. The time difference between t and t₁ was 0.008 s. Note this was not the time step for the controller to update the DBS frequency, and the controller updated the DBS frequency only after the former interpulse interval ended.

Results

Prediction Performance of the CAR Model

The performance of the CAR model during the model training process is shown in Figure 4A. The correlation coefficient between the actual and estimated data in the model training process was r(y, y_e) = 0.84. In addition, we generated different sequences of random stimulation frequencies, and delivered the corresponding stimulation signals to the network model to calculate the resulting sequences of beta power. The same sequences of stimulation frequencies were also delivered to the trained CAR model. The prediction performance of the trained CAR model on two example data sequences is shown in Figures 4B,C. In this testing phase, the correlation coefficient between the two outputs were r(y, y_e) = 0.82 and r(y, y_e) = 0.80. Thus, the prediction accuracy of the CAR model was ~80%.

FIGURE 4

Figure 4. The prediction performance of the CAR model during model training (A) and testing (B,C). The datasets used to train and test the CAR model were generated as described in section Computational Model of the Cortex-Basal Ganglia-Thalamus Network. The black line represented the beta power calculated from the original CTx-BG-Th network model, and the red line represented the beta power data predicted by the identified CAR model.

To create a quantitative comparator for the prediction accuracy of the identified CAR model, we delivered an identical test stimulation signal to the CTx-BG-Th network model five times. The mean correlation coefficient among any two output datasets was 0.95. Since the CTx-BG-Th was highly non-linear, while the structure of the CAR model presented here was linear, the difference between 0.95 and 0.8 may reflect the unmodeled non-linear dynamics between the stimulation frequency and the beta power. However, since our aim in identifying the CAR model was as a tool to design the PI controller, the 80% accuracy was deemed sufficient.

Tracking of Constant Beta Power

The relationship between the DBS pulse repetition frequency and the beta power of GPi model neuron spike times in the CTx-BG-Th model is shown in Figure 5. The beta band power in the healthy and PD states of the CTx-BG-Th model were 162 and 222.5, respectively. Similar to the effects of DBS frequency on motor symptoms (Birdno and Grill, 2008), reductions in beta band oscillatory activity were observed only for higher frequencies of DBS. The target beta power was selected to be 110, which was approximately the value generated by DBS at 115 Hz. When the stimulation frequency was larger than 100 Hz, the variations of beta power with changes in frequency were quite small, and the selection of a specific beta power target level had no particular impact on the results.

FIGURE 5

Figure 5. The relationship between DBS frequency and the beta band power of GPi model neuron spike times. Standard error bars are shown for 50 trials. The dotted line labels the 110 target beta power value.

The spectrograms of the spike times from model GPi, GPe, and STN neurons in the parkinsonian condition, during 115 Hz DBS, and during CL DBS are shown in Figure 6. Under the parkinsonian condition, the model neurons in these three nuclei exhibited oscillatory activity around 20 Hz. During the 115 Hz DBS and CL DBS cases, stimulation began at t = 2 s, after which the oscillatory activity rapidly diminished. CL DBS produced intermittent oscillatory activity in model STN neurons in the low frequency band (3–12 Hz), which was 6.2 times larger than the lower frequency power present during open loop DBS (OL DBS) at 115 Hz. The dynamic sequence of stimulation frequencies during CL DBS (Figure 7A) exhibited peaks in the power spectrum both around 115 Hz and between 3 and 12 Hz. DBS (Figure 7D). The stimulation signal power 3–12 Hz generated oscillatory activity in model STN neurons in the low frequency band that was larger than during 115 Hz OL DBS. Thus, although both stimulation methods reduced the power in the beta band, they may act through different mechanisms.

FIGURE 6

Figure 6. Spectrograms of the spike times from model GPi, GPe, and STN neurons in the normal (A), parkinsonian condition (B), during 115 Hz DBS (C), and during CL DBS (D). In the parkinsonian condition, all neurons exhibited excessive oscillatory activity compared with the normal condition. The 115 Hz DBS and CL DBS began at t = 2 s, and greatly reduced the beta band oscillatory activity.

FIGURE 7

Figure 7. Variations of DBS frequency (A) and beta power of spike times from model GPi neurons during CL DBS (B) and regular 115 Hz DBS (C). The stimulation signal began at t = 2 s, and the beta power converged to the target of 110 at t = 2.66 s in (B) and at t = 2.09 s in (C). (D) The corresponding spectral power of the stimulation sequence for this CL DBS example. (E) The relationship between the initial frequency of the CL controller and response time in the CL DBS system. The initial frequency of CL DBS was randomly selected from a uniform distribution between 5 and 200 Hz, and the corresponding response time was calculated across 200 trials.

The variations of DBS frequency and the corresponding changes in beta band power in model GPi neurons during CL DBS are shown in Figures 7A,B, respectively. The stimulation began at t = 2 s, the initial stimulation frequency was set to 5 Hz, and the CL DBS system calculated the subsequent frequencies automatically to drive the beta band power to the target of 110. The mean stimulation frequency from 2 to 30 s was 118.7 Hz, and the mean beta power from 2 to 30 s was 114.3, while the mean beta power during OL DBS from 2 to 30 s was 111.3 (Figure 7C). Compared to OL DBS at 115 Hz, the CL DBS controller generated a wider distribution of power in the stimulation frequency sequence (Figure 7D), and the power present in the low frequency band of the stimulation signal generated low frequency power in STN model neurons during CL DBS (Figure 6D). The response time was shorter for open loop 115 Hz DBS (0.09 s) than for CL DBS (0.66 s); however, the response time was strongly dependent on the initial value of frequency during CL DBS (Figure 7E). As the initial frequency was increased the response time decreased, and when the initial frequency was ≥60 Hz, the response time for CL DBS was <0.15 s.

To assess the robustness of the PI controller, we changed the target beta power while keeping the PI parameters unchanged (Figure 8A). Figures 8B–E illustrate the stimulation frequency and beta power variation when the desired beta power was 140 and 180, respectively. When the target beta power was 140, the response time was 0.89 s, and the mean stimulation frequency was 74 Hz. When the target beta power was 150, the response time was 1.15 s, and the mean stimulation frequency was 56 Hz. When the target beta power was larger than 160, the tracking performance declined. Thus, as the desired beta power was larger, the convergence time of GPi beta power became longer. When the target beta power was set to 60 (i.e., a value not achievable with OL DBS, Figure 5), the calculated stimulation frequency varied between 155 and 200 Hz (mean = 177.8 Hz), the mean beta power from 2 to 30 s was 82.3, and with OL DBS at 177.8 Hz, the mean beta power was 87.91.

FIGURE 8

Figure 8. Performance of the PI controller across different levels of target beta power. (A) The dotted line represented the value of desired beta power, and the solid line represented the value of controlled beta power; standard error bars are shown for 50 trials. The variation of DBS frequency and beta power of model GPi neuron spike times when the desired beta power was 140 (B,C) and 180 (D,E), respectively. The red solid line in (C,E) are the desired beta power value. The red dotted line in (C,E) indicate the time when the controlled beta power reached the desired beta power.

Tracking of Dynamic Changes in Target Beta Power

Beta power in the BG exhibits dynamic changes prior to and during voluntary movement and a fixed target beta power may not be appropriate for functional control of DBS.

Therefore, we tested the performance of the control system with time-varying beta power. According to Figure 5, when the stimulation frequency of regular DBS increased from 50 to 130 Hz, the GPi beta band power decreased gradually from 220 to 110, and the beta band power tended to saturate at DBS frequencies larger than 130 Hz. Therefore, the target values randomly selected from a uniform distribution between 110 and 220. The duration of the target value varied from 10 to 1 s (Figures 9a1–h1). The correlation coefficient between the target beta power and actual beta power between 3.5 and 30 s was calculated, as this mitigated the confounding effects of the initial stimulation frequency. The tracking performance of the CL DBS declined with the duration of the target value, and when the duration was 10, 5, 2, 1, and 0.5 s, the correlation coefficients were 0.83, 0.82, 0.71, 0.69, and 0.49, respectively. Sinusoidal trajectories of target beta power with frequencies ranging from 0.05 to 1 Hz were also tested. The BG model can generate beta power between 90 and 200 during regular DBS, and the minimum and maximum amplitude of the target sinusoidal trajectories were therefore set to be 90 and 200, respectively. The tracking performance and variation of stimulation frequency of the CL DBS system are shown in Figures 9a2–h2, and the correlation coefficient between the actual beta power and the target trajectory was used to quantify the tracking accuracy. The tracking performance declined with the increase in target sinusoidal frequency, and the correlation coefficient was 0.85, 0.65, 0.49, and 0.17 for sinusoidal frequencies of 0.05, 0.3, 0.5, and 1 Hz, respectively.

FIGURE 9

Figure 9. (a1–h1) The performance of the PI controller tracking dynamic changes in target beta power and the corresponding variation in stimulation frequency by the CL DBS system. The duration of the target beta power was 10 s (a1), 5 s (c1), 2 s (e1) and 1 s (g1), respectively. The black line represents the desired beta power, and the red line represents the actual beta power during CL DBS. (b1,d1,f1,h1) show the respective variations in stimulation frequency during CL DBS. When the duration of target beta power declined from 10 to 0.5 s, the correlation coefficients between the desired and actual beta power were 0.83, 0.82, 0.71, 0.69, 0.49, respectively. (a2–h2) The performance of the PI controller tracking sinusoidal trajectories at different frequencies and the associated DBS frequencies determined by the CL controller. The black line in (a2,c2,e2,g2) represented the desired beta power, and the red line represented the actual beta power. When the frequency of target sinusoidal trajectories increased from 0.05 to 1 Hz, the correlation coefficients between the desired and actual beta power were 0.85, 0.65, 0.49, and 0.17, respectively. Sinusoidal trajectories of target beta power with frequencies ranging from 0.05 to 1 Hz were also tested. The BG model can generate beta power between 90 and 200 during regular DBS, and the minimum and maximum amplitude of the target sinusoidal trajectories were therefore set to be 90 and 200, respectively. The tracking performance and variation of stimulation frequency of the CL DBS system are shown in (a2–h2), and the correlation coefficient between the actual beta power and the target trajectory was used to quantify the tracking accuracy. The tracking performance declined with the increase in target sinusoidal frequency, and the correlation coefficient was 0.85, 0.65, 0.49, and 0.17 for sinusoidal frequencies of 0.05, 0.3, 0.5, and 1 Hz, respectively.

Discussion

Beta band oscillatory activity in the BG is correlated with motor symptoms in PD and may be a suitable biomarker for CL DBS in PD (Little and Brown, 2012; Hoang et al., 2017). For example, Arlotti et al. (2016b) and Little et al. (2013) used the beta oscillation amplitude to control the on time of DBS. DBS was delivered only when the beta-band oscillation amplitude was larger than a pre-set threshold, which reduced energy consumption compared to continuous DBS, while increasing the therapeutic effects on motor symptoms. Subsequently, Dan et al. demonstrated that this approach was also effective in a PD patient with chronically implanted DBS (Piña-Fuentes et al., 2017). In complementary modeling studies, Grant and Lowery designed a CL DBS system to modulate the amplitude of DBS based on beta band oscillations of LFPs, where the coupling strength within the cortico-basal ganglia network was altered to illustrate the ability of CL DBS to respond to changes in network activity (Grant and Lowery, 2013).

However, beta oscillatory activity exhibits dynamic changes (desynchronization) during movement, and Johnson et al. found that a constant beta set point may not be suitable as CL DBS performed poorly during reaching behavior (Johnson et al., 2016). Therefore, if beta power is to be used as a feedback control signal, a constant reference value might not be appropriate. In more recent studies, DBS voltage was adjusted proportionally to the STN LFP beta power, and this adaptive DBS reduced side effects compared to traditional open-loop DBS (Rosa et al., 2015; Arlotti et al., 2018). In another alternative to simply reducing oscillatory activity below a fixed threshold, Santaniello et al. automatically adjusted the stimulation voltage in a mathematical model to match a desired profile of oscillatory neuronal activity (Santaniello et al., 2011). During go/no-go voluntary movements, dynamic changes in beta band power occur at 0.3–1Hz (Sanes and Donoghue, 1993; Zaepffel et al., 2013). The proposed controller could track dynamic changes slower than 1 Hz, and thus such an approach may account for the dynamic changes in beta oscillatory power that occur during movement. Instead of simply switching the stimulation on and off, or adjusting the stimulation amplitude, the controller regulated the stimulation frequency in real time. If the variation in beta band power during a wide range of movements was known a priori, such a closed-loop system that modulates stimulation frequency to track dynamic beta oscillatory activity may facilitate a wide range of individual patient motor behaviors.

The proposed closed-loop stimulation algorithm was simulated using a validated CTx-BG-Th model (Kumaravelu et al., 2016). There are several other potential models of the network effects of DBS, which might be used for development and evaluation of closed-loop controllers. Hahn and McIntyre developed a network model of the effects of DBS in the STN of the parkinsonian non-human primate, and demonstrated that effective DBS suppressed burst activity in the GPi (Hahn and McIntyre, 2010). Subsequently, Holt and Netoff implemented a mean field version of this model and analyzed the effects of different frequencies of DBS (Holt and Netoff, 2014). Similarly, Santaniello et al. (2015) implemented a network model of the effects of STN DBS in the parkinsonian non-human primate and demonstrated the importance of both antidromic and orthodromic activation. We selected the Kumaravelu et al. network model because it replicated a wide range of electrophysiological data from the unilateral 6-OHDA lesioned rat model of PD (Kumaravelu et al., 2016) thereby facilitating subsequent in vivo evaluation of the controller.

The proposed CL DBS controller was successful at regulating the beta oscillatory activity of spike times of model GPi neurons to track different beta reference values. The stimulation frequency was automatically calculated by the PI controller, and PI parameters were calculated using stability analysis of the system rather than trial-and-error adjustment (Gorzelic et al., 2013). However, there were several potential limitations of the proposed CL DBS method. The identified linear CAR model described only 80% of the relationship between the stimulation frequency and the beta power. Therefore, although the PI controller was robust to changes in the reference beta power, the dynamic changes in beta power could be tracked well only at frequencies of ≤ 1 Hz. When the target beta power changed faster than 1 Hz, the tracking error increased, likely as a result of the unmodeled dynamics. In subsequent trial-and-error tuning, it appeared that the best PI controller parameters were different for different beta power targets. Thus, adaptive controllers that modulate the PI controller parameters with the variation of target beta power may improve the tracking performance for dynamic reference signals. The CTx-BG-Th network was highly non-linear, and performance might also be improved using a non-linear controller. The beta oscillatory power was selected as the biomarker in this study, however, other biomarkers such as the spike time entropy (Dorval et al., 2008) and phase amplitude coupling (de Hemptinne et al., 2015) are also correlated with parkinsonian symptoms, and might be suitable feedback control signals. The application of other biomarkers or multiple biomarkers in the design of closed-loop stimulation for PD is worth exploring (Hoang et al., 2017). The controller regulated the stimulation frequency, but the effects of DBS are also dependent on the pulse amplitude, pulse duration, and stimulation pattern (Kuncel and Grill, 2004; Grill, 2018). Further, Holt et al. demonstrated that the effects of burst DBS in a network model of the basal ganglia (Hahn and McIntyre, 2010) were strongly dependent on timing relative to the phase of oscillatory activity (Holt et al., 2016).

We demonstrated successful tracking of different dynamic beta power reference signals, and the simulated dynamic targets could represent different movements of PD patients. Thus, an important challenge to implement the proposed CL DBS approach experimentally or clinically is to determine the relationship between reference beta oscillation power and the movement. In addition to real-time electrophysiological recording, movement sensors might also be useful to establish the dynamic reference signal.

Conclusion

CL DBS was proposed to reduce energy consumption and alleviate side effects compared to continuous fixed-parameter DBS. This requires design of a suitable closed-loop system that can account for dynamic changes in the feedback signal that occur during voluntary movement. We used the beta oscillatory power of GPi model neuron spike times as a biomarker of model state, and used a PI controller to calculate the DBS frequency according to dynamic variations in the beta power. This closed-loop adjustment of stimulation frequency approach was tested in a computational model of the CTx-BG-Th network and was able to track constant as well as dynamic beta oscillatory activity.

Data Availability

The raw data supporting the conclusions of this manuscript will be made available by the authors, without undue reservation, to any qualified researcher.

Author Contributions

All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.

Funding

This work was supported by the US National Institutes of Health under Grant R37 NS040894 and Grant UH3 NS103468 and by the National Natural Science Foundation of China under Grant 61801273 and Grant 61701336.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

References

Arlotti, M., Marceglia, S., Foffani, G., Volkmann, J., Lozano, A. M., Moro, E., et al. (2018). Eight-hours adaptive deep brain stimulation in patients with Parkinson disease. Neurology 90, e971–e6. doi: 10.1212/WNL.0000000000005121

PubMed Abstract | CrossRef Full Text | Google Scholar

Arlotti, M., Rosa, M., Marceglia, S., Barbieri, S., and Priori, A. (2016a). The adaptive deep brain stimulation challenge. Parkinson. Relat. Disord. 28, 12–17. doi: 10.1016/j.parkreldis.2016.03.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Arlotti, M., Rossi, L., Rosa, M., Marceglia, S., and Priori, A. (2016b). An external portable device for adaptive deep brain stimulation (aDBS) clinical research in advanced Parkinson's disease. Med. Eng. Phys. 38, 498–505. doi: 10.1016/j.medengphy.2016.02.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Aström, K. J., and Hägglund, T. (1995). PID Controllers: Theory Design and Tuning, 2nd Edn. Research Triangle Park, NC: Instrument Society of America.

Google Scholar

Birdno, M. J., and Grill, W. M. (2008). Mechanisms of deep brain stimulation in movement disorders as revealed by changes in stimulus frequency. Neurotherapeutics 5, 14–25. doi: 10.1016/j.nurt.2007.10.067

PubMed Abstract | CrossRef Full Text | Google Scholar

Brittain, J. S., and Brown, P. (2014). Oscillations and the basal ganglia: motor control and beyond. Neuroimage 85, 637–647. doi: 10.1016/j.neuroimage.2013.05.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Brocker, D. T., Swan, B. D., Turner, D. A., Gross, R. E., Tatter, S. B., Miller Koop, M. M., et al. (2013). Improved efficacy of temporally non-regular deep brain stimulation in Parkinson's disease. Exp. Neurol. 239, 60–67. doi: 10.1016/j.expneurol.2012.09.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Carron, R., Chaillet, A., Filipchuk, A., Pasillas-Lépine, W., and Hammond, C. (2013). Closing the loop of deep brain stimulation. Front. Syst. Neurosci. 7:112. doi: 10.3389/fnsys.2013.00112

PubMed Abstract | CrossRef Full Text | Google Scholar

Cyron, D. (2016). Mental side effects of deep brain stimulation (DBS) for movement disorders: the futility of denial. Front. Integr. Neurosci. 10:17. doi: 10.3389/fnint.2016.00017

PubMed Abstract | CrossRef Full Text | Google Scholar

Davie, C. A. (2008). A review of Parkinson's disease. Br. Med. Bull. 86, 109–127. doi: 10.1093/bmb/ldn013

PubMed Abstract | CrossRef Full Text | Google Scholar

de Hemptinne, C., Swann, N. C., Ostrem, J. L., Ryapolova-Webb, E. S., San Luciano, M., Galifianakis, N. B., et al. (2015). Therapeutic deep brain stimulation reduces cortical phase-amplitude coupling in Parkinson's disease. Nat. Neurosci. 18, 779–786. doi: 10.1038/nn.3997

PubMed Abstract | CrossRef Full Text | Google Scholar

Deuschl, G., Herzog, J., Kleiner-Fisman, G., Kubu, C., Lozano, A. M., Lyons, K. E., et al. (2006). Deep brain stimulation: postoperative issues. Mov. Disord. 21, S219–S237. doi: 10.1002/mds.20957

PubMed Abstract | CrossRef Full Text | Google Scholar

Dorval, A. D., Russo, G. S., Hashimoto, T., Xu, W., Grill, W. M., and Vitek, J. L. (2008). Deep brain stimulation reduces neuronal entropy in the MPTP-primate model of Parkinson's disease. J. Neurophysiol. 100, 2807–2818. doi: 10.1152/jn.90763.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Goldberg, J. A., Rokni, U., Boraud, T., Vaadia, E., and Bergman, H. (2004). Spike synchronization in the cortex/basal-ganglia networks of Parkinsonian primates reflects global dynamics of the local field potentials. J. Neurosci. 24, 6003–6010. doi: 10.1523/JNEUROSCI.4848-03.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Gopal, M. (2002). Control Systems: Principles and Design, 2nd Edn. New Delhi: Tata McGraw-Hill Education.

Google Scholar

Gorzelic, P., Schiff, S. J., and Sinha, A. (2013). Model-based rational feedback controller design for closed-loop deep brain stimulation of Parkinson's disease. J. Neural Eng. 10:26016. doi: 10.1088/1741-2560/10/2/026016

PubMed Abstract | CrossRef Full Text | Google Scholar

Grant, P. F., and Lowery, M. M. (2013). Simulation of cortico-basal ganglia oscillations and their suppression by closed loop deep brain stimulation. IEEE Trans. Neural Syst. Rehabil. Eng. 21, 584–594. doi: 10.1109/TNSRE.2012.2202403

PubMed Abstract | CrossRef Full Text | Google Scholar

Grill, W. M. (2018). Temporal pattern of electrical stimulation is a new dimension of therapeutic innovation. Curr. Opinion. Biomed. Eng. 8, 1–6. doi: 10.1016/j.cobme.2018.08.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Hahn, P. J., and McIntyre, C. C. (2010). Modeling shifts in the rate and pattern of subthalamopallidal network activity during deep brain stimulation. J. Comput. Neurosci. 28, 425–441. doi: 10.1007/s10827-010-0225-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Hebb, A. O., Zhang, J. J., Mahoor, M. H., Tsiokos, C., Matlack, C., Chizeck, H. J., et al. (2014). Creating the feedback loop: closed-loop neurostimulation. Neurosurg. Clin. N. Am. 25, 187–204. doi: 10.1016/j.nec.2013.08.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Hoang, K. B., Cassar, I. R., Grill, W. M., and Turner, D. A. (2017). Biomarkers and stimulation algorithms for adaptive brain stimulation. Front. Neurosci. 11:564. doi: 10.3389/fnins.2017.00564

PubMed Abstract | CrossRef Full Text | Google Scholar

Holt, A. B., and Netoff, T. I. (2014). Origins and suppression of oscillations in a computational model of Parkinson's disease. J. Comput. Neurosci. 37, 505–521. doi: 10.1007/s10827-014-0523-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Holt, A. B., Wilson, D., Shinn, M., Moehlis, J., and Netoff, T. I. (2016). Phasic burst stimulation: a closed-loop approach to tuning deep brain stimulation parameters for Parkinson's disease. PLoS Comput. Biol. 12:e1005011. doi: 10.1371/journal.pcbi.1005011

PubMed Abstract | CrossRef Full Text | Google Scholar

Jankovic, J. (2008). Parkinson's disease: clinical features and diagnosis. J. Neurol. Neurosurg. Psychiatry 79, 368–376. doi: 10.1136/jnnp.2007.131045

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, N., and Brown, P. (2011). New insights into the relationship between dopamine, beta oscillations and motor function. Trends Neurosci. 34, 611–618. doi: 10.1016/j.tins.2011.09.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, L. A., Nebeck, S. D., Muralidharan, A., Johnson, M. D., Baker, K. B., and Vitek, J. L. (2016). Closed-loop deep brain stimulation effects on Parkinsonian motor symptoms in a non-human primate-is beta enough? Brain Stimul. 9, 892–896. doi: 10.1016/j.brs.2016.06.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Kita, H., and Kita, T. (2011). Cortical stimulation evokes abnormal responses in the dopamine-depleted rat basal ganglia. J. Neurosci. 31, 10311–10322. doi: 10.1523/JNEUROSCI.0915-11.2011

PubMed Abstract | CrossRef Full Text | Google Scholar

Krack, P., Batir, A., Van Blercom, N., Chabardes, S., Fraix, V., Ardouin, C., et al. (2003). Five-year follow-up of bilateral stimulation of the subthalamic nucleus in advanced Parkinson's disease. N. Engl. J. Med. 349, 1925–1934. doi: 10.1056/NEJMoa035275

PubMed Abstract | CrossRef Full Text | Google Scholar

Kühn, A. A., Trottenberg, T., Kivi, A., Kupsch, A., Schneider, G. H., and Brown, P. (2005). The relationship between local field potential and neuronal discharge in the subthalamic nucleus of patients with Parkinson's disease. Exp. Neurol. 194, 212–220. doi: 10.1016/j.expneurol.2005.02.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumaravelu, K., Brocker, D. T., and Grill, W. M. (2016). A biophysical model of the cortex-basal ganglia-thalamus network in the 6-OHDA lesioned rat model of Parkinson's disease. J. Comput. Neurosci. 40, 207–229. doi: 10.1007/s10827-016-0593-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuncel, A. M., and Grill, W. M. (2004). Selection of stimulus parameters for deep brain stimulation. Clin. Neurophysiol. 115, 2431–2441. doi: 10.1016/j.clinph.2004.05.031

PubMed Abstract | CrossRef Full Text | Google Scholar

Leventhal, D. K., Gage, G. J., Schmidt, R., Pettibone, J. R., Case, A. C., and Berke, J. D. (2012). Basal ganglia beta oscillations accompany cue utilization. Neuron 73, 523–536. doi: 10.1016/j.neuron.2011.11.032

PubMed Abstract | CrossRef Full Text | Google Scholar

Levy, R., Ashby, P., Hutchison, W. D., Lang, A. E., Lozano, A. M., and Dostrovsky, J. O. (2002). Dependence of subthalamic nucleus oscillations on movement and dopamine in Parkinson's disease. Brain 125, 1196–1209. doi: 10.1093/brain/awf128

PubMed Abstract | CrossRef Full Text | Google Scholar

Little, S., and Brown, P. (2012). What brain signals are suitable for feedback control of deep brain stimulation in Parkinson's disease? Ann. N. Y. Acad. Sci. 1265, 9–24. doi: 10.1111/j.1749-6632.2012.06650.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Little, S., Pogosyan, A., Neal, S., Zavala, B., Zrinzo, L., Hariz, M., et al. (2013). Adaptive deep brain stimulation in advanced Parkinson disease. Ann. Neurol. 74, 449–457. doi: 10.1002/ana.23951

PubMed Abstract | CrossRef Full Text | Google Scholar

Ljung, L. (1999). System Identification Theory for User, 2nd Edn. Englewood Cliffs, NJ: Prentice Hall.

Google Scholar

Mallet, N., Pogosyan, A., Márton, L. F., Bolam, J. P., Brown, P., and Magill, P. J. (2008). Parkinsonian beta oscillations in the external globus pallidus and their relationship with subthalamic nucleus activity. J. Neurosci. 28, 14245–14258. doi: 10.1523/JNEUROSCI.4199-08.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Massano, J., and Garrett, C. (2012). Deep brain stimulation and cognitive decline in Parkinson's disease: a clinical review. Front. Neurol. 3:66. doi: 10.3389/fneur.2012.00066

PubMed Abstract | CrossRef Full Text | Google Scholar

McConnell, G. C., So, R. Q., Hilliard, J. D., Lopomo, P., and Grill, W. M. (2012). Effective deep brain stimulation suppresses low-frequency network oscillations in the basal ganglia by regularizing neural firing patterns. J. Neurosci. 32, 15657–15668. doi: 10.1523/JNEUROSCI.2824-12.2012

PubMed Abstract | CrossRef Full Text | Google Scholar

McQuarrie, A. D. R., and Chih-Ling, T. (1998). Regression and Time Series Model Selection. Hackensack, NJ: World Scientific. doi: 10.1142/3573

CrossRef Full Text | Google Scholar

Odekerken, V. J., Boel, J. A., Schmand, B. A., de Haan, R. J., Figee, M., van den Munckhof, P., et al. (2016). GPi vs. STN deep brain stimulation for Parkinson disease: three-year follow-up. Neurology 86, 755–761. doi: 10.1212/WNL.0000000000002401

PubMed Abstract | CrossRef Full Text | Google Scholar

Okun, M. S., and Foote, K. D. (2010). Parkinson's disease DBS: what, when, who and why? The time has come to tailor DBS targets. Expert. Rev. Neurother. 10, 1847–1857. doi: 10.1586/ern.10.156

PubMed Abstract | CrossRef Full Text | Google Scholar

Parastarfeizabadi, M., and Kouzani, A. Z. (2017). Advances in closed-loop deep brain stimulation devices. J. Neuroeng. Rehabil. 14:79. doi: 10.1186/s12984-017-0295-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Piña-Fuentes, D., Little, S., Oterdoom, M., Oterdoom, M., Neal, S., Pogosyan, A., et al. (2017). Adaptive DBS in a Parkinson's patient with chronically implanted DBS: a proof of principle. Mov. Disord. 32, 1253–1254. doi: 10.1002/mds.26959

PubMed Abstract | CrossRef Full Text | Google Scholar

Rodriguez-Oroz, M. C., Obeso, J. A., Lang, A. E., Houeto, J. L., Pollak, P., and Rehncrona, S. (2005). Bilateral deep brain stimulation in Parkinson's disease: a multicentre study with 4 years follow-up. Brain 128, 2240–2249. doi: 10.1093/brain/awh571

PubMed Abstract | CrossRef Full Text | Google Scholar

Rosa, M., Arlotti, M., Ardolino, G., Cogiamanian, F., Marceglia, S., Di Fonzo, A., et al. (2015). Adaptive deep brain stimulation in a freely moving parkinsonian patient. Mov. Disord. 30, 1003–1005. doi: 10.1002/mds.26241

PubMed Abstract | CrossRef Full Text | Google Scholar

Rosin, B., Slovik, M., Mitelman, R., Rivlin-Etzion, M., Haber, S. N., Israel, Z., et al. (2011). Closed-loop deep brain stimulation is superior in ameliorating Parkinsonism. Neuron 72, 370–384. doi: 10.1016/j.neuron.2011.08.023

PubMed Abstract | CrossRef Full Text | Google Scholar

Rossi, P. J., Gunduz, A., Judy, J., Wilson, L., Machado, A., and Giordano, J. J. (2016). Proceedings of the third annual deep brain stimulation think tank: a review of emerging issues and technologies. Front. Neurosci. 10:119. doi: 10.3389/fnins.2016.00119

PubMed Abstract | CrossRef Full Text | Google Scholar

Sanes, J. N., and Donoghue, J. P. (1993). Oscillations in local field potentials of the primate motor cortex during voluntary movement. Proc. Natl. Acad. Sci. U.S.A. 90, 4470–4474. doi: 10.1073/pnas.90.10.4470

PubMed Abstract | CrossRef Full Text | Google Scholar

Santaniello, S., Fiengo, G., Glielmo, L., and Grill, W. M. (2011). Closed-loop control of deep brain stimulation: a simulation study. IEEE Trans. Neural Syst. Rehabil. Eng. 19, 15–24. doi: 10.1109/TNSRE.2010.2081377

PubMed Abstract | CrossRef Full Text | Google Scholar

Santaniello, S., McCarthy, M. M., Montgomery, E. B. Jr., Gale, J. T., Kopell, N., and Sarma, S. V. (2015). Therapeutic mechanisms of high-frequency stimulation in Parkinson's disease and neural restoration via loop-based reinforcement. Proc. Natl. Acad. Sci. U.S.A. 112, E586–595. doi: 10.1073/pnas.1406549111

PubMed Abstract | CrossRef Full Text | Google Scholar

So, R. Q., McConnell, G. C., August, A. T., and Grill, W. M. (2012). Characterizing effects of subthalamic nucleus deep brain stimulation on methamphetamine-induced circling behavior in hemi-Parkinsonian rats. IEEE Trans. Neural Syst. Rehabil. Eng. 20, 626–635. doi: 10.1109/TNSRE.2012.2197761

PubMed Abstract | CrossRef Full Text | Google Scholar

Stein, E., and Bar-Gad, I. (2013). β oscillations in the cortico-basal ganglia loop during parkinsonism. Exp. Neurol. 245, 52–59. doi: 10.1016/j.expneurol.2012.07.023

PubMed Abstract | CrossRef Full Text | Google Scholar

Zaepffel, M., Trachel, R., Kilavik, B. E., and Brochier, T. (2013). Modulations of EEG beta power during planning and execution of grasping movements. PLoS ONE 8:e60060. doi: 10.1371/journal.pone.0060060

PubMed Abstract | CrossRef Full Text | Google Scholar

Zaidel, A., Spivak, A., Grieb, B., Bergman, H., and Israel, Z. (2010). Subthalamic span of oscillations predicts deep brain stimulation efficacy for patients with Parkinson's disease. Brain 133, 2007–2021. doi: 10.1093/brain/awq144

PubMed Abstract | CrossRef Full Text | Google Scholar