|Home | About | Journals | Submit | Contact Us | Français|
High-frequency (HF) deep brain stimulation (DBS) is the gold standard for the treatment of medically refractory movement disorders like Parkinson’s disease, essential tremor, and dystonia, with a significant potential for application to other neurological diseases. The standard setup of HF DBS utilizes an open-loop stimulation protocol, where a permanent HF electrical pulse train is administered to the brain target areas irrespectively of the ongoing neuronal dynamics. Recent experimental and clinical studies demonstrate that a closed-loop, adaptive DBS might be superior to the open-loop setup. We here combine the notion of the adaptive high-frequency stimulation approach, that aims at delivering stimulation adapted to the extent of appropriately detected biomarkers, with specifically desynchronizing stimulation protocols. To this end, we extend the delayed feedback stimulation methods, which are intrinsically closed-loop techniques and specifically designed to desynchronize abnormal neuronal synchronization, to pulsatile electrical brain stimulation. We show that permanent pulsatile high-frequency stimulation subjected to an amplitude modulation by linear or nonlinear delayed feedback methods can effectively and robustly desynchronize a STN-GPe network of model neurons and suggest this approach for desynchronizing closed-loop DBS.
Synchronization is a fundamental natural phenomenon in interacting networks [1–4]. Synchronization plays a crucial role in the human brain in, e.g., processing of sensory information , motor control , and cognitive function . However, excessive pathological neuronal synchrony may severely impair brain function and is a hallmark of several neurological disorders, such as Parkinson’s disease (PD) [8, 9], essential tremor , epilepsy , and tinnitus [12–14]. The standard therapy for the treatment of medically refractory PD is high frequency (HF) deep brain stimulation (DBS), where electrical HF pulse trains are administered at frequencies >100 Hz via depth electrodes chronically implanted in target areas such as the thalamic ventralis intermedius (VIM) nucleus, the subthalamic nucleus (STN), or the globus pallidus (GP) [15–18]. HF DBS has been developed empirically, and the clinical and electrophysiological mechanisms of the symptom suppression by HF DBS are still a matter of intensive research [17, 19, 20]. A large number of studies are devoted to an improvement of the therapeutic effects of HF DBS by appropriate calibration of the stimulation parameters such as stimulation frequency and intensity, the width and shape of the stimulation pulses, spatial spread and localization of the stimulation current in the neuronal tissue, as well as selection of appropriate stimulation targets, etc. [15–17, 21–25]. A key aspect for further improvement of DBS is the reduction of side effects. HF DBS may not only cause side effects by the spread of electrical current outside of the target region, but also by chronic stimulation of the target itself as well as due to functional disconnection of the stimulated structure [26–29]. Hence, it is crucial to reduce the integral current required.
In contrast to the standard open-loop HF DBS, the major goal of closed-loop, demand controlled DBS is to stimulate only when necessary and/or to adapt the strength of stimulation to the amount of abnormal neuronal synchrony. Demand-controlled DBS was initially introduced in computational studies with different types of specifically designed desynchronizing stimuli and different types of closed-loop control modes, e.g. demand-controlled timing of stimulus delivery or demand-controlled adaptation of stimulus strength during period stimulus delivery [30–32]. So far demand-controlled DBS was experimentally tested by means of conventional high-frequency stimulation and denoted as adaptive DBS (aDBS) [33–36]. In monkeys rendered parkinsonian with the neurotoxin 1-methyl-4-phenyl-1,2,3,6-tetrahydropyridine (MPTP), closed-loop DBS was tested under acute conditions , where a short stimulation pulse train (7 pulses at 130 Hz) was delivered through a pair of electrodes located in the globus pallidus internal (GPi) with an optimal time delay of 80 ms following the occurrence of an action potential recorded either from the GPi or the primary motor cortex (M1). This type of stimulation caused a strong decrease of the firing rate of pallidal neurons together with a pronounced decrease of the oscillatory neuronal activity at the tremor frequency (4−7 Hz) and at the double tremor frequency (9−15 Hz) along with an amelioration of the MPTP-induced akinesia . In contrast, standard continuous 130 Hz DBS caused a less pronounced decrease of the pallidal firing rate, the oscillatory neuronal activity and the amelioration of the akinesia .
Another study reported on a successful proof of principle of a closed-loop aDBS in PD patients, where the onsets and offsets of HF DBS were triggered by threshold crossings by local field potential (LFP) assessing beta-band STN activity . The onset of the stimulation was delayed by 30 to 40 ms with respect to the threshold crossing by LFP. The average improvement in clinical motor scores in the aDBS condition was significantly better by about 30% despite delivering less than 50% of the stimulation current as compared to the conventional continuous HF DBS (cDBS) condition . Clinical and electrophysiological (suppressing of beta-band LFP oscillations) effects of aDBS were also stronger compared to the intermittent random DBS, where random DBS bursts were not triggered by the LFP threshold crossings. This indicates that the intermittency itself is not a key determinate of the efficacy observed with aDBS, and the real-time feedback from the ongoing LFP oscillations is necessary .
The on-demand DBS was applied for suppression of essential tremor , where the electromyographic (EMG) signal was used to predict the onset of the tremor and initiation of the DBS epoch lasting for a few tens of seconds. Tremor was suppressed by such an intermittent stimulation, which indicated the feasibility of the EMG-based predictive on-off control of DBS in essential tremor patients. The on-demand closed-loop DBS was also tested in patients with intention tremor, where the HF DBS was switched on and off when the power of the recorded EMG activity exceeded or decreased below the EMG trigger threshold, respectively . The on-demand control system triggered the switching on/off of DBS accurately and controlled intention tremor completely.
Closed-loop aDBS also proved to be superior to conventional open-loop DBS (cDBS) in a freely moving patient with PD stimulated with aDBS up to 6 days after DBS electrode implantation . During aDBS the patient experienced a more stable condition with better control of symptoms and dyskinesia as well as improvement of bradykinesia than during cDBS. Instead of the on-off strategy of the papers [34, 37–39], the stimulation voltage was linearly adapted each second based on the LFP beta band power . When the beta-band LFP activity was reduced, the voltage got diminished, in this way avoiding unnecessary hyperstimulation, and the aDBS did not elicit side effects and was well tolerated.
Another direction of research is the model-based development of novel stimulation techniques aiming at specifically counteracting abnormal neuronal synchrony by desynchronization . Initially these techniques were designed to achieve a demand-controlled desynchronizing stimulation [30–32]. However, by taking into account spike timing-dependent plasticity (STDP) [40, 41] in the model neural networks, computationally it was shown that coordinated reset (CR) stimulation , a spatio-temporally patterned, desynchronizing stimulation technique developed in the framework of the model-based approach, reduces the rate of coincidences and, hence, down-regulates abnormal synaptic weights [42, 43]. In this way cumulative, long-lasting, sustained desynchronizing effects can be obtained [42, 43]. Based on numerous modeling studies, that formed a solid foundation for the application of CR stimulation with different stimulation modalities [32, 42, 44–46], these computational predictions were verified both pre-clinically and clinically.
Long-lasting CR-induced desynchronization was confirmed in vitro in rat hippocampal slices rendered epileptic by magnesium withdrawal . CR stimulation caused beneficial therapeutic long-lasting after-effects in parkinsonian MPTP-treated monkeys [48, 49] and in PD patients , whereas the standard HF DBS only has acute effects, i.e., neither clinical [17, 51] nor electrophysiological [52, 53] effects reliably persisted after cessation of stimulation. Based on the same principles of desynchronization-induced unlearning of abnormal synaptic connectivity and neuronal synchrony, i.e. anti-kindling , non-invasive, acoustic CR stimulation was successively applied in a clinical proof of concept study in tinnitus patients . In accordance with the results of modeling studies [44, 45], acoustic CR stimulation can significantly counteract both tinnitus symptoms and the underlying pathological neuronal synchronization by normalizing the effective connectivity and restoring the functional patterns of activity .
In line with the above promising model-based development of desynchronizing CR stimulation, other methods have been developed for the control of abnormal neuronal synchronization. They are based on feedback techniques, where the mean field of synchronized population is measured, preprocessed and fed back as stimulation signal [55–64], or on phase response properties of neurons, where the stimulation signal can be derived from the phase response curve (PRC) [65–67]. The feedback methods can be applied under a variety of conditions and possess an intrinsic demand-controlled character, where the stimulation signal is significantly reduced or even vanishes as soon as desynchronization is achieved. These methods have in common that the stimulation current is a smooth and slowly oscillating signal (although the PRC approach allows for pulsatile stimuli as well), whereas the standard HF DBS and CR stimulation utilize trains of short charge-balanced pulses [22, 68]. The experimental and clinical realization of these feedback methods is a challenging task, first of all, from the technical side, since stimulation signals have to fulfill all safety aspects like charge density limits [16, 69, 70]. These charge density limits actually crucially restrict the applicability of the slow feedback signals. As opposed to the timescales of HF DBS and CR pulses, the feedback signals are slow, so that, if the feedback signal is used directly for stimulation, during the comparably long feedback stimulation periods an irreversible charge deposit can exceed safety limits.
In this paper we resolve this problem and show that the desynchronizing impact and the demand-controlled character of feedback techniques from one side can be combined with the advantages of the charge-balanced property of the HF DBS signal from another side. For this, the amplitude of the HF train of charge-balanced pulses is modulated by the slow feedback signal, which represents a pulsatile feedback stimulation appropriate for electrical DBS. We consider two feedback stimulation methods, linear delayed feedback (LDF) [55, 56] and nonlinear delayed feedback (NDF) [57, 61], and illustrate their efficacy in desynchronizing the model network of STN-GPe neurons suggested in the papers [71, 72]. We administer the delayed feedback stimulations to strongly or weakly and intermittently synchronized neuronal populations and demonstrate that they have a robust desynchronizing impact with respect to parameter variation.
In this way we question the main conclusion of the paper by Dovzhenok et al.  considering the same model and claiming that the delayed feedback DBS is unlikely to be clinically successful because the delayed feedback can desynchronize strongly synchronized neurons, but the same feedback can boost synchronization for intermittently synchronized neurons. This claim is disorienting and requires clarification, because it depreciates the great amount of modeling [55–57, 61, 63, 74] and experimental [75, 76] studies of LDF and NDF as well as recent promising experimental and clinical attempts to advance standard HF DBS to a closed-loop stimulation setup, which attracts the attention of many researchers [34, 35, 37–39, 77–80]. We numerically demonstrate that the pulsatile LDF and NDF are equally effective in counteracting synchronization in strongly as well as weakly and intermittently synchronized neuronal populations. The feedback techniques (LDF, NDF) with smooth stimulation signal are also effectively causing desynchronization in our model. To this end, however, safety requirements for electrical stimulation of neuronal tissue have to be satisfied [16, 69, 70]. Hence, based on our computational study we suggest to use pulsatile LDF and NDF to achieve desynchronizing effects, thereby preserving mandatory safety requirements. In this paper we disprove the argumentation by Dovzhenok et al. . However, more importantly, we suggest a safe solution that paves the way for properly modified DBS feedback techniques to be applied in a first in man study. Our results contribute to the growing field of adaptive closed-loop DBS and show that the suggested approach of pulsatile delayed feedback stimulation could be a possible enhancement of the standard HF DBS, where an abnormal neuronal synchronization can reliably be suppressed by a minimal amount of stimulation.
We consider a network of two neuronal populations, which models the dynamics of STN and GPe neurons. Each cell is modeled by the following system :
where v is a membrane potential of the neuron, IL, IK, INa, IT, ICa, IAHP, Isyn, and Iapp are the corresponding leak, potassium, sodium, low threshold calcium, high threshold calcium, afterhyperpolarisation potassium, synaptic, and external current, respectively. The stimulation current Istim will be defined below. [Ca] is the intracellular concentration of Ca2+ ions, and X = n, h, r are the gating variables.
The following currents from Eq (1) attain the same form for both types of neurons:
whereas current IT is given by different expressions for the excitatory STN cells and for the inhibitory GPe cells:
where b∞(r) = 1/(1 + exp[(r − θb)/σb]) − 1/(1 + exp[−θb/σb]). The functions X∞(v) and τX(v) used in Eq (3) and in the above definition of the currents read
For GPe neurons τr(v) = τr is a constant parameter.
In our study we consider coupled populations of N = 200 STN and 200 GPe neurons. The STN and GPe neuronal ensembles and coupling among them are schematically illustrated in Fig 1. Each STN neuron excites a single GPe neuron, whereas each GPe neuron inhibits three neighboring STN neurons. We also consider periodic boundary conditions. Microscopic models of this type were introduced and investigated in a number of papers [71–73, 81], where STN neurons receive an inhibitory input from GPe neurons and, in turn, give an excitatory output to the GPe network.
The coupling among the neurons is realized via synaptic currents Isyn defined in the following way:
for STN and GPe cells, respectively. j is the index of neurons and summations are taken over all presynaptic neurons. The synaptic weights gS→G = 0.4 nS/μm2 and gG→S = 1.7 nS/μm2 reflect the strength of the coupling from STN neurons to GPe neurons, and in the opposite direction, respectively. We also consider a weaker coupling gG→S = 1.28 nS/μm2. The reversal potentials vS→G = 0 mV and vG→S = −100 mV reflect the excitatory coupling from STN to GPe neurons and inhibitory coupling from GPe to STN, respectively. The equation for the synaptic variables sj reads:
We suppose that the neurons in the STN and GPe ensembles are nonidentical. For this, the applied currents Iapp = Iapp, j for STN cells are Gaussian distributed with the mean 10 pA/μm2 and the standard deviation 0.015 pA/μm2. The parameter ε = εj for GPe neurons are also Gaussian distributed with the mean 0.0055 ms−1 and the standard deviation 2 10−5 ms−1. The values of the other parameters for the STN and GPe neurons are given in S1 Appendix.
In Fig 2 we illustrate the dynamics of STN neurons for the considered sets of parameters. For a strong coupling gG→S = 1.7 nS/μm2 the STN neurons synchronize and fire bursts nearly simultaneously [Fig 2A]. Such a synchronized dynamics of individual neurons results in a rhythmic activity of the STN as reflected, for example, by the collective firing rate [Fig 2B], which is the relative number of neurons firing a spike at a given time. The local field potential (LFP, see below for definition) of synchronized STN neurons also demonstrates well-pronounced oscillations of large amplitude [Fig 2B, red solid curve], which can serve as an indicator of a synchronized neuronal dynamics. For a weak coupling gG→S = 1.28 nS/μm2 the neurons are much less synchronized [Fig 2C], and STN does not produce a pronounced rhythmic output as illustrated by the firing rate in Fig 2D. The same applies to the dynamics of the LFP which exhibits low-amplitude oscillations [Fig 2D, red solid curve] indicating a desynchronized dynamics of the individual neurons.
In this study we focus on the control of the collective synchronized dynamics of the STN-GPe network. The extent of synchronization can be estimated either by the LFP amplitude [Fig 2B and 2D, red solid curves] or by the order parameter [2, 3, 82], where ψj(t) are the phases of individual neurons calculated from the neuronal bursting dynamics. The phase ψj(t) of the jth neuron attains the values ψj(tn) = 2πn, n = 0, 1,… at the time moments tn of the burst onsets, i.e., the first spikes in the bursts, and linearly increases between the neighboring bursts ψj(t) = 2π(t − tn)/(tn+1 − tn) + 2πn for t (tn, tn+1), n = 0, 1,… . The order parameter R(t) ranges from 0 to 1, where 1 corresponds to a perfect in-phase synchronization and 0 indicates a desynchronized state. For example, the time-averaged order parameter R(t) ≈ 0.69 and R(t) ≈ 0.21 for strongly and weakly synchronized regimes of STN neurons illustrated in Fig 2A, 2B, 2C and 2D, respectively.
We consider the two above cases of strongly and weakly synchronized regimes and show how the external stimulation based on delayed feedback can suppress the neuronal synchronization. We investigate the impact of two different feedback stimulation techniques counteracting neuronal synchronization. The first stimulation protocol is a linear delayed feedback (LDF) [55, 56]. To calculate the LDF stimulation signal, the LFP of synchronized STN neurons is measured. We model the LFP as an averaged synaptic activity of neurons , where sj(t) are the synaptic variables Eq (4) of STN neurons, see also papers [84, 85] for a more sophisticated approach. The measured LFP is on-line filtered by applying a linear damped oscillator
Parameter ω approximates the frequency of the LFP oscillations ω = 2π/T, where T is the mean period of LFP. For the strongly synchronized state illustrated in Fig 2A and 2B, T ≈ 110 ms. In a real application the mean period T for filtering can be selected, for example, based on the central frequency of the pronounced spectral peak in the β-band of the measured LFP observed in Parkinson’s disease [52, 53, 78, 86, 87]. As the output signal of Eq (5), that is the filtered LFP, we use the variable , which has a zero phase shift with respect to the original LFP signal . The damping and scaling coefficients in Eq (5) were chosen αd = kf = 0.008 which preserves the amplitude of the input raw LFP signal, see S1 Text for a parameter optimization approach for filtering. The filtered LFP signal is illustrated in Fig 2B and 2D by red curves for strongly and weakly synchronized regimes of STN neurons, respectively.
where K is the parameter of the stimulation intensity, and τ is the stimulation delay.
Another control method considered in this study is based on nonlinear delayed feedback (NDF) suggested in Refs. [57, 61] for the control of pathological neuronal synchronization. To construct the stimulation signal, we represent the measured mean field of the synchronized neuronal population in the form of an analytic complex signal Z(t) = x(t) + iy(t), where the variable x(t) is the filtered LFP signal obtained with the help of Eq (5) as for the case of LDF stimulation, and the corresponding y(t) signal can be calculated from x(t) by means of Hilbert transform . In a simple realization, which we use in this study, y(t) can be approximated by the time-shifted filtered LFP, y(t) = x(t − T/4), where T is the mean period of LFP. The stimulation signal of the NDF reads S(t) = KZ2(t)Z*(t − τ), where the asterisk denotes the complex conjugacy. In our case we consider only the real part of S(t) as the stimulation signal
where, as before, K is the stimulation intensity, and τ is the stimulation delay.
In this section we address the desynchronizing impact of the LDF and NDF stimulations with the stimulation signals Eqs (6) and (7), respectively. The stimulation signals are smooth and slowly oscillating since they are calculated from the filtered and, thus, smooth LFP of STN neurons. The stimulation is administered to STN neurons only, and the GPe neurons are not stimulated. We compare the effect of the stimulation for two conditions of strongly and weakly coupled neurons, where the neuronal population exhibits strongly and weakly synchronized intrinsic dynamics, respectively, as illustrated in Fig 2.
The desynchronizing impact of the differential LDF Eq (6) on the STN neurons is illustrated by the example in Fig 3. The order parameter R(t) of the STN without stimulation (stimulation intensity K = 0) saturates at R ≈ 0.69 for the case of strong coupling [Fig 3A, red circles] and fluctuates around R ≈ 0.21 for the case of weak coupling [Fig 3C, red circles]. In the latter case the STN neurons exhibit weak and intermittent synchronization, see Fig 2C and 2D. For the fixed stimulation delay τ = 95 ms and stimulation intensity K = 20, the differential LDF can significantly suppress the neuronal synchronization in the strongly coupled regime, where the values of the order parameter are reduced to R ≈ 0.02 (averaged during the last 20 s) [Fig 3A, black diamonds]. The stimulation-induced desynchronized dynamics is also reflected by the STN firing rate and behavior of LFP, where no rhythmic oscillations are observed any longer, see Fig 3B and compare to Fig 2B. Also for the case of weak coupling and weak and intermittent initial synchronization the differential LDF is still effective in counteracting neuronal synchronization. For the same stimulation parameters, the values of the order parameter are reduced to R ≈ 0.05 [Fig 3C, black diamonds], and the low-amplitude oscillations of the firing rate and LFP are further suppressed, see Fig 3D and compare to Fig 2D. The time course of desynchronization is monotonous in the case of strong coupling, where the feedback signal is regular [Fig 3A], while for weak coupling, where the feedback signal is contaminated by statistically relevant fluctuations, the time course of desynchronization displays fluctuations (and even epochs of slight increase of synchrony) as well [Fig 3C].
The LDF stimulation can have different impact on the collective dynamics of the stimulated neurons depending of the stimulation parameters. It is known that in the parameter plane of the stimulation delay and intensity (τ, K), there are regions of desynchronization complemented by domains, where synchronization is enhanced by stimulation [55, 56, 61]. This property is illustrated for the considered model in Fig 4. For example, for the strongly coupled regime, the LDF can induce a nearly perfect desynchronization in the blue parameter islands depicted in Fig 4A, where the order parameter practically vanishes, see also Fig 3A. In the red parameter domain in Fig 4A, on the other hand, the LDF stimulation forces the neurons to synchronize, and the order parameter closely approaches 1. Therefore, for a successful suppression of the neuronal synchronization the stimulation parameters have to be calibrated to the desynchronization regions. The situation is similar for the case, when the LDF stimulation is applied to weakly coupled and weakly synchronized STN neurons, which is illustrated in Fig 4B. The desynchronization regions are however larger for this regime, in particular, for a weak to moderate stimulation intensity.
It was also shown [55, 56] that the properties of the desynchronization regions, e.g., size, location, form, etc., depend on the properties and parameters of the stimulated system, such as oscillation frequency, extent of synchronization, coupling, etc. For example, changing the oscillation frequency of the stimulated neurons may shift the desynchronization regions in the (τ, K)-parameter plane with respect to the delay parameter τ. This effect is observed in Fig 4B for weakly coupled regime. Indeed, the period T of the LFP oscillations in the weakly coupled regime is approximately T ≈ 100 ms [Fig 2D] which is smaller than that for the strongly coupled neurons, where T ≈ 110 ms [Fig 2B]. This causes a shift of the desynchronization regions toward smaller values of the stimulation delay τ, as shown in Fig 4B (compare blue regions with those in Fig 4A bounded by white curves). This effect is more pronounced for large values of τ. It is also known that large values of stimulation delay and stimulation intensity are not favorable for a desynchronizing effect of LDF [55, 56, 61]. Therefore, for the optimal range of the stimulation parameters (small τ and K), the desynchronization induced by the LDF stimulation in the strongly synchronized neuronal population is preserved if the intrinsic neuronal dynamics changes to a weakly synchronized regime, except for a narrow boundary region, which has to be avoided also for the strongly synchronized case. Note, the authors of Ref.  considered the fixed stimulation delay τ = 50 ms which does not belong to the desynchronization regions [Fig 4A, blue domains] and, therefore, is not representative for the desynchronization effect of the LDF stimulation for the considered model.
For the case, when the stimulation parameters of the LDF are sub-optimally selected or shifted out of the desynchronization regions [Fig 4] because of the fluctuation of system properties (frequency, synchronization strength, etc.), they have to be recalibrated. This can be done automatically as suggested by Montaseri et al. , where the stimulation parameters can be adapted depending on the extent of the stimulation-induced desynchronization. The proposed calibration algorithm adjusts the stimulation parameters based on, for example, the amplitude of the LFP. Following Ref. , together with Eq (5) we consider an additional variable w(t):
where μ is some large coefficient such that μω 1, for instance, μ = 2000, and Eq (8) acts as an integrator. Since the variable is the filtered LFP, its amplitude can be found as . The rule of the automatic parameter calibration adapted from  reads
where parameters X1 = 0.005 and X2 = 0.028 define the cutoff thresholds for the growth of the stimulation parameters K and τ, respectively, a1 = 4000 and a2 = 500 determine the widths of the cutoffs, and b1 = 7.5 10−6 and b2 = 5 10−4 govern the speed of the parameter variation.
The performance of the parameter calibration algorithm Eqs (8) and (9) is illustrated in Fig 5. Starting at (τ(0), K(0)) = (0, 0), both parameters increase [Fig 5B] until the LFP amplitude [Fig 5A] falls below the cutoff threshold for the delay (parameter X2 in Eq (9), upper dashed line in Fig 5A). Then the growth of τ is slowed down, whereas the stimulation intensity K continues to increase and saturates when the LFP amplitude approaches the cutoff threshold for K (parameter X1 in Eq (9), lower dashed line in Fig 5A). At the end of the parameter calibration process the stimulation parameters reach the desynchronization region [Fig 5C], and the values of the order parameter are successfully lowered from R ≈ 0.69 at the initial state to R ≈ 0.16. The illustrated algorithm can also be used in the case when the parameters of the stimulated neuronal population are non-stationary, see Ref.  for details.
In this section we illustrate the desynchronizing effect of the NDF stimulation Eq (7) administered to the considered model of STN neurons Eqs (1)–(3). In Fig 6 the time-averaged order parameter R is depicted by color versus the parameters of the stimulation delay τ and stimulation intensity K for two cases of strong coupling [Fig 6A] and weak coupling [Fig 6B], which correspond to strongly and weakly synchronized STN neurons without stimulation, respectively, see Fig 2. Starting from the mentioned initial states of the STN, the NDF stimulation can effectively desynchronize the stimulated neurons with stimulation parameters selected from the desynchronization parameter regions showing up in blue color in Fig 6. The desynchronization regions for the two considered regimes are located at nearly the same places in the (τ, K)-parameter plane although the stimulation-free dynamics of the strongly and weakly coupled neurons differ from each other, e.g., by oscillation frequency and extent of synchronization. The NDF desynchronization method is thus less sensitive to the variation of these properties of the neuronal population: If strongly synchronized neurons can be desynchronized by the NDF stimulation, the desynchronization is also induced by stimulation with the same parameters in the case of weakly synchronized neurons (compare blue regions bounded by the white curves in Fig 6A and 6B).
The desynchronization mechanism of NDF differs from that of LDF. The theory of LDF says that the desynchronization regions are bounded by bifurcation curves, where the amplitude of the ensemble mean field vanishes, which corresponds to a state of perfect desynchronization, see Refs. [55, 56, 88]. The NDF stimulation, on the other hand, does not cause any bifurcation, but the amplitude of the mean field is gradually suppressed and decays ~K−1/2 as the stimulation intensity K increases [57, 61]. For a network of two interacting populations, where only one population is measured and stimulated, the mean field of the stimulated neurons decays ~K−1/3 as K grows . The latter case corresponds to the stimulation setup considered in this study, where only STN neurons from the STN-GPe network are recorded and stimulated.
The behavior of the time-averaged order parameter R, the absolute value of the filtered LFP |LFP|, and the absolute value of the stimulation signal |Istim| from Eq (7) of the NDF stimulation is shown in Fig 7A versus the stimulation intensity K. As expected, order parameter and LFP decay as K grows according to a power law ~Kγ, see the insert in Fig 7A. For the considered STN population of N = 200 neurons, the exponent γ slightly deviates from the predicted −1/3 , where the direct numerical fit gives γ ≈ −0.3, which may be related to the finite-size effect . We however evaluate the properties of the feedback stimulation for the considered finite-size neuronal population. Because of such a value of γ, the amplitude of the stimulation signal |Istim| slowly increases when K grows [Fig 7A, blue triangles], cf. Ref. . Nevertheless, even for large values of the stimulation intensity K, where the order parameter can reach R ≈ 0.035 indicating an about 20-fold suppression of synchronization, the amplitude of the stimulation signal |Istim| is still smaller than the amplitude of the LFP |LFP| ≈ 0.018 of the synchronized STN neurons without stimulation, see Fig 7A. Therefore, the NDF technique can effectively suppress the neuronal synchronization by a weak stimulation as compared, for instance, to HF DBS.
We evaluate the dynamics of the quantities R, |LFP|, and |Istim| also for the LDF stimulation Eq (6), which are depicted in Fig 7B versus K. The order parameter and the amplitude of LFP decay as K increases much faster than for the NDF stimulation, and the mentioned level of desynchronization R ≈ 0.035 can be reached already for K = 10 and ranges to R ≈ 0.0035 for K = 100. However, the amplitude of the stimulation signal |Istim| also increases as the stimulation intensity grows [Fig 7B, blue triangles]. To compare the efficacy of the LDF and NDF stimulations we plot in Fig 7C the amount of the administered stimulation as given by the amplitude of the stimulation signal |Istim| from Eqs (6) and (7), respectively, versus the reached extent of desynchronization as given by values of R in Fig 7A and 7B. For a range of the order parameter R > 0.055 the NDF stimulation is more efficient in inducing desynchronization than the LDF stimulation for a given ensemble size: To reach the same extent of desynchronization the NDF requires a smaller amount of the stimulation. However, if we need a very strong desynchronization, then the LDF stimulation is more appropriate.
Stronger desynchronization can be obtained by the LDF stimulation for a larger stimulation intensity K, which also results in a large amplitude of the stimulation signal |Istim|, see Fig 7B. Furthermore, for too large K the desynchronization regions [Fig 4A, blue domains] shrink, and the LDF stimulation becomes ineffective in inducing desynchronization, see Refs. [55, 56, 61]. For example, for τ = 95 ms the LDF stimulation can boost synchronization instead of suppressing it for K > 102 [Fig 7B]. In Fig 7D we plot the fraction of the parameter τ values, where the LDF and NDF stimulations can suppress synchronization with R < 0.35, i.e., at least twice as compared to the stimulation-free regime, where R ≈ 0.69. The desynchronization regions of the LDF stimulation can occupy up to 35% of the (τ, K)-parameter plane for small K, whereas they can be reduced to only 4% of the parameter volume for large K [Fig 7D, red diamonds]. For the NDF stimulation, the desynchronization regions occupy more than 65% of the (τ, K)-parameter plane, and this property seems to be independent of the parameter of the stimulation intensity K [Fig 7D, black circles]. Therefore, large values of parameter K are not desirable for the LDF stimulation, whereas they can be beneficial for the NDF stimulation.
The stimulation signals Eqs (6) and (7) of the LDF and NDF, respectively, are derived from the filtered LFP. The latter is a smooth signal and oscillates at the frequency ≈ 9−10 Hz for the considered model and parameters [Fig 2], which results in smooth stimulation signals Eqs (6) and (7) oscillating at these frequencies. When administered to the stimulated neurons as electrical stimulation, e.g., via a deeply implanted electrode, the neurons receive a stimulation current of the same polarity during a half of the oscillation period, which lasts several tens of milliseconds. Such a long stimulation phase of the same polarity may cause an irreversible and possibly impairing charge deposit in the vicinity of the electrode, which may violate safety requirements of the electrical stimulation of the neuronal tissue and lead to its damage [16, 69, 70].
The stimulation signal of the standard HF DBS consists of a pulse train of biphasic electrical pulses administered at high frequencies above 100 Hz, often at 130 Hz . The cathodic and anodic phases of each pulse should deliver the same charge of opposite polarity providing a charge-balanced stimulation with zero net charge injection into the stimulated tissue, which is a necessary requirement to avoid a tissue damage [16, 69, 70, 90].
Combining both approaches, we construct a stimulation signal which inherits the advantages of the charge-balanced property of the HF DBS signal and the desynchronizing impact of the delayed feedback stimulation. In such a way, the amplitude of the high-frequency pulse train is modulated by the slowly oscillating smooth signal of the LDF Eq (6) or NDF Eq (7), which is schematically shown in Fig 8. We consider four types of the asymmetric biphasic charge-balanced pulses illustrated in the inserts in Fig 8. The typical waveform of the biphasic charge-balanced stimulation pulses used for the standard HF DBS consists of a short first pulse (1st phase) of duration 60 to 450 μsec  followed by a charge-balancing 2nd phase of opposite polarity such that the total charge of the biphasic pulse is zero . We consider two cases of the anodic 1st phase [Fig 8A] and the cathodic 1st phase [Fig 8D]. The latter pulse shape is a standard pulse waveform widely used for HF DBS . In our numerical simulations we also consider additional pulse shapes illustrated in Fig 8B and 8C, where the 1st and the 2nd phases of the stimulation pulses in Fig 8A and 8D are exchanged with each other, respectively, and the charge-balancing phase advances the short stimulation pulse. We consider the frequency of the pulse train 130 Hz (the inter-pulse interval 1000/130 ≈ 7.69 ms). The width of the short pulse (the 1st phase in Fig 8A and 8D and the 2nd phase in Fig 8B and 8C) is taken PW = 0.5 ms or 0.2 ms, and relates to the duration of its long counterpart as 1: 10, which is found to be energy efficient . Potentially, charge-balanced waveforms could be incorporated in other feedback methods, such as PRC based ones .
For a given smooth signal Istim(t) of the delayed feedback calculated according to Eqs (6) or (7) [Fig 8, dashed curves], the amplitude of a stimulation pulse is calculated at the time t = tp of the pulse onset as Istim(tp). Then, depending on the pulse shape, see insets in Fig 8, this value is assigned either to the short phase of the biphasic pulse as in Fig 8A and 8B, or to the long phase of the pulse as in Fig 8C and 8D, as indicated by red dots in the inserts. The amplitude of the other counterpart of the pulse is obtained from the charge-balancing property such that the square delineated by the biphasic pulse is zero. Examples of the values tp (pulse onset) and Istim(tp) are indicated in Fig 8B by dotted lines, where the value Istim(tp) is assigned in this case to the short (second) phase of the pulse.
We illustrate the impact of the pulsatile LDF stimulation on the STN neurons Eqs (1)–(3) in Fig 9 for the case of strong coupling gG→S = 1.7 nS/μm2. Starting from the initially synchronized regime [Fig 2A and 2B] the stimulation of STN neurons with the pulse train modulated by the LDF signal Eq (6) can suppress the neuronal synchronization. The efficacy of the stimulation depends on the type of the charge-balanced biphasic pulses. For the pulses of type A or C [Fig 8A and 8C] the pulsatile LDF stimulation is hardly effective [Fig 9A and 9C]. On the other hand, for the pulses of type B or D [Fig 8B and 8D], the pulsatile LDF can pronouncedly desynchronize the stimulated neurons in large parameter regions [Fig 9B and 9D]. The size of the desynchronization regions approaches that for the NDF stimulation with smooth stimulation signal, compare the size of blue regions in Fig 9B and 9D with that in Fig 6A and in Fig 4A.
The impact of the pulsatile NDF stimulation on the STN neurons Eqs (1)–(3) is illustrated in Fig 10 for the case of strong coupling gG→S = 1.7 nS/μm2. When the pulse train is modulated by the NDF signal Eq (7), the neurons can effectively be desynchronized irrespectively of the considered type of the the charge-balanced biphasic pulses. The desynchronization regions are very similar in size (although they can be shifted with respect to parameter τ) in all four plots in Fig 10 (blue domains), which indicates that the pulsatile NDF stimulation is little sensitive to the waveform of the stimulation pulses.
To address the problem mentioned by Dovzhenok et al.  that the desynchronizing feedback stimulation may boost synchronization when administered to initially intermittently and weakly synchronized regimes, we apply the pulsatile LDF and NDF stimulations to the STN neurons in the weakly coupled regime for gG→S = 1.28 nS/μm2 [Fig 2C and 2D] for the same stimulation parameters as in Figs Figs99 and and10.10. We consider the pulses of type D [Fig 8D] widely used for the standard HF DBS . The results of the stimulation are shown in Fig 11A and 11B for the pulsatile LDF and NDF stimulations, respectively. Both pulsatile LDF and NDF stimulations are robust with respect to the variation of the extent of synchronization in the neuronal population. Indeed, for the parameters from the desynchronization regions delimited by the white curves in Fig 9D for LDF stimulation of initially strongly coupled and strongly synchronized neurons, the pulsatile LDF stimulation does not cause any enhancement of synchronization also for the ensembles of weakly coupled and weakly synchronized neurons, see Fig 11A, where the boundaries of desynchronization regions from Fig 9D are depicted by white curves for comparison. Also for the pulsatile NDF stimulation, the largest continuous desynchronization regions obtained for initially strongly coupled and strongly synchronized neurons [Fig 10D] nearly perfectly fit to the desynchronization regions for initially weakly coupled and weakly synchronized neurons, see Fig 11B (white curves). Our results cannot confirm the conclusions of Dovzhenok et al.  for the pulsatile LDF and NDF stimulations. For an appropriate selection of the stimulation parameters causing a pronounced desynchronization of initially strongly synchronized neurons, the stimulation by pulsatile LDF and NDF is not expected to cause any problem when the neuronal population runs into a regime of weak or intermittent synchronization.
The suggested pulsatile stimulation by LDF or NDF are also robust with respect to a variation of the pulse parameters. To illustrate this we consider a shorter pulse width PW = 0.2 ms. Also with such pulses of type D [Fig 8D] the pulsatile NDF stimulation can effectively desynchronize the STN neurons, where the quality of desynchronization is improved (the order parameter R and the amplitude of the LFP |LFP| decay) as the stimulation intensity K increases, see Fig 12A. As for the smooth NDF [Fig 7A], the decay of R and |LFP| obeys a power law ~Kγ, see the insert in Fig 12A. Since the considered neuronal model has a finite size, the numerically obtained exponent γ slightly deviates from the theoretically predicted γ = −1/3 , which leads to a slow increase of the amplitude |Istim| of the stimulation signal [Fig 12A, blue triangles]. Hence, for a better desynchronization the neurons have to be stimulated somewhat stronger. The same holds for the pulsatile LDF stimulation, and we illustrate the efficacy of both pulsatile feedback methods in Fig 12B, where the amount of the administered stimulation as given by |Istim| is plotted versus the reached extent of desynchronization as given by R. If compared to the case of smooth stimulation signals [Fig 7C], the pulsatile stimulation protocol diminishes the difference in the efficacy of NDF and LDF, where both methods require approximately the same amount of the stimulation to reach the given level of desynchronization. A similar conclusion can be drawn with respect to the size of desynchronization regions, except for the sensitivity of the pulsatile LDF and NDF to the shape of individual pulses, see Figs Figs99 and and1010.
To further verify the robustness of the considered stimulation methods and to model the effect of slow variations of system parameters, we calculate the order parameter R of the STN neurons when the stimulation delay τ slowly changes. The initial conditions of the model for the next value of τ correspond to the last state of the system for the previous value of τ. The results of such a continuation by parameter τ are illustrated in Fig 13A for pulsatile LDF and in Fig 13B for pulsatile NDF. The variation of τ, which can correspond to the changes of the oscillation frequency of STN neurons, does not cause any problem with respect to the stimulation-induced desynchronization. Indeed, the pulsatile LDF stimulation with slowly varying delay [Fig 13A, red circles and blue squares] demonstrates the same desynchronizing effects as for the case where the stimulation is administered to initially synchronized population [Fig 13A, black solid curve]: The size and location of desynchronization regions are preserved. For the pulsatile NDF stimulation the situation can even be improved, where the variation of parameters significantly extends the desynchronization regions [Fig 13B, red circles and blue squares]. With such an approach of slowly varying stimulation parameters, the NDF stimulation can desynchronize the stimulated neurons for any value of the stimulation delay as has also been reported for other models [57, 61].
In this paper we showed that the delayed feedback methods preserve their desynchronizing impact on synchronized neurons also in the case when the smooth feedback stimulation signal is replaced by a high-frequency pulse train of charge-balanced pulses with amplitude modulated according to the corresponding feedback algorithm. This may pave the way to pre-clinical and clinical tests of the delayed feedback methods for counteracting abnormal neuronal synchronization in the framework of closed-loop DBS. The latter paradigm is in focus of research nowadays, as discussed in the introduction. For instance, in a particular, successfully tested aDBS approach , the stimulation voltage was linearly adapted to the patient’s beta band power each second. From a technical standpoint this is close to our approach, where the amplitude of the HF DBS pulse train is directly modulated by the delayed feedback signal [Fig 8]. Our results suggest that an appropriately processed and time-delayed adaptation of the stimulation amplitude to the beta power amplitude might further increase the stimulation outcome. We presented a detailed investigation and comparison of the properties of smooth and pulsatile LDF and NDF administered to a physiologically based model of STN-GPe neurons  contributing to the abnormally synchronized neuronal dynamics characteristic for Parkinson’s disease [8, 9, 18].
We showed that both smooth and pulsatile LDF and NDF techniques can robustly desynchronize the stimulated neurons and revealed the corresponding desynchronization regions in the parameter plane of the stimulation delay τ and stimulation intensity K. To address the possible problem mentioned by Dovzhenok et al.  that the feedback stimulation may be ineffective because of a possible variation of the extent of synchronization in neuronal populations, we calculated the desynchronization parameter regions for initially strongly and weakly (and intermittently) synchronized neurons. We showed that the desynchronization regions for initially weakly synchronized neurons do not differ much in size and location from those for initially strongly synchronized neurons, whereas the former can even be larger for an optimal range of the stimulation parameters. For smooth LDF, the desynchronization regions can be displaced by the delay parameter τ (especially for suboptimally large delay) with respect to each other for the above two cases of strong and weak initial synchronization [Fig 4]. This however is caused by a shift in the oscillation frequency and not by the extent of synchronization among neurons itself since the optimal value of the stimulation delay τ is related to the oscillation period as known from several publications [55, 56, 61]. In order to observe the effects reported by Dovzhenok et al.  for the optimal range of the stimulation parameters, one has to particularly search for the corresponding parameter values at the border of the desynchronization region, which significantly devalues the generality of the claims of Dovzhenok et al. . For a practical application we hence need two types of closed-loop control processes: (i) On a fast time scale the measured signal is used for amplitude modulation of the HF pulse train according to the delayed feedback algorithm. (ii) On a slow time scale the delay of the feedback algorithm is adapted to the period of the abnormal oscillation in order to make sure the system is kept in a desynchronization region.
When parameters of the neuronal system significantly change to an extent considered by Dovzhenok et al. , the stimulation will in fact be administered to a completely different system with different dynamics, and the stimulation parameters have to be recalibrated. For the standard HF DBS, for example, selecting a clinically effective parameter set is an iterative process, and a systematic reprogramming of stimulation parameters is required in order to optimize clinical benefit [21–23, 78, 80, 91, 92]. Adaptive optimization of the stimulation parameters is a focus of research in the framework of the closed-loop DBS [78–80]. By the same token, for the delayed feedback stimulation characteristic stimulation parameters, especially the delay, should be adapted to measurable variables, such as the period of the abnormal oscillation. Using just one fixed set of stimulation parameters, irrespective of the system’s state, as suggested by Dovzhenok et al. , is neither an appropriate stimulation strategy nor does it take account of fundamental properties of delayed feedback stimulation mechanisms. The parameter adaptation for the delayed feedback stimulation can be done automatically as proposed in Ref. , which we successfully tested also for the considered STN-GPe model, see Fig 5. For smooth NDF as well as for pulsatile LDF and NDF, such a parameter recalibration seems to be unnecessary, where the stimulation-induced desynchronization obtained for initially strongly synchronized neurons is also preserved for initially weakly synchronized regime, see Figs Figs66 and and11.11. We therefore refute the conclusions of Dovzhenok et al.  for smooth NDF and do not confirm them for pulsatile LDF and NDF. Note, the authors of Ref.  do not distinguish between LDF and NDF which are clearly two different control techniques. Moreover, the desynchronizing impact of the pulsatile NDF is little sensitive to the stimulation delay when the latter slowly changes, which models the slow variations of the oscillation frequency [Fig 13B]. We also verified the efficacy of the considered approach for the case when neurons are perturbed by an independent noise modeling a noisy environment inevitable in real systems, see S2 Text. The investigated stimulation techniques are thus robust desynchronizing methods which can be promising candidates for closed-loop DBS.
For a successful utilization of any control method, the application conditions of the method have to be handled with care by precisely following all application requirements of a given method. Otherwise, any approach can easily be spoiled. The LDF and NDF techniques, for example, are essentially based on measuring of a mean field of a sufficiently representative, i.e., large enough oscillatory population such that the measured signal reliably reflects the properties of synchronized dynamics. The stimulation should also be administered to a large enough (sub)population with stimulation parameters appropriately related to the other parameters and properties of the stimulated system. All these fundamental requirements were not taken into account by the authors of Ref. . By simulating a model containing 10 STN neurons during maximally 5 sec (less than 50 oscillation periods), measuring 5 or 3 neurons as a mean field (LFP) contaminated by the stimulation signal and stimulating only a few neurons with nonlinearly modified stimulation signal and suboptimally selected stimulation parameters, it was concluded that the delayed feedback “is unlikely to be clinically successful”  which is inappropriate and might, in fact, prevent other researchers from further fruitful theoretical investigations and pre-clinical as well as clinical tests of delayed feedback control methods.
However, due to the substantial stimulation-induced artifact , it could be difficult to perform the pulsatile LDF and NDF with the same stimulation and recording electrode. To remove the stimulation-induced artifacts, a special filtering has to be applied , or the stimulation and recording can be arranged by adjacent contacts in a specific bipolar configuration . The pulsatile LDF and NDF can also be delivered by means of separate stimulation-registration setups [58, 63] or by an act-and-wait protocol .
The smooth LDF has been tested for many different models and stimulation setups demonstrating a pronounced desynchronizing effect [55, 56, 74]. The same has been demonstrated for smooth NDF [57, 61, 63]. The structure of the parameter space of LDF was also experimentally confirmed for arrays of coupled electrochemical oscillators . The desynchronizing impact of smooth NDF was also confirmed experimentally for suppression of synchronized alpha-rhythm in visual cortex by visual stimulation in healthy subjects . Several studies on closed-loop DBS revealed promising results [34, 35, 37–39, 77–79]. However, there are still relevant open questions [93, 95], and our approach, extending the LDF and NDF methods to physiologically relevant balanced pulsatile stimulation protocols, might help to resolve some of these issues. In particular, our approach requires an oscillatory biomarker sufficiently representing abnormal states or conditions. However, as yet, it is not clear whether such a biomarker is available . For instance, so far it is not clear whether low or high frequency beta might be more appropriate as feedback signal . Applying our method to different frequency (sub-)bands might help to probe and find possible oscillatory biomarkers. Actually, the beta oscillation alone might possibly not be an optimal feedback signal , e.g. because enhanced beta oscillations are not consistently observed in all Parkinson’s patients [52, 93]. For instance, characteristic changes of two distinct bands of high frequency oscillations (HFO), around 250 Hz and 350 Hz, were observed after levodopa administration, and the power ratio of these two bands was significantly correlated with the Unified Parkinson’s Disease Rating Scale hemibody akinesia/rigidity subscore, but not with beta power . As suggested by this study (prokinetic) HFO (>200 Hz) and beta oscillations might potentially react in an independent manner to levodopa administration. The potentially differential response characteristics of HFO and beta oscillations to different types of aDBS still remain to be investigated. Another relevant aspect is patient phenotype : In tremor dominant (as opposed to akinetic rigid) Parkinson’s patients resting state beta power might be reduced during epochs with tremor [53, 97]. In addition, low-frequency (theta and beta) and high-frequency oscillations interact under both physiological  and pathological  conditions. Physiological processes in motor circuits may be perturbed, since the HFO amplitude is coupled to the phase of beta oscillations [99, 100]. Hence, phase amplitude coupling (PAC) might be used as biomarker to trigger demand-controlled DBS . Obviously, PAC-based closed-loop DBS does not easily translate to the feedback approach presented in our paper, since the latter requires a feedback signal representing the synchronized activity of a large neuronal population. Johnson et al.  performed another study demonstrating that STN beta LFP power alone will likely not be an appropriate biomarker. They compared traditional continuous DBS (tDBS) with closed-loop DBS (CL-DBS) in MPTP monkeys. CL-DBS reduced rigidity to an amount comparable to tDBS, but with approximately half the stimulation ON time. However, only tDBS improved bradykinesia during reaching tasks, probably because beta power was reduced related to the reaching process, in this way reducing the biomarker .
It would be interesting to compare the efficacy in suppressing pathological neuronal oscillations of the methods considered in this paper to other control methods, for instance, to those relying on event-based or phase-locked stimulation [66, 67, 101], where the stimuli are administered at a particular phase of the oscillation cycle. Vulnerable tremor phases were studied with non-invasive approaches [102, 103] as well as with thalamic DBS . In the latter study  10 patients with essential tremor were stimulated with low-frequency stimulation with thalamic electrodes, while tremor amplitude and phase were recorded. Stimulation with a stimulation frequency close to the postural tremor frequency caused a tremor entrainment. Tremor amplitude was modulated by the tremor phase at which the stimulus was delivered. Stimuli delivered in the favorable half of the tremor cycle caused a reduction of the median tremor amplitude by approximately 10%, whereas stimuli delivered in the opposite half of the tremor cycle caused an increase of the tremor amplitude by a similar amount. Stimuli administered at an optimal tremor phase caused a tremor reduction of 27%. Improvement in tremor severity by the standard HF DBS was on average 70%, supporting satisfactory DBS electrode placement. This approach might be further improved by taking into account the dependence on the tremor amplitude itself. According to the theory of phase resetting in noise-free  and noisy oscillatory populations [3, 30] it is known that the vulnerable phase of a synchronized neuronal population strongly depends on the extent of synchrony. If not properly adapted to the amount of synchrony, the single pulse may, in fact, massively boost synchrony .
In summary, as yet there is no an established biomarker representing the patients’ abnormal states or conditions . In this context, we also have to take into account that increased beta power as well as PAC of beta phase to HFO amplitude may not be disease specific, but were also found in patients with isolated dystonia . This may well limit the diagnostic value of these LFP features , but still they might be useful biomarkers for closed-loop DBS. However, applying our approach to different frequency (sub-)bands might eventually help to probe possible oscillatory candidates and help to reduce side effects by inducing long-lasting, sustained therapeutic effects by a specifically desynchronizing effect, as demonstrated in pre-clinical studies [48, 49] and a clinical proof of concept study  on CR-DBS.
In this study we presented a novel method for demand-controlled desynchronization stimulation. To this end, we combined desynchronizing linear and non-linear delayed feedback stimulation techniques with cDBS by using the delayed feedback signal to modulate the amplitude of the cDBS pulse train. Such a stimulation signal satisfies mandatory safety requirements of charge deposit in neuronal tissue by electrical stimulation. In addition, it inherits the desynchronizing properties of the delayed feedback. We investigated the impact of smooth and pulsatile desynchronizing delayed feedback stimulation on a physiologically based model of a neuronal network of interacting STN-GPe neurons for linear and nonlinear delayed feedback techniques. We demonstrated that both feedback techniques can desynchronize both initially strongly and initially weakly synchronized STN neurons. When comparing smooth LDF and NDF stimulations, the latter appears to be less sensitive to parameter variations and stimulation conditions and has larger desynchronization parameter regions. The former, however, is more effective in inducing strong desynchronization, but less effective if desynchronization of a moderate extent is desirable. In the case of pulsatile stimulation, the difference between LDF and NDF is diminished. The pulsatile LDF is however hardly inducing synchronization for the considered range of parameters for initially weakly synchronized neurons, whereas the desynchronizing impact of pulsatile NDF is robust with respect to stimulation pulse shape and to slow parameter variations. Both stimulation techniques can be suggested for closed-loop desynchronizing DBS.
The study was funded by the Helmholtz Society (OVP, BL, PAT). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
All relevant data are within the paper and its Supporting Information files.