5 Nonlinear Integrate-and-Fire Models

5.4 Summary

The standard leaky integrate-and-fire model is rather limited in scope, since it has one universal voltage threshold. Nonlinear integrate-and-fire neurons, however, can account for the fact that in real neurons the effective voltage threshold for repetitive firing is different than the voltage threshold found with short current pulses. These two voltage thresholds, which are related to the minimum of the nonlinearity f(u)f(u) and to the unstable fixed point, respectively, are intrinsic features of nonlinear integrate-and-fire models. Once the membrane potential is above the intrinsic threshold, the upswing of the membrane potential starts. The integration is stopped at a numerical threshold θreset\theta_{\rm reset} which is much higher and conceptually very different than the intrinsic firing threshold of the model. In fact, the exact value of the numerical threshold does not matter, since, without such a threshold, the membrane potential would go to infinity in finite time.

In principle, many different forms of nonlinearity are imaginable. It turns out, however, that many neurons are well described by a linear term (the ‘leak’ term) combined with an exponential term (the ‘activation’ term); cf. Fig. 5.4. Therefore, the exponential integrate-and-fire model has the ‘correct’ nonlinearity, whereas the quadratic integrate-and-fire model is too nonlinear in the subthreshold regime and too slow in the superthreshold regime.

Both the exponential and the quadratic integrate-and-fire model show a frequency-current curve of type I. Indeed, close to the bifurcation point the two models become identical and can be mapped onto the canonical type I model.

Literature

Nonlinear generalizations of the one-dimensional equation of the leaky integrate-and-fire model can be found in the paper of Abbott and van Vreeswijk (1993). While the form of the nonlinearity was left open at that time, analytical methods from bifurcation theory suggested that there was a canonical one-dimensional model which describes the saddle-node-onto-limit-cycle bifurcation of type I models. This is called the canonical type I model or Theta model or Ermentrout-Kopell model (140; 141). While the quadratic nonlinearity in the voltage appeared in many early papers on bifurcations and type I models (141; 140; 499; 228) the active use of the quadratic integrate-and-fire model for simulations seems to have been started by Latham et al. in 2000 and later popularized by Izhikevich (239).

Based on biophysical and mathematical arguments, the exponential integrate-and-fire model was introduced in 2003 by Fourcaud-Trocme, Hansel, Vreeswijk and Brunel. While the quadratic integrate-and-fire model is the canonically correct type I model close to the bifurcation point, i.e. at the transition to repetitive firing, there is no fundamental reason that the quadratic nonlinearity should also be correct further away from the threshold, and it was unclear at that time what the nonlinearity of real neurons would look like. The issue was settled by the work of Badel et al. (35, 34) who designed an experimental method to measure the nonlinearity directly in experiments. The nonlinearity found in experiments is extremely well matched by the exponential integrate-and-fire model.

The link between ion channel activation or inactivation in conductance-based neuron models and the parameters of an exponential integrate-and-fire model with refractoriness (35; 34) is discussed in the overview paper of Platkiewicz and Brette (401).

We will see in the next chapter that the single-variable exponential integrate-and-fire model as it stands is not sufficient to account for the wide variety of neuronal firing patterns, but needs to be complemented by a second variable to account for slow processes such as adaptation - and this leads eventually to the Adaptive Exponential Integrate-and-Fire Model

Exercises

  1. 1.

    Quadratic vs. exponential integrate-and-fire For a comparison of the two models, take a look at Fig. 5.8B and answer the following questions:

    (a) Show for the exponential integrate-and-fire model that the minimum of the nonlinearity f(u)f(u) occurs at u=ϑrhu=\vartheta_{rh}. Calculate the curvature d2/du2f(u)d^{2}/du^{2}f(u) at u=ϑrhu=\vartheta_{rh}.

    (b) Find parameters of the quadratic integrate-and-fire model so that it matches the location and curvature of the exponential model in (a).

    (c) Suppose that the value of the numerical threshold is θreset=ϑrh+2ΔT\theta_{\rm reset}=\vartheta_{rh}+2\Delta_{T}. When the threshold is reached the membrane potential is reset to ur=ϑrh-2ΔTu_{r}=\vartheta_{rh}-2\Delta_{T}. Sketch qualitatively the trajectory u(t)u(t) for both neuron models in the repetitive firing regime. Pay particular attention to the following questions: (i) For u>ϑrhu>\vartheta_{rh} which of the two trajectories rises more rapidly towards the numerical threshold? (ii) For u=ϑrhu=\vartheta_{rh}, do both trajectories have the same speed of rise du/dtdu/dt or a different one? (iii) For u<ϑrhu<\vartheta_{rh}, which of the two trajectories rises more rapidly? (iv) Should your drawing of the trajectories follow any symmetry rules? Compare your results with Fig. 5.9.

  2. 2.

    Zero-curvature and rheobase threshold . Experimentalists sometimes determine the threshold voltage for spike initiation by searching for the zero-curvature point in the voltage trajectory, just before the upswing of the action potential. Show for a general nonlinear integrate-and-fire model that the zero-curvature point is equivalent to the rheobase threshold ϑrh\vartheta_{rh}. Hints: (i) Minimal curvature means that the voltage trajectory passes through a point d2u/dt2=0d^{2}u/dt^{2}=0. (ii) Note that the slope of the voltage trajectory du/dtdu/dt is given by the function f(u)f(u). (iii) Recall that the minimum of ff is taken at u=ϑrhu=\vartheta_{rh}.

  3. 3.

    Exponential integrate-and-fire and sodium activation. To derive the exponential integrate-and-fire model from a Hodgkin-Huxley model with multiple channels, follow the procedure indicated in Section 5.2.2 and perform the following steps.

    a) Show that NN linear currents can always be rewritten as one effective leak current kgk(u-Ek)=geff(u-Eeff)\sum_{k}g_{k}\,(u-E_{k})=g^{\rm eff}\,(u-E^{\rm eff}). Determine the effective conductance geffg^{\rm eff} and the effective reversal potential EeffE^{\rm eff}.

    b) Assume that the sodium activation is rapid and approaches an equilibrium value m0(u)=1/{1+exp[-β(u-θact)]}=0.5{1+tanh[-β(u-θact)/2]}m_{0}(u)=1/\{1+\exp[-\beta\,(u-\theta_{act})]\}=0.5\,\{1+\tanh[-\beta\,(u-% \theta_{act})/2]\} where θact\theta_{act} is the activation threshold of the sodium channel and β\beta the sharpness of the threshold. Using a Taylor expansion, show that for u<θact-β-1u<\theta_{act}-\beta^{-1} the activation function can be approximated by an exponential m0(u)=exp[β(u-θact)]m_{0}(u)=\exp[\beta\,(u-\theta_{act})].

    c) Assume that u<θact-β-1<ENau<\theta_{act}-\beta^{-1}<E_{\rm Na} and show that the sodium current

    INa=gNa[m0(u)]3hrest(u-ENa)I_{\rm Na}=g_{\rm Na}[m_{0}(u)]^{3}\,h_{\rm rest}\,(u-E_{\rm Na}) (5.21)

    gives rise to the exponential term.

    d) Using the steps (a) - (c), map the equation (5.15) to the exponential integrate-and-fire model

    dudt=f~(u)=-u-urestτ+ΔTτexp(u-ϑrhΔT){du\over dt}=\tilde{f}(u)=-{{u-u_{\rm rest}}\over\tau}+{\Delta_{T}\over\tau}\,% \exp\left({u-\vartheta_{rh}\over\Delta_{T}}\right) (5.22)

    Determine the parameters urestu_{\rm rest}, τ\tau and ΔT\Delta_{T}.

  4. 4.

    Refractory exponential integrate-and-fire model and sodium inactivation . In the previous exercise, we have assumed that h=hresth=h_{\rm rest} is constant throughout the spiking process. Repeat the steps of the previous calculation, but set h(t)=hrest+x(t)h(t)=h_{\rm rest}+x(t) and show that the sodium inactivation dynamics leads to an increase in the the firing threshold a few milliseconds after the spike.

    Hints. (i) Linearize the dynamics of gating variable hh in the Hodgkin-Huxley model around h=hresth=h_{\rm rest} so as to derive dx/dt=-x/τ0dx/dt=-x/\tau_{0}. (ii) Assume that during each spike, the inactivation variable hh increases by a fixed amount Δh\Delta h.

  5. 5.

    Quadratic integrate-and-fire model .

    Consider the dimensionless quadratic integrate-and-fire model

       𝑑      𝑑   tu=(u-urest)(u-uc)+I,{{\text{d}}\over{\text{d}}t}u=\,(u-u_{\rm rest})\,(u-u_{c})+I\,, (5.23)

    with urest=-1u_{\rm rest}=-1 and uc=1u_{c}=1 and I=0I=0.

    (a) Suppose that a trajectory starts at time t=0t=0 at --\infty. How long does it take to get close to the resting state and reach the value u=-1-ϵu=-1-\epsilon where ϵ1\epsilon\ll 1?

    (b) Suppose that a trajectory is initialized at u(0)=1+ϵu(0)=1+\epsilon. where ϵ1\epsilon\ll 1 how long does it take to reach u=?u=\infty?

    (c) Initialize as in (b) but stop the integration at u=θresetu=\theta_{\rm reset}. How long does the trajectory take to reach the reset threshold? Is the difference between (b) and (c) important?

    Hint: Use dx/[1-x2]=arcoth(x)=0.5[ln(x+1)-ln(x-1)]\int dx/[1-x^{2}]={\rm arcoth}(x)=0.5[{\rm ln}(x+1)-{\rm ln}(x-1)] and cothh(x)=[ex+e-x]/[ex-e-x]{\rm cothh}(x)=[e^{x}+e^{-x}]/[e^{x}-e^{-x}]

  6. 6.

    Gain function of the quadratic integrate-and-fire model

    We consider a quadratic integrate-and-fire model in the super-threshold regime

       𝑑      𝑑   tu=u2+I,{{\text{d}}\over{\text{d}}t}u=u^{2}+I\,, (5.24)

    with constant current I>0I>0.

    (a) Give a transformation of variables from Eq. (5.23) to (5.24).

    (b) Calculate the duration TT of a trajectory which starts at time t=0t=0 at u(0)=-u(0)=-\infty and ends at time TT at u(T)=+u(T)=+\infty.

    Hint: Use dx/[1+x2]=artan(x)\int dx/[1+x^{2}]={\rm artan}(x)

    (c) Plot the frequency ν=1/T\nu=1/T as a function of II.