Anti-kindling Induced by Two-Stage Coordinated Reset Stimulation with Weak Onset Intensity

Introduction

Several brain disorders are characterized by abnormally strong neuronal synchronized activity (Uhlhaas and Singer, 2006). Examples of such neuronal disorders are epilepsy (Wong et al., 1986; Schomer and Lopes da Silva, 2011), Parkinson's disease (Lenz et al., 1994; Nini et al., 1995; Hammond et al., 2007), and tinnitus (Ochi and Eggermont, 1997; Llinas et al., 1999; Weisz et al., 2005; Eggermont and Tass, 2015). Neuronal dynamics essentially depends on the neurons' connectivity patterns (Sporns, 2011). Furthermore, the mechanism by which the neurons interact crucially determines whether and how they synchronize. For instance, as shown computationally the synchronization behavior of neuronal networks may strongly differ depending on whether neurons are coupled via gap-junctions or synapses (Belykh et al., 2005; Belykh and Hasler, 2011). Accordingly, computational studies were devoted to the role of gap junctions in the emergence of epileptic seizures (Volman et al., 2011) or on the impact of the interplay between gap junctions and delayed inhibitory synaptic coupling on the emergence of different patterns of neuronal synchrony (Guo et al., 2012).

Connectivity and related function may undergo plastic changes that are not restricted to periods early in life (Hübener and Bonhoeffer, 2014). The timing pattern of neuronal activity may strongly shape the strength of neuronal connections (Hebb, 1949; Bliss and Lomo, 1973). Spike timing-dependent plasticity (STDP) is a fundamental mechanism by which neurons adapt the strength of their synapses to the relative timing of their action potentials (Gerstner et al., 1996; Markram et al., 1997; Bi and Poo, 1998; Feldman, 2000). A number of computational studies addressed the impact of activity dependent coupling/synaptic strength on the collective dynamics (Seliger et al., 2002; Tass and Majtanik, 2006; Maistrenko et al., 2007; Masuda and Kori, 2007; Ren and Zhao, 2007; Aoki and Aoyagi, 2011; Belykh and Hasler, 2011; Volman et al., 2011; Bayati and Valizadeh, 2012; Guo et al., 2012; Knoblauch et al., 2012; Popovych et al., 2013; Zhang et al., 2013). In the presence of STDP a variety of different dynamical phenomena and regimes emerge, such as complex clustering phenomena, see e.g., Maistrenko et al. (2007) and Belykh and Hasler (2011), and the emergence of traveling waves, see e.g., Zhang et al. (2013). Neuronal networks as well as oscillator networks with STDP typically display multistability. For instance, multistability was found in phase oscillator networks with symmetric as well as asymmetric phase difference-dependent plasticity (a time continuous approximation of STDP; Seliger et al., 2002; Maistrenko et al., 2007) as well as in phase oscillator networks with STDP (Tass and Majtanik, 2006) and, subsequently, in different types of neuronal networks with STDP (Tass and Hauptmann, 2006; Popovych and Tass, 2012).

Based on a computational approach aiming at the development of stimulation techniques that specifically counteract abnormal neuronal synchrony by desynchronization (Tass, 1999), finally Coordinated Reset (CR) stimulation was designed (Tass, 2003a,b). CR stimulation means that phase resetting signals are delivered to different sub-populations within the abnormally synchronized neural network in a spatiotemporal manner (Tass, 2003a,b), such that the phase of each subpopulation is reset once within a stimulation ON-cycle. Computational studies showed that CR stimulation in a network with STDP (Gerstner et al., 1996; Markram et al., 1997; Bi and Poo, 1998; Feldman, 2000) not only causes a desynchronization of the abnormally strong synchronized neuronal activity but due to STDP also decreases the average synaptic weight. This anti-kindling, i.e., unlearning of abnormal synaptic connectivity and of abnormal neuronal synchronized activity by repetitive stimulation, is a long-lasting sustained effect (Tass and Majtanik, 2006; Hauptmann and Tass, 2007, 2009; Popovych and Tass, 2012). Abnormal synaptic connectivity in this context means abnormally strong synaptic connectivity that induces abnormal neuronal synchrony. The anti-kindling approach is not restricted to epilepsy. Rather it was designed to be tested in a number of disease conditions characterized by abnormally strong neuronal synchrony. In addition, computational studies showed that anti-kindling can robustly be obtained in networks with plastic excitatory and inhibitory synapses, no matter whether CR stimulation is delivered directly to the soma or via synapses (Popovych and Tass, 2012; Tass and Popovych, 2012). In line with these computational findings, long-lasting CR-induced anti-kindling was achieved invasively as well as non-invasively: in pre-clinical and clinical studies with rat hippocampal slices rendered epileptic by magnesium withdrawal (Tass et al., 2009), with parkinsonian non-human primates (Tass et al., 2012b), and with parkinsonian patients (Adamchic et al., 2014a), as well as in a proof of concept-study with tinnitus patients (Tass et al., 2012a; Adamchic et al., 2012a,b, 2014b; Silchenko et al., 2013).

Under certain circumstances, like in the case of side effects, such as speech and gait deterioration in Parkinson's patients receiving deep brain stimulation (Mahlknecht et al., 2015) or in tinnitus patients suffering from hyperacusis (Sheldrake et al., 2015), it might be desirable or even unavoidable to use weak stimulation intensities for the CR stimulation. For instance, to date the goal of clinical programming recipes for traditional, permanent high-frequency deep brain stimulation is to select stimulation sites and parameters such that target regions are maximally covered, and stimulation of adjacent regions is minimized (Saint-Cyr et al., 2000; Dujardin et al., 2001; Schroeder et al., 2003; Rodriguez-Oroz et al., 2005; Deuschl et al., 2006a; Hershey et al., 2008; van Nuenen et al., 2008; McIntyre et al., 2011; Mikos et al., 2011; Xu et al., 2011; Homer et al., 2012; Tommasi et al., 2012; Pinsker et al., 2013; Jahanshahi et al., 2015). Accordingly, to adapt the spatial extent of current spread according to anatomical borders of stimulation targets new electrode designs were suggested (Martens et al., 2011; van Dijk et al., 2015). However, there are side effects that are at least partly due to stimulation of the sensorimotor target region itself (Moreau et al., 2008; Jahanshahi et al., 2015); Accordingly, in addition to spatially shaping stimulation currents according to anatomical borders of target regions, we pursue an approach that tends to minimize side effects by using qualitatively different stimulation patterns which induce therapeutic effects that outlast cessation of stimulation and require weaker intensities compared to standard deep brain stimulation, e.g., a third of the standard stimulation amplitude (Tass et al., 2012a). In this study, we want to drive our approach forth with the final goal to enable clinically effective stimulation at further and substantially reduced stimulation intensities. CR stimulation essentially requires to stimulate neuronal subpopulations in a sufficiently separated manner. Hence, for a given, fixed spatial alignment of the stimulation sites, like in the case of deep brain stimulation, an appropriate stimulation amplitude has to be chosen in order to avoid a global, spatially non-specific stimulation and rather manipulate the different subpopulations separately (Tass et al., 2012b). The situation is different and more complex if stimulation sites can easily be modified as, e.g., in the case of acoustic CR stimulation for the treatment of tinnitus. In that case, not only stimulation intensity, but also stimulation sites have to be calibrated properly. For instance, based on the tonotopic organization of the central auditory system, tones are audiologically calibrated to the tinnitus pitch (Tass et al., 2012b), and a lot of efforts were put into improving reliable audiological calibration procedures (Hauptmann et al., 2016).

Computational studies have shown that CR stimulation is most effective for intermediate stimulation intensities when applied with rapidly varying sequences (RVS; Lysyansky et al., 2011; Popovych and Tass, 2012; Ebert et al., 2014; Zeitler and Tass, 2015) as well as with slowly varying sequences (SVS; Zeitler and Tass, 2015). In this study a sequence denotes the sequence of activating the spatially evenly distributed stimulation sites within the neuronal population once during the stimulation CR-on cycle. For the RVS CR stimulation the sequence changes continuously, whereas for the SVS-n CR stimulation the sequence is repeated n-times before another sequence is applied. For intermediate stimulation intensities the SVS CR stimulation, with its many repetitions of each sequence is more effective and robust in inducing an anti-kindling than the completely random RVS CR stimulation. However, for weak stimulation intensities the probability that a pronounced anti-kindling is achieved is higher for the RVS than for the SVS CR stimulation. So, replacing the commonly used RVS CR by the SVS CR is rather a disadvantage at weak stimulation intensities. However, these numerical results led us to the idea of sequentially combining the two CR protocols by delivering a two-stage CR stimulation in which after a preparatory stimulation at weaker intensities we switch over to a different CR protocol at intermediate intensities.

In this study we set out to improve the CR stimulation protocol for weak stimulation intensities, so that the same amount of anti-kindling can be induced as for intermediate stimulation intensities. We do this by studying systematically the anti-kindling impact of the two-stage CR stimulation, which we compare with the impact of the single-stage CR stimulation. The latter contains only one CR stimulation type. In addition, we also study the anti-kindling effect of stepwise increase of the stimulation intensity during the stimulation period while retaining the same CR protocol (CR stimulation with weak onset). For comparison and illustration, as control condition we use sham stimulation, i.e., we deliver no stimulation at all. Sham stimulation means that a medical device in controlled clinical studies or a stimulator in a pre-clinical study act as a placebo or stim off control, see e.g., McCracken and Kiss (2014). We use this (in the context of our computational study) simple control condition, since already comparatively uncomplex stimulation protocols such as periodic stimulation delivered through one stimulation side may elicit complex dynamics (Popovych and Tass, 2011).

To investigate the effects of a two-stage CR stimulation, we performed numerical simulations of Hodgkin-Huxley neurons in a network which incorporates excitatory as well as inhibitory spike timing-dependent plasticity (Section Materials and Methods). First, we demonstrate the extent of anti-kindling as a function of the stimulation intensity for the single-stage RVS as well as for the single-stage SVS CR stimulation (Section Short Single-Stage CR-Stimulation). Subsequently, we compare the effects of single-stage CR stimulation with weak onset and two-stage CR stimulation with weak onset (section Single- and Two-Stage CR Stimulation with Weak Onset). In Section Single- and Two-Stage CR Stimulation with Weak Onset we also attempt to disentangle the roles played by the individual CR approaches for very weak and weak stimulation intensities. We compare, the anti-kindling induced by short and long single-stage CR stimulation as well as by the two-stage CR-stimulation for different stimulation intensities (Section Two-Stage CR Stimulation with Constant Stimulation Intensity). Finally, we discuss possible clinical applications of the two-stage CR approach presented here.

Materials and Methods

In this section, we describe the equations used to model the neuronal and synaptic dynamics, the characteristics of the RVS and the SVS CR stimulation signals including the equations used to model the CR-stimulation, additional simulation details as well as the methods used for data analysis.

Conductance-Based Hodgkin-Huxley Neuron Model

The conductance-based Hodgkin-Huxley neuron model (Hodgkin and Huxley, 1952) used in this work is that of Hansel et al. (1993). For clarity we briefly recapitulate part of their model here. The differential equation of the membrane potential is given by:

$\begin{array}{l}C\frac{dV}{dt}=I-g_{Na}{m}^{3}h\left(V-V_{Na}\right)-g_K{n}^4\left(V-V_K\right)\hfill \\ -g_l\left(V-V_l\right)+I(t).\hfill & (1)\end{array}$

The letters C, V, t, g, I, and I(t) denote membrane capacitance, membrane potential, time, maximum conductance for the ion or leak channel, the constant depolarizing current injected into the neuron which fixes the intrinsic firing rate of the neuron, and a time-varying current. The latter will be described in Equations (4), (10), and (11). The sodium and potassium reversal potentials are V_Na = 50 mV and V_K = −77 mV, and the leak reversal potential is V_l = −54.4 mV. The maximum conductances per unit area for the sodium, potassium, and leak currents are g_Na = 120 mS/cm², g_K = 36 mS/cm², and g_l = 0.3 mS/cm². Other values used in this article are C = 1 F/cm² and I ∈ [I₀ − ε_I, I₀ + ε_I] with I₀ = 11.0 A/cm² and ε_I = 0.45 A/cm². See Hansel et al. (1993) or Popovych and Tass (2012) for the equations of the time-varying gate variables m, h, and n. In this study the initial conditions of all N neurons were randomly drawn from uniform distributions (n_i, m_i, h_i, s_i ∈ [0, 1];V_i ∈ [−65, 5]mV; I_i ∈ [I₀ − σ_I, I₀ + σ_I]). We define the “spike times” of the neuron as times at which V = 0.0 mV and dV∕dt < 0 mV/ms.

Network

We use the network model described in Popovych and Tass (2012), which we briefly summarize here. The 1D neural network consists of N = 200 conductance-based Hodgkin-Huxley neurons and has periodic boundary conditions to minimize boundary effects. The neurons are initially all-to-all coupled with short-range strong excitatory and long-range weak inhibitory synaptic couplings according to a Mexican Hat (Wilson and Cowan, 1973; Dominguez et al., 2006; De la Rocha et al., 2008) like a weighting factor M_ij for the synaptic coupling from neuron j to i. The spatial profile of coupling is given by:

$\begin{array}{rc}{M}_{ij}=\left(1-d_{ij}^2∕{\sigma }_1^2\right)exp\left(-d_{ij}^2∕\left(2{\sigma }_2^2\right)\right),& (2)\end{array}$

where d_ij = d|i − j| denotes the distance between neurons i and j, and

$\begin{array}{rc}d=d_0∕\left(N-1\right),& (3)\end{array}$

the lattice distance between adjacent neurons within the ensemble. The length of the neuronal chain is given by d₀ (d₀ = 10), σ₁ = 3.5, and σ₂ = 2.0 as used in Popovych and Tass (2012).

Each neuron i receives post-synaptic currents (PSCs) from the other neurons in the network. The weighted ensemble average of all PSCs received by neuron i is (part of) the time-varying current I(t) (see Equation 1) and given by the coupling term S_i (Popovych and Tass, 2012)

$\begin{array}{rc}{S}_i=N^{-1}{\displaystyle \sum _{j=1}^N}\left(V_{r,j}-V_i\right)c_{ij}\left|M_{ij}\right|s_j,& (4)\end{array}$

where N is the number of neurons, V_{r, j} is the reversal potential of the synaptic coupling (20 mV for excitatory and - 40 mV for inhibitory synapses), V_i is the membrane potential of neuron i, and c_ij is the synaptic coupling strength from neuron j to neuron i. There are no self-connections within the network. In this study the initial synaptic weights c_ij between the neurons were drawn from a normal distribution ($c_{ij}~N(\mu =0.5\mu {\text{A}∕\text{cm}}^2,\sigma =0.01\mu {\text{A}∕\text{cm}}^2)$)

The synaptic variable s_j is given by

$\begin{array}{rc}\frac{d{s}_j}{dt}=\frac{0.5(1-s_j)}{1+exp\left[-\left(V_j+5\right)∕12\right]}-2s_j.& (5)\end{array}$

Spike Timing-Dependent Plasticity

In general, synapses are dynamic and change their connection strengths depending on the precise timing of pre-and post-synaptic spikes (Markram et al., 1997; Bi and Poo, 1998). All synapses within the network are sensitive to the time difference (△t_ij) between post- and pre-synaptic spike times t_i, respectively, t_j. According to the spike timing-dependent plasticity (STDP) learning rule (Bi and Poo, 1998) the changes in connection strength is given by

$\begin{array}{ll}△c_{ij}=\{\begin{array}{cc}{\beta }_1{e}^{\frac{-△t_{ij}}{{\gamma }_1\tau }},& △t_{ij}\ge 0\\ {\beta }_2\frac{△t_{ij}}{\tau }{e}^{\frac{△t_{ij}}{{\gamma }_2\tau }},& △t_{ij}<0\end{array}& (6)\end{array}$

(see Popovych and Tass, 2012) where △t_ij = t_i − t_j is the difference between post-synaptic spike time t_i and pre-synaptic spike time t_j. Synaptic weights c_ij are updated in an event-based manner by adding δ·△c_ij for excitatory connections and −δ·△c_ij for inhibitory connections with learning rate δ > 0 every time a neuron spikes. To avoid an unbounded strengthening or weakening, the synaptic weights are restricted to the interval c_ij ∈ [0, 1] mS/cm². In this study we use values for the STDP parameters as in Popovych and Tass (2012): β₁ = 1, β₂ = 16, γ₁ = 0.12, γ₂ = 0.15, τ = 14 ms, and δ = 0.002.

Coordinated Reset

The CR stimulation was delivered to the network in a spatiotemporal manner via N_s equidistantly spaced stimulation sites (Tass, 2003a). These stimulation sites are activated one after the other, such that after one stimulation ON-cycle of duration T_s each stimulation site is stimulated exactly once. This spatiotemporal activation of stimulation sites is represented by the indicator functions ρ_k(t) (kϵ {1, …, N}):

$\begin{array}{ll}{\rho }_k\left(t\right)=\{\begin{array}{cc}1,& k^{th}\text{stimulation}\text{site}\text{is}\text{active}\text{at}t\\ 0,& \text{otherwise}\end{array}.& (7)\end{array}$

These CR stimulation signals induce single brief excitatory PSCs, which evoke post-synaptic normalized conductances in the post-synaptic membrane. The conductances are represented by α-functions (Popovych and Tass, 2012 given by):

$\begin{array}{rc}{G}_{stim}\left(t\right)=\frac{t-t_k}{{\tau }_{stim}}{e}^{-(t-t_k)∕\tau },t_k\le t\le t_{k+1},& (8)\end{array}$

where τ_stim = T_s ∕ (6N_s) is the time-to-peak of G_stim, and t_k is the onset of the k^th activation of the stimulation site. The spatial spread of the induced excitatory PSCs in the network is defined by a quadratic spatial decay profile (see Popovych and Tass (2012) for motivation):

$\begin{array}{rc}D\left(i,x_k\right)=\frac{1}{1+d^2{(i-x_k)}^2∕{\sigma }_d^2},& (9)\end{array}$

were d denotes the lattice distance between two neighboring neurons as defined in Equation 3, i is the index of neuron i, x_k is the index of the neuron at stimulation site k, and σ_d = 0.08 d₀ denotes the spatial decay rate of the stimulation current.

The current evoked in neuron i by the CR stimulation with stimulation intensity K, is defined like in Popovych and Tass (2012):

$\begin{array}{rc}{F}_i=\left[V_r-V_i(t)\right]·K{\displaystyle \sum _{k=1}^{N_s}}D\left(i,x_k\right){\rho }_k\left(t\right)G_{stim}\left(t\right),& (10)\end{array}$

where V_r = 20 mV denotes the excitatory reverse potential, V_i the membrane potential of neuron I, and ρ, G, D are given by Equations (7), (8), and (9). This implies that during the period with CR stimulation (briefly CR-on period) the total time-dependent current I(t) of neuron i in Equation 1 is the sum of the coupling term and the stimulation current

$\begin{array}{rc}{I}_i\left(t\right)=S_i\left(t\right)+F_i\left(t\right).& (11)\end{array}$

During the period after the CR stimulation is switched off (briefly CR-off period) F(t) = 0. In this paper we use the RVS CR (Tass and Majtanik, 2006; Tass and Hauptmann, 2006; Popovych and Tass, 2012; Tass and Popovych, 2012), and the SVS CR (Zeitler and Tass, 2015). A stimulation site sequence, briefly sequence, is the spatiotemporal sequence by which all stimulation sites are stimulated exactly once. In the RVS CR approach, the sequences are drawn randomly (see Figure 1A). In the SVS-n CR approach a randomly drawn sequence is consecutively repeated n-times before the next sequence is randomly drawn (see Figure 1B for n = 4). In this paper we use n = 100. Since n is constant in this study we will write SVS instead of SVS-100 for the remainder of this study. The stimulation sites used in this study are located at neuron indices 25, 75, 125, and 175.

FIGURE 1

Figure 1. Spatiotemporal stimulation signals of CR stimulation. An example sequence order for the (A) RVS CR stimulation with a new randomly drawn sequence for each stimulation ON-cycle and for the (B) SVS CR stimulation with every randomly drawn sequence repeated during four consecutive stimulation ON-cycles (SVS-4) before the next sequence is drawn. Vertical dashed lines separate stimulation ON-cycles. Newly drawn sequences are represented by a different color. In this study the CR-on periods exclusively contain stimulation ON-cycles: no stimulation OFF-cycles. For the remainder of this study each sequence is repeated consecutively 100 times in the SVS CR approach.

Simulation Details

After an initial equilibration period of 2 s, STDP was included for the rest of the simulation. During the first 60 s with STDP the network rewires its connections. At the end of this STDP-only period the network activity is highly synchronized and at t = 0 ms the CR simulation is applied for the duration of the total CR-on period, T_{CR − on}. In this study this CR-on period consists only of stimulation ON-cycles, each lasting T_s = 16 ms. After T_{CR − on} the single-stage CR-stimulation is stopped permanently and the CR-off period starts. The short single-stage CR-stimulation lasts for T_{CR − on} = 64 s, the other single-stage CR-stimulation are twice as long (T_{CR − on} = 128 s). The CR-off period lasts as long as the CR-on period.

For comparison a sham stimulation is used as control stimulation: during sham stimulation no stimulus is delivered, so that this control condition is just a continuation of the STDP-only period since the stimulation intensity, K equals zero.

For the two-stage CR-stimulation the first CR-on period starts at t = 0 s and lasts for 64 s. At t = 64 s the first CR stimulation is switched off and the second CR-stimulation is applied. At t = 128 s the second CR-stimulation is then switched off, and the 128 s lasting CR-off period starts. The stimulation intensity of the CR-stimulation during the i-th CR-on period is denoted as K_i, the CR-approach during the i-th CR-on period of the two-stage CR as RVS_i if RVS CR was applied and SVS_i in case of SVS CR (with ∈ {1, 2}). For the single-stage CR it is also possible to have two CR-on periods with different stimulation intensities K₁ and K₂. However, by definition only one CR approach is used in the single-stage CR and therefore, the CR-approach is denoted without index i, e.g., simply as RVS or SVS CR.

The initial network conditions as well as the sequence order have an influence on the anti-kindling (Zeitler and Tass, 2015). Therefore, eleven simulations are executed for the same CR-protocol and stimulation intensity, K, but for different combinations of initial network conditions and sequence orders, hereafter referred to as “samples.” The obtained distributions for C_av and R_av are represented in boxplots (Tukey, 1977).

All simulations were done in Matlab R2007a. The differential equations were solved by the built-in function ODE45 with a relative tolerance of 10⁻⁵.

The present study is part of a series of studies addressing different topics and issues investigated with the same numerical model. As soon as the related papers will be in press, the Matlab Code will be uploaded to ModelDB (https://senselab.med.yale.edu/ModelDB/).

Data Analysis

In plastic neural networks CR stimulation can lead to anti-kindling, i.e., an unlearning of abnormal synaptic connectivity and abnormal synchronization. In this study the synaptic weights change constantly due to STDP and the different intrinsic firing rates of the neurons. Therefore, the dynamics of the synaptic connectivity is monitored on a population level by the synaptic weight averaged over the population:

$\begin{array}{rc}{C}_{av}\left(t\right)=N^{-2}{\displaystyle \sum _{i,j}}sgn\left(M_{ij}\right)c_{ij}\left(t\right),& (12)\end{array}$

where M_ij as defined in Equation 2, sgn denotes the sign-function and c_ij is the synaptic coupling strength from neuron j to neuron i.

The degree of synchronization of the neuronal activity is influenced by synaptic weights and can be measured by the order parameter R (Haken, 1983; Kuramoto, 1984)

$\begin{array}{rc}R\left(t\right)\text{exp}\left[i\Phi \left(t\right)\right]=N^{-1}{\displaystyle \sum _{j=1}^N}\text{exp}\left[i{\phi }_j\left(t\right)\right],& (13)\end{array}$

where φ_j(t) = 2π(t − t_{j, m})∕(t_{j, m + 1} − t_{j, m}) for t_{j, m} ≤ t<t_{j, m + 1} is a linear approximation of the phase of neuron j between its m^th and (m + 1)^th spikes at spike times t_{j, m} and t_{j, m + 1}. A phase in case of a rhythmically active neuron can be approximated along the lines of the phase oscillator concept (Kuramoto, 1984; Tass, 1999) or simply by a piece-wise linear function of time (Rosenblum et al., 2001). Φ(t) is the circular mean phase of the entire group of N neurons. An R-value of zero indicates a complete desynchronization and a value of one indicates perfect in-phase synchronization. For our data analysis the order parameter is averaged over the last 1.6 s of the CR-off period and is called the average order parameter R_av.

Statistical significances of differences between the obtained distributions for C_av and R_av (represented in boxplots) are determined by the one-sided Mann-Whitney-U- test with p < 0.05.

To investigate the mechanisms of RVS and SVS CR stimulation we study the stimulus-locked dynamics on a mesoscopic scale by considering subgroups of neurons as given by their proximity to the four different stimulation sites. The latter are located at neuron indices 25, 75, 125, and 175. Accordingly, we divide the whole population into 4 subgroups comprising N/N_s = 50 neurons each, where the number of stimulation sites is given by N_s = 4. The mean phase of subgroup sg, Φ_sg(t) is determined by

$\begin{array}{ll}{R}_{sg}\left(t\right)\text{exp}\left[i{\Phi }_{sg}\left(t\right)\right]=N_{sg}^{-1}{\sum }_{j=1}^{N_{sg}}\text{exp}\left[i{\phi }_j\left(t\right)\right],& (14)\end{array}$

where the summation runs only over the neurons within subgroup sg, and sg = 1, …, 4. We perform a cross-trial analysis of the stimulus responses of each subgroup separately. In particular, we focus on the stimulus-locked phase dynamics of the subgroups, to reveal e.g., phase resetting or entrainment processes that are tightly time locked to the repetitively administered stimuli. To study the phase dynamics of a subgroup in an ensemble of stimulus-locked responses for each time Δt in a time window W = [–32 ms, +32 ms] attached to each stimulus onset τ_k we determine the cross-trial distributions of the phases of the different subgroups

$\begin{array}{ll}{\left\{{\Phi }_{sg}\left({\tau }_k^{sg}+\Delta t\right)\left(\text{mod }2\pi \right)\right\}}_{k=1,\dots ,L}& (15)\end{array}$

where ${\tau }_k^{sg}$ denotes the onset of the k-th stimulus delivered to subgroup sg, and L is the total number of stimuli delivered to the subgroup sg (Tass, 2003d). These cross-trial distributions will be displayed in a color plot as a function of △t ∈ W. In our data we only found phase distributions that were either close to uniform or had (only) one pronounced peak. Accordingly, to quantify the amount of stimulus-locking of the phase dynamics, we use the resetting index E_sg(△ t) of subgroup sg defined by

$\begin{array}{ll}{E}_{sg}\left(△t\right)=\left|L^{-1}{\displaystyle {\sum }_{k=1}^L\text{exp}\left[i{\Phi }_{sg}\left({\tau }_k^{sg}+\Delta t\right)\right]}\right|,& (16)\end{array}$

where L is the number of administered stimuli (Tass, 2003d; see also Tallon-Baudry et al., 1996). If the phase dynamics is not stimulus-locked at time △t, the corresponding phase distribution is uniform. Otherwise, one observes a distribution that significantly differs from a uniform distribution and has, e.g., one dominant peak. For instance, in case of a stimulus-induced phase reset the phase distribution will turn from uniform in the pre-stim period to unimodal after stimulus onset (Tallon-Baudry et al., 1996; Tass, 2003d). In contrast, a permanent stimulus-locking of the phase dynamics in terms of a phase entrainment can show up as unimodal distribution throughout the whole time window W.

As a motivation for the subsequent analysis, let us consider possible differences of the mechanisms and effects of RVS CR and SVS CR:

SVS CR (during an epoch with a constant sequence): Each subgroup is stimulated periodically. Periodic stimulation of an isolated (i.e., not synaptically interacting) subgroup may cause an entrainment that may lead to an increase of the synchronization, provided stimulation is strong enough. SVS CR (with constant sequence) delivered to mutually isolated subgroups may, thus, cause a phase-shifted entrainment of the different subgroups at high enough stimulation intensities. The mutual synaptic interaction of the different subgroups may cause a perturbation of the respective entrainment processes. Note, the effect of periodic stimulation on an isolated subgroup may be complex, as shown e.g., in the context of periodic stimulation of neurons interacting through excitatory synapses (Popovych and Tass, 2011): Depending on the mismatch of intrinsic frequencies and stimulation frequency periodic stimulation may cause a desynchronization—not just an entrainment. Accordingly, it is not justified to simplify the mechanism of SVS CR as being a combination of exclusively synchronizing direct effects of stimulation and desynchronizing indirect, synaptically mediated effects.

RVS CR: The direct stimulation effect on each subgroup is fundamentally different compared to SVS CR, since it is not a periodic stimulation.

By design, SVS CR and RVS CR stimulation affect subgroups differently: Accordingly, we study synchronization processes not only on the macroscopic level of the entire neuronal population, but also on the mesoscopic level of the different subgroups. To assess the amount of synchronization of the j-th subgroup we introduce

$\begin{array}{l}{\overline{R}}_j^{pre}=<R_j{>}_{last10sbeforeCR-on}\\ {\overline{R}}_j^{on}=<R_j{>}_{last10sCR-on}\\ {\overline{R}}_j^{off}=<R_j{>}_{first10sCR-off}\end{array}$

where R_j denotes the synchronization order parameter of the j-th subgroup as determined by Equation 14. The acute CR effect on the j-th subgroup is given by

$\begin{array}{ll}1-\frac{{\overline{R}}_j^{on}}{{\overline{R}}_j^{pre}}& (17)\end{array}$

whereas the acute after-effect on the j-th subgroup reads

$\begin{array}{ll}1-\frac{{\overline{R}}_j^{off}}{{\overline{R}}_j^{pre}},& (18)\end{array}$

where an outcome of zero means no acute CR effect (after-effect, respectively), a positive outcome indicates a desynchronizing and a negative value a synchronizing effect on the j-th subgroup due to CR stimulation.

Since CR stimulation induces phase shifts between the different subgroups, we study synchronization processes between the different subgroups by calculating the phase difference between the j-th subgroup and the k-th subgroup defined by

$\begin{array}{ll}{\vartheta }_{jk}\left(t\right)={\Phi }_j\left(t\right)-{\Phi }_k\left(t\right),& (19)\end{array}$

with Φ_sg(t) the mean phase of subgroup sg as determined by Eq. 14. We will determine the distribution of ϑ_jk(t) during the last 10 s with CR-on as well as during the first 10 s with CR-off and show this distribution by an angular distribution plot with bins of π∕18 rad.

To detect whether a distribution of ϑ_jk has several peaks we introduce the indices

$\begin{array}{cc}{\lambda }_{jk}^{\left(\nu \right)}=\left|\frac{1}{P}{\displaystyle {\sum }_{l=1}^P\text{exp}\left[i\nu {\vartheta }_{jk}\left(t_l\right)\right]}\right|,& \left(20\right)\end{array}$

where P is the number of sampling points in the 10 s window (P = 10,000) and ν is the order of λ (Batschelet, 1981; Kuramoto, 1984; Tass, 1999, 2003d). Perfect phase synchronization between the j-th and k-th subgroup with ϑ_jk(t) = const is associated with ${\lambda }_{jk}^{\left(\nu \right)}=1$ for ν = 1, 2, 3, … . Conversely, a perfect uniform distribution of ϑ_jk is related to ${\lambda }_{jk}^{\left(\nu \right)}=0$ for ν = 1, 2, 3, … . A bimodal distribution of ϑ_jk with anti-phase peaks can be detected with ${\lambda }_{jk}^{\left(2\right)}-{\lambda }_{jk}^{\left(1\right)}$ (Tass, 2003d). Given the numerical results obtained in the present study (see Results), we seek to detect distributions with three rectangularly aligned peaks located at π∕2, π and3π∕2. To this end, we introduce the index

$\begin{array}{ll}{\alpha }_{jk}={\left[{\lambda }_{jk}^{\left(4\right)}-{\lambda }_{jk}^{\left(1\right)}\right]}_{+},& (21)\end{array}$

where [z]₊, stands for the half–wave rectification operation: [z]₊ = z if z > 0 and [z]₊ = 0 else.

To motivate the choice of α_jk and study its characteristics, let us consider two different families of synthetic phase difference distributions with (i) four equidistantly aligned peaks at 0, π∕2, π and 3π∕2 and (ii) three rectangularly aligned peaks located at π∕2, π and 3π∕2.

(i) 4 equidistantly aligned peaks: The distribution of the phase differences is given by ${\left\{{\psi }_l(m)\right\}}_{l=1}^{4Z}$, where m = 0, …, Z serves as parameter. The phase differences read ψ_l(m) = ϱ_l(m) + ξ_l, with random variable ξ_l uniformly distributed and fulfilling −ε ≤ ξ_l ≤ ε, where the parameter ε determines how narrow the distribution is, and l = 0, …, 4Z. The symmetrical 4-peak distribution is created by

$\begin{array}{ll}{\varrho }_l\left(m\right)=\left({\left[l\right]}_{\text{mod}4}-1\right)\frac{\pi }{2}{\mu }_l(m)& (22)\end{array}$

where

$\begin{array}{ll}{\mu }_l\left(m\right)=\{\begin{array}{ll}0:\hfill & l>4m\hfill \\ 1:\hfill & \text{else}\hfill \end{array}& (23)\end{array}$

determines which index groups of four actually contribute to a 4 peak distribution (in the case of μ_l(m) = 1) or to a 1 peak distribution instead (in the case of μ_l(m) = 0). The limiting cases are a perfect 1-peak (Dirac-type) distribution for m = 0 and ε = 0 and a perfectly symmetrical distribution with 4 rectangularly arranged peaks for m = Z and ε = 0. For each given m we determine the parameters

$\begin{array}{cc}{\lambda }^{\left(\nu \right)}(m)=\left|\frac{1}{4Z}{\displaystyle {\sum }_{l=1}^{4Z}\text{exp}\left[i\nu {\psi }_l\left(m\right)\right]}\right|& \left(24\right)\end{array}$

for ν = 1 and ν = 4 and calculate

$\begin{array}{ll}\alpha \left(m\right)={\left[{\lambda }^{(4)}\left(m\right)-{\lambda }^{(1)}(m)\right]}_{+}& (25)\end{array}$

For illustration we plot λ⁽¹⁾(m) vs. x = m∕Z and α(m) vs. x = m∕Z (Figure 2A).

FIGURE 2

Figure 2. Detection of rectangulary aligned peaks in a synthetic phase difference distribution. (A) The first order indices λ⁽¹⁾ (Equation 24) as well as the indices α (Equation 25) are shown as a function of x for the 4-peak phase difference distribution with different narrowness, ε, of the peaks. (B) Indices as in (A), but determined by Equations 28 and 29 for the 3-peak phase difference distributions. Two extreme cases are shown: (i) perfect Dirac-type peaks (ε = 0); (ii) a homogenous distribution without any pronounced peaks (ε = π∕4). We used Z = 100 for both panels. See section Data analysis for more details.

(ii) Asymmetric 3-peak distribution: For comparison, we consider the distribution given by ${\left\{{\psi }_l(m)\right\}}_{l=1}^{3Z}$, with parameter m = 0, …, Z. The phase differences are given by ψ_l(m) = ϱ_l(m) + ξ_l, with l = 0, …, 3Z. The uniformly distributed random variable ξ_l is defined as above: −ε ≤ ξ_l ≤ ε. The rectangularly aligned 3 peaks of the distribution are generated by

$\begin{array}{ll}{\varrho }_l\left(m\right)=\left({\left[l\right]}_{\text{mod }3}-1\right)\frac{\pi }{2}{\mu }_l(m)& (26)\end{array}$

where

$\begin{array}{ll}{\mu }_l\left(m\right)=\{\begin{array}{ll}0:\hfill & l>3m\hfill \\ 1:\hfill & \text{else}\hfill \end{array}& (27)\end{array}$

determines which index groups of three actually contribute to a 3-peak distribution [in the case of μ_l(m) = 1] or to a 1 peak distribution instead [in the case of μ_l(m) = 0]. Analogously to (i), the limiting cases are given by a perfect 1-peak (Dirac-type) distribution for m = 0 and ε = 0 and an asymmetric distribution with 3 peaks at π∕2, π, and 3π∕2 for m = Z and ε = 0. For each given m we determine the parameters

$\begin{array}{cc}{\lambda }^{\left(\nu \right)}(m)=\left|\frac{1}{3Z}{\displaystyle \sum _{l=1}^{3Z}\text{exp}\left[i\nu {\psi }_l\left(m\right)\right]}\right|& \left(28\right)\end{array}$

for ν = 1 and ν = 4. With this we obtain

$\begin{array}{ll}\alpha \left(m\right)={\left[{\lambda }^{(4)}\left(m\right)-{\lambda }^{(1)}(m)\right]}_{+}& (29)\end{array}$

For comparison, we plot λ⁽¹⁾(m) vs. x = m∕Z and α(m) vs. $x=\frac{m}{Z}$(Figure 2B).

From Equations (24) and (28) it follows that in the case of 4 equidistantly aligned peaks as well as in the case of an asymmetric 3-peak distribution for ε = 0 the maximum of λ⁽¹⁾ is attained for x = 0 and reads λ⁽¹⁾ = 1. For ε = 0 the minimum of λ⁽¹⁾ is reached at x = 1 and reads 0 in case of the symmetric 4-peak distribution and 1∕3 for the asymmetric 3-peak distribution. Hence, from Equations (25) and (29) we read off the maximum of α is reached for x = 1 for both distributions, but differs in size: α = 1 for the symmetric 4-peak distribution, but α = 2∕3 for the asymmetric 3-peak distribution (see Figure 2). Conversely, the minimum of α is 0 and attained for x = 0. Figure 2 illustrates how λ⁽¹⁾ and α vary with increasing ε. A limiting case is ε = π∕4, which renders both distributions uniform, so that α = 0 for 0 ≤ x ≤ 1.

Results

Short Single-Stage CR-Stimulation

Without stimulation (K = 0.0) the average synaptic weight, C_av, does not change, but it decreases by the RVS CR stimulation during the CR-on period (Figure 3A). After switching off the CR stimulation the network evolves spontaneously. During this CR-off period the average synaptic weight remains more or less the same for certain stimulation intensities, for others it increases. Since the same sequence order is applied to the same initial network, the differences observed are only caused by the different stimulation intensities, K.

FIGURE 3

Figure 3. Dynamics of the average synaptic weight, C_av, for different stimulation intensities, K, for the short RVS and short SVS CR-stimulations. (A) C_av as a function of time for the RVS CR with different stimulation intensities, K. (B) C_av dynamics for the SVS CR stimulation with different K. The red horizontal bar represents the CR-on period, which start at t = 0 s and is switched off at t = 64 s (dashed vertical line). No stimulation signals are delivered during the subsequent 64 s lasting CR-off period and C_av evolves spontaneously. The initial network is identical for all simulations. The sequence order of the RVS (SVS, respectively) CR stimulation is identical for all K-values.

SVS CR with 100 consecutive repetitions also decreases the average synaptic weight, C_av, for the stronger stimulation intensities, but for weak stimulation intensities it increases the average synaptic weight, C_av, during the CR-on period (Figure 3B). For K = 0.05 and K = 0.10 the C_av dynamics are even almost identical and still different compared to the case without stimulation (K = 0.0). Again, the differences in Figure 3B are due to the different values of K, since one SVS sequence order is applied to the same initial network at different stimulation intensities. Since the initial network used for Figures 3A,B is identical, these figures suggest that for this particular initial network RVS CR is more effectively decreasing C_av for weak stimulation intensities, whereas SVS CR is more effective at stronger stimulation intensities.

Repeating the simulations for eleven samples (different combinations of initial networks and sequence orders) confirms that for K = 0.10 the RVS CR stimulation induces smaller C_av values than the SVS CR stimulation, and that for K ≥ 0.20 the SVS CR induces smaller C_av values than the RVS CR stimulation. Both reductions of C_av are statistically significant (one-sided Mann-Whitney test, p < 0.05; Figure 4A). The average order parameter, R_av, is smaller for the RVS than for the SVS CR for K = 0.10, and again this difference is statistically significant (one-sided Mann-Whitney test, p < 0.05; Figure 4B). Combining Figures 4A,B shows that for a stimulation intensity of K = 0.10 only the RVS CR approach can decrease the average synaptic weight and reduce the degree of synchrony of the neuronal activity. For stimulation intensities K ≥ 0.20 the SVS approach is superior with respect to a reduction of the mean synaptic connectivity. However, this does not induce a significantly stronger reduction of the mean synaptic connectivity.

FIGURE 4

Figure 4. Comparison of the average synaptic weight, C_av, and the average order parameter, R_av, at the end of the CR-off period (t = 128 s) for the short RVS and the short SVS CR stimulation. (A) Boxplot of C_av(t = 128 s) for different stimulation intensities, K, for the RVS (red) and SVS CR stimulation. (B) As (A) for the order parameter averaged over the last 1.6 s of the CR-off period, R_av. For each condition (K value and CR approach) eleven samples are used. The distribution of the obtained C_av and R_av values are represented by a boxplot. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

To study the mechanism by which CR stimulation desynchronizes the neuronal activity we compare the spike times of the individual neurons at the end of the CR-on period for a stimulation intensity K = 0.25, for which RVS as well as SVS CR stimulation work well, with our control stimulation (K = 0.00). The control stimulation (K = 0.00) “applied” to the initially synchronized network does not change the form of the spike volleys of the neurons (see Figures 5A–D). This is to demonstrate that the model does not change its dynamic mode spontaneously. By applying CR stimulation at an intensity K > 0 the network is stimulated in a spatiotemporal way: stimulation at different positions at different times as shown by the red diamonds in Figure 5. Stimulation (with K > 0) induces post-synaptic currents in the neurons. The amplitude of this induced current depends not only on the stimulation intensity, but due to the spatial decay profile (see Equation 9) also on the distance between the neuron and the stimulation site. Therefore, the spiking behavior of neurons near a stimulation site is more strongly affected as opposed to distant neurons. For our control stimulation the spike times of the neurons show some variations. The neurons are not in perfect synchrony, but display a certain time jitter. Due to the CR stimulation the spiking jitter increases, in particular, in neurons with conflicting impact from neighboring stimulation sites (see Figure 5A for RVS and Figure 5B for SVS CR stimulation). Since the inter-spike-intervals of strongly synchronized neurons are (almost) identical, this jitter can only be reduced if the neurons are all reset to a certain internal phase. To which state a neuron is reset depends on the timing of the stimulation in relation to the internal cycle of the neuron's dynamics. An excitatory stimulation of a Hodgkin-Huxley model neuron can advance or delay a spike: an excitatory stimulation received in the first part of the internal cycle will delay the next spike, in the last part of the internal cycle it will advance the spike time. This implies that a spike volley with some jitter can be split in two groups, usually these two groups will merge under the influence of further stimulation (see e.g., the neurons near stimulation site 175 in Figure 5A). So, the main reason why CR desynchronizes a neuronal population is because stimulation delivered at different sites within the network at different times divides the population into several phase-shifted synchronized subpopulations (located close to the stimulation sites) and groups of desynchronized neurons (located in between stimulation sites; see Figure 5A for the RVS CR stimulation and Figure 5B for the SVS CR stimulation).

FIGURE 5

Figure 5. Effect of RVS and SVS CR stimulation on the spike times at the end of the 64 s lasting CR-on period.(A) Raster plot of the last 150 ms of the CR-on period with the RVS CR applied with stimulation intensity K = 0.25 (blue dots) and for the control stimulation with K = 0.00 (cyan dots). Each blue (cyan) dot represents the spike time of that particular neuron under the condition K = 0.25 (K = 0.00, respectively). (B) as (A) for SVS CR. (C) As (A) for K = 0.10. (D) as (B) for K = 0.10. The stimulation onset of the stimulation sites are represented by red diamonds in (A–D). (E) Spike counts for the last 150 ms of the RVS CR-on period show how many neurons fire within each time interval of 1 ms. (F) as (E) for the SVS CR stimulation. (G–I) Distributions of the mean phase of neurons 51–100, Φ_{sg = 2}(Δt), averaged across the 4000 activation times of subgroup 2 within the RVS CR-on period for different applied stimulation intensities. (J) Resetting index E(Δt) of subgroup 2 determined for the RVS CR-on period for different stimulation intensities. (K–M) as (G–I) for the SVS CR stimulation. (N) as (J) for the SVS CR stimulation. For all panels are the same initial network conditions and sequence orders as used in Figure 3.

By decreasing the stimulation intensity the induced synchronized subpopulations which are in general out-of-phase with each other, become smaller and will disappear first for the SVS and for weaker stimulation intensities also for the RVS CR stimulation. For K = 0.10 only the RVS CR is still able to divide the whole population into some desynchronized and some synchronized groups (Figure 5C). For the SVS CR stimulation a stimulation intensity of K = 0.10 is too weak to divide the whole population into distinct, phase-shifted groups (Figure 5D).

The spike count plots in Figures 5E,F show the number of spikes that are fired by the network within each time bin of size 1 ms. For the control stimulation there are high peaks, which can be reduced by the RVS CR stimulation applied at K = 0.10 (Figure 5E), but the SVS CR stimulation applied at the same intensity does not significantly reduce the overall synchronization (Figure 5F). However, for intermediate stimulation intensities like K = 0.25, both CR approaches desynchronize the neuronal network activity (Figures 5E,F)

The analysis with the cross-trial phase distributions, from Equation (15), reveals a qualitative difference of the effects induced by RVS as opposed to SVS CR stimulation: for sufficient stimulation strength the non-periodically delivered RVS stimuli cause phase resets, whereas the periodically delivered SVS stimuli induce an entrainment. The cross-trial phase distribution in a window from −32 to +32 ms relative to the stimulation onset of the site belonging to the neuronal subgroup comprising neurons 51–100 shows that for the sham stimulation (K = 0.0) the phase distribution is flat throughout the entire time window (Figure 5G). As the stimulation intensity of the RVS CR increases, the post-stimulus cross-trial phase distribution displays a feature characteristic of a phase reset: irrespective of the pre-stimulus phase at stimulation onset, the stimulus resets the phase to a preferred value (Figures 5H,I). The pre-stimulus resetting index E from Equation 14 is close to zero irrespective of the stimulation intensity, whereas the post-stimulus resetting index increases more for increasing stimulation intensities (Figure 5J). The SVS CR stimulation shows a different phenomenon. For sufficiently strong intensities the phase of the stimulated subgroup gets increasingly, and finally tightly stimulus-locked, i.e., entrained across trials, throughout the entire stimulus-locked time window (Figures 5L,M). Accordingly, the resetting index does not show a difference between the pre- and the post-stimulus period (Figure 5N). Obviously, the SVS CR stimulation results in an entrainment of the subgroup to the SVS CR stimulation signal. The other three neuronal subgroups show similar effects (results not shown).

Since RVS and SVS CR affect the subgroups differently, we analyze the desynchronization of all subgroups by CR also for different stimulation intensities. Figure 6A illustrates that the acute RVS CR effect on the subgroups of all eleven initial networks as well as the acute after-effect are present already for K = 0.10. The acute effect of RVS CR peaks for K = 0.15 and slightly decreases for higher intensities. The acute after-effect of RVS CR builds up and saturates with increasing intensity. In contrast, for the SVS CR stimulation an acute effect and acute after-effect shows up only at K = 0.15. With increasing intensity the acute effect of SVS CR further decreases, whereas the acute after-effect peaks at K = 0.20 and slightly decreases above (Figure 6B). Another difference between the RVS and SVS CR stimulation is that for K = 0.20 and 0.25 the differences between the acute CR effect and the acute after-effect are larger for the SVS CR (Figure 6B).

FIGURE 6

Figure 6. Typical examples of the desynchronizing effect of the RVS and SVS CR stimulation on the mesoscopic level. (A) Boxplots (n = 44 samples) show the distributions of the acute CR effect (Equation 17) in red and of the acute after-effect (Equation 18) in blue, induced in all subgroups of all eleven networks by the RVS CR stimulation at different stimulation intensities K. The CR-on period lasts for only 64 s. (B) as (A) for the SVS CR stimulation. (C–G) Angular histogram plots of the phase differences ϑ₂₃ between subgroups 2 and 3 for different stimulation intensities K. Each plot shows the distribution of the 10.000 phase differences ϑ₂₃, obtained during the last 10 s of the 64 s lasting RVS CR-on period (n = 10.000 phase difference samples). (H–L) as (C–G) for the SVS CR stimulation within the same initial network.

In order to study the effect of CR stimulation on the phase-shifts between the different subgroups within the networks, we determined the phase differences between subgroup 2 and 3 for each ms during the last 10 s of the RVS as well as the SVS CR stimulation at different K-values in one network. Without stimulation the average and standard deviation of all 10,000 values of the phase difference ϑ₂₃ is 6.2379±0.0008 rad (results not shown). Accordingly, the corresponding circular distribution of the phase difference ϑ₂₃ has one prominent peak (where all values are comprised in one bin when displayed as in Figure 6, data not shown). The angular phase difference histograms in Figures 6C–G illustrate that the RVS CR induces a certain jitter of the phase differences ϑ₂₃ already for K = 0.05, and that the distribution is not homogeneous for any shown K-value. The SVS CR induces also some jitter in ϑ₂₃ at K = 0.05 (Figure 6H) but less than the RVS CR does (Figure 6C). Increasing the strength of K for the SVS CR (Figures 6H–L) induces a different pattern for ϑ₂₃ than for the RVS CR stimulation: the majority of phase bins has no or just a very little number of entries, and for K = 0.20 and 0.25 a clear asymmetric 3-peak distribution shows up. Interestingly, the circular distribution does not peak at 0, which indicates that the SVS CR stimulation actually causes a non-vanishing phase shift, predominantly at π∕2, π, and 3π∕2 rad. Note during the last 10 s of the 64 lasting SVS CR-on period, seven different sequences are applied and not just one. Accordingly, the network is subjected to different phase-shifted entrainment stimulation patterns which, in turn, gives rise to the different rectangularly aligned peaks. For the other ten networks and the other combinations of subgroups we found similar asymmetric 3-peaks distributions for K = 0.20 and K = 0.25, although the height of the peaks vary, due to the applied sequences during the last 10 s, but the position of the peaks are always at π∕2, π, and 3π∕2 rad (results not shown). For the phase differences between other subgroups we obtain very similar results (not shown).

The indices ${\lambda }_{23}^{\left(1\right)}$ from Equation 28 and α₂₃ from Equation 29 shed further light on the stimulation-induced changes of the synchronization between subgroups 2 and 3. In accordance with the angular histogram plots from Figures 6C–G, RVS CR stimulation slightly reduces the phase synchronization between subgroups 2 and 3 and, hence, ${\lambda }_{23}^{\left(1\right)}\text{for }K=0.05$ (Figure 7A). For K = 0.1 and, in particular, for K > 0.1 due RVS CR the phase synchronization between subgroups 2 and 3 is significantly reduced (Figure 7A). Only for K > 0.1 the desynchronization between subgroups 2 and 3 persists during the acute OFF phase (Figure 7C). RVS CR stimulation does not cause 3-peak distributions (with peaks at π∕2, π and 3π∕2 rad) as detected with the index α₂₃ from Equation 29. In contrast to RVS CR and in accordance with the corresponding angular histogram plots from Figures 6H–L, the phase synchronization between subgroups 2 and 3 and, hence, ${\lambda }_{23}^{\left(1\right)}$ hardly decreases for stimulation intensities of SVS CR up to K = 0.1. For K = 0.15 SVS CR stimulation causes a more pronounced reduction of the phase synchronization between subgroups 2 and 3 as evidenced by a reduced ${\lambda }_{23}^{\left(1\right)}$ with α₂₃ from Equation 29 being still close to zero (Figure 7B). For higher intensities K 3-peak distributions (with peaks at π∕2, π and 3π∕2 rad) emerge and get more pronounced, so that α₂₃ increasingly exceeds the zero line (Figure 7B). For K > 0.15 we observe pronounced a desynchronization between subgroups 2 and 3 in the acute off phase following cessation of stimulation (Figure 7D). As indicated by the vanishing index α₂₃, 3-peak distributions (with peaks at π∕2, π and 3π∕2 rad) are not detected after stimulation is turned off, irrespective of the stimulation intensity K (Figure 7D). Very similar results were obtained for the phase difference between the other subgroups (data not shown).

FIGURE 7

Figure 7. Detection of rectangulary aligned peaks in the distribution of the phase differences between subgroups 2 and 3 during the last 10 s of the 64 s lasting CR-on period as well as for the first 10 s of the subsequent CR-off period as a function of the stimulation intensity K. (A) The first order indices ${\lambda }_{23}^{\left(1\right)}$ (Equation 20) as well as the indices α₂₃ (Equation 21) are determined for the last 10 s of the 64 s lasting RVS CR-on period as a function of the stimulation intensity K. (B) As (A) for the SVS CR stimulation. (C) As (A) for the first 10 s of the CR-off period. (D) As (C) for the SVS CR stimulation. The same initial network conditions are used for all panels.

Single- and Two-Stage CR Stimulation with Weak Onset

RVS CR stimulation is already effective at weak stimulation intensities, whereas SVS CR stimulation requires higher intensities to unfold its superior efficacy. We specifically hypothesize that RVS CR stimulation might render the network more susceptible to SVS CR stimulation by inducing a preparatory desynchronization. In this way SVS CR stimulation might be even more effective and require shorter stimulation durations. To study whether this approach is actually specific, from a more general perspective, we study whether CR stimulation can be empowered by a preparatory CR stimulation, of RVS or SVS type, at weaker intensities. To investigate this idea, the first CR stimulation is applied with a stimulation intensity K₁ = 0.10 and the second CR stimulation is applied with a stimulation intensity K₂ = 0.15. Since both CR approaches are applied to the same initial network as used in Figure 3 as well as with the same sequence orders as used in Figures 3A,B, the effect of the preceding CR approach with weaker onset can be determined. The obtained results of the different CR combinations show that all combinations do influence the average synaptic weight C_av and suggest that C_av is reduced in two cases (Figure 8A). The strongest reduction is obtained for the two-stage CR-stimulation with weak onset and RVS₁ and SVS₂. This results in a smaller final C_av-value compared to the short single-stage SVS CR stimulation with intensity K = 0.15 as shown in Figure 3B. Also the long single-stage RVS CR-stimulation with weak onset (Figure 8A) induces a smaller final C_av-value than the short single-stage RVS CR-stimulation with intensity K = 0.15 as shown in Figure 3A. A reduction of the long-lasting desynchronizing effect is observed for the long single-stage SVS CR-stimulation (compare Figure 8A with Figure 3B for K = 0.15). Obviously the different network conditions at the start of the 64 s lasting SVS CR-on period with K = K₂ = 0.15 have a relevant effect on the stimulation outcome. The less advantageous network conditions at the start of the second SVS CR-on period (K₂ = 0.15) were caused by the preceding SVS CR-on period with K₁ = 0.10. This shows that a too weak stimulation intensity of the first SVS CR stimulation can, in fact, be unfavorable for the second SVS CR stimulation with a stronger stimulation.

FIGURE 8

Figure 8. Dynamics of C_av for the two-stage CR stimulation with weak onset. (A) Results for the same initial network conditions and sequence orders as used in Figure 3. (B) Results obtained for other initial network conditions and sequence orders than used in (A). The horizontal bar represents the two consecutive CR-on periods. During the first CR on period (0–64 s) the stimulation intensity equals K₁ = 0.10, during the second CR-on period (64–128 s) K₂ = 0.15, except for the no stimulation condition for which K₁ = K₂ = 0.00. The dashed vertical lines indicate the end of a CR-on period. No stimulation signals are delivered during the CR-off period (128–256 s) and C_av evolves spontaneously.

For the initial network and sequence orders as used in Figures 3, 8A, the two-stage CR-stimulation (RVS₁, SVS₂) with weak onset causes the strongest reduction of the C_av-value. The second strongest reduction is achieved by the single-stage RVS CR-stimulation with weak onset (Figure 8A). However, different initial network conditions or other sequence orders can lead to different results (Figure 8B). Therefore, all simulations are executed for eleven samples. All two-stage CR-stimulation (RVS₁, SVS₂) trials with weak onset revealed the most favorable results (one-sided Mann-Whitney test, p < 0.05; Figure 9). The two-stage CR-stimulation (RVS₁, SVS₂) with weak onset causes the strongest reduction of C_av-values for (K₁, K₂) = (0.10, 0.15), and the difference with the other combinations is statistically significant (one-sided Mann-Whitney test, p < 0.05), (Figure 9A). The resulting final R_av distributions are statistically not different for the two-stage CR stimulation and single-stage RVS CR-stimulation with weak onset. However, R_av caused by the two-stage CR stimulation with RVS₁ and SVS₂ is smaller than the other two combinations of the CR-approaches and the difference is statistically significant (one-sided Mann-Whitney test, p < 0.05; Figure 9B).

FIGURE 9

Figure 9. Comparison of anti-kindling effects of the CR-stimulations with weak onset. (A) Boxplots of the average synaptic weights, C_av, at t = 256 s for the four possible combinations of the two CR approaches (K₁ = 0.10, K₂ = 0.15) and the control stimulation (K₁ = K₂ = 0.00). (B) As (A) for the average order parameter, R_av, at t = 256 s. For each combination of CR approaches eleven samples were used. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

Figure 10 illustrates that the RVS CR is not so sensitive to different initial conditions as the SVS CR stimulation (K₁ = 0.10 and K₂ = 0.15). The second CR-on periods of Figures 10A,B started with spike timing patterns as shown in Figure 5C, those from Figures 10C,D with spike timing patterns as in Figure 5D. Despite these differences applying a RVS CR during the second CR-on period results in similar raster plots and spike counts (Figures 10A,C). If a SVS CR stimulation is applied during the second CR-on period, the different conditions at the beginning of this CR-on period will result in different raster plots and spike counts (Figures 10B,D). It seems that a slightly reduced synchronization, induced by the RVS₁ CR, makes it possible for the SVS₂ CR stimulation to desynchronize the neuronal population, where the phase of the (sub)group of neurons gets time locked to the stimulation onsets. The distributions of the mean phase of neuronal subgroup 2, Φ_{sg = 2}(△t) during the second CR-on period display patterns characteristic of phase resets for the RVS₂ and of phase entrainment for the SVS₂ CR stimulation (Figures 10E–H). The final phase reset is stronger if first RVS₁ is applied compared to SVS₁ CR stimulation (Figure 10I). The entrainment structure of the cross-trial phase distribution during the SVS₂ CR-on period is much more pronounced if RVS is applied during the first CR-on period compared to a SVS₁ CR stimulation (Figure 10I). So, the SVS₂ CR stimulation is much more sensitive to the (change in) initial conditions at the beginning of the second CR-on period than the RVS₂ CR stimulation and can obtain a resetting index comparable with the short (64 s lasting) single-stage SVS CR stimulations with a stronger stimulation intensity (compare Figures 5N at K = 0.20 with Figure 10I at K₂ = 0.15).

FIGURE 10

Figure 10. Effect of single-and two-stage RVS and SVS CR stimulations with weak onset on the spike times at the end of the second 64 s lasting CR-on period. (A) Raster plot and spike counts of the last 150 ms of the second CR-on period for the single-stage RVS CR stimulation (K₁ = 0.10, K₂ = 0.15; blue dots) and for the control stimulation with K = 0.00 (cyan dots). The stimulation onsets are represented in all panels by red diamonds. The spike counts are determined for each time bin of 1 ms. (B) as (A) for the two-stage (RVS₁, SVS₂), CR (K₁ = 0.10, K₂ = 0.15; blue dots). (C) As in (A) with blue dots for the spike times of the two-stage (SVS₁, RVS₂), CR (K₁ = 0.10, K₂ = 0.15). (D) As in (B) with blue dots for the spike times of the single-stage SVS CR (K₁ = 0.10, K₂ = 0.15) (E, G) Distributions of the mean phase of neurons 51–100, Φ_{sg = 2}(Δt), averaged across the 4000 activation times of stimulation site 2 (= neuron 75) within the RVS₂ CR-on period under the condition that first RVS₁ (E) or SVS₁ (G) was applied. (F, H) as (E, G) within the SVS2 CR-on period for RVS₁ (F) and SVS₁ (H). (I) Resetting index E(Δt) determined for the second CR-on period for the different combinations of RVS and SVS CR stimulations with K₁ = 0.10 and K₂ = 0.15. For all panels are the same initial network conditions and sequence orders as used in Figure 3.

Figure 11 sheds some light on the effects of the two-stage CR with weak onset on the average synaptic weight C_av. For a very weak stimulation intensity (K₁ = 0.10) during the first CR-on period, the SVS₁ CR-stimulation is not able to decrease the average synaptic strength, C_av (Figure 11A). But a RVS₁ CR-stimulation with K₁ = 0.10 decreases the initial C_av value (Figure 11A). The effectiveness of reducing C_av during the second CR-on period with K₂ = 0.15 depends on the C_av -value at the onset of this second CR-on period, especially the SVS₂ CR-stimulation can result in a large range of C_av -values for initially large C_av-values (C_av > 0.1; Figure 11B). However, for initially small C_av -values it induces a smaller C_av than the RVS₂ CR-stimulation (Figure 11B) and also has a superior aftereffect on the average synaptic weight during the CR-off period (Figure 11C). Except for large C_av-values at the onset of the SVS₂ CR-stimulation, the SVS₂ CR-stimulation with K₂ = 0.15 induces a very robust and long-lasting reduction of the C_av-values (Figure 11D). The equivalent for the RVS₂ CR stimulation is less effective and less robust (Figure 11D).

FIGURE 11

Figure 11. The effect of the RVS and SVS CR-stimulations for weak stimulation intensities, K, on the average synaptic connectivity. (A) C_av at the end of the first CR-on period vs. C_av at the beginning of the first CR-on period for RVS₁ (red asterisks) and SVS₁ (blue circles) CR-stimulation signals with K₁ = 0.10. (B) C_av at the end of the second CR-on period vs. C_av at the beginning of the second CR-on period under the condition that RVS₂ (red asterisks) or SVS₂ (blue circles) CR-stimulation signals are delivered during the second CR-on period with stimulation intensity K₂ = 0.15. (C) C_av at the end of the CR-off period vs. C_av at the end of the second CR-on period under the condition that RVS₂ (red asterisks) or SVS₂ (blue circles) CR-stimulation signals are delivered during the second CR-on period with K₂ = 0.15. (D) C_av at the end of the CR-off period vs. C_av at the beginning of the second CR-on period under the condition that RVS₂ (red asterisks) or SVS₂ (blue circles) CR-stimulation signals with K₂ = 0.15 were delivered during the second CR-on period.

In summary, for the two-stage CR-stimulation with weak onset (K₁ = 0.10; K₂ = 0.15) it is best to first apply the RVS CR and then the SVS CR-stimulation because the RVS₁ CR-stimulation with K₁ = 0.10 reduces the amount of synchronization of the initially very strongly synchronized neuronal activity (Figure 5C) as well as the very strong C_av (Figure 11A), and the SVS₂ CR stimulation is robustly and effectively performing for C_av ≤ 0.1 (Figure 11D) as well as entraining the mean phases of the neuronal subpopulations to the stimulation onsets (Figures 10F,I).

This two-stage CR-stimulation (RVS₁ with K₁ = 0.10 and SVS₂ with K₂ = 0.15) performs better on C_av than the short as well as the long single-stage RVS CR-stimulation with K = 0.15 and also better than the short single-stage SVS CR-stimulation with K = 0.15 (Figure 12A). The differences are statistically significant (one-sided Mann-Whitney-U-Test, p < 0.05). Although this optimal two-stage CR-stimulation is statistically not different from the long single-stage SVS CR-stimulation with K = 0.15, it has the advantage that it uses a smaller stimulation intensity during the first CR-on period. The average order parameter, R_av, also suggests that the two-stage CR-stimulation (RVS₁ with K₁ = 0.10 and SVS₂ with K₂ = 0.15) has a stronger desynchronizing effect in comparison with the other CR approaches (Figure 12B). Increasing the stimulation intensities for the two-stage CR-stimulation (RVS₁ with K₁ = 0.15 and SVS₂ with K₂ = 0.20) shows similar effects Figures 12C,D).

FIGURE 12

Figure 12. Comparison of anti-kindling effects of four single-stage CR-stimulations with one two-stage CR-stimulation with weaker onset stimulation intensity (K₁ < K₂). (A) Boxplots of C_av(t = end of CR-off period) for different CR-stimulation approaches with K = 0.15 for the single-stage CR-stimulations. The CR-off period lasts as long as the total CR-on period. From the left to the right: 64 s lasting RVS CR, 64 s lasting SVS CR, two-stage CR stimulation with RVS₁ CR with K₁ = 0.10 and SVS₂ CR with K₂ = 0.15 and finally a 128 s lasting CR-off period, 128 s lasting RVS CR, 128 s lasting SVS CR. (B) as (A) for R_av(t = end of CR-off period). (C) as (A) for K = 0.20 and the middle boxplot represent the results of the two-stage CR-stimulation with K₁ = 0.15 and K₂ = 0.20.(D) as (B) for K = 0.20 and the middle boxplot represent the results of the two-stage CR-stimulation with K₁ = 0.15 and K₂ = 0.20. For each condition (K value, CR approach and duration) eleven samples are used. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

The advantage of the two-stage approach becomes even more pronounced for particularly short stimulation durations. Decreasing the duration of each CR-on period to 32 s or 48 s results in more reduced C_av as well as R_av values at the end of the CR-off period for the two-stage CR-stimulation (RVS₁ with K₁ = 0.15 and SVS₂ with K₂ = 0.20) than for the single-stage RVS CR-stimulation with K₁ = 0.10 and K₂ = 0.15 with a total CR-on duration of 64 s or 96 s, respectively (Figure 13). This is because for K = 0.15 the SVS CR stimulation decreases C_av much faster than the RVS CR stimulation does. The mentioned reductions caused by decreasing the CR-on period are statistically significant (one-sided Mann-Whitney-U-Test, p < 0.05).

FIGURE 13

Figure 13. Comparison of anti-kindling effects induced by the two-stage CR and the single-stage RVS CR stimulation, both with weak onset intensity (K₁ = 0.10, K₂ = 0.15), as a function of the duration of the total CR-on period, T_{CR − on}. (A) Boxplots of C_av(t = end of CR-off period) for different lasting CR-on periods. (B) Boxplots of R_av(t = end of CR-off period) different lasting CR-on periods. Results obtained by the two-stage (RVS₁, SVS₂) CR stimulation are shown in black, for the single-stage CR-stimulation in red. K₁ = 0.10 and K₂ = 0.15. The CR-off period lasts as long as the total CR-on period T_{CR − on}. For each condition (CR approach and duration) eleven samples are used. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

Two-Stage CR Stimulation with Constant Stimulation Intensity

So far the two-stage CR-stimulation with weak onset seems to be a favorable choice compared to the single-stage CR-stimulation in the range of K = K₂ ∈ [0.15, 0.20]. In this section we will investigate the effect of the two-stage CR-stimulation with identical stimulation intensities. This is to study whether the two-stage CR approach might be applicable beyond the initial concept of a two-stage CR stimulation with weak onset.

Figure 14 shows that except for K = 0.10 the two-stage CR-stimulation with RVS₁ and SVS₂ CR causes a stronger reduction of C_av-values than the two-stage CR with SVS₁ and RVS₂ CR, where the differences are statistically significant (one-sided Mann-Whitney-U-test, p < 0.05). However, for the very weak stimulation intensity of K = 0.10 the results are not very pronounced, and statistically there is no difference between the results obtained by the two different two-stage CR-stimulations. In what follows we will only compare the best two-stage CR stimulation variant (RVS₁and SVS₂) with the single-stage CR stimulations.

FIGURE 14

Figure 14. Comparison of anti-kindling effects induced by the two different two-stage CR-stimulations as a function of constant stimulation intensity K₁ = K₂ = K. (A) shows the boxplots of the C_av results in black for the two-stage CR-stimulation with RVS₁ and SVS₂, and in cyan for the two-stage CR-stimulation with SVS₁ and RVS₂, as a function of the constant stimulation intensity, K.(B) As (A) for the R_av results. Eleven samples are used for each condition (K value and order of the CR approaches). The distribution of the obtained C_av and R_av values are represented by a boxplot. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

Figure 4 showed that for K = 0.15 the average synaptic weight, C_av, at the end of the CR-off period is not very robust against initial network conditions and sequence orders. Extending the total CR-on period from 64 to 128 s increases the robustness (compare the two boxplots at the left to the two boxplots at the right side of Figure 15A). The C_av-value resulting from the most effective two-stage CR-stimulation for K = 0.15 (RVS₁and SVS₂) is smaller and more robust against initial network conditions and sequence orders than the short and long single-stage RVS and the short single-stage SVS CR stimulation (Figure 15A). These differences are statistically significant (one-sided Mann-Whitney-U-test, p < 0.05). The median order parameter, averaged over the last 1.6 s of the CR-off period, R_av, is smaller for the two-stage CR stimulation than for the single stage CR stimulations, however, the difference is only statistically significant for the single-stage RVS CR stimulations, independent of the duration of the single-stage CR stimulation (Figure 15B). For stronger stimulation intensities (K = 0.20 and K = 0.25) the advantage of using the two-stage CR-stimulation (RVS₁and SVS₂) is less pronounced, but still statistically not worse than the single-stage SVS CR-stimulation (for K = 0.20 see Figures 15C,D; for K = 0.25 results are similar and therefore not shown), although the short (64 s) single-stage SVS CR-stimulation should be preferable because of its shorter stimulation duration.

FIGURE 15

Figure 15. Comparison of anti-kindling effects of four single-stage CR-stimulation with one two-stage CR-stimulation with constant stimulation intensity (K₁ = K₂). (A) Boxplots of C_av(t = end of CR-off period) for different CR-stimulation approaches with K = 0.15 for the single-stage CR-stimulations. The CR-off period lasts as long as the total CR-on period. From the left to the right: 64 s lasting RVS CR, 64 s lasting SVS CR, two-stage CR-stimulation with RVS₁ and SVS₂ CR-on period and finally a 128 s lasting CR-off period, 128 s lasting RVS CR, 128 s lasting SVS CR. (B) as (A) for R_av(t = end of CR-off period). (C) as (A) for K = 0.20. (D) as (B) for K = 0.20. Eleven samples are used for each boxplot. The gray horizontal line within the box represents the median, the box itself the middle 50% and the whiskers below (above) the box the first (last, respectively) 25%. Outliers are defined as 1.5 times the length of the box below or above the box and represented by open circles. An asterisk indicates a statistically significant lower C_av- or R_av-value (one-sided Mann-Whitney test, p < 0.05).

Discussion

In this study we investigated if a two-stage CR-stimulation induces better anti-kindling effects than a single-stage CR-stimulation for weak to intermediate stimulation intensities (K ≤ 0.25). The most optimal CR-approach would be the one that induces the best long-lasting anti-kindling results for the shortest and weakest CR-stimulation. The main result of this study is that the two-stage CR-stimulation with first a RVS CR-stimulation with intensity K₁ = 0.10, followed by a SVS CR-stimulation with intensity K₂ = 0.15 is the optimal CR-approach for weak stimulation intensities (Figure 12), even for short CR-stimulations (Figure 13). This is because it uses the advantage of a very weak RVS CR-stimulation to reduce the average synaptic strength, C_av, and the advantage of the weak SVS CR-stimulation of inducing stronger and more robust anti-kindling effects. For weaker K₂ values the SVS CR-stimulation cannot induce stronger anti-kindling effects (Figure 14) than the short single-stage RVS CR-stimulation (Figure 4). Therefore, the short single-stage RVS CR-stimulation is preferred if no stronger stimulation intensity can be used than K ≥ 0.10. For stronger stimulation intensities (K = K₁ = 0.20) the two-stage CR-stimulation induces statistically similar anti-kindling effects as the short single-stage SVS CR-stimulation (Figures 15C,D), but because of the shorter duration the single-stage SVS CR-stimulation is preferred in case no weaker stimulation intensities can be used than K = 0.20.

From a previous study (Zeitler and Tass, 2015) with single-stage CR stimulation it is known that for very weak stimulation intensity (K = 0.10) only the RVS CR-stimulation could induce small anti-kindling effects. It showed further that for intermediate stimulation intensities (0.20 ≤ K ≤ 0.50) the SVS CR-stimulation is the most optimal CR-approach. In that study we applied after three stimulation ON-cycles two stimulation OFF-cycles during the CR-on period. In the present study we used no stimulation OFF-cycles during the CR on-period(s), because preliminary results showed that without stimulation OFF-cycles the anti-kindling effects are increased, mainly for the very weak and weak RVS CR-stimulation, but also for the weak SVS CR stimulation (results not shown). The reason for this might be that the CR stimulation, without stimulation OFF-cycles as used in the current study, contains about 67% more ON-cycles. To investigate this we performed additional simulations in which we applied two stimulation OFF-cycles after each group of three simulation-ON cycles until in total 6000 stimulation ON-cycles were applied. The resulting CR-on period lasts 120 s and corresponds to a CR-on period of 96 s in case no stimulation OFF-cycles were applied. These additional simulations showed that if the two-stage CR-stimulation with weak onset contains the same amount of ON-cycles, similar anti-kindling effects are induced by the stimulation signals with ON:OFF ratios of 3:2 and 5:0, but for a ratio of 3:2 slightly stronger stimulation intensities (K₁ = 0.13 and K₂ = 0.16) are needed than for a ratio of 5:0 (K₁ = 0.10 and K₂ = 0.15; results not shown).

Lysyansky and coworkers showed in their computational model with phase oscillators that a permanent CR stimulation with high frequency bursts as stimulation signal induces desynchronization for weak stimulation intensities (Lysyansky et al., 2011). This is in good agreement with our finding that for a permanent CR stimulation with weak intensity a good desynchronization can be obtained.

As shown in numerous studies in neuroscience, psychology and education, learning effects can be reinforced by means of the spacing principle, i.e., by administering repeated stimuli spaced by pauses as opposed to administering a massed stimulus in a single long stimulation session (Ebbinghaus, 1913; Itoh et al., 1995; Frey and Morris, 1997; Menzel et al., 2001; Scharf et al., 2002; Cepeda et al., 2006, 2009; Pavlik and Anderson, 2008; Xue et al., 2011; Naqib et al., 2012). By the same token, as demonstrated computationally the spacing principle may also boost the CR-induced unlearning of abnormal synaptic connectivity and, hence, of abnormal neuronal synchrony (Popovych et al., 2015). For this demonstration, CR stimulation was applied at particularly weak intensities, so that neither acute effects (i.e., desynchronization under stimulation) nor long-lasting effects (i.e., anti-kindling) could be elicited with permanently delivered CR. Intriguingly, spaced CR stimulation at these particularly weak intensities can effectively induce anti-kindling. The dynamical mechanism behind this stimulation protocol is that the spaced stimulation causes the neuronal population to consecutively bounce from one attractor to another one, ultimately approaching desired attractors characterized by down-regulated synaptic connectivity and synchrony. The downside to this approach is the immense time it takes to build up anti-kindling effects (Popovych et al., 2015). For implanted stimulation treatments this might be tolerable. However, for non-invasive treatments, requiring compliance, e.g., requiring the patient to wear a stimulation device, it is indispensable to minimize treatment duration. In contrast, the two-stage CR stimulation is effective even at particularly short stimulation durations (Figure 13). Accordingly, one might use the spacing principle for single-stage CR stimulation at higher stimulation intensities. In addition, one could also apply the spacing principle to the two-stage CR stimulation in an attempt to further increase its efficacy.

The two-stage CR stimulation comprises two stimulation types, RVS and SVS CR stimulation, that both build up their effects on time scales, that are large compared to the period of the synchronized oscillation, and both act, in principle, in the same direction, i.e., they have desynchronizing and anti-kindling effects, at least for appropriate parameter ranges. Combining the two CR variants enables to induce anti-kindling at particularly weak intensities. In contrast, in the field of desynchronizing stimulation protocols preparatory stimuli have been used in computational studies (in the absence of STDP) that act on the short time scale of the period of the synchronized oscillation and enable an acute desynchronization when repetitively delivered as part of a double pulse stimulation (Tass, 2001a,b, 2002a,b), but typically have synchronizing effects if solely applied. The first, preparatory stimulus of a double pulse is stronger and causes a reset (and often a synchronization) of the collective dynamics within one (Tass, 2001a,b) or a few cycles of the collective oscillation, irrespective of the network's initial dynamic state (Tass, 2002a,b). The first stimulus might be a strong pulse (Tass, 2001a, 2002b) or a low-frequency pulse train (Tass, 2002a,b). The weaker, second stimulus is applied after a fixed time delay and hits the network in its vulnerable state in a stereotyped way, in this way causing a desynchronization. Although the two-stage CR stimulation and the double pulse stimulation share some conceptual similarities (empowered stimulation by combined stimulation protocols), both approaches are, in fact, quite different. In fact, it adds to the value of the two-stage CR approach that the preparatory type of stimulus, the RVS CR stimulation causes an acute desynchronization rather than a synchronization as compared to the first, strong and globally resetting stimulus of a double pulse.

Our analysis of the dynamics induced by RVS CR stimulation and SVS CR stimulation revealed clear differences (Figures 5–7). At sufficient intensities K RVS CR stimulation causes repetitive phase resets of the different neuronal subgroups (located close to the respective stimulation sites; Figure 5) along with a desynchronization between different subgroups (Figures 6, 7) as well as an overall desynchronization of the whole population (Figure 4) combined with a desynchronization on the mesoscopic level, i.e., of the different subgroups (Figure 6). In contrast, applied at sufficient intensity K, SVS CR stimulation induces a phase entrainment between subgroups and corresponding stimuli (Figure 5). Changing the sequence of the SVS CR stimulation gives rise to different phase shifted entrainment patterns and, hence, characteristic 3-peak distributions with peaks at π∕2, π and 3π∕2 rad (Figure 6) as detected by an appropriately designed index α₂₃ from Equation 29 (Figure 7). The subgroups' mutually phase shifted entrainment is combined with a reduction of the overall synchronization (Figure 4) as well as of the synchronization on the mesoscopic level of subgroups (Figure 6). Combining a preparatory RVS CR stimulation at weak intensities (first stage) with a SVS CR stimulation at higher intensity (second stage) renders the latter more effective (Figures 9, 12) and furthers the SVS characteristic phase entrainment mechanism (Figure 10).

For our analysis of the stimulus-locking of the subgroups' phases we used the resetting index E_sg(Δt) for each subgroup sg, see Equation (16). This index is nothing but a circular mean of the phases at time Δt relative to the stimulus onsets. In case of a unimodal distribution of phase (i.e., periodic) variables, a circular mean is a standard measure (Batschelet, 1981). A quantity comparable to E_sg(Δt), based on a wavelet transformation and denoted as “phase-locking factor,” has been used to study phase resetting in EEG signals registered during sensory stimulation experiments (Tallon-Baudry et al., 1996; Makeig et al., 2002). However, since in all our simulations we only found phase distributions that were either close to uniform or unimodal, we used the circular mean for quantification. In contrast, multimodal distributions of the phase of an oscillatory signal may occur in the case of stimulus-locked response clustering (Tass, 2003d). The detection of the corresponding bi- or multimodal cross-trial distributions of the phases require different statistical methods (Tass, 2004). Applying a circular mean-based analysis to the phase difference of two oscillatory signals across trials led to quantities like the phase locking value (Lachaux et al., 1999) and the stimulus locked-phase synchronization index (Tass, 2003d).

A global pattern of intermittent bursting of synchronization was observed in the context of Parkinson's (Park et al., 2011). Accordingly, our stimulation approach could be tested in computational models showing that type of dynamics, or straight away in Parkinsonian monkeys or Parkinson's patients in comparison to standard CR approaches (see Tass et al., 2012a; Adamchic et al., 2014a). However, our approach is not limited to Parkinson's (see Introduction).

The dynamics of our model is quite complex due to STDP, and it is difficult to predict what the exact long-lasting anti-kindling effect is for different individual neuronal initial conditions as well as different sequence orders. Individually randomly varying background noise might even make it more complicated. Nevertheless, the impact of noise on stimulus responses is a relevant topic. For instance, it was shown that cross-trial response clustering in two coupled oscillators subjected to CR stimulation obeys a stochastic-like behavior (Tass, 2003c). In addition, in the presence of STDP noise may have unexpected, e.g., coupling- stabilizing effects (Popovych et al., 2013). Hence, further studies should address the impact of different noise processes on the outcome of different desynchronizing stimulation approaches.

In this paper, we investigated how the CR stimulation algorithm for weak stimulation intensities could be optimized by combining the RVS and the SVS CR approaches and by using permanent CR stimulation signals during the CR-on period. The results predict for weak onset better anti-kindling by the two-stage CR stimulation than by the single-stage RVS CR stimulation. Stimulation protocols of this kind might be tested for invasive as well as non-invasive CR approaches. In the case of acoustic CR stimulation for the treatment of tinnitus (Tass et al., 2012b) stimulation at particularly weak intensities might be beneficial for several reasons: In general, cognitive functions, communication and social interaction are less hindered by low-intensity sound therapy. In tinnitus with pronounced hearing loss two-stage CR might enable to reduce the stimulation intensity and, hence, reduce the risk of sound-induced hearing impairment. About 30–40% of the tinnitus patients additionally suffer from hyperacusis, an over-sensitivity to certain volume and frequency ranges of sound (Sheldrake et al., 2015). In such patients low-intensity sound therapy may prevent from unpleasant or painful perceptions. For invasive CR approaches the two-stage protocol might be advantageous, too. Standard high-frequency permanent deep brain stimulation for the treatment of Parkinson's disease (Benabid et al., 1991; Deuschl et al., 2006b) may cause side effects such as speech and gait deterioration (Mahlknecht et al., 2015). Some of these side effects might be caused by particular anatomical relations between targets and neighboring areas and/or fiber tracts passing by, others by sub-optimally positioned leads. In any case, the reduction of the amplitude of the stimulation current may help to reduce side effects. This computational study is intended to provide results that can be tested in pre-clinical (animal) studies and phase 1 (first in man) and subsequently phase 2 (proof of concept and dose finding) studies. Our model is qualitatively rather than biophysically realistic. Hence, our experimentally testable predictions are qualitatively and refer to, e.g., the strategy of two-stage stimulation and the related ratio between the stimulation intensity of the first vs. second stage—as opposed to absolute intensity values. A comparable qualitatively prediction of the ratio between the optimal intensity of standard high-frequency deep brain stimulation and an optimal intensity for CR deep brain stimulation was verified in a pre-clinical study in parkinsonian non-human primates and enabled to reveal pronounced long-lasting after-effects of CR stimulation (Tass et al., 2012b). By the same token, a dose finding study might reveal intensity ranges for RVS and SVS CR stimulation, separately. One could then select a RVS intensity close to the lower border of the RVS range and a SVS intensity that exceeds the RVS intensity, say, by 30% or is selected from the lower middle range of SVS intensities. These intensities could then be tested in a first man study of two-stage CR stimulation. Of course, based on the computational results obtained here, we cannot make any claims concerning safety, tolerability and efficacy in potential clinical applications.

Author Contributions

PT came up with the initial ideas for this work. MZ designed and performed the simulations. MZ and PT analyzed the data, discussed findings and interpretations, and wrote the manuscript.

Conflict of Interest Statement

PT works at Jülich Research Center and as consulting professor at Stanford University. Several patents protect invasive and non-invasive CR neuromodulation. The main inventor of the CR patent portfolio is PT the assignee is Jülich Research Center. MZ works at Jülich Research Center. PT and MZ are co-inventors of recently filed two-stage CR patents.

References

Adamchic, I., Hauptmann, C., Barnikol, U. B., Pawelcyk, N., Popovych, O. V., Barnikol, T., et al. (2014a). Coordinated Reset has lasting aftereffects in patients with Parkinson's Disease. Mov. Disord. 29, 1679–1684. doi: 10.1002/mds.25923

CrossRef Full Text

Adamchic, I., Langguth, B., Hauptmann, C., and Tass, P. A. (2012a). Psychometric evaluation of Visual Analog Scale for the assessment of chronic tinnitus. Am. J. Audiol. 2, 215–225. doi: 10.1044/1059-0889(2012/12-0010)

CrossRef Full Text | Google Scholar

Adamchic, I., Tass, P. A., Langguth, B., Hauptmann, C., Koller, M., Schecklmann, M., et al. (2012b). Linking the tinnitus questionnaire and the clinical global impression: which differences are clinically important? Health Qual. Life Outcomes 10:79. doi: 10.1186/1477-7525-10-79

PubMed Abstract | CrossRef Full Text | Google Scholar

Adamchic, I., Toth, T., Hauptmann, C., and Tass, P. A. (2014b). Reversing pathologically increased EEG power by acoustic CR neuromodulation. Hum. Brain Mapp. 35, 2099–2118. doi: 10.1002/hbm.22314

CrossRef Full Text | Google Scholar

Aoki, T., and Aoyagi, T. (2011). Self-organized network of phase oscillators coupled by activity-dependent interactions. Phys. Rev. E 84:066109. doi: 10.1103/PhysRevE.84.066109

PubMed Abstract | CrossRef Full Text | Google Scholar

Batschelet, E. (1981). Circular Statistics in Biology. London: Academic Press.

Google Scholar

Bayati, M., and Valizadeh, A. (2012). Effect of synaptic plasticity on the structure and dynamics of disordered networks of coupled neurons. Phys. Rev. E 86:011925. doi: 10.1103/PhysRevE.86.011925

PubMed Abstract | CrossRef Full Text | Google Scholar

Belykh, I., de Lange, E., and Hasler, M. (2005). Synchronization of bursting neurons: what matters in the Network Topology. Phys. Rev. Lett. 94:188101. doi: 10.1103/PhysRevLett.94.188101

PubMed Abstract | CrossRef Full Text | Google Scholar

Belykh, I., and Hasler, M. (2011). Mesoscale and clusters of synchrony in networks of bursting neurons. Chaos 21:016106. doi: 10.1063/1.3563581

PubMed Abstract | CrossRef Full Text | Google Scholar

Benabid, A. L., Pollak, P., Gervason, C., Hoffmann, D., Gao, D. M., Hommel, M., et al. (1991). Long-term suppression of tremor by chronic stimulation of the ventral intermediate thalamic nucleus. Lancet 337, 403–406. doi: 10.1016/0140-6736(91)91175-T

PubMed Abstract | CrossRef Full Text | Google Scholar

Bi, G. Q., and Poo, M. M. (1998). Synaptic modifications in cultured hippocampal neurons: dependence on spike timing, synaptic strength, and postsynaptic cell type. J. Neurosci. 18, 10464–10472.

PubMed Abstract | Google Scholar

Bliss, T. V., and Lomo, T. (1973). Long-lasting potentiation of synaptic transmission in the dentate area of the anaesthetized rabbit following stimulation of the performant path. J. Physiol. 232, 331–356. doi: 10.1113/jphysiol.1973.sp010273

PubMed Abstract | CrossRef Full Text | Google Scholar

Cepeda, N. J., Coburn, N., Rohrer, D., Wixted, J. T., Mozer, M. C., and Pashler, H. (2009). Optimizing distributed practice: theoretical analysis and practical implications. Exp. Psychol. 56, 236–246. doi: 10.1027/1618-3169.56.4.236

PubMed Abstract | CrossRef Full Text | Google Scholar

Cepeda, N. J., Pashler, H., Vul, E., Wixted, J. T., and Rohrer, D. (2006). Distributed practice in verbal recall tasks: a review and quantitative synthesis. Psychol. Bull. 132, 354–380. doi: 10.1037/0033-2909.132.3.354

PubMed Abstract | CrossRef Full Text | Google Scholar

De la Rocha, J., Marchetti, C., Schiff, M., and Reyes, A. D. (2008). Linking the response properties of cells in auditory cortex with network architecture: cotuning versus lateral inhibition. J. Neurosci. 28, 9151–9163. doi: 10.1523/JNEUROSCI.1789-08.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

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

PubMed Abstract | CrossRef Full Text | Google Scholar

Deuschl, G., Schade-Brittinger, C., Krack, P., Volkmann, J., Schäfer, H., Bötzel, K., et al. (2006b). A randomized trial of deep-brain stimulation for Parkinson's disease. N. Engl. J. Med. 355, 896–908. doi: 10.1056/NEJMoa060281

PubMed Abstract | CrossRef Full Text | Google Scholar

Dominguez, M., Becker, S., Bruce, I., and Read, H. (2006). A spiking neuron model of cortical correlates of sensorineural hearing loss: spontaneous firing, synchrony, and tinnitus. N. Comp. 18, 2942–2958. doi: 10.1162/neco.2006.18.12.2942

PubMed Abstract | CrossRef Full Text | Google Scholar

Dujardin, K., Defebvre, L., Krystkowiak, P., Blond, S., and Destée, A. (2001). Influence of chronic bilateral stimulation of the subthalamic nucleas on cognitive funcion in Parkinson's disease. J. Neurol. 248, 603–611. doi: 10.1007/s004150170139

PubMed Abstract | CrossRef Full Text | Google Scholar

Ebbinghaus, H. (1913). Memory: A Contribution to Experimental Psychology. Columbia University. Teachers College. Educational reprints. no. 3. New York, NY : Teachers College, Columbia University.

Ebert, M., Hauptmann, C., and Tass, P. A. (2014). Coordinated reset stimulation in a large-scale model of the STN-GPe circuit. Front. Comput. Neurosci. 8:154. doi: 10.3389/fncom.2014.00154

PubMed Abstract | CrossRef Full Text | Google Scholar

Eggermont, J. J., and Tass, P. A. (2015). Maladaptive neural synchrony in tinnitus: origin and restoration. Front. Neurol. 6:29. doi: 10.3389/fneur.2015.00029

PubMed Abstract | CrossRef Full Text | Google Scholar

Feldman, D. E. (2000). Timing-based LTP and LTD at vertical inputs to layer II/III pyramidal cells in rat barrel cortex. Neuron 27, 45–56. doi: 10.1016/S0896-6273(00)00008-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Frey, U., and Morris, R. G. M. (1997). Synaptic tagging and long-term potentiation. Nature 385, 533–536. doi: 10.1038/385533a0

PubMed Abstract | CrossRef Full Text | Google Scholar

Gerstner, W., Kempter, R., Van Hemmen, J. L., and Wagner, H. (1996). A neuronal learning rule for sub-millisecond temporal coding. Nature 383, 76–78. doi: 10.1038/383076a0

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, D., Wang, Q., and Perc, M. (2012). Complex synchronous behavior in interneuronal networks with delayed inhibitory and fast electrical synapses. Phys. Rev. E 85:061905. doi: 10.1103/PhysRevE.85.061905

PubMed Abstract | CrossRef Full Text | Google Scholar

Haken, H. (1983). Advanced Synergetices. Berlin: Springer.

Hammond, C., Bergman, H., and Brown, P. (2007). Pathological synchronization in parkinson's disease: networks, models and treatments. Trends Neurosci. 30, 357–364. doi: 10.1016/j.tins.2007.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansel, D., Mato, G., and Meunier, C. (1993). Phase dynamics of weakly coupled Hodgkin-Huxley neurons. Europhys. Lett. 23, 367–372. doi: 10.1209/0295-5075/23/5/011

CrossRef Full Text | Google Scholar

Hauptmann, C., and Tass, P. A. (2007). Therapeutic rewiring by means of desynchronizing brain stimulation. Biosystems 89, 173–181. doi: 10.1016/j.biosystems.2006.04.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Hauptmann, C., and Tass, P. A. (2009). Cumulative and after-effects of short and weak coordinated reset stimulation - a modeling study. J. Neural Eng. 6:016004. doi: 10.1088/1741-2560/6/1/016004

PubMed Abstract | CrossRef Full Text | Google Scholar

Hauptmann, C., Wegener, A., Poppe, H., Williams, M., Popelka, G., and Tass, P. A. (2016). Validation of a mobile device for acoustic coordinated reset neuromodulation tinnitus therapy. J. Am. Acad. Audiol. doi: 10.3766/jaaa.15082

CrossRef Full Text | Google Scholar

Hebb, D. O. (1949). The Organization of Behavior: A Neurophysiological Approach. New York, NY: John Wiley and Sons.

Google Scholar

Hershey, T., Wu, J., Weaver, P. M., Perantie, D. C., Karimi, M., Tabbal, S. D., et al. (2008). Unilateral vs. bilateral STN DBS effects on working memory and motor function in Parkinson disease. Exp. Neurol. 210, 402–408. doi: 10.1016/j.expneurol.2007.11.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Hodgkin, A., and Huxley, A. F. (1952). A quantitative description of membrane current and application to conduction and excitation. J. Physiol. 117, 500–544. doi: 10.1113/jphysiol.1952.sp004764

PubMed Abstract | CrossRef Full Text | Google Scholar

Homer, M. A., Rubin, S. S., Horowitz, T. D., and Richter, E. (2012). Linguistic testing during ON/OFF states of electrical stimulation in the associative portion of the subthalamic nucleus. Neuromodulation 15, 238–245. doi: 10.1111/j.1525-1403.2012.00441.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Hübener, M., and Bonhoeffer, T. (2014). Neuronal plasticity: beyond the critical period. Cell 159, 727–737. doi: 10.1016/j.cell.2014.10.035

PubMed Abstract | CrossRef Full Text | Google Scholar

Itoh, K., Stevens, B., Schachner, M., and Fields, R. D. (1995). Regulated expression of the neural cell adhesion molecule L1 by specific patterns of neural impulses. Science 270, 1369–1372. doi: 10.1126/science.270.5240.1369

PubMed Abstract | CrossRef Full Text | Google Scholar

Jahanshahi, M., Obeso, I., Baunez, C., Alegre, M., and Krack, P. (2015). Parkinson's disease, the subthalamic nucleus, inhibition, and impulsivity. Mov. Disord. 30, 128–140. doi: 10.1002/mds.26049

PubMed Abstract | CrossRef Full Text | Google Scholar

Knoblauch, A., Hauser, F., Gewaltig, M.-O., Körner, E., and Palm, G. (2012). Does spike-timing-dependent synaptic plasticity couple or decouple neurons firing in synchrony? Front. Comput. Neurosci. 6:55. doi: 10.3389/fncom.2012.00055

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuramoto, Y. (1984). Chemical Oscillations, Waves and Turbulence. Berlin; Heidelberg; New York, NY:Spinger.

Google Scholar

Lachaux, J.-P., Rodriguez, E., Martinerie, J., and Varela, F. J. (1999). Measuring phase synchrony in brain signals. Hum. Brain Mapp. 8, 194–208.

PubMed Abstract | Google Scholar

Lenz, F. A., Kwan, H. C., Martin, R. L., Tasker, R. R., Dostrovsky, J. O., and Lenz, Y. E. (1994). Single unit analysis of the human ventral thalamic nuclear group. Brain 117, 531–543. doi: 10.1093/brain/117.3.531

PubMed Abstract | CrossRef Full Text | Google Scholar

Llinas, R. R., Ribary, U., Jeanmonod, D., Kronberg, E., and Mitra, P. P. (1999). Thalamocortical dysrhythmia: a neurological and neuropsychiatric syndrome characterized by magnetoencephalography. Proc. Natl. Acad. Sci. U.S.A. 96, 15222–15227. doi: 10.1073/pnas.96.26.15222

PubMed Abstract | CrossRef Full Text | Google Scholar

Lysyansky, B., Popovych, O. V., and Tass, P. A. (2011). Desynchronizing anti-resonance effect of m: n ON-OFF coordinated reset stimulation. J. Neural Eng. 8:036019. doi: 10.1088/1741-2560/8/3/036019

PubMed Abstract | CrossRef Full Text | Google Scholar

Mahlknecht, P., Limousin, P., and Foltynie, T. (2015). Deep brain stimulation for movement disorders: update on recent discoveries and outlook on future developments. J. Neurol. 262, 2583–2595. doi: 10.1007/s00415-015-7790-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Maistrenko, Y. L., Lysyansky, B., Hauptmann, C., Burylko, O., and Tass, P. A. (2007). Multistability in the Kuramoto model with synaptic plasticity. Phys. Rev. E 75:066207. doi: 10.1103/PhysRevE.75.066207

PubMed Abstract | CrossRef Full Text | Google Scholar

Makeig, S., Westerfield, M., Jung, T.-P., Enghoff, S., Townsend, J., Courchesne, E., et al. (2002). Dynamic brain sources of visual evoked responses. Science 29, 690–694. doi: 10.1126/science.1066168

CrossRef Full Text | Google Scholar

Markram, H., Lübke, J., Frotscher, M., and Sakmann, B. (1997). Regulation of synaptic efficacy by coincidence of postsynaptic Aps and EPSPs. Science 275, 213–215. doi: 10.1126/science.275.5297.213

PubMed Abstract | CrossRef Full Text | Google Scholar

Martens, H. C., Toader, E., Decré, M. M. J., Anderson, D. J., Vetter, R., Kipke, D. R., et al. (2011). Spatial steering of deep brain stimulation volumes using a novel lead design. Clin. Neurophysiol. 122, 558–566. doi: 10.1016/j.clinph.2010.07.026

PubMed Abstract | CrossRef Full Text | Google Scholar

Masuda, N., and Kori, H. (2007). Formation of feedforward networks and frequency synchrony by spike-timing-dependent plasticity. J. Comput. Neurosci. 22, 327–345. doi: 10.1007/s10827-007-0022-1

PubMed Abstract | CrossRef Full Text | Google Scholar

McCracken, C. B., and Kiss, Z. H. T. (2014). Time and frequency-dependent modulation of local field potential synchronization by deep brain stimulation. PLoS One 9:e102576. doi: 10.1371/journal.pone.0102576

PubMed Abstract | CrossRef Full Text | Google Scholar

McIntyre, C. C., Richardson, S. J., Frankemolle, A. M., Varga, G., Noecker, A. M., and Alberts, J. (2011). Improving postural stability via computational modeling appraoch to deep brain stimulation programming. Conf. Proc. IEEE Eng. Med. Biol. Soc. 2011, 675–676. doi: 10.1109/IEMBS.2011.6090151

PubMed Abstract | CrossRef Full Text | Google Scholar

Menzel, R., Manz, G., Menzel, R., and Greggers, U. (2001). Massed and spaced learning in honeybees: the role of cs, us, the intertrial interval, and the test interval. Learn. Mem. 8, 198–208. doi: 10.1101/lm.40001

PubMed Abstract | CrossRef Full Text | Google Scholar

Mikos, A., Bowers, D., Noecker, A. M., McIntyre, C. C., Won, M., Chaturvedi, A., et al. (2011). Patient-specific analysis of the relationship between the volume of tissue activated during DBS and verbal fluency. Neuroimage 54, S238–S246. doi: 10.1016/j.neuroimage.2010.03.068

PubMed Abstract | CrossRef Full Text | Google Scholar

Moreau, C., Defebvre, L., Destée, A., Bleuse, S., Clement, F., Blatt, J. L., et al. (2008). STN-DBS frequency effects on freezing of gait in advanced Parkinson disease. Neurology 71, 80–84. doi: 10.1212/01.wnl.0000303972.16279.46

PubMed Abstract | CrossRef Full Text | Google Scholar

Naqib, F., Sossin, W. S., and Farah, C. A. (2012). Molecular determinants of the spacing effect. Neural Plast. 2012:581291. doi: 10.1155/2012/581291

PubMed Abstract | CrossRef Full Text | Google Scholar

Nini, A., Feingold, A., Slovin, H., and Bergman, I. (1995). Neurons in the globus pallidus do not show correlated activity in the normal monkey, but phase-locked oscillations appear in the MPTP model of Parkinsonsism. J. Neurophysiol. 74, 1800–1805.

PubMed Abstract | Google Scholar

Ochi, K., and Eggermont, J. J. (1997). Effects of quinine on neural activity in cat primary auditory cortex. Hear. Res. 105, 105–118. doi: 10.1016/S0378-5955(96)00201-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, C., Worth, R. M., and Rubchinsky, L. L. (2011). Neural dynamics in Parkinsonian brain: the boundary between synchronized and nonsynchronized dynamics. Phys. Rev. E 83:042901. doi: 10.1103/PhysRevE.83.042901

PubMed Abstract | CrossRef Full Text | Google Scholar

Pavlik, P. I., and Anderson, J. R. (2008). Using a model to compute the optimal schedule of practice. J. Exp. Psychol. Appl. 14, 101–117. doi: 10.1037/1076-898X.14.2.101

PubMed Abstract | CrossRef Full Text | Google Scholar

Pinsker, M., Amtage, F., Berger, M., Nikkhah, G., and van Els, L. T. (2013). Psychiatric side-effects of bilateral deep brain stimulation for movement disorders. Acta Neurochirgica. Suppl. 117, 47–51. doi: 10.1007/978-3-7091-1482-7_8

PubMed Abstract | CrossRef Full Text | Google Scholar

Popovych, O., and Tass, P. A. (2011). Macroscopic entrainment of periodically forced oscillatory ensembles. Prog. Biophys. Mol. Biol. 105, 98–108. doi: 10.1016/j.pbiomolbio.2010.09.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Popovych, O., and Tass, P. A. (2012). Desynchronizing electrical and sensory coordinated reset modulation. Front. Hum. Neurosci. 6:58. doi: 10.3389/fnhum.2012.00058

CrossRef Full Text | Google Scholar

Popovych, O. V., Xenakis, M. N., and Tass, P. A. (2015). The spacing principle for unlearning abnormal neuronal synchrony. PLoS One 10:e0117205. doi: 10.1371/journal.pone.0117205

PubMed Abstract | CrossRef Full Text | Google Scholar

Popovych, O. V., Yanchuk, S., and Tass, P. A. (2013). Self-organized noise resistance of oscillatory neural networks with spike timing-dependent plasticity. Sci. Rep. 3:2926. doi: 10.1038/srep02926

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, Q., and Zhao, J. (2007). Adaptive coupling and enhanced synchronization in coupled phase oscillators. Phys. Rev. E 76:016207. doi: 10.1103/PhysRevE.76.016207

PubMed Abstract | CrossRef Full Text

Rodriguez-Oroz, M. C., Obeso, J. A., Lang, A. E., Houeto, J.-L., Pollak, P., Rehncrona, S., et al. (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

Rosenblum, M. G., Pikovsky, A. S., Schäfer, C., Kurths, J., and Tass, P. A. (2001). “Phase synchronization: from theory to data analysis,” in Handbook of Biological Physics, ed F. Moss (Amsterdam: Elsevier), 279–321.

Google Scholar

Saint-Cyr, J. A., Trépanier, L. L., Kumar, R., Lozano, A. M., and Lang, A. E. (2000). Neurophysiological consequences of chronic bilateral stimulation of the subthalamic nucleus in Parkinson's disease. Brain 123, 2091–2108. doi: 10.1093/brain/123.10.2091

CrossRef Full Text | Google Scholar

Scharf, M. T., Woo, N. H., Lattal, K. M., Young, J. Z., Nguyen, P. V., and Abel, T. (2002). Protein synthesis is required for the enhancement of long-term potentiation and long-term memory by spaced training. J. Neurophysiol. 87, 2770–2777. doi: 10.1152/jn.00792.2001

PubMed Abstract | CrossRef Full Text | Google Scholar

Schomer, D. L., and Lopes da Silva, F. H. (2011). Niedermeyer's Electroencephalography, 6th rev. Edn. Philadelphia, PA: Lippincott Williams & Wilkins.

Schroeder, U., Kuehler, A., Lange, K. W., Haslinger, B., Tronnier, V. M., Krause, M., et al. (2003). Subthalamic nucleus stimulation affects a frontotemporal network: a PET study. Ann. Neurol. 54, 445–450. doi: 10.1002/ana.10683

PubMed Abstract | CrossRef Full Text | Google Scholar

Seliger, P., Young, S. C., and Tsimring, L. S. (2002). Plasticity and learning in a network of coupled phase oscillators. Phys. Rev. E 65:041906. doi: 10.1103/PhysRevE.65.041906

PubMed Abstract | CrossRef Full Text | Google Scholar

Sheldrake, J., Diehl, P. U., and Schaette, R. (2015). Audiometric characteristics of hyperacusis patients. Front. Neurol. 6:105. doi: 10.3389/fneur.2015.00105

PubMed Abstract | CrossRef Full Text | Google Scholar

Silchenko, A. N., Adamchic, I., Hauptmann, C., and Tass, P. A. (2013). Impact of acoustic coordinated reset neuromodulation on effective connectivity in a neural network of phantom sound. Neuroimage 77, 133–147. doi: 10.1016/j.neuroimage.2013.03.013

PubMed Abstract | CrossRef Full Text

Sporns, O. (2011). Networks of the Brain. Cambridge, MIT Press.

Google Scholar

Tallon-Baudry, C., Bertrand, O., Delpuech, C., and Pernier, J. (1996). Stimulus specificity of phase-locked and non-phase-locked 40 Hz visual responses in Human. J. Neurosci. 16, 4240–4249.

PubMed Abstract | Google Scholar

Tass, P. A. (1999). Phase Resetting in Medicine and Biology: Stochastic Modelling and Data Analysis. Berlin: Springer.

Google Scholar

Tass, P. A. (2001a). Effective desynchronization with a resetting pulse train followed by a single pulse. Europhys. Lett. 55, 171–177. doi: 10.1209/epl/i2001-00397-8

CrossRef Full Text | Google Scholar

Tass, P. A. (2001b). Desynchronizing double-pulse phase resetting and application to deep brain stimulation. Biol. Cybern. 85, 343–354. doi: 10.1007/s004220100268

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A. (2002a). Effective desynchronization with a stimulation technique based on soft phase-resetting. Europhys. Lett. 57, 164–170. doi: 10.1209/epl/i2002-00557-x

CrossRef Full Text | Google Scholar

Tass, P. A. (2002b). Desynchronization of brain rhythms with soft phase-resetting techniques. Biol. Cybern. 87, 102–115. doi: 10.1007/s00422-002-0322-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A. (2003a). A model of desynchronizning deep brain stimulation with a demand-controlled coordinated reset of neural subpopulations. Biol. Cybern. 89, 81–88. doi: 10.1007/s00422-003-0425-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A. (2003b). Desynchronization by means of a Coordinated Reset of neural sub-populations. Prog. Theor. Phys. Suppl. 150, 281–296. doi: 10.1143/PTPS.150.281

CrossRef Full Text | Google Scholar

Tass, P. A. (2003c). Stochastic phase resetting of stimulus-locked responses of two coupled oscillators: transient response clustering, synchronization, and desynchronization. Chaos 13, 364–376. doi: 10.1063/1.1505813

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A. (2003d). Stochastic phase resetting of two coupled phase oscillators stimulated at different times. Phys. Rev. E 67:051902. doi: 10.1103/PhysRevE.67.051902

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A. (2004). Transmission of stimulus-locked responses in two coupled phase oscillators. Phys. Rev. E 69:051909. doi: 10.1103/PhysRevE.69.051909

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A., Adamchic, I., Freund, H. J., von Stackelberg, T., and Hauptmann, C. (2012b). Counteracting tinnitus by acoustic coordinated reset neuromodulation. Restor. Neurol. Neurosci. 30, 137–159. doi: 10.3233/RNN-2012-110218

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A., and Hauptmann, C. (2006). Long-term anti-kindling effects induced by short-term, weak desynchronizing stimulation. Nonl. Phen. Compl. Syst. 9, 298–312.

Google Scholar

Tass, P. A., and Majtanik, M. (2006). Long-term anti-kindling effects of desynchronization brain stimulation: a theoretical study. Biol. Cybern. 94, 58–66. doi: 10.1007/s00422-005-0028-6

CrossRef Full Text | Google Scholar

Tass, P. A., and Popovych, O. V. (2012). Unlearning tinnitus-related cerebral synchrony with acoustic coordinated reset stimulation – theoretical concept and modelling. Biol. Cybern. 106, 27–36. doi: 10.1007/s00422-012-0479-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Tass, P. A., Qin, L., Hauptmann, C., Doveros, S., Bezard, E., Boraud, T., et al. (2012a). Coordinated reset neuromodulation has sustained effects in parkinsonian non-human primates. Ann. Neurol. 72, 816–820. doi: 10.1002/ana.23663

CrossRef Full Text | Google Scholar

Tass, P. A., Silchenko, A. N., Barnikol, U., Hauptmann, C., and Speckmann, E. J. (2009). Long-lasting desynchronization in rat hippocampal slice induced by coordinated reset stimulation. Phys. Rev. E 80:011902. doi: 10.1103/PhysRevE.80.011902

PubMed Abstract | CrossRef Full Text | Google Scholar

Tommasi, G., Krack, P., Fraix, V., and Pollak, P. (2012). Effects of varying subthalamic nucleus stimulation on apraxia of lid opening in Parkinson's disease. J. Neurol. 259, 1944–1950. doi: 10.1007/s00415-012-6447-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Tukey, J. W. (1977). Exploratory Data Analysis. Reading; London: Addison-Wesley Publishing Company, Inc.

Google Scholar

Uhlhaas, P. J., and Singer, W. (2006). Neural synchrony in brain disorders: relevance for cognitive dysfunctions and pathophysiology. Neuron 52, 155–168. doi: 10.1016/j.neuron.2006.09.020

PubMed Abstract | CrossRef Full Text | Google Scholar

van Dijk, K. J., Verhagen, R., Chaturvedi, A., McIntyre, C., Bour, L. J., Heida, C., et al. (2015). A novel lead design enables selective deep brain stimulation of neural populations in the subthalamic region. J. Neural Eng. 12:046003. doi: 10.1088/1741-2560/12/4/046003

PubMed Abstract | CrossRef Full Text | Google Scholar

van Nuenen, B. F. L., Esselink, R. A., Munneke, M., Speelman, J. D., van Laar, T., and Bloem, B. R. (2008). Postoperative gait deterioration after bilateral subthalamic nucleus stimulation in Parkinson's disease. Mov. Disord. 23, 2404–2406. doi: 10.1002/mds.21986

PubMed Abstract | CrossRef Full Text | Google Scholar

Volman, V., Perc, M., and Bazhenov, M. (2011). Gap junctions and epileptic seizures – Two sides of the same coin? PLoS One 6:e20572. doi: 10.1371/journal.pone.0020572

PubMed Abstract | CrossRef Full Text | Google Scholar

Weisz, N., Moratti, S., Meinzer, M., Dohrmann, K., and Elbert, T. (2005). Tinnitus perception and distress is related to abnormal spontaneous brain activity as measured by magnetoencephalography. PLoS Med. 2, 546–553. doi: 10.1371/journal.pmed.0020153

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilson, H. R., and Cowan, J. D. (1973). A mathematical theory of the functional dynamics of cortical and thalamic nervous tissue. Kybernetik 13, 55–80. doi: 10.1007/BF00288786

PubMed Abstract | CrossRef Full Text | Google Scholar

Wong, R. K., Traub, R. D., and Miles, R. (1986). Cellular basis of neuronal synchrony in epilepsy. Adv. Neurol. 44, 583–592.

PubMed Abstract | Google Scholar

Xu, W., Miocinovic, S., Zhang, J., Baker, K., McIntyre, C. C., and Vitek, J. L. (2011). Dissociation of motor symptoms during deep brain stimulation of the subthalamic nucleus in the region of the internal capsule. Exp. Neurol. 228, 294–297. doi: 10.1016/j.expneurol.2010.08.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Xue, G., Mei, L., Chen, C., Lu, Z. L., Poldrack, R., and Dong, Q. (2011). Spaced learning enhances subsequent recognition memory by reducing neural repetition suppression. J. Cogn. Neurosci. 23, 1624–1633. doi: 10.1162/jocn.2010.21532

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeitler, M., and Tass, P. A. (2015). Augmented brain function by coordinated reset stimulation with slowly varying sequences. Front. Syst. Neurosci. 9:49. doi: 10.3389/fnsys.2015.00049

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, H., Wang, Q., Perc, M., and Chen, G. (2013). Synaptic plasticity induced transition of spike propagation in neuronal networks. Commun. Nonlinear Sci. Numer. Simulat. 18, 601–615. doi: 10.1016/j.cnsns.2012.08.009

CrossRef Full Text | Google Scholar