The Spacing Principle for Unlearning Abnormal Neuronal Synchrony

Desynchronizing stimulation techniques were developed to specifically counteract abnormal neuronal synchronization relevant to several neurological and psychiatric disorders. The goal of our approach is to achieve an anti-kindling, where the affected neural networks unlearn abnormal synaptic connectivity and, hence, abnormal neuronal synchrony, by means of desynchronizing stimulation, in particular, Coordinated Reset (CR) stimulation. As known from neuroscience, psychology and education, learning effects can be enhanced by means of the spacing principle, i.e. by delivering repeated stimuli spaced by pauses as opposed to delivering a massed stimulus (in a single long stimulation session). To illustrate that the spacing principle may boost the anti-kindling effect of CR neuromodulation, in this computational study we carry this approach to extremes. To this end, we deliver spaced CR neuromodulation at particularly weak intensities which render permanently delivered CR neuromodulation ineffective. Intriguingly, spaced CR neuromodulation at these particularly weak intensities effectively induces an anti-kindling. In fact, the spacing principle enables the neuronal population to successively hop from one attractor to another one, finally approaching attractors characterized by down-regulated synaptic connectivity and synchrony. Our computational results might open up novel opportunities to effectively induce sustained desynchronization at particularly weak stimulation intensities, thereby avoiding side effects, e.g., in the case of deep brain stimulation.


Introduction
Efficient and safe treatment requires proper dosage (see, e.g., [1][2][3]). Typically, an optimal dosage enables high efficacy at a minimal or at least tolerable rate of side effects. Different types of interventions require different dosage patterns. For instance, periodic regimes for non-invasive (peroral, topical etc.) drug delivery take into account washout characteristics to enable appropriate plasma levels [1][2][3]. Fractionation plays a major role in radiation therapy to both reduce side effects and increase therapeutic effects [4].
Deep brain stimulation (DBS) is the golden standard for the treatment of medically refractory movement disorders [5][6][7]. DBS protocols were empirically developed and basically advise to permanently stimulate at high frequencies (>100 Hz) at intensities close to sensory thresholds [8,9]. The mechanism of standard high-frequency (HF) DBS is still a matter of debate [10][11][12]. However, HF DBS only has acute effects, i.e. neither clinical [13] nor electrophysiological [14,15] effects persist after cessation of stimulation.
Long-lasting CR-induced desynchronization was confirmed in vitro in rat hippocampal slice rendered epileptic by magnesium withdrawal [26]. In the context of Parkinson's disease (PD), beneficial therapeutic long-lasting after-effects of electrical CR stimulation were observed in 1-methyl-4-phenyl-1,2,3,6-tetrahydropyridine (MPTP)-treated macaque monkeys, in contrast to standard HF DBS [27]. It was shown that unilateral CR stimulation delivered to the subthalamic nucleus (STN) of parkinsonian MPTP monkeys for only 2 h per day during 5 consecutive days leads to significant and sustained bilateral therapeutic after-effects for at least 30 days, while standard 130 Hz DBS had no after-effects [27]. Lasting after-effects of electrical CR stimulation of the STN were also confirmed in PD patients [28].
To develop CR neuromodulation further, we take advantage of a principle playing an important role in different research traditions including psychology, neuroscience and education. As known from pioneering studies on retrieving memories [29], experimental studies in several species addressing various aspects ranging from behavior to the underlying biochemical machinery [30][31][32][33][34], functional magnetic resonance imaging (fMRI) studies on face recognition [35], studies on recall of word pairs and other tasks [36][37][38] as well as studies in educational science [39], learning effects can be enhanced by means of the spacing principle, i.e. by delivering repeated stimuli spaced by pauses as opposed to delivering a massed stimulus (in a single long stimulation session). Another and similar principle is the pausing principle from eduction [40].
Given these findings, we hypothesize that appropriate spacing might enhance the effects of CR neuromodulation. To illustrate that the spacing principle may, in fact, boost the antikindling effect of CR neuromodulation, in this computational study we carry this approach to extremes. To this end, we use CR neuromodulation at particularly weak intensities which make CR ineffective if permanently delivered. We compare very weak CR neuromodulation delivered either permanently or delivered according to a spaced protocol. Intriguingly, spaced CR neuromodulation at these particularly weak intensities causes an effective anti-kindling while massed CR neuromodulation is ineffective. The spacing principle utilizes the pronounced multistability arising due to STDP. During the delivery of spaced CR neuromodulation the neuronal population successively bounces from one attractor to another one, and finally approaches an attractor characterized by down-regulated synaptic connectivity and synchrony. In this paper we first analyze the STDP-induced multistability and then continue with the specific spacing-induced effects. Finally, we discuss our findings in a broader context of DBS, noninvasive CR neuromodulation and further stimulation interventions.
The strength of the coupling of neuron j to neuron i is controlled by the time-dependent synaptic weight s ij (t) which is updated according to s ij ! s ij + δ Δs ij , i 6 ¼ j, with δ = 0.1, when the pre-synaptic neuron j or post-synaptic neuron i produces a burst. The increment Δs ij is governed by the STDP rule which depends on the time difference Δt ij of the nearest burst onsets of the pre-and post-synaptic neurons. The burst onsets are detected by the upward zero crossing by the membrane potential ν i following an inter-burst hyperpolarization. Parameters c p = 0.038, c d = 0.02, τ p = 10 ms, and τ d = 25 ms. The synaptic weights s ij are confined to the interval [s min , s max ] = [0.05, 1.5] by setting s ij to s min as soon as s ij tends to go below s min via STDP or, respectively, to s max if s ij exceeds this value. The amount of the coupling in the ensemble and the extent of synchronization among neurons are measured by the mean synaptic weight (mean coupling) S ¼ 1 [46], respectively. ψ j is the phase of neuron j and defined as c j ðtÞ ¼ 2p tÀt j;n t j;nþ1 Àt j;n ; t j,n t < t j,n+1 , where t j,n , n = 1, 2, . . ., are the times of the burst onsets [47].

CR stimulation
The neuronal ensemble with STDP (1) is stimulated by coordinated reset (CR) stimulation [17,18] via four stimulation sites. The N = 256 neurons are arranged in a 16 × 16 2Dim square lattice, and the considered four stimulation sites are placed at the center of each quadrant, as illustrated in Fig. 1A. The stimulation current administered via a stimulation site m decays as c n exp(-2d im ) with c n = 0.588 [44] and distance d im = dkm-ik from the stimulation site m to the stimulated neuron i, where d = 1/16 is the lattice distance. An example of the spatial stimulation profile of the active stimulation site 1 is depicted in color in Fig. 1A. The stimulation signals consist of short high-frequency (130 Hz) trains of charge-balanced electrical pulses, as schematically illustrated in Fig. 1B. Each positive rectangular pulse of unit amplitude and length 0.2 ms is followed by a negative counterpart of length 3 ms and amplitude-1/15 with inter-pulse interval 1/130 s.
All four stimulation sites are consequently activated in a randomized way once per stimulation period of length T s . Activation of a stimulation site means to deliver a burst, i.e., a stimulation train of six pulses. A CR stimulus is a sequence of bursts delivered at different sites with a delay of T s / 4 [ Fig. 1B]. The considered stimulation period T s = 200 ms approximates the oscillation period of the synchronized neurons. CR stimulation is administered intermittently, where a few periods of activated stimulation sites (ON-periods) are followed by a few periods of deactivated stimulation sites (OFF-periods), see

STDP-induced multistability
STDP has a great impact on the intrinsic dynamics of the neuronal ensemble, in particular, by inducing a pronounced multistability. The coexisting stable states are characterized by different coupling topologies and strengths as well as by different amounts of synchronization. We illustrate the STDP-induced multistability by simulating system (1) without stimulation (I stim i ¼ 0) for different initial conditions: The individual synaptic weights s ij are initially uniformly and randomly distributed in the interval [S 0 -σ, S 0 + σ], and the initial values of the other variables are uniformly and randomly distributed in the interval [-0.1, 0.1]. For long enough simulation we observe that the mean synaptic weight S(t) saturates at some value and stays close to it hereafter. In Fig. 1A we plot the time courses of S(t) for 12 selected initial conditions such that 6 different regimes are realized in ensemble (1). Each of these states can attract trajectories starting from different initial conditions, which supports the stability of the established coupling regimes as revealed by the mean coupling.
To explore possible coupling regimes that can be established in the ensemble for the mentioned initial conditions, we perform a large number of simulations by varying parameters S 0 and σ and collect the values of the time-averaged mean synaptic weight hSi after skipping long enough transients. The obtained mean coupling hSi is depicted in color in Fig. 2B versus parameters S 0 and σ, which illustrates a multitude of coexisting states with mean coupling ranging from small to large values. This diagram also provides an estimate of the structure of the basins of attraction of the coexisting stable states. Two main regions in other. Nevertheless, an important tendency can be observed: Larger initial mean coupling S 0 is expected to lead to a stronger coupling established in the ensemble due to STDP.
The volume of the basins of attraction from Fig. 2B can be quantified by calculating the frequency of occurrence of a given state. For this, we calculate the distribution density of hS i obtained for the considered initial conditions and plot its histogram in Fig. 2C. The distribution clearly demonstrates three well-pronounced peaks at weak, strong, and intermediate coupling.
Examples of the time dynamics of the mean coupling S(t) for such states are depicted in Fig. 2A with mean values and standard deviations, for instance, S = 0.123 ± 0.001 (filled green squares), S = 0.546 ± 0.0009 (filled violet up-triangles), S = 0.822 ± 0.0002 (filled black down-triangles), S = 0.933 ± 0.0002 (filled blue diamonds), S = 1.195 ± 0.0004 (magenta asterisks), and S = 1.499 ± 0.00003 (filled red circles). As follows from the presented values, the mean synaptic weight S(t) indeed exhibits very small time fluctuations, so that S(t) is nearly constant, and the established coupling matrices are suspected to be nearly constant as well, in particular, for large coupling.
The extent of synchronization in the neuronal ensemble depends on the amount of coupling as illustrated in Fig. 2D. The time-averaged order parameters hRi tends to grow as the mean coupling of the STDP-induced states increases. In other words, the neurons get more synchronized in the regime of larger mean coupling. Depending on the coupling strength the order parameter R(t) exhibits time fluctuations of different magnitude as illustrated by error bars in Fig. 2D. These fluctuations are very small for strongly coupled regimes (the order parameter is nearly constant), somewhat larger for weakly coupled regimes, and well pronounced for intermediately coupled regimes. Strongly coupled neurons are phase locked to each other, and the distribution of their phases is narrow and stationary. For the weakly coupled regimes the   variation of the order parameters can be accounted for by a finite-size effect leading to the fluctuation of the order 1= ffiffiffiffi N p of magnitude [48]. On the other hand, the neuronal dynamics for intermediately coupled regimes, where the order parameter is relatively large, but, nevertheless, exhibits large-amplitude oscillations, requires additional investigations. For this we consider the structure of the coupling matrices and the neuronal phase dynamics in more detail below.
In Fig. 3 we present distributions of the individual synaptic weights s ij of the coupling regimes illustrated in Fig ]. Accordingly, the coupling matrices have four indices (I, j, i 0 , j 0 ) identifying a synaptic weight from the pre-synaptic neuron (i 0 , j 0 ) to the post-synaptic neuron (i, j). It thus becomes difficult to illustrate such matrices in 2Dim plots. We simplify this problem by labeling the neurons in the lattice by a single index k according to an algorithm adopted for the sake of illustration. For example, a simple labeling can be k = 16(i-1) + j, i, j = 1,. . .,16, where the neurons are ordered "row-by-row". Such a labeling is used in Fig. 4A, B, where the coupling matrices of weakly coupled and strongly coupled regimes are presented, respectively. The coupling matrices for these cases have a relatively simple structure, and the neurons are either weakly or strongly coupled, see also To illustrate such a cluster structure of synaptic weights we assign strongly coupled neurons to separate groups, where neurons within the same group (cluster) are coupled with the strength above a certain threshold, e.g., when s ij > 1. 49. In such a way several clusters of strongly coupled neurons can be detected. We depict these clusters in the 2 Dim lattice by different colors in Fig. 4C, E, G, and I. Scanning the neuronal lattice row-by-row as explained above, we label the neurons "cluster-by-cluster" such that the neurons labeled by index k of the first cluster are followed by the neurons of the second cluster and so on. Coupling matrices for the neurons labeled in such a way are shown in Fig. 4D, F, H, and J for the coupling regimes of intermediate strength illustrated in Fig. 3B-E, respectively.
From the shape of the coupling matrices shown in Fig. 4 it follows that neurons can form a few strongly coupled clusters with or without weak interaction between clusters [ Fig. 4D, F]. The inter-cluster coupling may get somewhat stronger, but still remains much smaller than the intra-cluster coupling [ Fig. 4H]. An interesting case of two-cluster coupling is illustrated in Fig. 4J, where a group of neurons from the first cluster strongly interact with the second cluster (note the range of synaptic weights in the attached color bar). In the corresponding plot in Fig. 4I of the clusters in the neuronal lattice, this particular subgroup of neurons of the first cluster is colored in green, whereas the reminder of the first cluster is shown in dark blue.
The reported cluster structure of the coupling topology strongly influences the synchronization dynamics of the neurons in the ensemble. To illustrate this effect, we calculate the phases ψ j of the neurons and explore their dynamics for coupling matrices shown in Fig Fig. 5A. One can see that the extent of phase synchronization in the neuronal ensemble, as addressed by the values of the order parameter, strongly oscillates with relatively large periods. These oscillations originate from the time-varying distribution of the neuronal phases, where large values of the order parameter R reflect a narrow phase distribution, whereas a symmetric (or close to that) distribution of the phases on the unit circle leads to small values of R. We found that the cluster structure of  Fig. 2D]. In this paper we suggest an approach for shifting the neuronal ensemble from the strongly coupled and synchronized regime to the weakly coupled and desynchronized regime by a remarkably weak external simulation. Single-epoch CR stimulation The impact of CR stimulation on the neuronal ensemble with STDP is illustrated in Fig. 6. We are particularly interested in the behavior of synaptic coupling and the extent of synchronization during and after CR stimulation. It is known that CR stimulation has a desynchronizing effect on the stimulated oscillators [17,18,49,50]. CR stimulation can induce an anti-kindling [23], i.e., it can initiate a process, where neuronal networks "unlearn" pathologically strong interactions and synchronization. During CR stimulation the stimulated neurons get rewired in such a way that the connectivity gets on average substantially weakened. This process is illustrated in Fig. 6A where the time course of the mean synaptic weight S(t) is shown for different stimulation intensities c. The population is initially strongly coupled and synchronized with hSi = 1.499 [compare Figs. 2, 3F, and 4B]. We show that there is a range of the CR stimulation intensity, within which a robust and long-lasting suppression of coupling and synchronization among neurons can be achieved [Fig. 6A, green curves].
We distinguish between on-stimulation and post-stimulation regimes of the ensemble during long enough CR stimulation and after long enough transients after stimulation cessation, respectively, see Fig. 6A. Accordingly, the time-averaged mean synaptic weight hSi and order parameter hRi are calculated during the last 500 s of the CR stimulation epoch and the post-stimulation transients to characterize the on-and post-stimulation regimes, respectively. The on-stimulation states just at the end of CR stimulation epoch can be considered as initial conditions for the corresponding post-stimulation states, and the on-and post-stimulation states demonstrate to a certain extent similar behavior when stimulation parameters are varied, see Fig. 6B, C. For the considered range of the stimulation intensity c, one can distinguish three intervals of c with respect to qualitatively different effects of CR stimulation on the neuronal ensemble.
• For c < 0.5, CR stimulation is weak and has no significant post-stimulation impact on the mean coupling and synchronization [Figs. 6A, red curves, and 6B, C, blue squares]. In this case, however, the synchronization is somewhat perturbed during the stimulation, but it is restored after a post-stimulation transient [ Fig. 6C]. To explain the observed effects of strong CR stimulation, we illustrate the coupling topology established in the neuronal ensemble during and after stimulation in Fig. 7. Already at optimal stimulation intensities, e.g., for c = 0.95, four coupling clusters are formed during the stimulation in the vicinity of the stimulation sites [ Fig. 7A]. These clusters are visualized in Fig. 7B, where the individual synaptic weights s mk from pre-synaptic to post-synaptic neurons in the vicinity of the stimulation sites are encoded in color. To obtain such a plot we considered a 1Dim section of the neuronal lattice along a square curve running through all stimulation sites, see  The coupling clusters induced by moderately strong, i.e., optimal CR stimulation are small and comprise neurons within the closest vicinity of the stimulation sites [ Fig. 7A, B]. Moreover, these clusters do not interact in such a way that the neurons assigned to different stimulation sites do not develop connections to each other. Rather such a structure of the on-stimulation coupling topology results in a weak mean coupling [Fig. 6A, green curves], which is preserved if CR stimulation is completely switched off. The coupling matrix of the post-stimulation state after the optimal CR stimulation does not demonstrate any pronounced clustered structure whatsoever [ Fig. 4A]. In such a way, the on-stimulation coupling clusters induced by optimal CR stimulation get dissolved during the post-stimulation transient, and the neuronal ensemble converges to a weakly coupled and desynchronized state. The anti-kindling effect of CR stimulation depends on the length of the stimulation epoch as illustrated in Fig. 8. The time-averaged mean synaptic weight hSi of the post-stimulation regimes is calculated after skipping long enough post-stimulation transients for varying length L CR of single-epoch CR stimulation. For an optimal stimulation strength, the stimulated neuronal ensemble converges to a weakly coupled and desynchronized state if CR stimulation lasts long enough [ Fig. 8C]. However, even for optimal stimulation intensity, e.g., for c = 0.7, CR stimulation has to be administered at least for 2000 s in order to shift the system to a weakly coupled regime. Therefore, to come up with parameters for effective CR stimulation one should find both the optimal stimulation intensity and the corresponding length of the stimulation epoch.
If CR stimulation is weak, even very long single-epoch CR stimulation does not cause any significant change of the coupling and synchronization of strongly coupled neurons, see Fig. 8A, B. In fact, the post-stimulation regime is apparently independent of the stimulation duration. Only very occasionally the observed post-stimulation regimes exhibit somewhat smaller coupling, but in most cases the neuronal ensemble returns to the initial strongly coupled state [ Fig. 8A, B]. We conclude that weak CR stimulation is ineffective in inducing antikindling if administered to strongly coupled neurons. Note, CR stimulation of the same intensity c = 0.35 as in Fig. 8B can effectively decouple and desynchronize the neuronal ensemble if administered to less strongly coupled neurons, see Fig. 8D for L C > 3000 s. The considered pre-stimulation state (the state of coupling and synchrony before the onset of CR stimulation) is characterized by hSi = 0.944 and has a coupling matrix similar to that shown in Fig. 4H, where three coupling clusters lead to two clusters of 213 and 43 neurons synchronized at frequencies 4.81Hz and 5.17Hz, respectively.
The presented results suggest that the notion of weak stimulation strongly relates to the initial conditions at stimulation onset. For example, CR stimulation of intensity c = 0.35 is weak and ineffective for strongly coupled neurons [ Fig. 8B], whereas it is very effective in inducing an anti-kindling for intermediately coupled neurons [ Fig. 8D]. It is thus reasonable to consider the notion of weak stimulation in the above broad sense, when the stimulation is ineffective in inducing the desirable effect either due to small stimulation intensity or stability properties of the stimulated system and its resistance to the external stimulation. We consider weak CR stimulation for the case when it is administered with intensity c = 0.35 to strongly coupled and synchronized neurons [ Fig. 8B] and suggest a novel approach to improve its efficacy.

Spaced CR stimulation
We show that weak CR stimulation can induce anti-kindling if it is appropriately applied. Taking into account the great multistability induced in the neuronal ensemble by STDP [ Fig. 2] we suggest to administer the stimulation in a spaced manner, where stimulation-epochs of CR stimulation with the stimulation signals as in Fig. 1B are recurrently followed by pause-epochs, where the stimulation is completely switched off. Time courses of the mean synaptic weight S and the order parameter R during spaced CR stimulation are illustrated in Fig. 9A and B (red curves), respectively. For long enough stimulation-and pause-epochs the neuronal ensemble converges to a weakly coupled and desynchronized regime after a few spaced stimulation cycles. The latter consist of a stimulation-epoch and the subsequent pause-epoch. In Fig. 9 the duration of the stimulation-and pause-epochs are chosen large enough, L stim = L pause = 3600 s, to explore and illustrate the dynamics of the coupling during each epoch. For stimulation-and pause-epochs of this duration the mean coupling saturates during many of these epochs, and the neuronal ensemble approaches the target state of weak coupling in a stepwise manner. Moreover, during each pause-epoch the system is attracted to one of the multiple coexisting stable states which may have a weaker connectivity than that of the previous pause-epoch, see pause-epochs 6, 8, 10, 12, and 14. Accordingly, as the number of stimulation-epochs increases, stimulation is delivered to a neuronal population with progressively weaker synaptic connectivity. Eventually, the spaced procedure decouples and desynchronizes neurons [ Fig. 9].
We address the mechanism of the spaced stimulation by inspecting the dynamics of the coupling in the ensemble at each transition from stimulation-to pause-epochs. The pre-stimulation state (before the first stimulation-epoch) is characterized by strong coupling such that all neurons form a single cluster with all-to-all coupling with s ij = s max , see Fig. 4B. Weak CR stimulation slightly perturbs and weakens the coupling in the ensemble during stimulation-epochs. During the following pause-epochs the coupling further evolves, and some of the neurons may split (decouple) from the main (largest) coupling cluster. We illustrate such transitions in  Fig. 8D.
We explore the properties of the suggested spaced protocol of CR stimulation by performing a large number of numerical simulations for varying lengths L stim and L pause of the stimulationand pause-epochs, respectively. In each run the neuronal ensemble with STDP (1) starts in the strongly coupled regime and undergoes 100 cycles of length L stim + L pause , where CR stimulation of intensity c = 0.35 is delivered during the stimulation-epochs only. The time-averaged mean synaptic weight hSi of the post-stimulation states is collected and analyzed. We found that for short pause-epochs the stimulated neurons remain strongly coupled, see Fig. 11A for L pause < 40 s. On the other hand, if L pause increases, the coupling in the neuronal ensemble gets weaker, see Fig. 11A for L pause > 40 s. Hence, effective spaced CR stimulation requires pauses (i.e., pause-epochs) of sufficient length.
The coupling in the neuronal ensemble can further be suppressed for larger L pause [Fig. 11B]. However, the stimulation-epochs also have to be of sufficient length in order to reach the states with smaller coupling [Fig. 11C]. For example, for L stim = 80 s, spaced CR stimulation can shift the connectivities close to 0.8-0.9 only, although the length of the pauseepochs was varied up to L pause = 960 s [ Fig. 11A-C]. On the other hand, for L stim = 560 s, the   Fig. 11C]. Further prolongation of the stimulation-epochs can lead the system to the desired weakly coupled and desynchronized state, for instance, with hSi % 0.13 for L stim = 1440 s and L pause = 160 s [ Fig. 11D]. However, if the stimulation-epochs are too long with respect to the pause-epochs, the spaced CR stimulation does no longer induce an anti-kindling, see Fig. 11D. Accordingly, a certain balance between the lengths of the stimulation-and pause-epochs is required for a pronounced decoupling and desynchronization.
In Fig. 12 we summarize our findings obtained from the set of 540 simulations for different pairs (L pause ,L stim ) mentioned above. The distribution of hSi of the post-stimulation states reached by the neuronal ensemble during these simulations is shown in Fig. 12A. Based on the shape of the distribution, one can divide the entire range of the mean coupling into 5 domains: hSi 2 (1.4, 1.5)-D0, hSi 2 (1.0, 1.4)-D1, hSi 2 (0.6, 1.0)-D2, hSi 2 (0.2, 0.6)-D3, and hSi 2 (0.05, 0.2)-D4. We ordered all values (L pause , L stim ) with respect to these domains and plot them in Fig. 12B. In particular, anti-kindling does not occur for insufficiently long stimulationand pause-epochs [ Fig. 12B, red circles]. Too long stimulation-epochs with respect to pauseepochs will also lead to failure of the stimulation in inducing anti-kindling, and the neuronal population stays in domain D0, i.e., it remains strongly coupled and synchronized. Even for very long pause-epochs, but short stimulation-epochs only a weak suppression of the coupling can be achieved, i.e., only transitions to D1 and, occasionally, to D2 are possible [Fig. 12B, blue squares and violet diamonds]. The tendency to reach domains D1, D2, and D3 grows as both L stim and L pause increase. For L pause ! 160 s the transition to D4 becomes possible for appropriately large L stim . We actually found a large region of values (L pause , L stim ), where weak spaced CR stimulation at intensity c = 0.35 decouples and desynchronizes strongly coupled neurons and reliably induces an anti-kindling [ Fig. 12B, magenta asterisks and green triangles]. Note, The stimulation outcome clusters in Fig. 12B, and, hence, for spaced CR stimulation we find optimal values of L stim and L pause for the transitions to the different domains D1-D4, respectively. An important quantity in this context is the total amount of stimulation current required. We estimate it by measuring the total stimulation time, i.e., we calculate the total length of the stimulation-epochs T stim required before the neuronal population converged to the final state D1-D4, respectively. We thus follow the mean synaptic weight S(t), detect the time moment when it enters to and remains in a small vicinity of one of the different final states D1-D4, and calculate the total stimulation time T stim . The different final domains D1-D4 are related to different characteristic properties of the stimulation time T stim [Fig. 12C, D]. In particular, the distribution of T stim is skewed for D1 and D2 [ Fig. 12C], while it is more symmetric for D3 and D4 [ Fig. 12D]. Statistical properties of the stimulation time T stim and the corresponding parameters L pause and L stim are given in Table 1.
Stronger anti-kindling requires more stimulation time T stim . In fact, the mean stimulation time hT stim i necessary to reach the domains D1-D4 nearly doubles from one domain to the next with subsequently weaker coupling [ Table 1]. There is a relatively broad spread of the stimulation time as given by its standard deviation and large differences between the mean and minimal values of T stim [Fig. 12C, D]. Small T stim provides optimal combinations of parameters The Spacing Principle for Unlearning Abnormal Neuronal Synchrony L stim and L pause , so that a given domain can be reached with a minimal amount of stimulation. Interestingly, the regimes D1 and D2 can be reached already by very short stimulation-epochs and for relatively short stimulation time T stim . Accordingly, to further optimize spaced CR stimulation one could increase the duration of the stimulation-and pause-epochs in the course of the stimulation. This might help to clear the hurdles D1-D4 at lower cost T stim . Note, effective spaced CR neuromodulation may require more total stimulation as opposed to optimal conventional (massed) CR neuromodulation. For instance, the stimulation currents administered in these two cases may relate as c spaced T stim : c opt L CR = 4822: 1400, where c spaced = 0.35 and c opt = 0.7 are the stimulation intensities, and T stim = 13777 s [ Table 1] and L CR = 2000 s [ Fig. 8C] are the total stimulation times required to induce a pronounced anti-kindling by the considered spaced CR neuromodulation and optimal conventional CR neuromodulation, respectively.

Discussion
In this paper we have computationally shown that the spacing principle may significantly enhance the desynchronizing effect of CR neuromodulation. In particular, we have demonstrated that the spacing principle may enable an anti-kindling [23], i.e. an unlearning of abnormal synaptic connectivity and abnormal neuronal synchrony, by using CR neuromodulation at very weak intensities which render permanently delivered CR neuromodulation ineffective. Put otherwise, even ineffective CR neuromodulation may become effective by using a spacing protocol. The CR approach fundamentally relies on the multistability of neuronal networks with STDP, which is a common phenomenon in such networks and has been observed for a variety of models and STDP rules [23][24][25][51][52][53][54]. In fact, multistable dynamics is relevant to a variety of learning and memory processes in the brain and is essentially involved in the emergence of different functional and pathological regimes [55][56][57][58][59][60][61][62][63].
Multistability may lead to complex stimulus responses. Initially, anti-kindling was studied by using massed CR neuromodulation of sufficient intensity and duration [23]. Stimuli of this kind shift the neuronal population close to an attractor characterized by down-regulated synaptic connectivity and neuronal synchrony. Subsequent studies employed slightly modified CR protocols to probe fundamental dose-response features of CR neuromodulation: Transitory short-term rebound of synchrony [51,52]: Shorter and/or weaker CR stimuli may still be able to shift the neuronal population into the basin of attraction of a favorable attractor (with down-regulated synaptic connectivity and neuronal synchrony), but not close to the attractor. Accordingly, after cessation of stimulation the neuronal population spontaneously relaxes towards that attractor. However, on its trajectory it may undergo different types of transients, e.g., even a transitory rebound of synchrony [51,52]. This illustrates that acute after-effects and sustained after-effects of CR neuromodulation may differ substantially (see also [64]). Accordingly, acute and sustained after-effects were separately analyzed in parkinsonian monkeys, showing that after-effects of standard high-frequency STN DBS were acute and vanished within 30 min after offset of stimulation, whereas STN CR DBS applied for only 2 h per day on 5 consecutive days led to sustained after-effects lasting for up to 30 days [27].
Cumulative effects [44]: In studies on CR DBS in parkinsonian monkeys [27] and PD patients [28] CR neuromodulation had to be delivered via externalized leads, since at that time no implantable pulse generator (IPG) dedicated for CR was available, so that CR stimulation had to be delivered through a portable stimulator specifically developed for this application [65]. For practical purposes stimulation via externalized leads could not be performed 24 h/day. Hence, the question arose as to whether an anti-kindling or at least long-term desynchronizing effects could be achieved by repetitively delivered CR stimuli of insufficient duration. For instance, imagine the simplified picture of a bistable system, in particular, with the dynamics being governed by a double-well potential, massed CR stimulation of sufficient duration is optimal, since during premature pauses the neuronal population relaxes back into the initial, pathological attractor state [44]. In this clinically motivated context it was computationally shown that CR neuromodulation has cumulative effects [44]: Two or more subsequently delivered short CR stimuli with appropriately bounded pauses have an effect greater than the sum of the effects of the individual CR stimuli. Cumulative effects of this kind were observed in PD patients treated with unilateral CR DBS of the STN [28]. Note, the motivation to repetitively deliver short stimuli was to come up with sort of a workaround, since for practical reasons stimuli of sufficient duration could not be delivered. In contrast, here we have shown that pauses need not degrade the desynchronizing effect of CR neuromodulation. In fact, the spacing principle may significantly enhance the efficacy of CR neuromodulation. This intriguing effect essentially relies on the STDP-induced multistability rather than bistability.
In the present study we used a network of excitatorily coupled Morris-Lecar bursting neurons with a symmetric STDP rule (Equation (2), see Ref. [44] for motivation). It would be interesting to investigate whether the spacing principle is also beneficial in networks composed of other types of model neurons and/or an asymmetric STDP rule. In our computational studies [25,54] we have shown that CR neuromodulation can have very similar desynchronizing effects for both direct electrical stimulation and excitatory or inhibitory synaptically mediated stimulation, which was tested in neuronal ensembles with excitatory and inhibitory interactions. We therefore suspect that the presented effects of spaced stimulation could be observed for such complicated coupling and stimulation modalities, too.
Spaced CR neuromodulation might enable to further decrease the stimulation intensity of the pulses used for CR DBS and might, hence, presumably further decrease the rate of side effects caused by, e.g., current spread. Note, the total amount of stimulation current required to achieve an anti-kindling with spaced CR neuromodulation might exceed the total amount of stimulation current required for the case of optimal conventional (massed) CR neuromodulation. Put otherwise, spaced CR DBS is a strategy which might be favorable in case of, e.g., sub-optimal electrode placement or target constrains required very weak stimulation intensities, rather than a strategy to reduce the total amount of stimulation current and, hence, the battery demand.
As shown computationally, CR neuromodulation is also effective if exclusively delivered via synapses along the lines of a non-invasive approach [25,54]. Accordingly, acoustic CR neuromodulation was developed for the treatment of chronic subjective tinnitus [54,66]. In a proof of concept-study it was shown that acoustic CR neuromodulation causes a significant and sustained reduction of tinnitus symptoms along with a concomitant reduction of abnormal neuronal synchrony [66,67] and abnormal effective connectivity within a tinnitus-related network of brain areas [68]. The spacing principle might also be beneficial to acoustic CR neuromodulation, for instance, by enabling to deliver this treatment at minimal loudness levels or by enabling optimized dosage patterns.
The spacing principle might also be used to improve other computationally developed desynchronizing stimulation techniques utilizing phase resetting [64,[69][70][71][72], delayed feedback [73][74][75] or proportional-integro-differential feedback [76] principles. However, so far, these techniques have only been applied to neuronal networks with fixed synaptic weights, and it still remains to be shown that they cause an anti-kindling in neuronal networks with STDP. In fact, uncorrelated noise delivered to a neuronal network with fixed synaptic weights is a perfect means for desynchronization, while the same stimulation pattern delivered to the same network, but equipped with STDP causes an increase of the synaptic weights and, hence, a stabilization of the neuronal synchrony [77].
In a computational study in a neuronal network without STDP it was shown that standard HF DBS may have desynchronizing effects as opposed to blocking or inhibitory effects, provided stimulation affects (nearly) exclusively excitatory targets [24]. Let us assume this holds for neuronal networks with STDP, too. In that case one might speculate whether spaced HF DBS delivered through an electrode implanted, e.g., in a fiber tract projecting onto a target population via excitatory synapses might cause an anti-kindling. In a study on closed-loop adaptive DBS (aDBS) in PD patients onsets and offsets of HF DBS were triggered by threshold crossings of the local field potential assessing beta band STN activity [78]. Effects of different stimulation protocols were studied during a few minutes of stimulus delivery. Clinical and electrophysiological (beta suppressing) effects of aDBS were significantly stronger compared to standard HF DBS as well as intermittent random DBS, where random DBS bursts were not triggered by beta LFP threshold crossings [78]. Interestingly, during aDBS on the time scale of a few minutes there was a tendency of the duration of the stimulation epochs to decrease, possibly indicative of plastic changes [78]. It might be interesting to investigate whether spaced HF DBS might be superior to intermittent random DBS (which was used as specific control for aDBS) and compare it with aDBS.
Spaced desynchronizing stimulation might, in general, serve as a control for corresponding studies or possibly even an alternative to adaptive [78] or responsive [79,80] or demand-controlled [17,23] brain stimulation techniques where threshold crossings are used to trigger stimulation, in particular, if the temporal specificity of closed-loop stimulation is actually not required to induce long-term plastic changes. By the same token, the spacing principle might be applied to design dosing protocols for other types of intermittent neuromodulation, e.g., vagus nerve stimulation for the treatment of epilepsies [81].
For the sake of illustration, in this study we used a simple spacing protocol with fixed duration of stimulation (L stim ) and pause (L pause ) epochs [ Fig. 9]. However, our results [ Fig. 12 and Table 1] suggest that for the neuronal network under consideration temporally non-uniform spacing protocols might be feasible and possibly even superior. In fact, we have shown that the neuronal network can be shifted from the strongly coupled and synchronized regime to intermediately coupled states [ Fig. 12A, regions D1 and D2] already by means of relatively short stimulation-epochs with sufficiently long pauses, requiring only small amounts of the administered stimulation as measured by the stimulation time T stim [ Fig. 12 and Table 1]. However, a more pronounced anti-kindling [ Fig. 12A, regions D3 and D4] requires much longer stimulation-epochs and pauses together with greater amounts of stimulation time T stim . Put otherwise, one might possibly use a temporally non-uniform spacing protocol where in the course of the intervention the duration of the stimulation-epochs increases. Further improvement might be obtained by means of closed-loop protocols by incorporating feedback signals appropriately reflecting the state of the network. This might eventually lead to time-dependent thresholds governing the stimulation (L stim ) and pause (L pause ) epochs.