Bifurcations in the modified Hodgkin-Huxley neuron model

BIFURCATIONS IN THE MODIFIED HODGKIN-HUXLEY NEURON MODEL

The MHH model has been simulated numerically first in the absence of noise. For an easy comparison with electrophysiological experiments, we use the interspike intervals, $\tau_n$ , as a state variable. The interspike intervals are registered by the membrane potential crossing a threshold (at -20mV) with positive derivative (Poincaré surface of section). Note that for oscillatory neuron models they can be associated with the instantaneous period of the action potential, while $\tau_{n+1}$ versus $\tau_n$ represents the first return map of the system.

Fig. 1. Bifurcation diagram of the MHH system. The deterministic case.

The bifurcation diagram of the MHH is shown in Fig. 1 and represents a classical route to chaos through a period-doubling cascade³⁰ up to a temperature $T \approx 10^oC$ . The first period doubling bifurcation occurs at $T_1\approx 6.765^oC$ , the second at $T_2\approx 7.195^oC$ . Finally for $T>T_{ch}\approx 7.31^oC$ the system becomes chaotic. Inside the chaotic region we observe several periodic windows opened by a saddle-node bifurcation and closed by a global bifurcation, namely an interior crisis. This picture is very similar to that of the logistic map, for example, the return map $\tau _{n+1}=f(\tau _n)$ looks like a logistic map for T<8.5^oC, that is, a parabola with a single maximum. However, with further increase of the temperature the shape of the return map changes from parabolic to a curve with multiple maxima.

Fig. 2. The return map for two parameters values: (a) T=10.6 ^oC and (b) T=10.7 ^oC.

The focus of our interest is the abrupt increase of the duration of the interspike intervals in the region 10.658^oC<T<10.9^oC (see Fig. 1.). For T<T_cr=10.6589^oC the duration of the interspike intervals is not larger than 1600 ms. At T=T_cr an abrupt explosion of the interspike intervals occurs. In Fig. 2. we show the return map before and after this transition. We see that after the transition an additional, well-separated part of the return map with longer time intervals appears. This transition is well pronounced in the return map, but it is difficult to recognize in the phase space of the system: before and after the transition the phase portraits of the system have nearly the same structure. Thus, the signature of this transition in phase space is the appearance of a new type of orbits within the attractor.

Let us study the explosion in the interspike intervals in detail. As mentioned above, long time intervals between consecutive returns to a Poincaré section can be considered as an indicator for the approach of a homoclinic orbit. Such a homoclinic orbit is associated with an equilibrium state of saddle type and its stable and unstable manifolds. In our case this saddle equilibrium is embedded in the already existing chaotic attractor. It has two stable and two unstable directions. The eigenvalues, corresponding to the stable directions, are real ( $\lambda_1<0, \lambda_2<0$ ) while the eigenvalues, corresponding to the unstable directions, are complex ( $\gamma_{1,2}=\rho \pm i\omega, \rho>0$ ). Thus we have an equilibrium state of saddle-focus type. If we suppose that there exists a chaotic set with an embedded homoclinic orbit then the eigenvalues of the equilibrium state should satisfy the Shilnikov condition:¹⁸ $-\lambda>\rho$ . Indeed, the eigenvalues at a temperature value T=10.7456^oC which is, as we will show, close to the homoclinic bifurcation are: $\lambda_1=-0.182$ , $\lambda_2=-0.146$ , $\gamma_{1,2}=0.327\times10^{-2}\pm {\rm i}0.282\times10^{-2}$ (the units of the eigenvalues are in (ms)^-1). From those values one can see, that the two real eigenvalues are comparable in magnitude. For this reason one cannot easily reduce the problem to the known cases of homoclinic bifurcations in R³.

Fig. 3. Periouds of different limit cycles versus temperature: (1, solid line) the period-1 cycle, (2, dashed line) the period-2 cycle,and (3, dotted line) the period-3 cycle.

In this model we have to deal with an equilibrium state of saddle-focus type having a two-dimensional stable and a two-dimensional unstable manifold. The return time for the homoclinic orbit itself is infinity. A separatrix loop is built up, where the outgoing path is tangent to the unstable manifold while the returning path is tangent to the stable manifold of the equilibrium state. Nearby trajectories on the attractor approach the equilibrium state along its stable manifold. The motion slows down so that the trajectory spends a long time in the neighborhood of the equilibrium state before it is ejected along the unstable manifold. The same happens to saddle limit cycles embedded in the attractor. If they are close to the homoclinic orbit then at least one point of the saddle limit cycle is close to the saddle equilibrium state. As a consequence the Poincaré return time of these saddle limit cycles, that is their period, becomes also very large. The motion on the saddle limit cycle slows down approaching the vicinity of the equilibrium state along its stable manifold. In Fig. 3. we show the dependencies of the periods of different limit cycles versus the control parameter calculated using the continuation software CONTENT.³¹ In particular, we show the period-1 cycle (which is born from the equilibrium point in a Hopf bifurcation), the period-2 cycle (which is born from the period-1 cycle trough a period-doubling bifurcation at $T \approx 6.765^oC$ ), and the period-3 cycle (which is born via saddle-node bifurcation at $T\approx 7.58^oC$ (see Fig. 1.). Although these cycles possess turning points at different temperature values, their periods diverge at the same temperature value $T \approx 10.7456^oC$ . Extensive calculations have shown qualitatively the same behavior for other limit cycles of higher periods: they are accumulated and their periods diverge at $T \approx 10.7456^oC$ .

Fig. 4. Two dimensional projection of the phase portrait of the system at T=10.65oC (a) and T=10.7456oC (b). A neighborhood of the saddle-focus equilibrium point is also shown. The gray dots represent the chaotic attractor; the black filled circle is the saddle-focus equilibrium state; the unstable limit cycles of periods 1, 2 and 3 are shown by blue, red and green bold lines, respectively.

The computation of the homoclinic orbit itself using direct numerical methods or perturbation methods is rather difficult since the saddle-focus has a two-dimensional stable and two-dimensional unstable manifold with two stable real eigenvalues very close to each other. Therefore, we can show the evidence of the homoclinic bifurcation only by the discussed indirect indicators.

Fig. 5. (a): Minimal distance between the saddle-focus equilibrium state and the trajectories on the chaotic attractor versus control parameter. (b): Maximal residence time of the trajectories in a small neighborhood ( $\delta =0.01$ )of the saddle-focus.

The phase portrait of the system reflects the behavior discussed above. In Fig. 4. we show a two dimensional projection of the chaotic attractor for two parameter values: one before the suspected homoclinic bifurcation and the other at the temperature value where we find the divergent behavior discussed above. We also plot the saddle cycles and the saddle-focus equilibrium point. At T=10.65^oC, that is before the suspected bifurcation, the saddle-focus is not included into the saddle cycles. The saddle cycles are located far from each other. However, at the temperature T=10.7456^oC the saddle cycles are extremely close to the saddle-focus. Moreover, they are extremely close to each other, that is, the period-1 cycle is a part of the period-2 cycle, the period-2 cycle is a part of the period-3 cycle and so on. At the same time, as it can be seen from the figure, a large fraction of the phase trajectories passes the vicinity of the saddle focus. Unfortunately we were not able to compute the homoclinic orbit itself, but each of the presented saddle orbits can be considered as a candidate for it.

We note, that as soon as the phase trajectory approaches the vicinity of the equilibrium state the motion of the system slows down, so that the phase trajectory spends a long time in the neighborhood of the equilibrium state. This gives rise to very long interspike intervals. To support this argument we present in Fig. 5. the dependence of the minimal distance between the saddle-focus and the trajectory on the chaotic attractor versus control parameter and the maximal residence time of the phase trajectory spent in a small neighborhood of the equilibrium state. As can be seen, the distance becomes very small at a certain parameter value and the phase trajectories of the system spend a maximal time near the equilibrium state. These plots together with the discussed bifurcation diagram suggest that the homoclinic bifurcation occurs at $T \approx 10.7456^oC$ .