5.

​ Simplified Neuron Models

In chapter 3, we have introduced the HH model which is the earliest model of action
potential (AP) generation in an axon. It shows how voltage-sensitive dynamics of
sodium and potassium channels results in generation of APs. But real neurons have a
much large variety of ion channels. There are dozens of voltage- and Ca2+-gated
channels known today. Combinations of these channels can give rise to an
astronomically large number of neuron models. Even the HH model, the simplest
model of AP generation that we encountered so far, has 4 differential equations, three
of them being nonlinear. More realistic neuron models with larger number of ion
channels can easily become mathematically untractable. It would be desirable to
construct simplified models, with fewer variables, and milder nonlinearities, in such a
way that the reduced model preserves the essential dynamics of their more complex
versions. One of the first reduced model of that kind is the Fitzhugh-Nagumo neuron
model.

5.1​​ FitzHugh-Nagumo model

FitzHugh-Nagumo model is a two-variable neuron model, constructed by reducing

the 4-variable HH model, by applying suitable assumptions.

Hodgkin- Huxley model:

Assumptions:-
The time scales for m, h and n variables are not all of the same order. These
disparities provide a basis for eliminating some of the gating variables.

1)​ Since the time scale for m is much smaller than that of the other two, we assume
that m relaxes faster than the other two gating variables. Therefore, we let,
 i.e.,

0=

​ ​ ​ ​ ​ ​ ​ ​ (5.1.1)

2)​ h varies too slowly. Therefore, we let h to be a constant, h = h0

The resulting system has only two variables (v, n). After transformation to dimensionless
variables, and some approximations, the resulting FN model may be defined as,

​ ​ ​ ​ ​ ​ ​ (5.1.2)

where​ ​ ​ ​ ​ ​ (5.1.3)
​

​ ​ ​ ​ ​ ​ ​ ​ (5.1.4)

In the above system, v is analogous to the membrane voltage, and w represents all
the three gating variables.
 Figure 5.1.1:​ The nullclines V and w​ ​ ​ ​ ​ ​

​ ​ ​ ​ ​

The two nullclines of the above system are,

​ ​ ​ ​ ​ ​

​ (F-nullcline)
​

​ ​ (g-nullcline)

To examine the stability of a stationary point, we calculate the Jacobian, A, of the

system at that point.

​ ​ ​ ​ ​ ​ ​ ​ (5.1.5)
​ ​ ​ ​ ​ ​ ​ ​ (5.1.6)

The type of the stationary point can be expressed in terms of determinant, Δ, and trace,
τ,of the Jacobian, using the following rules:
-​ ifΔ< 0, the stationary point is a saddle irrespective of the value of τ .
-​ if Δ> 0, τ < 0, stable point
 -​ ifΔ> 0, τ > 0, unstable point.
​ ​ ​ ​ ​

Δ> 0

​ ​
​ ​ ​ ​ ​ ​ ​ (5.1.7)

That is, Δ> 0, when the slope of the F-nullcline is lesser than slope of w-nullcline.

​ ​ ​ ​ ​ ​ ​ (5.1.8)
Now let us consider the behavior of FN model as external current Ia is gradually
increased.

5.1.1)​ Ia=0, Excitability

The phase-plane shown below depicts the situation when Ia=0. There is only one
stationary point at the origin.

​ Figure [Link]:​ Phase plane analysis: Excitability at a=0.5; b=0.1; r=0.1;

Ia=0​
 Stability of origin:

a)​ Slope of F-nullcline> slope of w-nullcline

b)​

​ ​ ​ Origin is stable.

​ ​ Figure [Link]:​ The points a,p1,p2 and their trajectories

Consider the evolution of the variable v (membrane voltage) when the initial condition
is at points p1 or p2 (fig. [Link]).

This behavior will first be anticipated using loose arguments, and then confirmed using
simulations.
The F- and w-nullclines in fig. ([Link],[Link]) above intersect only at one point (the
origin) and therefore divide the plane into four regions, numbered from 1 to 4 (fig.
[Link]). The flow patterns in the four regions can be seen to be as follows:

Region 1:
Region 2:
Region 3:
Region 4:
 We have just talked about the signs of , but it must be noted that, far from the null-clines,
the magnitude of is much smaller than that of , since b,r<<1.
We divide the F-null-cline into three segments: Segment 1 (to the left of point M), Segment 2
(between points M and N) and Segment 3 (to the right of point N). We will refer to these
segments in the following discussion.

Let us now consider the two initial conditions:

Initial condition at p1: This point is inside region 1. Therefore the flow is leftwards, with a small
upward component. The system state approaches the origin and settles there, confirming our
earlier result that the origin is the only stable point.

Initial condition at p2: This point is inside region 4, where the flow is rightwards with a small
upward component. Therefore the system state moves rightwards until it hits the F-nullcline
inside Segment 3. Since the flow has no horizontal component on the F-nullcline, the small
upward component pushes the state upwards. In Segment 3 of F-nullcine there is a tendency for
the state stay on the F-nullcline, since the flow is leftwards on its right, and rightwards on its left.
Therefore, the state creeps along the F-nullcline in the upward direction until it reaches the
topmost point, N, on the F-nullcline. Beyond this point the flow is still upwards and leftwards
while the F-nullcline bends downwards. Therefore, the state leaves the F-nullcline and drifts
leftwards until it reaches Segment 1. The situation now similar to what we encountered on
Segment 3, but with a downward flow component. Therefore the state creeps downward along
Segment 3 until reaches the origin where it finally settles down.
​ Therefore, when p2 is the initial condition, the system exhibits this large excursion by
which the membrane voltage reaches a maximum before it returns to the origin. Such an
excursion of membrane voltage resembles an action potentials.
​ Therefore the FitzHugh-Nagumo neuron model exhibits excitability. For the initial
membrane less than a threshold value (=a), the voltage quickly returns to zero (fig. [Link]).
When the initial voltage exceeds the threshold (=a), the voltage exhibits an action potential (fig.
[Link]).
 Initial voltages

​ ​ ​

Figure [Link]:​ a) The nullclines vand w, and 'v' simulation from b) p1 and c) p2

5.1.2)​ Limit Cycles (Ia> 0):

As increases, for a range of values of , the w-nullcline intersections the

F-nullcline in the “middle branch” where the F-nullcline has a positive slope. In this case
too there only a single intersection.
 Figure [Link]:​ Phase plane analysis: Oscillations at a=0.5; b=0.1; r=0.1; Ia=0.6​

Figure [Link]:​ a) The nullclines v and w and 'v' simulation with initial b)
v(1)=0.4

Therefore, the stationary point is unstable.

The rough ‘arrowplot’ in fig.([Link]) above shows that there is a ‘circulating field’
around the stationary point, which is unstable. Thus it can be expected that there is a limit
cycle enclosing the stationary point, which is actually true. Fig. [Link] shows the
oscillations in membrane voltage (v) produced by a MATLAB program.

5.1.3)​ Depolarization (higher ):

 As increases further, the two nullclines intersect in the “right branch” of the
F-nullcline where the slope of F-nullcline is negative.

Figure [Link]:​ Phase plane analysis: Depolarisation at a=0.5; b=0.1; r=0.6; Ia=0​
.3

Figure [Link]:​ a) The nullclines v and w and 'v' simulation with initial b)
v(1)=0.4

Since,
 We know that the stationary point is [Link] this case, the membrane voltage remains
stable at a high value. This corresponds to regime 4 (in fig [Link]) in the HH model
where for a sufficiently high current, the neuron does not fire but remains tonically
depolarized.

5.1.4)​ ​ Bistability

Some real neurons exhibit bistable behavior – their membrane voltage can remain
at tonically high (“UP” state) or tonically low (“DOWN” state) values. These
UP/DOWN neurons are found in for example in medium spiny neurons of Basal
ganglia striatum.

The FN model exhibits bistability for a certain range of model parameters.

Fig. [Link] below shows a configuration in which the null-clines intersect at three
points (p1, p2 and p3). It can easily be shown that p1 and p3 are stable, and p2 is
a saddle.

​ Figure [Link]:​ Phase plane analysis: Neuronal on-off bi-stable behavior at

a=0.5; b=0.01; r=0.8; Ia=0.02​
 ​ Figure [Link]:​ a) The nullclines v and w, and 'v' simulation with initial b)
p1: v(1)=0.4 c) p3: v(1) = 0.8

@ p1: (V(1) = 0.4)

@p2: (V(1) = 0.5)

@p3: (V(1) = 0.8)

​ ​

FN model in this case, can remain stable at either p1(low v value, DOWN state) or at
p3 (high v value, UP state).

5.2​ Morris-Lecar Model

The Morris-Lecar (ML) model (Morris and Lecar 1981) describes the membrane voltage
dynamics of the barnacle muscle fiber. It consists of three channel currents:
a)​ A fast activating Ca2+ current (activating variable – m)
 b)​ A delayed rectifying K+ current (activating variable – w)
c)​ A leakage current

Simplifying Assumption:
As in the case of FN model, we assume that the dynamics of m-variable is fast.
Therefore, m may be substituted by .

Equations of the ML model are given as,

​ ​ ​ ​ (5.2.1)

​ ​ ​ ​ ​ ​ ​ ​ ​ (5.2.2)

​ ​ ​ ​ ​ ​ ​ ​
(5.2.3)

​ ​ ​ ​ ​ ​ ​ ​
(5.2.4)

​ ​ ​ ​ ​ ​ ​ ​ (5.2.5)
v-nullcline:

​ ​ ​ ​ ​ ​ (5.2.6)

w-nullcine:

​ ​ ​ ​ ​ ​ ​ (5.2.7)

Note that the v-nullcline is ‘inverted N’-shaped, similar to the F-nullcline of FN model. The
w-nullcline has a sigmoidal shape. The v-nullcine can be divided into ‘left’, ‘middle’ and ‘right’
branches. Depending on the branch in which the intersection of V- and w-nullclines, we have
different dynamics. The v-nullcine rises with increasing Ia (fig. 5.2.2).

Figure 5.2.1:​ Bifurication plot for the ML model

For Ia in the range (0-68 pA): the intersection point is in the left branch. The point can be shown
to be stable.
Therefore, the ML model does not produce spikes (fig. 5.2.2).
 Figure 5.2.2:​ Nullclines and Voltage simulation for I = 60 pA

For Iain the range (69-217 pA): the intersection point is in the middle [Link], the ML
model produces limit cycle oscillations (fig. 5.2.3).

Figure 5.2.3:​ Nullclines and Voltage simulation for I = 150 pA

ForIain the range (218 pA ): the intersection point is in the right [Link], the membrane
voltage exhibits depolarization (fig. 5.2.4).

Figure 5.2.4:​ Nullclines and Voltage simulation for I = 300 pA

5.3​ INAPK Model

So far we have visited two simplified neuron models that are capable of exhibiting some
of the basic dynamic properties of a neuron: resting state, oscillations and depolarization.
Both the models are 2-variable systems which are preferred since they can be studied
conveniently using phase-plane techniques. The FitzHugh-Nagumo model was reduced
systematically and derived from the HH model, whereas the Morris-Lecas model is
slightly different from the HH model.

Let us now consider another model that is constructed by systematic reduction from the
HH model. This reduced model, known as the Ina,p + IK model has a persistent sodium
channel and a potassium channel. The “persistent” sodium channel (if you recall sodium
channel model from chapter 3), has only an activation variable; there is no inactivation
variable. The potassium channel also has a single activation variable.

Membrane voltage dynamics of such a model can be written as usual in the following
form:

​ ​ ​ ​ (5.3.1)

​ ​ ​ ​ ​ ​ ​ ​ ​ (5.3.2)

​ ​ ​ ​ ​ ​ ​ ​ ​ (5.3.3)​
If we assume that ‘m’ dynamics is fast, we replace m in eqn. (5.3.1) above with m∞ and
eliminate the ‘m’ dynamics (eqn. (5.3.2)). We are now left with the following two
equations:

​ ​ ​ ​ (5.3.4)

​ ​ ​ ​ ​ ​ ​ ​ ​ (5.3.5)

m∞ and τm are modeled on the lines of eqns. (5.3.6,5.3.7).

​​ ​ ​ ​ ​ ​ ​ (5.3.6)

​​ ​ ​ ​ ​ (5.3.7)

Thus, the parameters for m∞ are:

V1/2 = -20 mV and λ = 15,

And the parameters for n∞ are,
V1/2 = -25 mV and λ = 5. τn(V) = 1 for all V.
​
Rewriting eqns. (5.3.4, 5.3.5) above, we may express the V- and n- null-clines as follows:

V-nullcine:​

​ ​ ​ ​ ​ ​ (5.3.8)
and
n-nullcline:

​ ​ ​ ​ ​ ​ ​ ​ ​ ​ (5.3.9)

The V- and n-nullclines obtained from the above equations are shown in fig. 5.3.1, 5.3.2 for Ia =
60 pA and 500 pA.

Figure 5.3.1:​ Nullclines Ia = 60 pA

 Figure 5.3.2:​ Nullclines Ia = 500 pA

5.4​ A Simplified two-dimensional model:

The three models visited above in this chapter have certain common features. They
differ in the ion channel composition; they are derived from different cell models by
different steps of simplification. But all three share certain features of dynamics. All the
three models have:
-​ Resting state for I = 0
-​ Spikes/oscillations/limit cycles for a sufficiently large I.
-​ An N-shaped V-nullcline. The shape of the second nullcine is linear in the FN
model, and sigmoidal in ML and INAPK model.
-​ The intersection between the V-nullcine and the second nullcline near the U-shaped,
lower arch of the V-nullcline plays a crucial role in neuron dynamics.

​ When the second nullcline intersects the V-nullcline in the 1st branch there is a
resting state; when the intersection occurs in the middle branch there are
oscillations.
 Consider the role of the remaining part of the state space. In the FN model for
example, when the external current I = 0. When the initial voltage (V(0)) is less
than ‘a’, the system returns to the resting potential, 0. But if V(0) > a, the neuron
state (V, w) increases until it touches the 3rd branch of the V-nullcine, climbs up
towards the maximum of the V-nullcline, turns leftwards at the maximum, proceeds
up to the 1st branch of the V-nullcline, before it slides down the 3rd branch to the
resting state. Thus the remaining part of the phase-plane determines the downstroke
and the peak of the action potential.

When the intersection point shifts slightly from the left of the minimum of
V-nullcline (where the system exhibits excitability) to the right of the minimum
(where the system oscillates), the shape and size of the action potential are about the
same. Thus in the oscillatory regime too, the part of the phase-space other than the
shaded portion shown in fig. 5.4.1 merely contributes to downstroke and the peak
value of the action potential.

Figure 5.4.1 : Area of interest for the upstroke: Phase plane analysis- Oscillations
(FN model) at a=0.5; b = 0.1; r=0.1; Ia=0.6

Similar comments can be said about the dynamics in the other two models also.

​
Thus the above three systems can be reduced to a more general, simpler system as
follows. The U-shaped, lower arch of the V-nullcline is approximated by a quadratic
function. Since only the shape of the second null-cline at its intersection with V-nullcline
matters, and not its shape elsewhere, the second nullcline is approximated by a
straightline. Dynamics far from the intersection point are implemented by a simple
resetting once the membrane voltage hits a peak value or a minimum.

The quadratic approximation of the V-nullcline at its minimum is given as,

 ​ ​ ​ ​ ​ ​ ​ (5.4.1)

The linear approximation of the second nullcline may be expressed as,

​​ ​ ​ ​ ​ ​ ​ ​ (5.4.2)

The dynamics of the reduced system then becomes,

​ ​ ​ ​ ​ ​ (5.4.3)

​ ​ ​ ​ ​ ​ ​ (5.4.4)

Where τV and τu denote the time-scales of V- and u-dynamics. Since V and u denote the
fast and slow variables respectively, τV << τu.

If the V-dynamics of Eqn. (5.4.3) above is implemented without any auxiliary conditions,
V blows up to infinity in a finite time 't'. This can be shown easily on integrating Eqn.
(5.4.3) as:

​ ​ ​ ​ ​ ​ ​ ​ (5.4.5)

Therefore, the downstroke of the action potential is modeled by resetting V(t) when it
reaches the peak value as follows.

, when .​ ​ ​ ​ ​ (5.4.6)

When the membrane voltage exceeds , both Vand u are reset instantaneously as
specified by the above equation.
Eqns. (5.4.3, 5.4.4) may be rewritten in a simpler form as,

​ ​ ​ ​ ​ ​ ​ ​ (5.4.7)
 ​ ​ ​ ​ ​ ​ ​ ​ (5.4.8)
If , ​ ​ ​ ​ ​ ​ ​ (5.4.9)
a,b,c,d = constants;

The last set of eqns. (5.4.7- 5.4.9) is generally referred to as the Izhikevich neuron model in
current literature. Its merit lies in the low computational cost, and the ability to reproduce firing
patterns of a large variety of neurons (E.M. Izhikevich et al, 2004, 2004).

5.4.1​ Quadratic integrate and fire neuron:

There is a simpler, one-dimensional version of the Izhikevich model of eqns.

(5.4.6-5.4.8). This model, known as the quadratic integrate and fire neuron model (
Latham et al., 2000 ) consists of only the membrane voltage dynamics with quadratic
nonlinearity:

​ ​ ​ ​ ​ ​ ​ ​ ​ ([Link])

If , ​ ​ ​ ​ ​ ​ ​ ([Link])

Note that for any non-zero value of I or V(0), eqns. ([Link]) above blows up in a finite time, as
already shown above on Eqn. 5.4.5. The resetting condition of eqn. (5.4.6) prevents the blowing
up.
By appropriate choice of the model parameters - c- the above model can be made to express a
variety of neurodynamic behaviors.

If I < 0, at two values, . Let us call these roots, and

. These names can be easily justified.
​
For , . Therefore, V increases indefinitely.
Similarly, , and therefore V decreases towards . Therefore

is an unstable point.

For , . Therefore, V increases towards . Hence is a stable

state.

Now consider the dynamics of the neurons in the following 3 cases:

Case i:

For I < 0, and are real numbers.

If , V approaches and remain there forever.

If , V grows indefinitely and, unless clamped, reaches infinity in a finite

time.

The reset condition of eqn. ([Link]) prevents the runaway of V, and resets it to as
soon as V reaches .

Since , V now tends to and settles there. This latter behavior is

described as ‘excitability.’ It is analogous to the case of FN model when I = 0.

Excitable behavior may be produced not just by the initial condition, but also by giving a
series of pulses which progressively push the membrane voltage towards and
beyond, causing excitation.
 Figure [Link]: a) Reset condition b) Voltage simulation at I = -5; vo = 3; vreset = -3;
vthreshold = 4;

Case (ii):

In this case, the model behavior for V(0) < is similar to that of the previous case.

But when V(0) > , V quickly reaches and gets reset to . But instead of a slow

return to , it rises again and again to , exhibiting continuous, periodic firing. Thus the
model in case (ii) exhibits two stable states: one corresponding to the resting state, and the other
the state of continuous firing. The model is therefore said to have bistability.

Figure [Link]: a) Reset condition b) Voltage simulation at I = -5; vo = 3; vreset = 3;

vthreshold = 20;
Case (iii): I > 0

When I > 0 and and are unreal (Fig. 5.4.1.3a). Therefore V rapidly

approaches and gets reset to , again and again.

This neuron shows tonic firing, even without excitation.

Figure [Link]: a) Reset condition b) Voltage simulation at I = 5; vo = 3;

5.4.2​ Leaky Integrate and fire neuron model:

This is the simplest model of spike generation. An electric circuit implementation of it consists of
a capacitance charged by a current, and discharged whenever the voltage across the capacitance
exceeds a limit.

Figure [Link]: Leaky Integrate and Fire circuit diagram

Applying the Kirchoff’s current law of electric circuits, the external current, Iext, going into the
circuit may be expressed as,

​ ​ ​ ​ ​ ​ ​ ([Link])

Capacitance is discharged whenever, V >θ, a threshold value. It is assumed that whenever the

capacitor is charged, the “neuron” emits a spike. Note that in this model, there is no nonlinear,
explosive build of excitation reaching a peak producing an action potential. The spike generated
in this neuron model is more notional. This has always been one of the points of criticism about
the leaky integrate and fire model.

If the capacitance starts off at 0 voltage, let us consider the time taken by the capacitance to reach

the threshold, θ.

If Iext is a constant current, I0, voltage variation while the capacitance is charging may be
expressed as,

whereτ the time-constant of the circuit equals RC. Since the charging stops at V =θ, the time

taken, T, to reach this threshold is given by setting V(t) =θ, or,

Or

Since there is a spike every time the capacitance discharges, spike frequency, f, is the reciprocal
of T.

​ ​ ​ ​ ​ ​ ([Link])

A key property of a real neuron reproduced by the above model is thresholding effect. Note that

the model exhibits firing only when R I0>θ. But as I0 increases beyond R/θ, f increases

indefinitely (Fig. [Link]), instead of saturating as it happens in a real neuron.

Figure [Link]: frequency vs Io without inclusion of absolute refractory period

In order to restore the saturation property, the condition of absolute refractory period is
introduced into the above model. Accordingly the capacitance can begin to be charged only a
little while, Tref, after the external current is applied. Therefore, the time taken to charge is
incremented as,
And the new firing rate is,

​ ​ ​ ​ ​ ([Link])

With the inclusion of absolute refractory period, the plot of I0vs f, shows saturating behavior (Fig.
[Link]).

Figure [Link]: frequency vs Io with inclusion of absolute refractory period

5.4.3​ Binary neuron models:

The simplified neuron models visited so far in this chapter are models of action potential or spike
generation. If we wish to simply the single neuron model further, we may consider rate-coded
models which represent the neuron state in terms of the spike rate. These models typically
describe neurons as binary elements, with a high (excited state) and a low (resting) state.

[Link]​ Dynamic Binary neuron model:

A dynamic version of a binary neuron is a bistableneuron whose dynamics is described as,

​ ​ ​ ​ ​ ​ ​ ([Link].1)
​​ ​ ​ ​ ​ ​ ​ ([Link].2)

Figure [Link].1: u,v nullclines

V is an abstract quantity that denotes if the neuron is excited (V is close to 1) or in the resting
state (V is close to -1).

For I = 0 and λ> 1, the above system has three fixed points (Fig. [Link].1). The fixed point at the

origin can be shown to be unstable, while the ones on the flanks can be shown to be stable. When
released from a random initial condition, the neuron state, V, settles near +1 or -1.
For I > b, where b = 2.63, and λ> 1, Eqn ([Link].1) above has only a single, stable fixed point at

V close to +1.

For I < -b, andλ> 1, eqn ([Link].1) above has only a single, stable fixed point at V close to -1.

For –b < I < b, there are again three fixed points, two close to + 1 and -1 respectively, and the
third somewhere in the middle, not necessarily at the origin.

The above neuron model also shows hysteresis effect. If I is reduced from a large positive value
to a large negative value and back, the points at which V makes transitions from (+1 to -1 and
vice versa) are different in the forward pass and the reverse pass (Fig. [Link].2).

Figure [Link].2 : Hysteresis curve for a dynamic binary neuron

[Link]​ Static Binary Neuron Model:

The McCulloch and Pitts neuron model is a good example of a static binary neuron. It combines
the inputs, xi, that it receives from other neurons and computes its output, y. The relationship
between the inputs xi and the output y is given as,

​ ​ ​ ​ ​ ​ ​
Where wi s are the synaptic strengths, or the “weights” of the connections from the neurons that
send inputs to the neuron of interest; q is the threshold for excitation; g(.) is known as a transfer
function, which typically has a sigmoidal shape.

Four types of transfer functions are usually considered depending on the range of values that y is
permitted.

Hardlimiting nonlinearity:
y = 1/0:

​ ​ ​ Figure [Link].1: Hardlimiting nonlinearity- y = 1/0:

y = +1.

​ ​ ​ Figure [Link].2: Hardlimiting nonlinearity- y = +1

Smooth sigmoid nonlinearity:

Logistic function:

​ ​ ​ Figure [Link].3: Logistic function-

Tanh(.) function:

​ ​ ​ Figure [Link].4: Tanh(.) function -

References:
E. M. Izhikevich, J. A. Gally, and G. M. Edelman, “Spike-timing dynamics of neuronal
groups,” Cerebral Cortex, vol. 14, pp. 933–944,2004.

E. M. Izhikevich, N. S. Desai, E. C. Walcott, and F. C. Hoppensteadt, “Bursts as a unit of

neural information: Selective communication via resonance ,” Trends Neurosci., vol. 26,
pp. 161–167, 2003.
E. M. Izhikevich, “Simple model of spiking neurons,” IEEE Trans. Neural Networks, vol.
14, pp. 1569–1572, Nov. 2003.

Latham, P. E., B. J. Richmond, P. G. Nelson, and S. Nirenberg. (2000). "Intrinsic

Dynamics in Neuronal Networks. I. Theory.". Journal of Neurophysiology 88 (2):
808–27.