=Paper= {{Paper |id=Vol-1662/mod9 |storemode=property |title=Analysis of noise-induced transitions between spiking and bursting regimes in Hindmarsh-Rose neuron model |pdfUrl=https://ceur-ws.org/Vol-1662/mod9.pdf |volume=Vol-1662 |authors=Lev B. Ryashko,Evdokia S. Slepukhina }} ==Analysis of noise-induced transitions between spiking and bursting regimes in Hindmarsh-Rose neuron model== https://ceur-ws.org/Vol-1662/mod9.pdf

Analysis of noise-induced transitions between spiking
and bursting regimes in Hindmarsh–Rose neuron model

Lev B. Ryashko Evdokia S. Slepukhina
lev.ryashko@urfu.ru evdokia.slepukhina@urfu.ru
Ural Federal University (Yekaterinburg, Russia)

Abstract
The stochastic dynamics of the Hindmarsh–Rose model of neuronal ac-
tivity is studied. For the parametric zone of tonic spiking oscillations,
it is shown that random disturbances transform the spiking dynamic
regime into the bursting one. For a quantitative analysis of the noise-
induced bursting, we suggest a constructive approach based on the
stochastic sensitivity function technique and the method of confidence
domains. It allows us to give a geometric description for a distribu-
tion of random states around the deterministic attractors and estimate
critical values for the noise intensity corresponding to the qualitative
changes in stochastic dynamics. We show that the obtained estimations
are in a good agreement with direct numerical simulations.

1 Introduction
Random perturbations can considerably affect the behavior of dynamical systems. Neuronal models are very
sensitive to noise, and even small stochastic fluctuations can lead to significant qualitative changes in properties
of such systems.
One of the important mechanisms of neuronal behavior is excitability. Under the external impulse, a neuron
can generate either individual spikes (i.e. abrupt changes in the electrical potential across a cell’s membrane)
or sequences of periodic spikes. A special type of neuronal activity is bursting: a mode when periodic spiking
alternates with intervals of quiescence. [1].
The three-dimensional Hindmarsh–Rose (HR) [2] system is one of the simplest models representing the bursting
neural activity.
The deterministic HR model is fairly well studied. The detailed bifurcation analysis in different parametrical
zones can be found in [3–5]. As for the stochastic case, the effects of noise on the system with periodic external
impulse, such as stochastic resonance [6, 7] and coherence resonance [8] were studied.
The original deterministic HR system describes a wide range of dynamic regimes, such as periodic oscillations
of various types, oscillations zones with period doubling and adding, coexistence of several attractors, chaos.
This article studies the effect of random disturbances on the dynamics of the HR model in the parametric zone
of tonic spiking oscillations. We show that under noise the spiking dynamic regime transforms into the bursting
one. This phenomenon is confirmed by changes of distribution of random trajectories in the phase space.
The most detailed description of the stochastic attractors in terms of probability density function is given by
Kolmogorov-Fokker-Planck (KFP) equation. However, a direct usage of this equation is very difficult technically

Copyright c by the paper’s authors. Copying permitted for private and academic purposes.
In: A.A. Makhnev, S.F. Pravdin (eds.): Proceedings of the 47th International Youth School-conference “Modern Problems in
Mathematics and its Applications”, Yekaterinburg, Russia, 02-Feb-2016, published at http://ceur-ws.org

306
even in simple cases, therefore various asymptotics and approximations were developed. For an approximation of
KFP solutions, a well-known quasipotential method [9,10] and a stochastic sensitivity function technique [14–16]
can be used.
The stochastic sensitivity function (SSF) technique allows us to construct confidence domains that are simple
and evident geometrical models for a spatial description of a configurational arrangement of random states
around the deterministic attractors. In this article, for a quantitative analysis of the noise-induced phenomena
in the HR model, we suggest a constructive semi-analytical approach based on the stochastic sensitivity function
technique and the method of confidence domains.

2 Stochastic sensitivity analysis: theoretical background
Consider a nonlinear system of stochastic differential equations:

dx = f (x) dt + εσ(x) dw(t). (1)
Here, x is n-vector function, w(t) is n-dimensional standard Wiener process, ε is a scalar parameter of noise
intensity. Assume that a corresponding deterministic system (ε = 0) has a T -periodic solution x = ξ(t) with an
exponentially stable phase curve Γ (limit cycle).
Under stochastic disturbances, random trajectories of the system (1) leave the deterministic cycle Γ and form
some bundle around it. A detailed probabilistic description of the stochastic trajectories in this bundle is given
by Kolmogorov-Fokker-Planck (KFP) equation. In a steady regime, one can consider stationary probability
density function ρ(x, ε) governed by the stationary KFP equation.
To avoid well-known technical difficulties with the direct use of this equation, various asymptotics and ap-
proximations are developed [11–13]. For the approximation of KFP solutions, a well-known quasipotential
method [9, 10] and a stochastic sensitivity function (SSF) technique [14, 15] can be applied.
Let Πt be a hyperplane that is orthogonal to the limit cycle at the point ξ(t) (0 ≤ t < T ). For this plane, in
the neighborhood of the point ξ(t), a Gaussian approximation of the stationary probabilistic distribution can be
written [14] as:
(x − ξ(t))> W + (t)(x − ξ(t))

ρt (x, ε) = K exp −
2ε2
with the mean value mt = ξ(t) and the covariance matrix D(t, ε) = ε2 W (t). The stochastic sensitivity matrix
W (t) is a unique solution of the boundary problem

Ẇ = F (t)W + W F > (t) + P (t)S(t)P (t) (2)

with conditions
W (T ) = W (0), W (t)r(t) = 0.
Here
∂f
F (t) = (ξ(t)), S(t) = G(t)G> (t), G(t) = σ(ξ(t)),
∂x
rr>
r(t) = f (ξ(t)), P (t) = Pr(t) , Pr = I − .
r> r
The eigenvalues λi (t) and the eigenvectors vi (t) of the SSF matrix characterize the distribution of random
states in the Poincare section Πt near the point ξ(t) of the cycle.
The value M = max λ1 (t) is a useful characteristic of the cycle as a whole. We consider M as a stochastic
[0,T ]
sensitivity factor of the limit cycle Γ.
SSF matrix allows to construct the confidence ellipse with the center in point ξ(t). The equation of this ellipse
in plane Πt looks like
(x − ξ(t))> W + (t)(x − ξ(t)) = 2k 2 ,
where the parameter k determines a fiducial probability P = 1 − e−k . A set of these ellipses for t ∈ [0; T ) specify
some confidence torus around a deterministic cycle. This torus is a confidence domain in a phase space for the
stochastic cycle as a whole [16].

307
3 Deterministic Hindmarsh–Rose system
Consider the 3D Hindmarsh–Rose (HR) model [2]:
ẋ = y − x3 + 3x2 + I − z

ẏ = 1 − 5x2 − y (3)

ż = r(s(x − x0 ) − z),
where x is a membrane potential, variables y, z describe ionic currents, I is an external current; 0 < r 1 is a
time scale parameter; s, x0 are other parameters.
In this paper, we fix r = 0.002, s = 4, x0 = −1.6. The dynamics of the system (3) is considered under
variation of the parameter I.
Fig. 1 shows a bifurcation diagram of the deterministic system on the interval I ∈ (1.2, 4.0). Here, z-
coordinates of equilibrium points (green solid for stable equilibrium, red dashed for unstable one), and z-
coordinates of points of intersection of limit cycles and chaotic attractors with Poincare section plane x = 0
(blue solid) are plotted.
The system (3) has a single equilibrium which is stable for I < I1 ≈ 1.288. It loses stability due to the
subcritical Andronov–Hopf bifurcation when the parameter I passes from left to right through the point I1 .
For I < I0 ≈ 1.2677 the stable equilibrium is a unique attractor of the system. As the parameter passes
the point I0 , a stable limit cycle appears as a result of a saddle-node bifurcation. Such type of cycle is termed
”burst”. The bursting activity is observed in the system for I0 < I < I2 ≈ 3.292. For I0 < I < I1 the
system exhibits a coexistence of the stable limit cycle and the stable equilibrium. Initially, near I0 , a limit
cycle emerges with one spike in a burst. Then, with the increase of I, a sequence of period-adding bifurcations
occurs: bursts with 2,3,4, ..., 13 spikes appear. A transition from burst with 12 spikes to burst with 13 spikes is
accompanied with the transition to chaos. For I > I2 ≈ 3.292, a transition from the bursting chaos to the spiking
chaos occurs. With further increase of I, the intermittency of chaotic and regular spiking is observed, then, at
I ≈ 3.316, the transition from the chaotic spiking regime to regular one occurs. A sequence of backward period
doubling bifurcations takes place: 4-cycle transforms to 2-cycle at I ≈ 3.328, then a transition from 2-cycle to
1-cycle occurs at I ≈ 3.364. Finally, the limit cycle transforms to the stable equilibrium due to the supercritical
Andronov–Hopf bifurcation at I ≈ 25.261.
Fig. 2 shows the examples of typical regular bursting (Fig. 2a) and tonic spiking (Fig. 2b) attractors.
In this paper, we focus on the zone of transition from the bursting regime to the spiking one: I ∈ (3.0, 4.0).

Figure 1: Bifurcation diagram: z-coordinates of stable (green solid) and unstable (red dashed) equilibrium points,
and z-coordinates of points of intersection of limit cycles and chaotic attractors with x = 0 plane (blue solid).

308
a)
x z
z 1.5 3.5
3.4
3.4
1
3.3
3.3
3.2
0.5
3.1 3.2
3 0
3.1
2.9
−0.5
3
y −1
−2 x 2.9
−4 1
0.5
−6 0
−8 −0.5 0 200 400 600 800 t 0 200 400 600 800 t
−1

b)
x z
z 1.5 3.5

3.4
1
3.5
3.3
0.5
3.2
3.4 0
3.1
−0.5
3
y −1
−2 x 2.9
−4 1
0.5
−6 0
−8 −0.5 0 200 400 600 800 t 0 200 400 600 800 t
−1

Figure 2: Deterministic attractors and time series x(t), z(t) for a) I = 3.24, b) I = 3.4.
4 Stochastic Hindmarsh–Rose system
Consider now the stochastically forced HR model:

ẋ = y − x3 + 3x2 + I − z + εẇ,

ẏ = 1 − 5x2 − y (4)

ż = r(s(x − x0 ) − z),
where w is a standard Wiener process with E(w(t) − w(s)) = 0, E(w(t) − w(s))2 = |t − s| and ε is a noise
intensity.
Let us study the effect of random disturbances on the system (4) in the zone of tonic spiking oscillations.
Consider the value I = 3.7. Here, the spiking cycle with one spike per period is the attractor of the determin-
istic system. Fig. 3 shows random trajectories starting from this deterministic cycle (in the projection on xOz
plane) and corresponding x(t) plots for two values of noise intensity. For a sufficiently small noise (ε = 0.01),
random trajectories are localized in some small vicinity of the deterministic limit cycle. The type of oscillations
remains spiking, with one spike per period (Fig. 3a). This can be considered as noisy tonic spiking oscillations. If
the noise intensity is greater than some threshold value (ε = 0.1), random trajectories go far from the determin-
istic cycle, with a sharp decrease of z-coordinates, and then an oscillatory transition process of approaching to
the deterministic cycle is observed (Fig. 3b). On x(t) plot, one can observe clearly defined intervals of quiescence
and spiking phases. This indicates that under noise, the type of oscillations changed from spiking to bursting.
Fig. 4 shows how noise changes the dispersion of random trajectories for I = 3.7. Consider Poincare section
plane x = x̄, where x̄ is x-coordinate of the equilibrium. In Fig. 4a, z-coordinates of points of intersection of
random trajectories with this Poincare section are plotted. For small noise intensities, random states concentrate
near the deterministic limit cycle and have a sufficiently small dispersion. With an increase of noise, a significant
deviation of z-coordinates in the direction of smaller values is observed. This corresponds to the emergence of
bursting oscillations.
Fig. 4b demonstrates the details of the distribution of random trajectories. The probability density function
P (z) for z-coordinates of points of intersection of random trajectories with x = x̄ is plotted for I = 3.7 for

309
a)
z
x z
3.8

3.7 1
3.7

3.6
0
3.5 3.5

3.4
−1

3.3 3.3

−1 0 1 x 0 200 400 600 800 t 0 200 400 600 800 t
b)
z
x z
3.8

3.7 1
3.7

3.6
0
3.5 3.5

3.4
−1

3.3 3.3

−1 0 1 x 0 200 400 600 800 t 0 200 400 600 800 t

Figure 3: Stochastic generation of bursting oscillations for I = 3.7: a) ε = 0.01, b) ε = 0.1.

various noise intensities. For a weak noise, peaks of P (z) are located above the limit cycle. For greater noise
intensity values, one can observe a shift of the peak to the zone with smaller z-coordinates.
Note that for bursting oscillations, random states with x < −1 are observed (the quiescence zones in Fig.
3b). The value x = −1 can be used as a threshold that separates in the phase space the quiescence from the
spiking phase. For a quantification of the weight of spiking time in the total time of the observation, consider
T
the numerical characteristic ζ = Tq , where Tq is the time spending by the system in the region x < −1 and
T is the total time. In Fig. 5, functions ζ(ε) for different I are plotted. For a small noise, random states are
concentrated near the deterministic cycle, so ζ = 0. For greater noise intensity values, the bursting oscillations
are observed, and ζ increases.
Using the function ζ(ε), we can make empirical estimations for critical values of noise intensity corresponding
to the transition from the spiking regime to the bursting one. For I = 3.5 we estimate ε∗ ≈ 0.006, for I = 3.7
we have ε∗ ≈ 0.02 and for I = 3.9 we get ε∗ ≈ 0.04.
a) b)
z P ε = 0.01
3.8
ε = 0.02
3.7 ε = 0.03
ε = 0.05
3.6 ε = 0.1

3.5

3.4

3.3

−3
10 10
−2
ε 0
3.3 3.4 3.5 3.6 3.7 3.8 z

Figure 4: Distribution for points of intersection of random trajectories with Poincare section x = x̄ in dependence
of noise intensity for I = 3.7: a) z-coordinates of points, b) probability density functions.

310
ζ I = 3.9
I = 3.7
0.4 I = 3.5
I = 3.4
0.3

0.2

0.1

10
−3
10
−2
10
−1
ε
Figure 5: Function ζ.

The emergence of noise-induced bursting oscillations can be explained by peculiarities of the phase portrait
of the deterministic system. Trajectories starting close to the deterministic cycle, tend to it, but a character of
a transient process depends on the value of the initial deviation. Indeed, if the deviations are relatively small,
the trajectories tend to the cycle monotonically. If we take the initial deviations larger than some threshold, the
trajectory goes sufficiently far from cycle and approaches to it making several loops that correspond to spikes
(Fig. 6). Thus, considering different initial points, one can specify a border between these two different transient
regimes in the phase space. Let us define this border by term pseudo-separatrix.
To analyze the mechanism of stochastic generation of bursting oscillations, we apply the stochastic sensitivity
function (SSF) technique [14, 15]. SSF allows us to approximate a dispersion of random states around the
deterministic attractors.
Eigenvalues and eigenvectors of the SSF define a geometric configuration of the confidence domains, i.e., the
variance of random states around the stable attractor. For a stable limit cycle, eigenvalues and eigenvectors of
the SSF depend on time along the limit cycle. They form a family of confidence ellipses around the stable limit
cycle. The eigenvalues of SSF may vary nonuniformly along the limit cycle. This is shown on Fig. 7a, which
displays the non-zero eigenvalues of the SSF matrix of the limit cycle for I = 3.7. One can observe a significant
overfall of values along the cycle. The largest values of SSF correspond to the spiking phase of the limit cycle.
Let us consider a point ξ of the limit cycle from ”transition region”, i.e. a part of the cycle from which the
stochastic trajectories go off to the bursting zone in the phase space (see Fig. 7b). Note that the transition
region is not a part of the cycle with maximal stochastic sensitivity, but it is a zone in which large sensitivity
combines with the proximity to the pseudo-separatrix.

z z

3.5 3.59

3.4
3.58

3.3

3.57
3.2

−1 0 1 x −0.93 −0.92 x
Figure 6: Deterministic phase trajectories for I = 3.5.

311
a) b)
λ1,2
z
1
10
3.78

0 3.77
10 x
3.76
1
10
−1 3.75

0
−2
y 0
10 −2
0 5 10 15 20 25 t −4 −1

Figure 7: Stochastic sensitivity of limit cycle (a) and transition region (b) for I = 3.7.
Consider a hyperplane Π that is orthogonal to the cycle at the point ξ. For points of intersection of random
trajectories with the plane Π, we can construct confidence ellipses corresponding to different values of noise
intensity. For small noise, confidence ellipses give a good approximation for the dispersion of random states
around the deterministic cycle (see Fig. 8a).
Examining different initial points in the plane Π, one can define a pseudo-separatrix, specifying a border
between two different transient regimes in the phase space. To determine the type of the transient regime
(spiking of bursting), we use the threshold value x = −1: if a deterministic trajectory has points with x < −1,
we assume that the type of the transient regime is bursting.
Let us suggest a method allowing to estimate the critical value of the noise intensity ε∗ corresponding to the
onset of noise-induced bursting. Fig. 8b shows a point of cycle in the transition region for I = 3.7, pseudo-
separatrix, and confidence ellipses for two values of the noise intensity (a fiducial probability is P = 0.99). For
the noise intensity less than critical value (ε = 0.01), a confidence ellipse is close to the deterministic cycle and
do not intersect the pseudo-separatrix. With an increase of the noise intensity (ε = 0.02), the ellipse expands
and intersects the pseudo-separatrix (see Fig. 8b). This means that with the high probability (P = 0.99),
stochastic trajectories can go to the bursting zone of the phase space. The noise intensity that corresponds to
the intersection of the confidence ellipse with pseudo-separatrix can be used as an estimation of the threshold
value ε∗ . For I = 3.7, we get ε∗ ≈ 0.02. The obtained value is in a good agreement with the results of the
a) b)

ε = 5×10
−5 ε = 0.01
ε = 1×10−4 ε = 0.02

Figure 8: Confidence ellipses method for estimation of critical noise intensity for I = 3.7, P = 0.99: confidence
ellipses (solid, dashed) and pseudo-separatrix (dash-dotted).

312
numerical simulations.

Conclusion
The effect of random disturbances on the Hindmarsh–Rose model in the parametric zone of tonic spiking os-
cillations was studied. We showed that under noise, the spiking dynamic regime transforms into the bursting
one. For a quantitative analysis of the noise-induced bursting, we suggested a constructive approach based on
the stochastic sensitivity function technique and the method of confidence domains. It allows us to describe
geometrically a distribution of random states around the deterministic attractors and estimate critical values for
the noise intensity corresponding to the qualitative changes in stochastic dynamics. The obtained estimations
are in a good agreement with direct numerical simulations.

Acknowledgements
The work was supported by the Government of the Russian Federation (Act 211, contract No. 02.A03.21.0006)
and the Russian Foundation for Basic Research (project No. 16-31-00317 mol a).

References
[1] E. M. Izhikevich. Dynamical Systems in Neuroscience: The Geometry of Excitability and Bursting. Cam-
bridge: MIT Press, 2007.

[2] J. L. Hindmarsh, R. M. Rose. A model of neuronal bursting using three coupled first order differential
equations. Proc R Soc Lond B Biol Sci, 221:87–102, 1984.

[3] M. Desroches, T. Kaper, M. Krupa. Mixed-mode bursting oscillations: Dynamics created by a slow passage
through spike-adding canard explosion in a square-wave burster. Chaos, 23:046106, 2013.

[4] G. Innocenti, A. Morelli, R. Genesio, A. Torcini. Dynamical phases of the Hindmarsh–Rose neuronal model:
Studies of the transition from bursting to spiking chaos. Chaos, 17:043128, 2007.

[5] A. Shilnikov, M. Kolomiets. Methods of the qualitative theory for the Hindmarsh–Rose model: A case study
– a tutorial. Int J Bifurcation Chaos, 18:2141–2168, 2008.

[6] A. Longtin. Autonomous stochastic resonance in bursting neurons. Phys Rev E, 55:868–876, 1997.

[7] S. Reinker, E. Puil, R. M. Miura. Resonances and noise in a stochastic Hindmarsh–Rose model of thalamic
neurons. Bull Math Biol, 65:641–663, 2003.

[8] X. Shi, Q.-S. Lu. Coherence resonance and synchronization of Hindmarsh–Rose neurons with noise. Chinese
Physics, 14:1088–1094, 2005.

[9] M. Dembo, O. Zeitouni. Large deviations techniques and applications. Boston: Jones and Bartlett Publishers,
1995.

[10] M. I. Freidlin, A. D. Wentzell. Random Perturbations of Dynamical Systems. New York: Springer, 1984.

[11] C. Kurrer, K. Schulten. Effect of noise and perturbations on limit cycle systems. Physica D, 50:311–320,
1991.

[12] B. Lindner, L. Schimansky-Geier. Analytical approach to the stochastic FitzHugh–Nagumo system and
coherence resonance. Phys Rev E, 60:7270–7276, 1999.

[13] G. N. Milshtein, L. B. Ryashko. A first approximation of the quasipotential in problems of the stability of
systems with random non-degenerate perturbations. J Appl Math Mech, 59:47–56, 1995.

[14] I. A. Bashkirtseva, L. B. Ryashko. Stochastic sensitivity of 3D-cycles. Mathematics and Computers in
Simulation, 66:55–67, 2004.

[15] I. Bashkirtseva, L. Ryashko. Sensitivity analysis of stochastic attractors and noise-induced transitions for
population model with Allee effect. Chaos, 21:047514, 2011.

313
[16] L. Ryashko, I. Bashkirtseva, A. Gubkin, P. Stikhin. Confidence tori in the analysis of stochastic 3D-cycles.
Mathematics and Computers in Simulation, 80:256269, 2009.

314