Royal Society Publishing

Multi-frequency activation of neuronal networks by coordinated reset stimulation

Borys Lysyansky, Oleksandr V. Popovych, Peter A. Tass


We computationally study whether it is possible to stimulate a neuronal population in such a way that its mean firing rate increases without an increase of the population's net synchronization. For this, we use coordinated reset (CR) stimulation, which has previously been developed to desynchronize populations of oscillatory neurons. Intriguingly, delivered to a population of predominantly silent FitzHugh–Nagumo or Hindmarsh–Rose neurons at sufficient stimulation amplitudes, CR robustly causes a multi-frequency activation: different Arnold tongues such as 1 : 1 or n : m entrained neuronal clusters emerge, which consist of phase-shifted sub clusters. Owing to the clustering pattern the neurons' timing is well balanced, so that in total there is no synchronization. Our findings may contribute to the development of novel and safe stimulation treatments that specifically counteract cerebral hypo-activity without promoting pathological synchronization or inducing epileptic seizures.

1. Introduction

Normal brain function requires an appropriate and well-coordinated action of related neuronal populations. The level of neuronal activity, e.g. in terms of the mean neuronal discharge rate in a neuronal population, has to be appropriate. Neurons should neither spike or burst at abnormally high frequencies [1], nor be hypo-active or even silent [2,3]. Hypo-active brain regions can be found in a number of diseases, e.g. in dementia [2,3], schizophrenia [46] and bipolar affective disorders [5,7]. In addition, interactions within and between neuronal populations are crucial for normal brain function. Apparently, uncorrelated [8] or well-dosed, e.g. task-related, moderate and transient synchronization processes [9] are significant for neuronal information processing. In contrast, excessive, pathologically strong neuronal synchronization severely impairs brain function. This is well known from a number of brain diseases, e.g. from Parkinson's disease [1,8] and tinnitus [10,11]. By the same token, boosting pathological synchronization processes, e.g. by periodic stimulation, may cause symptoms to deteriorate [1214] and induce an increase in synchronization that outlasts stimulation offset [15]. Analogously, kindling stimulation protocols may evoke epileptiform activity [16].

This computational study is devoted to the question of whether it is possible to specifically counteract neuronal hypo-activity in a safe manner, i.e. without running the risk of inducing abnormal synchronization. More specifically, we computationally show that it is possible to stimulate predominantly silent neuronal populations in such a way that their mean firing rate increases without a net increase in their synchronization. Note, this is not a trivial task. Of course, a sufficiently strong stimulus periodically delivered to a predominantly silent neuronal network is able to increase the mean firing rate (e.g. [13]). However, the periodic stimulation protocols tested so far typically cause a significant increase in the neuronal synchronization. This is because periodic stimulation of sufficient strength causes neuronal discharges that are time-locked to the stimulus and, hence, mutually synchronized [17,18].

To overcome this issue, we use coordinated reset (CR) stimulation [19,20]. This stimulation technique was initially developed for a different purpose. CR stimulation specifically counteracts abnormally strong synchronization in populations of oscillatory neurons by desynchronization. The desynchronizing effect of CR stimulation has been shown computationally in networks of phase oscillators [19,20] and bursting neurons [2124] as well as experimentally in bursting electroreceptors of the paddle fish [25] and in rat hippocampal slice rendered epileptic by magnesium withdrawal [15]. We here consider two qualitatively different types of neuronal networks, a population of predominantly silent FitzHugh–Nagumo (FHN) neurons [26] and a population of predominantly silent Hindmarsh–Rose (HR) neurons [27,28]. Active FHN neurons are spiking neurons, whereas active HR neurons are bursters. Intriguingly, applied to a population of predominantly silent FHN or HR neurons at sufficient stimulation amplitudes, CR causes a multi-frequency activation, i.e. the emergence of clusters of different resonant frequencies that split into different phase-shifted sub-clusters, hence reliably giving rise to an overall desynchronization.

2. Activation of a single neuron

We start by considering the impact of a periodic stimulation on an intrinsically silent single neuron. For this, in the absence of an external input or stimulation the neuron does not generate any spike or burst. In this section, we consider the FHN model [26] of a spiking neuronEmbedded Image 2.1The variable vj models the dynamics of the membrane potential of the individual neuron and the synaptic current Isyn,j models the input from the other neurons of the network. As mentioned above, in this section, we focus on the dynamics of a single neuron and, thus, set Isyn,j = 0 (no interactions with other neurons) and omit the index j. The natural dynamics of model (2.1) without stimulation (Istim = 0) is determined by the parameter I0. For example, system (2.1) has a unique stable steady state for I0 < IH1 and for I0 > IH2, where the values IH1 ≈ 0.3313 and IH2 ≈ 1.4187 are the bifurcation points of a subcritical Hopf bifurcation. The parameter range I0 ∈ (IH1, IH2) is associated with a spiking behaviour of the model given by a stable limit cycle.

The above model neuron receives an external periodic stimulation given by the stimulation current Istim. To illustrate the impact of such a forcing, we first consider the neuron in an intrinsically oscillating regime, where system (2.1) with Istim = 0 exhibits spiking behaviour. We thus fix the parameter I0 = 0.34, where the natural spiking frequency (number of spikes per second) is fs ≈ 21.4 Hz (figure 1a).

Figure 1.

(a) Time course of the membrane potential v of the intrinsically oscillating FHN neuron (2.1) spiking at the frequency fs ≈ 21.4 Hz without coupling (Isyn = 0) and stimulation (Istim = 0). (be) Frequency-locking Arnold tongues versus stimulation frequency F and stimulation strength I for (b,c) a single intrinsically spiking and (d,e) a single silent FHN neuron (2.1). The neuron is stimulated with (b,d) a harmonic periodic signal and (c,e) a periodic pulse-train signal (see text for details). Four Arnold tongues are shown in colour corresponding to 1 : 1 (grey), 1 : 2 (red), 1 : 3 (green) and 1 : 4 (blue) entrainment. Magenta domains in plots (d) and (e) correspond to small sub threshold fluctuations. The inset in plot (b) is in linear-log scale, where the Arnold tongues are magnified for small values of the stimulation intensity and are shown versus the ratio F/fs of the forcing frequency F and the natural spiking frequency fs of the neuron. Vertical dashed lines in plot (e) indicate the stimulation frequencies F = 62.5 Hz and F = 83.3 Hz used below. Parameter I0 = 0.34 in (ac), which leads to the natural spiking frequency fs ≈ 21.4 Hz of the intrinsically oscillating FHN neuron (2.1), and I0 = 0.32 in (d) and (e).

It is well known that a periodic oscillator can be entrained by a sufficiently strong external periodic force [17,18]. Moreover, the entrainment can be of higher order such that the frequency of the entrained oscillator f rationally relates to the frequency of the forcing F, i.e. mF = nf for some positive integers m and n, which constitutes an m : n entrainment or frequency locking. The frequency-locking regions in the parameter space of the forcing frequency F and the forcing strength I attain the well-known form of the Arnold tongues [17,18]. Arnold tongues were originally found for the circle map and were intensively used in the context of, for example, cardiac activity [2932] and laser dynamics [33,34].

For the harmonic force Istim = I sin(2π Ft/1000), where the forcing frequency F is measured in hertz (Hz) and the time units in equation (2.1) are given in milliseconds (ms), the frequency-locking regions in the (F,I)-parameter plane are illustrated in figure 1b for the FHN neuron. For the sake of clarity only four Arnold tongues are shown, corresponding to the 1 : 1, 1 : 2, 1 : 3 and 1 : 4 entrainment. As expected, the Arnold tongues originate from the points on the horizontal axis with zero stimulation intensity I = 0, where the forcing frequency F equals the corresponding integer multiple of the oscillator's natural frequency (figure 1b inset). It is interesting to observe that the frequency-locking parameter regions have a form of bent stripes aligning one beside the other with respect to parameter I for large stimulation frequency F. This means that, if the stimulation intensity of the periodic force at a relatively high stimulation frequency is gradually increased, the forced neuron will first be entrained by the periodic force at an m : n frequency-locking ratio of rather high order. The order of the entrainment will successively decrease for larger stimulation strength, and, finally, the stimulated neuron will be 1 : 1 entrained by the periodic force and fire at the forcing frequency if the stimulation intensity is large enough.

The above property is common for different periodic stimulation forces. For example, the (F,I)-parameter plane of the forced neuron (2.1) has a similar structure for the harmonic forcing and for a periodic pulse-train stimulation, respectively (figure 1b,c). The stimulation current Istim in equation (2.1) of the periodic pulse-train stimulation consists of short high-frequency (HF) pulse trains of several rectangular pulses of the amplitude I and inter-pulse intervals Tp with pulse duration Tp/2, where Tp = 1 ms (1000 pulses s−1) for the FHN model. During every stimulation period of length T = 1/F a single pulse train of length T/4 is administered to a neuron such that the overall stimulation frequency of such a patterned periodic stimulation is F. A modelling and clinical study with this type of patterned low-frequency stimulation can be found in Barnikol et al. [13].

The discussed structure of the parameter space is preserved if the stimulation is delivered to a silent neuron (figure 1d,e). The main difference, however, is that the stimulation activates and entrains the neuron, which otherwise does not fire. Therefore, if the stimulation is weak, it cannot initiate a spiking oscillatory dynamics, and only small sub threshold fluctuations can be observed (figure 1d,e, magenta domains). Except for this, the stimulated silent neuron exhibits a high-order m : n entrainment which successfully progresses towards the 1 : 1 entrainment as the stimulation strength increases, similar to the forced intrinsically oscillating neuron (compare figure 1be). The property of the m : n activating entrainment is essentially used by the method of multi-frequency activation of neuronal networks by CR stimulation.

3. Neuronal network activation

In this section, we investigate the impact of the CR stimulation on weakly coupled neuronal networks that mostly consist of silent neurons. For the neuronal model (2.1) under consideration, the coupling is realized via the synaptic input current Isyn,j. To model a densely coupled network, where each neuron receives an input from many other neurons, we use a global coupling. Then the synaptic current Isyn,j can be related to the population field potential by taking an ensemble average of post-synaptic potentials sj,Embedded Image 3.1where j = 1, 2, … , N is the index of a neuron in the ensemble of N neurons, and each variable in equation (2.1) is equipped with such an index. K is the coupling strength and Vsyn is a reversal potential determining the type of synaptic connection [35], which is taken as Vsyn = 2 for excitatory coupling. The post-synaptic potential sj generated by neuron j is modelled in the standard way by an additional equation for sj(t) [35],Embedded Image 3.2

The neuronal ensemble (2.1)–(3.2) will be controlled by CR stimulation [19,20]. For this, the neurons are assumed to be arranged on a one-dimensional lattice, e.g. equally spaced on the linear segment [0,L] with coordinates xj = (j − 1) L/(N − 1), Embedded Image, with L = 10. The stimulation signals are delivered via Nc stimulation sites, which are equidistantly spaced within the above segment and located at the points ck = (k − 1/2) L/Nc, k = 1, 2, … , Nc. Below we consider four stimulation sites, i.e. Nc = 4. Following Tass [19,20], the stimulation currents Istim,j in equation (2.1) readEmbedded Image 3.3where I is the stimulation strength. The coefficients D(xj, k) define the impact of the stimulation signal delivered via the kth stimulation site on the jth oscillator with the spatial coordinate xj and are distance dependent, as found for brain tissue [36],Embedded Image 3.4where σ defines the spatial decay rate of the administered current with the distance from the stimulation contact (figure 2a).

Figure 2.

(a) Spatial profile functions D(x, k), Embedded Image of the current decay for four stimulation contacts equally spaced within the oscillator population arranged on a linear segment of length L = 10. Filled circles depict the positions of the stimulation sites. (b) Three stimulation cycles with four stimulation signals (brief HF pulse trains) of the CR stimulation protocol for Nc = 4 stimulation sites. Different cycles of the CR stimulation are separated by dashed lines. The encircled numbers indicate the stimulation sites in (a) and the corresponding stimulation signals in (b). The length of CR stimulation cycle T = 2 s and σ = 1.0 in equation (3.4).

CR stimulation signals consist of short HF pulse trains sequentially administered via different stimulation contacts [19,20]. The HF pulse train is given by function P(t) in equation (3.3). This function provides a periodic train of rectangular pulses of unit amplitude with inter-pulse interval Tp = 1 ms for the FHN model (2.1) with pulse duration Tp/2. The activation of the kth stimulation site is governed by the indicator function ρk (t), which attains 1 if this site is active at time t and zero otherwise. Within one stimulation cycle of length T all stimulation sites are sequentially activated, and the length of each short pulse train equals T/Nc. Exemplary time courses of the stimulation signals administered during three cycles of CR stimulation are shown in figure 2b.

We consider an ensemble of N = 400 FHN neurons (2.1), where the parameters I0,j are Gaussian distributed around the mean value I0 = 0.315 with s.d. 0.006. For such parameter values, the uncoupled neuronal ensemble (K = 0 in equation (3.1)) contains approximately 5–10% of intrinsically spiking neurons at the frequency fs ≈ 21 Hz (see §2), whereas the other neurons are silent. The fraction of the oscillating neurons in the ensemble also depends on the coupling strength, such that a strong coupling (e.g. K = 0.2) leads to an entirely oscillating network. We therefore consider a weak coupling K = 0.01, where the fraction of oscillating neurons in the ensemble remains below 10 per cent.

If a periodic pulse-train stimulation is administered to the neuronal population via only a single stimulation site, we observe a multi-frequency activation of the neuronal network (figure 3a,b). Such a stimulation corresponds to the stimulation protocol considered in §2 for a single neuron and can be obtained from the CR stimulation set-up, when three out of four stimulation sites are permanently deactivated, i.e. the corresponding indicator functions ρk(t) in equation (3.3) are permanently set to 0. Let only the second stimulation site c2 placed over the oscillator j = 150 be active. Then, depending on the distance from the stimulation site, the stimulated neurons can exhibit different dynamics. In particular, the neurons in the nearest proximity to the electrode are most strongly affected by the stimulation and, hence, can most easily follow the frequency of the external stimulation and be 1 : 1 activated and entrained by the periodic stimulation (figure 3a,b, dark red regions). In contrast, over a wide parameter range the more distant neurons exhibit a higher order m : n activation and entrainment (figure 3a,b, green and blue regions). In such a way the stimulation causes (i) an activation of the stimulated neurons and (ii) the emergence of a variety of neuronal clusters that fire at different frequencies and are m : n entrained by the stimulation.

Figure 3.

Multi-frequency activation of the network of FHN neurons (2.1) mostly consisting of silent neurons by (a,b) a periodic pulse-train stimulation delivered via a single stimulation site c2 placed over the oscillator j = 150 and (cf) by CR stimulation with Nc = 4 stimulation sites. The ratio μ = fs/F of the stimulation-induced spiking frequency fs of the neurons to the stimulation frequency F is encoded in colour for all neurons in the ensemble versus stimulation strength I ranging from dark red (1 : 1 activation and entrainment) to dark blue (no activation). Stimulation frequencies F ≈ 83.3 Hz corresponding to a stimulation period T = 12 ms (figure 1e, right vertical dashed line) in plots (a), (c) and (e) and F = 62.5 Hz corresponding to a stimulation period T = 16 ms (figure 1e, left vertical dashed line) in plots (b), (d) and (f). Coupling K = 0.01, spatial stimulation current decay rate σ = 1.0 in plots (a)–(d) and σ = 0.5 in plots (e) and (f). Parameters I0,j in equation (2.1) are Gaussian distributed around the mean value I0 = 0.315 with s.d. 0.006.

Since the coupling is weak, the collective dynamics is mainly dominated by the stimulation. Accordingly, the phenomenon of the multi-frequency activation of a neuronal network essentially relies on the higher order activation and entrainment of a single neuron by a periodic stimulation discussed in §2 (figure 1e). The distance between the stimulation site and the stimulated neuron within the network can be translated to the amount of stimulation received by the neuron according to the function D(x, k) of equation (3.4), such that the closer the neuron is located to the electrode, the more strongly it is stimulated. Therefore, if a stimulation frequency F is fixed, the displacement along the stimulated network away from the stimulation site corresponds to a decreasing stimulation strength in figure 1e (vertical dashed lines). Accordingly, a much broader spectrum of possible higher order m : n activation and entrainment can be observed for a larger stimulation frequency F. Indeed, for the considered F = 62.5 Hz (figure 3b) four active frequency clusters can easily be distinguished with the entrainment order 1 : 1, 1 : 2, 2 : 5 and 1 : 3. For the stimulation frequency F = 83.3 Hz (figure 3a) at least five neuronal frequency clusters can be found with the entrainment order 1 : 1, 1 : 2, 1 : 3, 2 : 7 and 1 : 4. The diagrams in figure 3a,b are thus tightly linked to figure 1e, and even small peculiarities, such as the thin dark red stripe of 1 : 1 entrainment within the green region of 1 : 2 entrainment in figure 3a, can be related to the form of the Arnold tongues for a periodically stimulated single FHN neuron (figure 1e, right vertical dashed line).

The activating CR stimulation administered via several stimulation sites induces a rich neuronal dynamics in the stimulated network. Around each of the stimulation sites, which are equidistantly spaced in the network, a structure similar to that in figure 3a,b is established (figure 3cf). The CR stimulation thus activates many spatially distinct neuronal clusters firing at different frequencies, which are m : n entrained by the corresponding stimulation signals. This takes place at moderate stimulation intensities, whereas a strong stimulation forces the entire neuronal network to fire at the stimulation frequency (figure 3c,d, dark red region for large I).

Obviously, the spatial decay rate of the stimulation current, given by parameter σ in equation (3.4), plays an important role concerning the spatial patterns of the multi-frequency activation. A more localized stimulation for small σ still leads to a multi-frequency neuronal dynamics in the network even for a relatively large stimulation intensity (compare figure 3c,d with 3e,f). Parameter σ and the stimulation intensity I show an interdependent impact on the activating effect of CR stimulation such that a similar cluster structure of the stimulated network for smaller σ can be obtained for stronger stimulation. However, the activation and entrainment threshold of CR stimulation with respect to the stimulation strength only weakly depends on values of σ for the neurons in the nearest proximity to the stimulation sites. Rather it depends on the structure of the Arnold tongues shown in figure 1e. From the experimental and clinical studies, it is known that a safe extracellular electrical stimulation can be several orders of magnitude larger (ranging up to a few volts) than the natural voltage (some tens of millivolts) of the neurons [27,37]. The stimulation current of an effective CR stimulation or patterned pulse-train stimulation is in the above range [13,15]. We therefore expect that the considered stimulation strength will not violate the safety requirements for clinical applications.

A larger stimulation frequency F also contributes to the number of the stimulation-induced frequency clusters. For example, for a fixed stimulation strength I = 5 (figure 3c,d, horizontal dashed line), CR stimulation with F = 83.3 Hz creates three large-frequency clusters firing at the frequencies 83.3, 46.65 and 27.77 Hz (figure 3c), whereas the neurons stimulated at F = 62.5 Hz mostly fire either at 62.5 or at 31.25 Hz (figure 3d). Figure 3cf shows a slight asymmetry in the structure of the frequency-locking domains, which is caused by the temporal activation order of the stimulation sites and which becomes less pronounced for a more localized stimulation, i.e. for smaller σ in equation (3.4) (compare figure 3c,e).

Another important characteristic of the multi-site activating CR stimulation is the phase dynamics of the stimulation-induced frequency clusters. The stimulation signals administered via different stimulation sites have a predefined phase shift of 2π/Nc caused by the activation times shifted by T/Nc from one signal to the next (see figure 2b). The activated neuronal clusters become m : n phase locked to the corresponding stimulation signals and, thus, also acquire the corresponding phase shifts. For this reason, the neurons activated by different stimulation sites do not fire simultaneously, but rather with predefined phase shifts, even if they oscillate at the same frequencies. This can clearly be observed in figure 4a, where an exemplary raster plot of the neuronal firing is shown before and during activating CR stimulation. The network is mostly silent in the stimulation-free regime for t < 1000 ms. The activating CR stimulation is switched on at t = 1000 ms and induces a regular neuronal firing, where three different spiking frequencies can be distinguished (figure 4a, red, green and blue horizontal stripes), which correspond to 1 : 1, 1 : 2 and 1 : 3 activating entrainment. Furthermore, the neurons exhibiting the same order of the entrainment, i.e. oscillating at the same frequencies, fire at different times, which are systematically shifted from one activated neuronal cluster to the next, compare, for example, the neuronal firing in different red stripes in figure 4a. Therefore, the activating CR stimulation has a twofold impact on the stimulated neuronal network. First, it induces a variety of spatially distinct frequency clusters, and, second, it induces a variety of phase clusters within the neuronal sub populations firing at the same frequency. On the population level, activating CR stimulation hence causes no net synchronization.

Figure 4.

Frequency and phase clustering induced by the activating CR stimulation in the network of FHN neurons (2.1) mostly consisting of silent neurons. (a) The raster plot of the neuronal firing in the stimulation-free regime (t < 1000 ms) and during CR stimulation (switched on at t = 1000 ms, vertical red line), where the spike onsets are indicated by black dots. The colour stripes indicate the neuronal clusters which are 1 : 1 (red), 1 : 2 (green) and 1 : 3 (blue) activated and entrained by the corresponding stimulation signals. Stimulation strength I = 5, which corresponds to the horizontal dashed line in figure 3c. (b,c) The time-averaged order parameter 〈R〉 given by equation (3.5) of the neuronal network (2.1) stimulated by (b) CR stimulation with Nc = 4 sites (as in figure 3cf) and (c) by a single-site pulse-train stimulation (as in figure 3a,b) versus the stimulation intensity I. Parameters σ = 1, F = 83.3 Hz and the other parameters as in figure 3.

The extent of synchronization in the stimulated neuronal ensemble can be estimated by the order parameter [38]Embedded Image 3.5where ψj(t) are the phases of the individual neurons, which are calculated from the oscillatory spiking dynamics. For this, the phase ψj(t) of neuron j attains the values ψj(tn) = 2πn, n = 0, 1, … , at the time moments tn of the spike onsets and linearly increases over 2π from one spike onset to the next [17]. The synchronization order parameter R ranges from 1 to 0, which corresponds to a perfect in-phase synchronization of the neurons (where all neurons fire precisely at the same time) or its complete absence, respectively.

The desynchronizing impact of the activating CR stimulation, as revealed by the order parameter R, is illustrated in figure 4b versus the stimulation strength I. As the latter increases, the stimulated neurons get activated at multiple frequencies and the order parameter decays (figure 4b for I > 2.5, compare with figure 3c), which indicates an onset of desynchronization in the activated neuronal ensemble. The desynchronization is preserved even in the uni-frequency activated regime, where, for a relatively large I, all neurons fire at the same frequency (figures 4b and 3c for I > 12). This is caused by the stimulation-induced phase clustering of CR stimulation, illustrated in figure 4a. A single-site pulse-train stimulation, on the other hand, does not have such an effect (figure 4c). For moderate stimulation strength it induces a frequency clustering (figure 3a) and slightly desynchronized dynamics, as reflected by moderate values of the order parameter (figure 4c). However, for large stimulation strength the single-site stimulation forces most of the neurons to fire at the same frequency, and the order parameter increases, which is indicative of an onset of synchronization in the neuronal network. Note that the ranges of the stimulation intensity I given in figure 4b,c for CR and single-site stimulation, respectively, correspond to the same range of the effective amount of stimulation current received on average by an individual neuron in the network per time unit. In such a way, we may directly compare the two graphs of the order parameter between the two stimulation protocols. Accordingly, single-site stimulation does not cause a pronounced desynchronization. Nevertheless, within a narrow range at weak stimulation amplitudes single-site stimulation is somewhat superior to CR stimulation with respect to the boost of synchrony. However, for sufficiently large stimulation amplitudes CR stimulation robustly activates a well-desynchronized neuronal firing within a wide parameter range. Similar results were obtained for the network of HR neurons.

4. Hindmarsh–Rose ensemble

We show the CR-induced multi-frequency activation in another, qualitatively different neuronal network. For this, we consider a network of N = 200 excitatory coupled and stimulated HR bursting neurons [27,28]Embedded Image 4.1The HR model (4.1) has a unique stable steady state, in particular, for z0 > 5.13, which corresponds to a silent regime of the neuron. The parameter value z0 ≈ 5.13 is a point of a subcritical Hopf bifurcation. For the parameter range z0 ∈ (3.16, 5.13), system (4.1) shows a bursting dynamics of the square-wave type [39] (figure 5a).

Figure 5.

(a) Time course of the membrane potential v of the intrinsically bursting HR neuron (4.1) at the frequency fb ≈ 1.6 Hz without coupling and stimulation (Isyn = Istim = 0). (b,c) Frequency-locking Arnold tongues versus stimulation frequency F and stimulation strength I for (b) a single intrinsically bursting and (c) a single silent HR neuron (4.1) stimulated with a periodic pulse-train signal (colouring as in figure 1). Stimulation frequencies F = 5 Hz and F = 6.7 Hz used in plots (d)–(f) are indicated by vertical dashed lines in (c). (d,e) Multi-frequency activation of the network of HR neurons (4.1), (3.1)–(3.3) mostly consisting of silent neurons by CR stimulation. The ratio μ = fb/F of the stimulation-induced bursting frequency fb to the stimulation frequency F is encoded in colour as in figure 3. (f) The time-averaged order parameter 〈R〉 given by equation (3.5) of the HR network stimulated by CR stimulation (solid green curve, lower horizontal axis) and by a single-site pulse-train stimulation (dashed red curve, upper horizontal axis) versus the stimulation intensity I. Parameter z0 = 5.1 in (a) and (b) and z0 = 5.2 in (c). In (d)–(f) coupling K = 0.05, and parameters z0,j are Gaussian distributed around the mean value z0 = 5.15 with s.d. 0.007. Number of the stimulation sites Nc = 4, stimulation frequency (d,f) F = 6.7 Hz and (e) F = 5 Hz, spatial stimulation current decay rate σ = 1, and inter-pulse interval Tp = 10 ms.

A single intrinsically bursting or silent HR neuron (4.1) stimulated with a periodic pulse-train stimulation with frequency F and stimulation intensity I shows the same structure of the Arnold tongues in the (F,I)-parameter plane as the FHN neuron (2.1) (compare figure 5b,c with figure 1c,e). Analogously to the ensemble of FHN neurons investigated in §3, the population of HR neurons (4.1), (3.1)–(3.3) demonstrates a very similar response to the activating CR stimulation (figure 5d,e). For parameters z0,j and coupling strength K (as in figure 5df) only about 6 per cent of the neurons in the network demonstrate an intrinsically bursting dynamics at the bursting frequency fb ≈ 1.6 Hz (figure 5a), whereas the other neurons are silent. CR stimulation activates a variety of neuronal frequency clusters, which are in agreement with the Arnold tongues found for a single HR neuron receiving a periodic pulse-train stimulation (see figure 5c, vertical dashed lines). A large stimulation frequency F amplifies the spectrum of the activating and clustering effects of CR stimulation similar to that for the FHN neurons (compare figure 5d,e). Moreover, CR stimulation has a well-pronounced desynchronizing impact on the activated neuronal population, as revealed by the synchronization order parameter R given by equation (3.5) (figure 5f, solid green curve). The single-site pulse-train stimulation, on the other hand, synchronizes the stimulated neurons as the stimulation strength increases (figure 5f, dashed red curve). Therefore, the conclusions drawn for the ensemble of FHN spiking neurons straightforwardly apply to the population of HR bursting neurons. The activating CR stimulation can thus be applied to neuronal networks of very different nature and can effectively induce an active neuronal firing on multiple frequencies based on the high-order m : n activating entrainment without net overall synchronization.

5. Discussion

Our computational results show that CR stimulation enables us to activate intrinsically silent neuronal populations without inducing synchronization. In particular, we have shown that CR stimulation at sufficient amplitudes robustly activates networks of predominantly silent FHN or HR neurons. Activation without synchronization is achieved, because CR causes a multi-frequency activation. In an Arnold tongue-like manner, different 1 : 1 or n : m entrained clusters of neurons emerge. These different resonant neuronal frequency clusters consist of phase-shifted sub clusters, respectively. Owing to the CR-induced equally spaced discharge pattern of the sub clusters, the overall timing of the neurons is well balanced, which reliably prevents it from synchronization. The multi-frequency activating effect of CR stimulation is based on the inclined structure of Arnold tongues observed in this paper for the models of periodically driven spiking and bursting neurons (figures 1 and 5). Such a parameter structure has received strong experimental support and has been found for periodically stimulated squid axon membrane and the Hodgkin–Huxley model [40] as well as for the respiratory system of cats and the corresponding model similar to the van der Pol oscillator [41]. This contributes to the generality of the method and the obtained results.

In contrast to CR stimulation, periodic stimulation is associated with a variety of dynamical phenomena that are typically related to a boost of synchronization. Systems of coupled self-sustained oscillators with global periodic forcing have been shown to generate a variety of complex dynamical effects on the way towards the full entrainment. Coupled and forced phase oscillators introduced in Sakaguchi [42] can demonstrate the emergence of so-called resonance clusters, which have been confirmed experimentally for coupled and forced electrochemical oscillators as well as numerically for the corresponding modelling system [43]. The dynamics of the individual oscillators or oscillatory neurons may significantly disagree with the dynamics of the large-amplitude field potential generated by this population [44]. On the other hand, strongly coupled and forced oscillators behave like a single giant oscillator, and, hence, the dynamics of the mean field effectively becomes low dimensional [17,18,4547]. The transition to the forced entrained regime is mediated by a complicated sequence of bifurcations of the mean field [45,46,48]. Nevertheless, this tranisition is characterized by relatively large values of the order parameter and, thus, by synchronized collective dynamics of the forced oscillators [43,44,49].

A population of intrinsically silent and periodically stimulated neurons has been studied in Barnikol et al. [13], where single-site periodic burst stimulation, as in figure 3a,b, was used both numerically and for target point localization for deep brain stimulation (DBS). In analogy to the above-mentioned studies on periodically forced active neurons, the periodically stimulated and intrinsically silent neurons are activated and synchronize as the forcing strength increases. More precisely, high levels of collective synchronization in the neuronal population require higher stimulus intensities, whereas the portion of synchronously active neurons nevertheless is strongly phase locked to the periodic force already at weak stimulus intensities. This robust entrainment effect was successfully applied as the sole functional method for DBS electrode localization to a patient with spinocerebellar ataxia type 2 (SCA2) with pronounced tremor that disappeared intraoperatively under general anaesthesia [13]. The stimulation-induced increase in synchrony corresponded to an increase in the amplitude of the evoked tremor, whereas the phase entrainment of the weak tremor at weak stimulation intensities (corresponding to only a local activation of brain tissue) displayed an optimal location of the stimulation contact.

The medical motivation behind the present study is to contribute to the design of a safe and efficient stimulation technique that enables us to activate hypo-active or even silent neuronal populations without inducing synchronized neuronal activity. For safety reasons an induction of synchronization should be avoided for the following reason. In a number of brain diseases the amount of pathological synchronization and the extent of the corresponding symptoms are closely correlated [5053]. Consequently, periodic stimulation, in particular, at resonant frequencies causes a deterioration of the symptoms [1214]. In addition, owing to spike timing-dependent plasticity (STDP) [54,55] neuronal networks may actually learn pathological synchronization: periodic stimulation causes an increase in the rate of coincidences and, in turn, an increase in the mean synaptic weights, which may shift the network from a stable state with weak synchrony and weak mean synaptic connectivity to a stable state with pathological synchrony and pathologically up regulated mean synaptic weight [2125,56]. Also, excessive synchronization may severely impair physiological information processing [57].

Activating CR stimulation might, for instance, be used for DBS. As yet, DBS has successfully been applied to patients suffering from neurological diseases like pain [58], Parkinson's disease [37,59,60], essential tremor [61], dystonia [62] and Huntington disease [63] as well as to patients with treatment-resistant depression [64,65]. Another new field of early applications of DBS is the treatment of different types of dementia, e.g. Alzheimer's disease [66] or Parkinson disease dementia (PDD) [67] and PDD-related apraxia [68]. In a patient with PDD, in addition to the bilateral subthalamic nucleus (STN) stimulation at 130 Hz, which is the typical DBS strategy for patients with Parkinson's disease, bilateral DBS of the nucleus basalis of Meynert (NBM) at 20 Hz was performed [67,68]. The STN stimulation at 130 Hz presumably exerted an inhibitory effect, while the low-frequency NBM stimulation at 20 Hz was performed to activate, for example, via cortical projections. In that patient STN DBS improved motor symptoms, but cognitive performance was almost unchanged. The additional 20 Hz NBM DBS caused a pronounced improvement in cognitive functions [67] and apraxia [68]. For safety reasons and in order to avoid excessive synchronization that might impair physiological information processing (e.g. [57]), CR stimulation at 20 Hz could be a beneficial alternative to the simple periodic stimulation at 20 Hz used so far. However, given the topology of the currently available depth electrodes for the STN or for the thalamus, delivering electrical CR stimulation to the NBM still requires the development of depth electrodes with suitable topology of the stimulation contacts. In particular, owing to the size of the NBM the contacts' spacing has to be sufficiently narrow in order to guarantee that an appropriate number of contacts, say four, are located within the target area. According to computational studies, CR stimulation with ineffective electrodes is sub optimal or even inefficient [69]. An extreme limiting case of CR stimulation, with only one contact being placed in the target area, is the single-site stimulation from figures 3a,b, 4c and 5f. As shown above, such a setting leads to an activation that is typically associated with a significant increase in synchrony.


  • Received October 10, 2010.
  • Accepted November 8, 2010.


View Abstract