- Email: [email protected]
A minimal model for a slow pacemaking neuron
Chaos, Solitons & Fractals 45 (2012) 640–644
Contents lists available at SciVerse ScienceDirect
Chaos, Solitons & Fractals Nonlinear Science, and Nonequilibrium and Complex Phenomena journal homepage: www.elsevier.com/locate/chaos
A minimal model for a slow pacemaking neuron D.G. Zakharov a,⇑,1, A. Kuznetsov b,2 a b
Nonlinear dynamics department, Institute of Applied Physics, RAS, 46 Ulyanov St., Nizhny Novgorod 603950, Russia Department of Mathematical Sciences and Center for Mathematical Biosciences, IUPUI, 402 N. Blackford St., Indianapolis, IN 46202, USA
a r t i c l e
i n f o
Article history: Available online 24 February 2012
a b s t r a c t We have constructed a phenomenological model for slow pacemaking neurons. These are neurons that generate very regular periodic oscillations of the membrane potential. Many of these neurons also differentially respond to various types of stimulation. The model is based on FitzHugh–Nagumo (FHN) oscillator and implements a nonlinearity introduced by a current that depends on an ion concentration. The comparison with the original FHN oscillator has shown that the new nonlinear dependence allows for differentiating responses to various stimuli. We discuss implications of our results for a broad class of neurons. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction There are several key groups of neurons in the central nervous system that generate very regular periodic oscillations of the membrane potential. The examples of these neurons are serotonin-containing neurons from the raphe nuclei [1], noradrenergic neurons located in the pontine nucleus locus coeruleus (LC) [2] and dopaminergic neurons from the substantia nigra pars compacta [3]. These slow pacemaking neurons mostly release modulatory neurotransmitters, such as dopamine or noradrenalin. The lowfrequency (1–4 Hz) metronomic activity is important for maintaining constant background levels of the neuromodulators. By contrast, pauses or burst-like frequency increases (above 20 Hz) cause sudden pulses in the transmitter release. This occurs in response to external stimuli. In the pacemaking neurons, such bursts are found to be very hard to elicit, reflecting the importance of the regularity of their oscillations. In particular, stimulation ⇑ Corresponding author. E-mail address: [email protected] (D.G. Zakharov). The first author was supported by the Ministry of Education and Science (Contracts No. 02.740.11.5188, 14.740.11.0348, P942) and the RFBR (grants No. 09-02-00719, 09-02-91061, 10-02-00643). 2 The second author was supported by NSF grant DMS-0817717 and International Development Grant of Indiana University. We thank Dr. Rubchinsky for useful discussions. 1
0960-0779/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.chaos.2012.01.007
by injecting a current into the cell body (applied somatic depolarization) is expected to elicit bursting similar to stimulation by inputs from other neurons (synaptic currents), but it does not [4,5]. Taken together, there is a subset of several key neuron types that differentially respond to applied depolarization and synaptic stimulations. The aim of this work is to construct a minimal model that reproduces the common properties of the slow pacemaking neurons. No simple model that mimics differential responses to various stimuli is available to our knowledge. The stimuli are represented by addition of corresponding terms to the model (see below). Depending on the stimulus type, the term is constant, linear or nonlinear function of the variables. The classical FHN equations [6] not only respond identically, but also hardly elevate the frequency in response to the stimuli. We reproduced the differential responses to the stimuli in a biophysical model for the dopaminergic (DA) neuron [7]. The difference in responses was based on the nonlinearity of the stimulus term. In particular, the activation of N-methyl-D-aspartate receptors (NMDAR) was mimicked by addition of a nonlinear function of voltage and elevated the frequency in the model. By contrast, addition of a constant term, which mimics tonic applied depolarization, did not elevate the frequency significantly. We reduce the model to focus on the general properties of the slow pacemaking neurons. The properties are (1) periodic oscillations in isolation; (2) more than
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
641
5-fold elevation in frequency during NMDAR stimulation; (3) less than 2-fold frequency elevation during applied depolarization.
This introduces new features into the dynamics, as described below. System (1) includes two external stimuli – an applied current and NMDAR synaptic stimulation [8]:
2. Model
jstim ¼ japp þ g N ðuÞðEN uÞ; g N ðuÞ ¼ g N =ð1 þ 0:1½Mgexp6u Þ;
We study the dimensionless neuron model
u_ ¼ f ðuÞ þ jKCa ðu; v Þ þ jstim ; v_ ¼ gðu; v Þ v ;
ð1Þ
where variable u describes membrane potential, and v corresponds to calcium ion concentration. is a small parameter reflecting the slowness of v relative to u. In contrast to the simplest FHN oscillator [6], the nonlinear term jKCa(u, v) replaces a linear function of v in system (1). Importantly, the function has a steep sigmoidal dependence on v, and we borrow its particular form from our model [7], where it represents an SK-type Ca2+-dependent potassium current: 4
jKCa ðu; v Þ ¼ g KCa ðEk uÞv =ðv þ k Þ: 4
4
ð2Þ +
Here, gKCa is its maximal current density, EK is the K reversal potential, k is the half activation of the current. Functions f(u) and g(u, v) are set as
f ðuÞ ¼ a1 ðu3 þ a2 u2 þ a3 u þ a4 Þ; u c; u P c; gðu; v Þ ¼ 0:01ðu cÞ v ; u < c:
ð3Þ
To explain our results, we use nullcline analysis. The nullclines are two curves on which either variable equilibrates (u_ ¼ 0 or v_ ¼ 0, respectively). As in FHN, the u-nullcline is folded to qualitatively mimic excitability of the neuron. The nonlinearity introduced by the SK current (2) does not unfold the u-nullcline, but changes its shape. First, another folded branch of the nullcline (symmetric about the u axis) emerges below the u axis (Fig. 1). Intersections with the v-nullcline in that region may become stable equilibrium states and attract trajectories, which is physiologically implausible. To avoid that, as well as the trajectories crossing into the region of negative v we changed the v nullcline (v = g(u, v)) from linear to piece-wise linear (3). Second, a minimum or maximum of the u nullcline may be stretched vertically by the transformation.
v =0
1.5 1
u =0
ð4Þ
where japp represents the applied current, gN is the maximal conductance of NMDAR, EN is its reversal potential. The conductance of NMDAR, gN(u), was borrowed from [9,7]. We fixed all parameters (a1 = 1, a2 = 1.35, a3 = 0.54, a4 = 0.0539, c = 0.585, [Mg] = 2, E N = 0, g KCa = 0.5, E KCa = 1, k = 10) except for japp and gN which vary as the control parameters of the system (1). For this parameter set, the model (1) without an external stimulus demonstrates robust low-frequency tonic (periodic) oscillations. 3. Responses of the model to stimulations The v equation in (1) depends on the morphology of the neuron. More precisely, parameter includes the ratio of the volume to the surface area of a neuron section. This is a common property of any equation written for an ion concentration as oppose to a gating variable or the voltage. The ratio is very different for the soma and a fine distal dendrite. As a result, a neuron is represented as a chain of electrically coupled compartments with different values of parameter . In this section, we do not consider coupled compartments, but examine the dynamics for two values of the small parameter ( = 0.01 and = 0.1) to understand the difference between a thin and a thick segment of the neuron. Because the small parameter closely controls the period of oscillations, the thin segment corresponds to a fast compartment, and the thick segment corresponds to a slow compartment. 3.1. Applied depolarization blocks oscillations In simulations, applied depolarization sharply reduces the amplitude of oscillations (solid) in the compartments with both values of (Fig. 2a). This occurs in the range of the applied current much lower compared to the amplitude of any other current in the model. Simultaneously, the frequency increases, but most of the frequency increase occurs at small amplitudes, especially in the slow compartment. Thus, applied depolarization causes very pronounced blockade of oscillations without a significant frequency elevation.
0.5 3.2. NMDAR stimulation evokes robust high-frequency firing
v 0 -0.5 -1 -1.5 -0.8
-0.6
u
-0.4
Fig. 1. The nullclines of model (1).
-0.2
In contrast to the previous type of stimulation, the NMDAR stimulation elicits robust firing in a wide interval of its conductance density (Fig. 2b). The amplitude of oscillations slightly grows with increasing NMDAR conductance. The growth is stronger in the fast compartment, and its amplitude becomes as high as that in the slow one. The frequencies grow first as the NMDAR current is introduced. When the current density exceeds a value
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
10
=0.01
0.4
=0.1
6
0.2
=0.1
2
=0.01 0
0.002
0
japp
0.004
Amplitude
=0.01
0.4
6
0.2 0
=0.1 0
0.5
1
gN
1.5
2
2
10
=0.1 0.4
6
0.2
2
=0.01 0
0.5
2.5
=0.01 and
0.6
1
japp
1.5
2
2.5
d
=0.1
14 10
=0.1
0.4
6
0.2 0
14
2
=0.01 0
0.5
1
gN
1.5
2
Frequency, 10 -3
10
=0.1
c
=0.1
0.8
Frequency, 10 -3
0.6
14
=0.01
=0.01 and
0.6
0
0.006
0.8
b
0.8
Amplitude
0.6
14
Frequency, 10 -3
a
Frequency, 10 -3
Amplitude
0.8
Amplitude
642
2.5
Fig. 2. The dependences of the amplitude (solid) and the frequency (dashed) on different stimuli in model (1) (a, b) and in the FHN equations (c, d).
around 1, its further elevation causes decrease in the frequencies. In both compartments, the frequency grows more than 7-fold. There is a wide interval of NMDAR conductance where the frequency is significantly elevated. Therefore, there is a range of optimal NMDAR conductance where robust high-frequency firing is achieved in the model. 3.3. The difference with FHN oscillator Here, we compare the dynamics of model (1) and the classical FHN oscillator, in which the dependence of the SK current (2) is replaced by a linear function (a1 = 50.5, a2 = 1.35, a3 = 0.54, a4 = 0.0472). Applied depolarization introduces no changes either to the amplitude or to the frequency of the compartments (Fig. 2d). It does not lead to blockade of oscillations either. NMDAR stimulation does not block the oscillations (Fig. 2e), similar to model (1). However, it does not elevate the frequency. Both the frequency and amplitude change little, and the frequency only decreases, in contrast to model (1). Therefore, the responses of the classical FHN oscillator to all three types of stimulation are qualitatively different from those in our model (1). In particular, the responses to stimulations are not significantly different from one another. 3.4. Influence of the SK current The difference between the models can be explained by the nullcline analysis. Oscillations circumscribe folding of u-nullcline as shown in Fig. 3. Without external stimulus, both our model (1) and the FHN oscillator displays robust oscillations. In FHN oscillator, applied depolarization only shifts the u-nullcline up without changing its shape (v = f(u) + japp, Fig. 3c). This does not change dynamics
because the v-nullcline is vertical (u = c) and the vertical shift of the u-nullcline does not change their relative positions (Fig. 3c). The amplitude and frequency of oscillations remain constant for all values of the applied current (Fig. 3d). Only the v range spanned by the oscillation shifts together with the folded region of the u-nullcline. In model (1), the u-nullcline not only shifts up, but also unfolds with growing applied depolarization (Fig. 3a). Its minimum shifts to the right. As a result, the intersection between the nullclines shifts across the minimum and becomes a stable equilibrium. It interrupts the oscillatory cycle, and the oscillation collapse. This occurs at very low applied depolarization and causes the blockade of oscillations. This analysis confirms the similarity of model (1) to the conductance-base model of the DA neuron [7]. The NMDAR current also shifts the u-nullcline up in the FHN oscillator (Fig. 3d). However, it folds the nullcline even stronger and, most importantly, shifts its extrema away from the intersection with the other nullcline. As a result, NMDAR stimulation does not block oscillations in the FHN oscillator. Although, it does not elevate the frequency either. The comparison with the transformation of the nullcline in model (1) (Fig. 3b) explains this result – in the latter, the nullcline flattens when NMDAR current is introduced (Fig. 3b). This transformation is slightly different from unfolding discussed above - the extrema become closer in vertical direction (in Ca2+ concentration), but remain well separated in the horizontal (in the voltage). The separation of extrema ensures robustness of oscillations. The smaller difference in Ca2+ concentration results in a frequency elevation because Ca2+ concentration is a slow variable and its amplitude mostly determines the period. Thus, the NMDAR-induced frequency elevation is based on shaping the voltage nullcline by the SK current in the model.
643
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
2.5
2.5
2
japp =0.03
1.5
v
1
a
v =0
japp =0.01
v japp =0
-0.8
-0.6
u
-0.4
v
-0.2
u =0
1
0
japp =0 -0.8
2.5
b
gN=0.5
1.5 1
japp =1
0.5
2.5 2
v
u
-0.4
-0.2
v =0
0
d
gN=5
1.5 1
gN=0
0.5 -0.8
-0.6
2
u =0
v =0
u =0
c
1.5
u =0
0.5
v =0
2
-0.6
-0.4
u
gN=0
0.5 -0.2
0
-0.8
-0.6
u
-0.4
u =0 -0.2
0
Fig. 3. The nullclines of model (1) (A, B) and the FHN equations (c, d) for different external stimuli – (a, c) applied current, (b, d) NMDAR. The thick curves are trajectories. The thin curves are nullclines.
4. Responses of the DA neuron As mentioned above, a neuron is represented as a chain of electrically coupled compartments with different parameters . The complexity of such coupled oscillator model and regimes generated by it is very high. In particular, the mismatch in the parameter excludes the possibility of exact synchronization. On the other hand, the oscillations in the whole neuron must stay very close to periodic. Following [7], we investigate the simplest coupled oscillator model – a pair of compartments. The compartments stand for the soma (the subscript s) and a distal dendrite (d):
u_ s;d ¼ f ðus;d Þ þ jKCa þ jstim þ Dðus;d ud;s Þ; v_ s;d ¼ s;d gðus;d ; v s;d Þ v s;d :
ð5Þ
The small parameters take the values introduced before – s = 0.01 and d = 0.1. The only new term in the equations is the coupling. We choose the coupling strength to achieve almost synchronous oscillations (D = 0.5). The model displays periodic oscillations in this regime, and the compartments oscillate with the same frequency (frequency synchronization). The coupled oscillator model (5) replicated the distinct responses to the stimuli observed in the single compartments. Fig. 4a shows the voltage time-series in response to a pulse of applied depolarization localized in the soma. The activation is assumed instantaneous. Depolarizing applied current is elevated to a value much lower than any other current in the model. The oscillations quickly fail during the stimulation, and a high frequency is not displayed even in the transient at the onset. After terminating the stimulation, low-frequency activity of the neuron
resumes. Fig. 4b shows the voltage time-series for the case of dendritic NMDAR stimulation. The model switches to the high frequency with the onset of the stimulation and quickly returns to the low frequency after the stimulation ends. These result reproduce those for the biophysical model of the DA neuron [7] and match experimental data [3]. For instance, one can compare them with experimental voltage time-series of a DA neuron in response to a pulse of applied depolarization (Fig. 5B from [5]) and NMDAR synaptic current (Fig. 7A from [10]). 5. Discussion In this paper, we have developed a minimal phenomenological model that reproduces responses of slow pacemaking neurons to various stimuli. During applied depolarization the neurons only slightly elevate the frequency because of an early onset of the depolarization block. During NMDAR stimulation, the neurons robustly fire and elevate the frequency to much greater values. The model is dimensionless, and there is no quantitative correspondence with the values of the frequency or conductance measured in experiments. Thus, we target the relative frequency elevation observed in experiments. During NMDAR stimulation, the model shows 7-fold increase in the frequency. Therefore, the model qualitatively mimics the wide range of frequencies and the specificity of responses to various stimuli. We have separated the most important property of the neuron that allows for the differential responses. In the model, the nonlinearity of a current that depends on an ion concentration is crucial for such differentiation. We have explicitly shown that when the nonlinear activation function is replaced by a linear dependence, the model
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
0 -0.2 -0.4 -0.6 -0.8 0
a us
us
644
0.04 1000
2000
japp 3000
4000
0 -0.2 -0.4 -0.6 -0.8 0
time
b 1
1000
2000
gNd 3000
4000
time
Fig. 4. The voltage time-series are shown in response to square pulses of the applied depolarization (a) and NMDAR activation (b).
responds similarly to all stimuli and shows very little frequency variation. These are the properties of the classical FHN oscillator. Along with introducing the nonlinear current, we also changed the linear v nullcline of the FHN oscillator to a piece-wise linear (3). This was done to avoid additional intersections between the nullclines that are physiologically implausible. Any other choice of the v nullcine that intersects the u nullcline only once as shown in Fig. 1 will qualitatively preserve the frequency responses. On the other hand, changing this nullcline alone in FHN oscillator does not result in differentiating responses to stimuli. This modification allows for a higher frequency growth in the oscillator, but the responses to the different types of stimuli are very similar to each other [11]. In the DA neuron, an SK-type Ca2+-dependent K+ current provides the necessary nonlinearity. Presumably, the same current works in serotonergic neurons [15]. No data is available for other neurons, and the SK current may not be the only current that differentiates the responses. The mechanism works for other currents with a sigmoidal activation function. In fact, any function that starts flat and then sharply increases its slope works. The saturation part of the sigmoidal dependence is not necessary for the frequency responses. The current may depend on ion concentration, on the voltage, or on both variables. This further expands the applicability of our results to neurons expressing various currents. Therefore, in a wide class of neurons, the sigmoidal nonlinearity of a conductance will cause distinct responses to stimuli. Currents that depend on slowly-changing activation variables also are thought to play a central role in shaping the low-frequency activity of slow pacemaking neurons. For instance, the SK current provides a low frequency robust periodic oscillations by keeping the neuron hyperpolarized in response to the elevation in Ca2+ concentration after each spike (see e.g. [16]). But how the long-lasting hyperpolarization is overcome during high-frequency firing? We have shown that even low-amplitude deviations in the slow activation variable are enough to support robust oscillations. Decreasing amplitude of oscillations in the activation variable strongly increases the frequency because the variable changes slowly. Taken together, a current activated by a slow variable supports robust oscillations in a wide range of frequencies. In the DA neuron, the slow activation variable is the Ca2+ concentration, and the current is the SK-type potassium current. The regularization of firing is an established role of the SK current [12], although it may support chaotic dynamics when the handling of intracellular Ca2+ is complex [13,14]. I any case,the
SK current is thought of as slowing the firing or providing the low frequency modulation. By contrast, our result specified for the DA neuron says that the SK current supports oscillations in a wide range of frequencies. This is a novel role of the SK-type Ca2+-dependent K+ current, and we suggest that this is a common feature of currents that depend on a slowly-changing variable. References [1] Hajos M, Gartside SE, Villa AE, Sharp T. Evidence for a repetitive (burst) firing pattern in a sub-population of 5-hydroxytryptamine neurons in the dorsal and median raphe nuclei of the rat. Neuroscience 1995;69:189–97. [2] Muntoni AL, Pillolla G, Melis M, Perra S, Gessa GL, Pistis M. Cannabinoids modulate spontaneous neuronal activity and evoked inhibition of locus coeruleus noradrenergic neurons. Eur J Neurosci 2006;23:2385–94. [3] Wilson CJ, Callaway JC. A coupled oscillator model of the dopaminergic neuron of the substantia nigra. J Neurophysiol 2000;83:3084–100. [4] Hajos M, Sharp T, Newberry NR. Intracellular recordings from burstfiring presumed serotonergic neurones in the rat dorsal raphe nucleus in vivo. Brain Res 1996;737:308–12. [5] Richards CD, Shiroyama T, Kitai ST. Electrophysiological and immunocytochemical characteristics of GABA and dopamine neurons in the substantia nigra of the rat. Neuroscience 1997;80:545–57. [6] FitzHugh R. Impulses and physiological states in theoretical models of nerve membrane. Biophys J 1961;1:445–66. [7] Kuznetsov AS, Kopell NJ, Wilson CJ. Transient high-frequency firing in a coupled – oscillator model of the mesencephalic dopaminergic neuron. J Neurophysiol 2006;95:932–47. [8] Destexhe A, Mainen ZF, Sejnowsky TJ. Kinetic models of synaptic transmission. In: Koch C, Segev I, editors. Methods in neuronal modeling: from ions to networks. Cambridge, MA: MIT Press; 1998. p. 1–27. [9] Li Y-X, Bertram R, Rinzel J. Modeling N-methyl-D-aspartate-induced bursting in dopamine neurons. Neuroscience 1996;71:397410. [10] Deister CA, Teagarden MA, Wilson CJ, Paladini CA. An intrinsic neuronal oscillator underlies dopaminergic neuron bursting. J Neurosci 2009;29:1588815897. [11] Zakharov DG, Kuznetsov AS, Nekorkin VI. A two-Compartment phenomenological model of a dopaminergic neuron. Biophysics 2010;55:241247. [12] Surmeier DJ, Mercer JN, Chan CS. Autonomous pacemakers in the basal ganglia: who needs excitatory synapses anyway? Curr Opin Neurobiol 2005;15:312–8. [13] Falcke M, Huerta R, Rabinovich MI, Abarbanel HDI, Elson RC, Selverston AI. Modeling observed chaotic oscillations in bursting neurons: the role of calcium dynamics and IP3. Biol Cybern 2000;82:517–27. [14] Varona P, Torres JJ, Huerta R, Abarbanel HDI, Rabinovich MI. Regularization mechanisms of spiking-bursting neurons. Neural Netw 2001;14:865–75. [15] Rouchet N, Waroux O, Lamy C, Massotte L, Scuvee-Moreau J, Leigeois J-F, et al. SK channel blockade promotes burst firing in dorsal raphe serotonergic neurons. Eur J Neurosci 2008;28:1108–15. [16] Wolfart J, Neuhoff H, Franz O, Roeper J. Differential expression of the small-conductance, calcium-activated potassium channel SK3 is critical for pacemaker control in dopaminergic midbrain neurons. J Neurosci 2001;21:344356.
Contents lists available at SciVerse ScienceDirect
Chaos, Solitons & Fractals Nonlinear Science, and Nonequilibrium and Complex Phenomena journal homepage: www.elsevier.com/locate/chaos
A minimal model for a slow pacemaking neuron D.G. Zakharov a,⇑,1, A. Kuznetsov b,2 a b
Nonlinear dynamics department, Institute of Applied Physics, RAS, 46 Ulyanov St., Nizhny Novgorod 603950, Russia Department of Mathematical Sciences and Center for Mathematical Biosciences, IUPUI, 402 N. Blackford St., Indianapolis, IN 46202, USA
a r t i c l e
i n f o
Article history: Available online 24 February 2012
a b s t r a c t We have constructed a phenomenological model for slow pacemaking neurons. These are neurons that generate very regular periodic oscillations of the membrane potential. Many of these neurons also differentially respond to various types of stimulation. The model is based on FitzHugh–Nagumo (FHN) oscillator and implements a nonlinearity introduced by a current that depends on an ion concentration. The comparison with the original FHN oscillator has shown that the new nonlinear dependence allows for differentiating responses to various stimuli. We discuss implications of our results for a broad class of neurons. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction There are several key groups of neurons in the central nervous system that generate very regular periodic oscillations of the membrane potential. The examples of these neurons are serotonin-containing neurons from the raphe nuclei [1], noradrenergic neurons located in the pontine nucleus locus coeruleus (LC) [2] and dopaminergic neurons from the substantia nigra pars compacta [3]. These slow pacemaking neurons mostly release modulatory neurotransmitters, such as dopamine or noradrenalin. The lowfrequency (1–4 Hz) metronomic activity is important for maintaining constant background levels of the neuromodulators. By contrast, pauses or burst-like frequency increases (above 20 Hz) cause sudden pulses in the transmitter release. This occurs in response to external stimuli. In the pacemaking neurons, such bursts are found to be very hard to elicit, reflecting the importance of the regularity of their oscillations. In particular, stimulation ⇑ Corresponding author. E-mail address: [email protected] (D.G. Zakharov). The first author was supported by the Ministry of Education and Science (Contracts No. 02.740.11.5188, 14.740.11.0348, P942) and the RFBR (grants No. 09-02-00719, 09-02-91061, 10-02-00643). 2 The second author was supported by NSF grant DMS-0817717 and International Development Grant of Indiana University. We thank Dr. Rubchinsky for useful discussions. 1
0960-0779/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.chaos.2012.01.007
by injecting a current into the cell body (applied somatic depolarization) is expected to elicit bursting similar to stimulation by inputs from other neurons (synaptic currents), but it does not [4,5]. Taken together, there is a subset of several key neuron types that differentially respond to applied depolarization and synaptic stimulations. The aim of this work is to construct a minimal model that reproduces the common properties of the slow pacemaking neurons. No simple model that mimics differential responses to various stimuli is available to our knowledge. The stimuli are represented by addition of corresponding terms to the model (see below). Depending on the stimulus type, the term is constant, linear or nonlinear function of the variables. The classical FHN equations [6] not only respond identically, but also hardly elevate the frequency in response to the stimuli. We reproduced the differential responses to the stimuli in a biophysical model for the dopaminergic (DA) neuron [7]. The difference in responses was based on the nonlinearity of the stimulus term. In particular, the activation of N-methyl-D-aspartate receptors (NMDAR) was mimicked by addition of a nonlinear function of voltage and elevated the frequency in the model. By contrast, addition of a constant term, which mimics tonic applied depolarization, did not elevate the frequency significantly. We reduce the model to focus on the general properties of the slow pacemaking neurons. The properties are (1) periodic oscillations in isolation; (2) more than
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
641
5-fold elevation in frequency during NMDAR stimulation; (3) less than 2-fold frequency elevation during applied depolarization.
This introduces new features into the dynamics, as described below. System (1) includes two external stimuli – an applied current and NMDAR synaptic stimulation [8]:
2. Model
jstim ¼ japp þ g N ðuÞðEN uÞ; g N ðuÞ ¼ g N =ð1 þ 0:1½Mgexp6u Þ;
We study the dimensionless neuron model
u_ ¼ f ðuÞ þ jKCa ðu; v Þ þ jstim ; v_ ¼ gðu; v Þ v ;
ð1Þ
where variable u describes membrane potential, and v corresponds to calcium ion concentration. is a small parameter reflecting the slowness of v relative to u. In contrast to the simplest FHN oscillator [6], the nonlinear term jKCa(u, v) replaces a linear function of v in system (1). Importantly, the function has a steep sigmoidal dependence on v, and we borrow its particular form from our model [7], where it represents an SK-type Ca2+-dependent potassium current: 4
jKCa ðu; v Þ ¼ g KCa ðEk uÞv =ðv þ k Þ: 4
4
ð2Þ +
Here, gKCa is its maximal current density, EK is the K reversal potential, k is the half activation of the current. Functions f(u) and g(u, v) are set as
f ðuÞ ¼ a1 ðu3 þ a2 u2 þ a3 u þ a4 Þ; u c; u P c; gðu; v Þ ¼ 0:01ðu cÞ v ; u < c:
ð3Þ
To explain our results, we use nullcline analysis. The nullclines are two curves on which either variable equilibrates (u_ ¼ 0 or v_ ¼ 0, respectively). As in FHN, the u-nullcline is folded to qualitatively mimic excitability of the neuron. The nonlinearity introduced by the SK current (2) does not unfold the u-nullcline, but changes its shape. First, another folded branch of the nullcline (symmetric about the u axis) emerges below the u axis (Fig. 1). Intersections with the v-nullcline in that region may become stable equilibrium states and attract trajectories, which is physiologically implausible. To avoid that, as well as the trajectories crossing into the region of negative v we changed the v nullcline (v = g(u, v)) from linear to piece-wise linear (3). Second, a minimum or maximum of the u nullcline may be stretched vertically by the transformation.
v =0
1.5 1
u =0
ð4Þ
where japp represents the applied current, gN is the maximal conductance of NMDAR, EN is its reversal potential. The conductance of NMDAR, gN(u), was borrowed from [9,7]. We fixed all parameters (a1 = 1, a2 = 1.35, a3 = 0.54, a4 = 0.0539, c = 0.585, [Mg] = 2, E N = 0, g KCa = 0.5, E KCa = 1, k = 10) except for japp and gN which vary as the control parameters of the system (1). For this parameter set, the model (1) without an external stimulus demonstrates robust low-frequency tonic (periodic) oscillations. 3. Responses of the model to stimulations The v equation in (1) depends on the morphology of the neuron. More precisely, parameter includes the ratio of the volume to the surface area of a neuron section. This is a common property of any equation written for an ion concentration as oppose to a gating variable or the voltage. The ratio is very different for the soma and a fine distal dendrite. As a result, a neuron is represented as a chain of electrically coupled compartments with different values of parameter . In this section, we do not consider coupled compartments, but examine the dynamics for two values of the small parameter ( = 0.01 and = 0.1) to understand the difference between a thin and a thick segment of the neuron. Because the small parameter closely controls the period of oscillations, the thin segment corresponds to a fast compartment, and the thick segment corresponds to a slow compartment. 3.1. Applied depolarization blocks oscillations In simulations, applied depolarization sharply reduces the amplitude of oscillations (solid) in the compartments with both values of (Fig. 2a). This occurs in the range of the applied current much lower compared to the amplitude of any other current in the model. Simultaneously, the frequency increases, but most of the frequency increase occurs at small amplitudes, especially in the slow compartment. Thus, applied depolarization causes very pronounced blockade of oscillations without a significant frequency elevation.
0.5 3.2. NMDAR stimulation evokes robust high-frequency firing
v 0 -0.5 -1 -1.5 -0.8
-0.6
u
-0.4
Fig. 1. The nullclines of model (1).
-0.2
In contrast to the previous type of stimulation, the NMDAR stimulation elicits robust firing in a wide interval of its conductance density (Fig. 2b). The amplitude of oscillations slightly grows with increasing NMDAR conductance. The growth is stronger in the fast compartment, and its amplitude becomes as high as that in the slow one. The frequencies grow first as the NMDAR current is introduced. When the current density exceeds a value
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
10
=0.01
0.4
=0.1
6
0.2
=0.1
2
=0.01 0
0.002
0
japp
0.004
Amplitude
=0.01
0.4
6
0.2 0
=0.1 0
0.5
1
gN
1.5
2
2
10
=0.1 0.4
6
0.2
2
=0.01 0
0.5
2.5
=0.01 and
0.6
1
japp
1.5
2
2.5
d
=0.1
14 10
=0.1
0.4
6
0.2 0
14
2
=0.01 0
0.5
1
gN
1.5
2
Frequency, 10 -3
10
=0.1
c
=0.1
0.8
Frequency, 10 -3
0.6
14
=0.01
=0.01 and
0.6
0
0.006
0.8
b
0.8
Amplitude
0.6
14
Frequency, 10 -3
a
Frequency, 10 -3
Amplitude
0.8
Amplitude
642
2.5
Fig. 2. The dependences of the amplitude (solid) and the frequency (dashed) on different stimuli in model (1) (a, b) and in the FHN equations (c, d).
around 1, its further elevation causes decrease in the frequencies. In both compartments, the frequency grows more than 7-fold. There is a wide interval of NMDAR conductance where the frequency is significantly elevated. Therefore, there is a range of optimal NMDAR conductance where robust high-frequency firing is achieved in the model. 3.3. The difference with FHN oscillator Here, we compare the dynamics of model (1) and the classical FHN oscillator, in which the dependence of the SK current (2) is replaced by a linear function (a1 = 50.5, a2 = 1.35, a3 = 0.54, a4 = 0.0472). Applied depolarization introduces no changes either to the amplitude or to the frequency of the compartments (Fig. 2d). It does not lead to blockade of oscillations either. NMDAR stimulation does not block the oscillations (Fig. 2e), similar to model (1). However, it does not elevate the frequency. Both the frequency and amplitude change little, and the frequency only decreases, in contrast to model (1). Therefore, the responses of the classical FHN oscillator to all three types of stimulation are qualitatively different from those in our model (1). In particular, the responses to stimulations are not significantly different from one another. 3.4. Influence of the SK current The difference between the models can be explained by the nullcline analysis. Oscillations circumscribe folding of u-nullcline as shown in Fig. 3. Without external stimulus, both our model (1) and the FHN oscillator displays robust oscillations. In FHN oscillator, applied depolarization only shifts the u-nullcline up without changing its shape (v = f(u) + japp, Fig. 3c). This does not change dynamics
because the v-nullcline is vertical (u = c) and the vertical shift of the u-nullcline does not change their relative positions (Fig. 3c). The amplitude and frequency of oscillations remain constant for all values of the applied current (Fig. 3d). Only the v range spanned by the oscillation shifts together with the folded region of the u-nullcline. In model (1), the u-nullcline not only shifts up, but also unfolds with growing applied depolarization (Fig. 3a). Its minimum shifts to the right. As a result, the intersection between the nullclines shifts across the minimum and becomes a stable equilibrium. It interrupts the oscillatory cycle, and the oscillation collapse. This occurs at very low applied depolarization and causes the blockade of oscillations. This analysis confirms the similarity of model (1) to the conductance-base model of the DA neuron [7]. The NMDAR current also shifts the u-nullcline up in the FHN oscillator (Fig. 3d). However, it folds the nullcline even stronger and, most importantly, shifts its extrema away from the intersection with the other nullcline. As a result, NMDAR stimulation does not block oscillations in the FHN oscillator. Although, it does not elevate the frequency either. The comparison with the transformation of the nullcline in model (1) (Fig. 3b) explains this result – in the latter, the nullcline flattens when NMDAR current is introduced (Fig. 3b). This transformation is slightly different from unfolding discussed above - the extrema become closer in vertical direction (in Ca2+ concentration), but remain well separated in the horizontal (in the voltage). The separation of extrema ensures robustness of oscillations. The smaller difference in Ca2+ concentration results in a frequency elevation because Ca2+ concentration is a slow variable and its amplitude mostly determines the period. Thus, the NMDAR-induced frequency elevation is based on shaping the voltage nullcline by the SK current in the model.
643
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
2.5
2.5
2
japp =0.03
1.5
v
1
a
v =0
japp =0.01
v japp =0
-0.8
-0.6
u
-0.4
v
-0.2
u =0
1
0
japp =0 -0.8
2.5
b
gN=0.5
1.5 1
japp =1
0.5
2.5 2
v
u
-0.4
-0.2
v =0
0
d
gN=5
1.5 1
gN=0
0.5 -0.8
-0.6
2
u =0
v =0
u =0
c
1.5
u =0
0.5
v =0
2
-0.6
-0.4
u
gN=0
0.5 -0.2
0
-0.8
-0.6
u
-0.4
u =0 -0.2
0
Fig. 3. The nullclines of model (1) (A, B) and the FHN equations (c, d) for different external stimuli – (a, c) applied current, (b, d) NMDAR. The thick curves are trajectories. The thin curves are nullclines.
4. Responses of the DA neuron As mentioned above, a neuron is represented as a chain of electrically coupled compartments with different parameters . The complexity of such coupled oscillator model and regimes generated by it is very high. In particular, the mismatch in the parameter excludes the possibility of exact synchronization. On the other hand, the oscillations in the whole neuron must stay very close to periodic. Following [7], we investigate the simplest coupled oscillator model – a pair of compartments. The compartments stand for the soma (the subscript s) and a distal dendrite (d):
u_ s;d ¼ f ðus;d Þ þ jKCa þ jstim þ Dðus;d ud;s Þ; v_ s;d ¼ s;d gðus;d ; v s;d Þ v s;d :
ð5Þ
The small parameters take the values introduced before – s = 0.01 and d = 0.1. The only new term in the equations is the coupling. We choose the coupling strength to achieve almost synchronous oscillations (D = 0.5). The model displays periodic oscillations in this regime, and the compartments oscillate with the same frequency (frequency synchronization). The coupled oscillator model (5) replicated the distinct responses to the stimuli observed in the single compartments. Fig. 4a shows the voltage time-series in response to a pulse of applied depolarization localized in the soma. The activation is assumed instantaneous. Depolarizing applied current is elevated to a value much lower than any other current in the model. The oscillations quickly fail during the stimulation, and a high frequency is not displayed even in the transient at the onset. After terminating the stimulation, low-frequency activity of the neuron
resumes. Fig. 4b shows the voltage time-series for the case of dendritic NMDAR stimulation. The model switches to the high frequency with the onset of the stimulation and quickly returns to the low frequency after the stimulation ends. These result reproduce those for the biophysical model of the DA neuron [7] and match experimental data [3]. For instance, one can compare them with experimental voltage time-series of a DA neuron in response to a pulse of applied depolarization (Fig. 5B from [5]) and NMDAR synaptic current (Fig. 7A from [10]). 5. Discussion In this paper, we have developed a minimal phenomenological model that reproduces responses of slow pacemaking neurons to various stimuli. During applied depolarization the neurons only slightly elevate the frequency because of an early onset of the depolarization block. During NMDAR stimulation, the neurons robustly fire and elevate the frequency to much greater values. The model is dimensionless, and there is no quantitative correspondence with the values of the frequency or conductance measured in experiments. Thus, we target the relative frequency elevation observed in experiments. During NMDAR stimulation, the model shows 7-fold increase in the frequency. Therefore, the model qualitatively mimics the wide range of frequencies and the specificity of responses to various stimuli. We have separated the most important property of the neuron that allows for the differential responses. In the model, the nonlinearity of a current that depends on an ion concentration is crucial for such differentiation. We have explicitly shown that when the nonlinear activation function is replaced by a linear dependence, the model
D.G. Zakharov, A. Kuznetsov / Chaos, Solitons & Fractals 45 (2012) 640–644
0 -0.2 -0.4 -0.6 -0.8 0
a us
us
644
0.04 1000
2000
japp 3000
4000
0 -0.2 -0.4 -0.6 -0.8 0
time
b 1
1000
2000
gNd 3000
4000
time
Fig. 4. The voltage time-series are shown in response to square pulses of the applied depolarization (a) and NMDAR activation (b).
responds similarly to all stimuli and shows very little frequency variation. These are the properties of the classical FHN oscillator. Along with introducing the nonlinear current, we also changed the linear v nullcline of the FHN oscillator to a piece-wise linear (3). This was done to avoid additional intersections between the nullclines that are physiologically implausible. Any other choice of the v nullcine that intersects the u nullcline only once as shown in Fig. 1 will qualitatively preserve the frequency responses. On the other hand, changing this nullcline alone in FHN oscillator does not result in differentiating responses to stimuli. This modification allows for a higher frequency growth in the oscillator, but the responses to the different types of stimuli are very similar to each other [11]. In the DA neuron, an SK-type Ca2+-dependent K+ current provides the necessary nonlinearity. Presumably, the same current works in serotonergic neurons [15]. No data is available for other neurons, and the SK current may not be the only current that differentiates the responses. The mechanism works for other currents with a sigmoidal activation function. In fact, any function that starts flat and then sharply increases its slope works. The saturation part of the sigmoidal dependence is not necessary for the frequency responses. The current may depend on ion concentration, on the voltage, or on both variables. This further expands the applicability of our results to neurons expressing various currents. Therefore, in a wide class of neurons, the sigmoidal nonlinearity of a conductance will cause distinct responses to stimuli. Currents that depend on slowly-changing activation variables also are thought to play a central role in shaping the low-frequency activity of slow pacemaking neurons. For instance, the SK current provides a low frequency robust periodic oscillations by keeping the neuron hyperpolarized in response to the elevation in Ca2+ concentration after each spike (see e.g. [16]). But how the long-lasting hyperpolarization is overcome during high-frequency firing? We have shown that even low-amplitude deviations in the slow activation variable are enough to support robust oscillations. Decreasing amplitude of oscillations in the activation variable strongly increases the frequency because the variable changes slowly. Taken together, a current activated by a slow variable supports robust oscillations in a wide range of frequencies. In the DA neuron, the slow activation variable is the Ca2+ concentration, and the current is the SK-type potassium current. The regularization of firing is an established role of the SK current [12], although it may support chaotic dynamics when the handling of intracellular Ca2+ is complex [13,14]. I any case,the
SK current is thought of as slowing the firing or providing the low frequency modulation. By contrast, our result specified for the DA neuron says that the SK current supports oscillations in a wide range of frequencies. This is a novel role of the SK-type Ca2+-dependent K+ current, and we suggest that this is a common feature of currents that depend on a slowly-changing variable. References [1] Hajos M, Gartside SE, Villa AE, Sharp T. Evidence for a repetitive (burst) firing pattern in a sub-population of 5-hydroxytryptamine neurons in the dorsal and median raphe nuclei of the rat. Neuroscience 1995;69:189–97. [2] Muntoni AL, Pillolla G, Melis M, Perra S, Gessa GL, Pistis M. Cannabinoids modulate spontaneous neuronal activity and evoked inhibition of locus coeruleus noradrenergic neurons. Eur J Neurosci 2006;23:2385–94. [3] Wilson CJ, Callaway JC. A coupled oscillator model of the dopaminergic neuron of the substantia nigra. J Neurophysiol 2000;83:3084–100. [4] Hajos M, Sharp T, Newberry NR. Intracellular recordings from burstfiring presumed serotonergic neurones in the rat dorsal raphe nucleus in vivo. Brain Res 1996;737:308–12. [5] Richards CD, Shiroyama T, Kitai ST. Electrophysiological and immunocytochemical characteristics of GABA and dopamine neurons in the substantia nigra of the rat. Neuroscience 1997;80:545–57. [6] FitzHugh R. Impulses and physiological states in theoretical models of nerve membrane. Biophys J 1961;1:445–66. [7] Kuznetsov AS, Kopell NJ, Wilson CJ. Transient high-frequency firing in a coupled – oscillator model of the mesencephalic dopaminergic neuron. J Neurophysiol 2006;95:932–47. [8] Destexhe A, Mainen ZF, Sejnowsky TJ. Kinetic models of synaptic transmission. In: Koch C, Segev I, editors. Methods in neuronal modeling: from ions to networks. Cambridge, MA: MIT Press; 1998. p. 1–27. [9] Li Y-X, Bertram R, Rinzel J. Modeling N-methyl-D-aspartate-induced bursting in dopamine neurons. Neuroscience 1996;71:397410. [10] Deister CA, Teagarden MA, Wilson CJ, Paladini CA. An intrinsic neuronal oscillator underlies dopaminergic neuron bursting. J Neurosci 2009;29:1588815897. [11] Zakharov DG, Kuznetsov AS, Nekorkin VI. A two-Compartment phenomenological model of a dopaminergic neuron. Biophysics 2010;55:241247. [12] Surmeier DJ, Mercer JN, Chan CS. Autonomous pacemakers in the basal ganglia: who needs excitatory synapses anyway? Curr Opin Neurobiol 2005;15:312–8. [13] Falcke M, Huerta R, Rabinovich MI, Abarbanel HDI, Elson RC, Selverston AI. Modeling observed chaotic oscillations in bursting neurons: the role of calcium dynamics and IP3. Biol Cybern 2000;82:517–27. [14] Varona P, Torres JJ, Huerta R, Abarbanel HDI, Rabinovich MI. Regularization mechanisms of spiking-bursting neurons. Neural Netw 2001;14:865–75. [15] Rouchet N, Waroux O, Lamy C, Massotte L, Scuvee-Moreau J, Leigeois J-F, et al. SK channel blockade promotes burst firing in dorsal raphe serotonergic neurons. Eur J Neurosci 2008;28:1108–15. [16] Wolfart J, Neuhoff H, Franz O, Roeper J. Differential expression of the small-conductance, calcium-activated potassium channel SK3 is critical for pacemaker control in dopaminergic midbrain neurons. J Neurosci 2001;21:344356.