Abstract

We introduce a simple two-dimensional model that extends the Poincaré oscillator so that the attracting limit cycle undergoes a saddle node bifurcation on an invariant circle (SNIC) for certain parameter values. Arbitrarily close to this bifurcation, the phase-resetting curve (PRC) continuously depends on parameters, where its shape can be not only primarily positive or primarily negative but also nearly sinusoidal. This example system shows that one must be careful inferring anything about the bifurcation structure of the oscillator from the shape of its PRC.

1.  Introduction

Stimuli delivered to biological oscillators affect the timing of subsequent cycles. Because of the functional importance of this phase resetting, theoretical and experimental studies have probed resetting in diverse settings (Hoppensteadt & Keener, 1982; Winfree, 2000; Glass & Mackey, 1988; Izhikevich, 2007; Canavier, 2006; Rinzel & Ermentrout, 1998; Ermentrout & Terman, 2010). The results of such experiments are typically represented by phase-resetting curves (PRCs), which give the change of the phase of an oscillation as a function of the phase of a stimulus, and phase transition curves (PTCs), which give the new phase of an oscillation as a function of the phase of a stimulus (precise definitions are given below). If the original oscillation rapidly reestablishes itself following a pulsatile stimulus, then these curves can be used to analyze the effects of inputs to biological oscillations as well as the effects of coupling and noise on systems of biological oscillators. For example, phase-resetting curves can be used to predict the effects of periodic stimuli to oscillating heart cells (Guevara, Glass, & Shrier, 1981) and can also be used to analyze synchronization in systems of oscillators with pulsatile coupling (say, neurons with short-lasting synapses or fireflies) (Strogatz, 1994). The shape of the PRC determines the manner in which a pair of mutually connected oscillators will synchronize and how sensitive the locking will be to noise and heterogeneities (Hansel, Mato, & Meunier, 1995; Ermentrout, 1996). The shape also plays a crucial role in how oscillators react to correlated noise and other signals (Ermentrout, Galán, & Urban, 2008; Marella & Ermentrout, 2008). Thus, there is a good deal of interest in what aspects of an oscillator determine the shape of the PRC.

A common rule of thumb about PRCs in neuroscience is that oscillators that are near a saddle node on an invariant circle (SNIC) are of one sign (usually positive) (Ermentrout, 1996; Brown, Moehlis, & Holmes, 2004; Canavier, 2006), while those in which the oscillation arises from a Hopf bifurcation are bimodal (Ermentrout & Kopell, 1984; Brown et al., 2004). Recently Achuthan, Butera, and Canavier (2011) showed that as one moves away from the SNIC, a substantial bimodality developed in the PRC of the Morris-Lecar neural model. The main goal of this letter is to extend the results of a recent model (Krogh-Madsen, Butera, Ermentrout, & Glass, 2012), which demonstrates that PRCs that are sinusoidal in shape can also occur in systems that are arbitrarily close to a SNIC bifurcation. This analytically solvable system allows one to continuously transition from PRCs that are of one sign to PRCs that are sinusoidal.

2.  Terminology

The PRC can be defined as follows for a given oscillator. Suppose the period of the oscillator is T and that we define a phase, , which is a time-like variable that increases at a constant rate of 1. At some phase , a brief perturbation is applied to one or more of the variables in the oscillator. This leads to a sequence of times, at which the oscillator crosses the zero phase line, which must be defined someplace. For excitable systems such as nerve cells or heart cells, it is common to take the point of maximal slope of the upstroke of an action potential as the point of zero phase, in which case it is called the fiducial point. Provided the state point following the stimulus is in the basin of attraction of the limit cycle, there will be an asymptotic return to the limit cycle in the limit . In the absence of the perturbation, Tj=jT, where j is an integer. From the timing of the subsequent action potential upstrokes, we define the PTC mapping to the new phase, :
2.1
For the situation in which the perturbation is small, we define the PRC, which is defined by asymptotic phase shift as
In some situations that arise in experiments (e.g., transients following a stimulus in which several beats are skipped before an oscillation is reestablished), this definition of the PRC curve will not be applicable.

In this letter, we scale the phase so that , , and , but in other contexts, the phase here can be multiplied by or . The degree to which the value of TjTj−1, j>1 differs from T depends on the rate at which the rhythm returns to the limit cycle oscillation.

In the case of a small perturbation and a rapid return to the limit cycle, the relationship between the PRC and PTC is
2.2
To help understand the geometric features of resetting curves, it is useful to define the concept of isochron (Guckenheimer, 1975; Winfree, 2000). Each point in the basin of attraction of a stable limit cycle will asymptotically approach the limit cycle as . Consider two trajectories in phase space xi(t) and xj(t) starting from points xi(0) and xj(0), where xi(0) is a point lying on the limit cycle and xj(0) is a point in the basin of attraction of the limit cycle. The points xi(0) and xj(0) are said to lie on the same isochron if the distance as . The isochrons foliate the basin of attraction of the limit cycle, so that each point in the basin of attraction of the limit cycle lies on an isochron.

3.  A Simple Model of a Nonlinear Oscillator

There are several mechanisms by which a system of differential equations with a stable equilibrium point can lose the stability of this equilibrium and start to oscillate. Two of the best-known bifurcations are the supercritical Hopf bifurcation and the saddle node on an invariant circle (SNIC, also known as SNIPER for infinite-period) bifurcation. In the former case, the equilibrium exists but loses stability, and a small-amplitude limit cycle emerges. There appears to be widespread belief that for systems with a SNIC bifurcation, the PRC is nonnegative so that all stimuli will lead to an advance in the timing of subsequent action potentials, while near a Hopf bifurcation, the PRC will have a positive and negative lobe (Canavier, 2006). Our main point in this letter is to explain the shape of the PRC for a system with a SNIC bifurcation.

Near a supercritical Hopf bifurcation, it is possible to reduce the dynamics of the full set of differential equations to a simple two-dimensional model, which we can write in polar coordinates as:
3.1
where a and q are positive. This equation (or similar equations with a circular limit cycle attractor) has been called the Poincaré oscillator (Glass & Mackey, 1988). In the limit when q=0 or , the PRC and PTC can be readily computed using simple trigonometry, and the resulting equation has been called the radial isochron clock (Winfree, 2000; Hoppensteadt & Keener, 1982). This model has been useful for understanding the effects of periodic stimulation of the heart (Guevara & Glass, 1982).
We describe a simple model that has a SNIC bifurcation but is such that the PRC can be analytically computed in the limit of small amplitude stimuli. Consider the following simple extension of the Poincaré oscillator:
3.2
where we will take without loss of generality. In rectangular coordinates, we can rewrite this as
This model has a SNIC bifurcation when b=1 and stable attracting limit cycle oscillations for . For b>1, the limit cycle is lost, but there are two additional fixed points centered around . Thus, the parameter determines where the dynamics slows down. Figure 1 shows the geometry of this model. The limit cycle of this model is given by r=1 and satisfying
which can be solved by quadrature
3.3
From these formulas, we find , , −T/2<t<T/2, and
defines the phase.
Figure 1:

Geometry of the simple model and illustration of phase resetting. The angle determines where on the limit cycle the slowing occurs. The thin black curves without arrowheads are isochrons and perturbations of magnitude s in the u-direction, occurring when the trajectory is at the dark gray circles, reset the phase to the corresponding light gray circles.

Figure 1:

Geometry of the simple model and illustration of phase resetting. The angle determines where on the limit cycle the slowing occurs. The thin black curves without arrowheads are isochrons and perturbations of magnitude s in the u-direction, occurring when the trajectory is at the dark gray circles, reset the phase to the corresponding light gray circles.

4.  Phase Resetting

Understanding the shape of the PRCs in models is facilitated by understanding the geometry of the isochrons. As is well known (Glass & Mackey, 1988), for the Poincaré oscillator with q=0, the isochrons are equally spaced around the cycle and are aligned along the radii of the circle defined by the limit cycle oscillation. For the modified equation defined by equation 3.2 with 1>b, the isochrons accumulate in the neighborhood of the incipient saddle node bifurcation at .

Figure 2 shows six panels with 10 computed isochrons each at equidistant time intervals for and a=1, 10, b=0.99, 0.999, 0.9999. These isochrons were computed using AUTO (Doedel & Oldeman, 2012) via the algorithm described in Osinga and Moehlis (2010). This algorithm needs the Floquet eigenfunction corresponding to the nontrivial (in this case, stable) Floquet multiplier of the periodic orbit. Specifically, let dX/dt=F(X) be a differential equation in , and let U(t) be an exponentially stable limit cycle solution with period T. Let A(t)≔DXF(U(t)) be the linearization of F evaluated along the limit cycle. Then the Floquet eigenfunction V(t) can be defined as the unique T-periodic solution to
4.1
where V(0) is the corresponding eigenvector for the time-T map evaluated at t=0. Moreover, V(0) is the tangent vector to the isochron at that point. The algorithm collects starting points of time-nT maps that, using boundary value methods, are forced to end at , where varies between 0 and and .
Figure 2:

Panels with 10 computed isochrons, each at equidistant time intervals for . For other values of , the isochrons rotate so they are always clustered in the region of .

Figure 2:

Panels with 10 computed isochrons, each at equidistant time intervals for . For other values of , the isochrons rotate so they are always clustered in the region of .

The Floquet multiplier becomes extremely small for large values of a and b, for example, for a=10, b=0.9999, which caused numerical difficulties when computed directly using the techniques described in Doedel, Kooi, van Voorn, and Kuznetsov (2008). We avoided these difficulties by first computing the Floquet eigenfunction for smaller parameter values: a=0.5 and b=0. Then the phase of the extended system consisting of the differential equations for the periodic orbit U(t) and its eigenfunction V(t) was shifted so that the desired isochron was at t=0. Finally, the extended system was continued in a, , and T and then in b, , and T, to reach the desired values of a and b. This homotopy approach, in combination with a high number of mesh points (the AUTO-constant NTST was set to 700), allowed us to reliably find the Floquet eigenfunction and then run the algorithm from Osinga and Moehlis (2010) that computes the actual isochrons.

From Figure 2, we note that as the value of b gets closer to 1, the isochrons are more tightly packed around , reflecting the slower angular velocity at that point, and as the value of a increases, the curvature of the isochrons decreases, as expected. For a different value of , the geometry can be found by an appropriate rotation of the images in Figure 2 so that the isochrons are always clustered in the region of .

As we have recently discussed in Krogh-Madsen et al. (2012), the shape of the PRC will depend on the location of the state point following a perturbation relative to and the fiducial point of the limit cycle oscillator. Thus, depending on the parameters, the PRC might be mostly positive, mostly negative, or bimodal. These situations are illustrated in Figure 3a, which shows the PRC for a stimulus applied to u with magnitude s=0.1 and several values of .

Figure 3:

(a) The PRC for an amplitude stimulus of s=0.1 (see Figure 1) and several values of . The right-hand scale of this panel rescales the PRC by dividing by s. (b) The adjoint zu for several values of . Away from , the adjoint loses its bimodal shape. Since the magnitude of the perturbation is small, the adjoint is a good approximation to the PRC/s. In this figure, we take a=10 and b=0.99.

Figure 3:

(a) The PRC for an amplitude stimulus of s=0.1 (see Figure 1) and several values of . The right-hand scale of this panel rescales the PRC by dividing by s. (b) The adjoint zu for several values of . Away from , the adjoint loses its bimodal shape. Since the magnitude of the perturbation is small, the adjoint is a good approximation to the PRC/s. In this figure, we take a=10 and b=0.99.

When , the PRC is almost entirely negative; when , the PRC is almost entirely positive; and when , the PRC is bimodal. The reason for this behavior is clear, and it results from the accumulation of the isochrons in the neighborhood of the incipient fixed point at . Thus, despite being very close to the saddle node bifurcation, the PRC can look very sinusoidal, and in particular, it is quite far from being strictly of one sign. To give the reader some sense of what the time traces of the modified oscillator look like, we plot U(t) for three values of in Figure 4 corresponding to the traces in Figure 3. The trace of is similar to the shape of an action potential.

Figure 4:

A plot of U(t) along the limit cycle for three different values of corresponding to the top three curves in Figures 3a and 3b. Curves have been shifted for clarity, but all have the same period.

Figure 4:

A plot of U(t) along the limit cycle for three different values of corresponding to the top three curves in Figures 3a and 3b. Curves have been shifted for clarity, but all have the same period.

The trajectories in Figure 4 suggest that the shape of the action potential might also be an indicator of whether the PRC has a negative lobe near the SNIC. As noted, for example, in Rinzel and Ermentrout (1998), near a SNIC, there is sometimes a small negative region after the peak of the action potential (compare the curve for in Figure 3). Geometrically, we would expect that if the slow part of the trajectory is at a place in the phase plane that is nearly opposite the peak of the potential, as, for example, is seen in the Hindmarsh-Rose model (Hindmarsh & Rose, 1984), then it may be possible to greatly increase the negative region of the PRC near the SNIC (compare the curve for in Figure 3).

For many biological applications in which oscillators weakly interact, it is adequate to consider small-amplitude perturbations of the limit cycle. From a mathematical perspective, this limit is of interest because it is possible to carry out analytical study of the PRCs.

If the stimulation magnitude s is small, we expect that the PRC will be of the form , that is, it depends linearly on the amplitude s. Here the function is the infinitesimal response to a stimulus applied to the nth component of the limit cycle. Several authors (Izhikevich, 2007; Ermentrout & Terman, 2010) have shown that is a component in a certain vector quantity, the adjoint , which is closely related to the limit cycle. Indeed, is the gradient of the isochron map evaluated on the limit cycle (Izhikevich, 2007). The adjoint is easy to compute numerically for stable limit cycles. Specifically, Z(t) is the unique T-periodic solution to
where U(t) and A(t) are defined above at the definition of the Floquet eigenfunction in equation 4.1.
It is sometimes possible to compute Z(t) explicitly because of its close relationship with the original limit cycle. For example, consider the scalar differential equation on the circle,
where f(x+1)=f(x) and f(x)>0 and continuous. Then there exists a solution U(t) defined implicitly by
where u(0)=0, u(t+T)=u(t)+1 for all t, and . The adjoint , and since f(u)>0, the adjoint is strictly positive. For with I>0, the oscillator has a period and , which is strictly positive. Adjoints and PRCs that are nonnegative are called class I PRCs.1 PRCs that have bimodal shapes (substantial regions where there is phase advance and phase delay) are called class II PRCs.
For the Poincaré oscillator, equation 3.1, we can compute the adjoint (it is also possible to analytically compute the PRC for finite stimuli), which will be proportional to
where depends on a and q. Thus, for the simple Poincaré oscillator, the PRC is bimodal.
We now proceed to analyze the adjoint for the modified Poincaré oscillator with a SNIC bifurcation defined by equations 3.2. The adjoint in rectangular coordinates can be found in terms of the adjoint in polar coordinates using
4.2
where r=1 and are the solutions to the adjoint equation in polar form:
4.3
We can readily solve for . Then we could solve the scalar equation for zr,
4.4
where and h(t+T)=h(t), and apply the appropriate normalization criteria. To begin, we use some simple approximations. We seek periodic solutions zr to equation 4.4. Let , and . Then we can rewrite equation 4.4 as

Now it follows from standard singular perturbation theory that (indeed, set , then Y(t)=−h(t) and , which is bounded as as it is independent of ). Thus, .

Suppose that the attraction to the limit cycle is strong so that . Then , so that we obtain the PRC for the u component of the limit cycle as approximately
If we choose , then with some trigonometric identities, we get
This says that for a strongly attracting limit cycle, even if it is very close to the SNIC, the adjoint, and thus the PRC, can have a nearly sinusoidal shape. The validity of the approximation can be seen by looking at the relative sizes of the two terms that comprise zu in equation 4.2. The term zr multiplying is O(1/a), so that for a large, the approximation will hold. However, as a decreases, the term involving zr will become important. In fact, for any fixed a (big or small), as , both zr and become large, but zr is much larger (zr is O(1/(1−b2)) while ). Thus, the term with zr will eventually dominate no matter what value of a is used. Indeed, as ,
and except near . Thus, independent of a, as , the adjoint will look like , as the theory (Ermentrout, 1996; Brown et al., 2004; Canavier, 2006) predicts. The point is that the magnitude of zr is proportional to 1/[a(1−b2)], and the magnitude of is proportional to so that for a large, the zr term will dominate only for . In the remainder, we numerically explore the adjoint.
In Figure 3b, we show the adjoint as we move the location of the saddle node bifurcation at around the circle. This figure also shows the relationship between the PRC and the adjoint as
where s is the amplitude of the perturbation in the u direction. Comparison of Figure 3a, which shows the PRC for s=0.1, and Figure 3b shows that the adjoint is still a fairly good approximation to the PRC even for a large perturbation of s=0.1.

Having established that the adjoint is a reasonable approximation of the PRC, we plot the adjoints for different values of and a, for b near the SNIC value of 1. We have chosen b fairly close to the value of the saddle node. Despite this, we see in Figure 5, where , that for the strongly attracting limit cycle (a=100), the adjoint or PRC is nearly sinusoidal, in agreement with our asymptotic approximation. The approximation is predicated on the idea that the contribution of zr is small for large a. As our calculations above show, as a decreases, this contribution increases and, depending on the value of , is seen through the appearance of a nearly nonnegative PRC or adjoint as predicted by the theory for systems near a saddle node.

Figure 5:

(a) The adjoint zu as the attraction to the limit cycle decreases, where and b=0.99. As the attraction decreases, the adjoint loses its strong negative region. (b) The adjoint zu as b approaches the saddle node case, b=1, where and a=10. The nearly sinusoidal adjoint when b=0.99 nearly disappears when b=0.9999.

Figure 5:

(a) The adjoint zu as the attraction to the limit cycle decreases, where and b=0.99. As the attraction decreases, the adjoint loses its strong negative region. (b) The adjoint zu as b approaches the saddle node case, b=1, where and a=10. The nearly sinusoidal adjoint when b=0.99 nearly disappears when b=0.9999.

With a=1, the negative region is a small part near 0 phase. The shape of the PRC curves in this case can be appreciated by the geometry of the isochrons in Figure 2. When a=1, the curvature of the isochrons in the neighborhood of leads to a shift to an isochron associated with a greater phase, and hence a positive PRC curve. However, as a increases, the isochrons have less curvature, and the PRC, as approximated by the adjoint, becomes more bimodal.

As with the case for a=1, for and strong attraction, a=10, the negative lobe of the adjoint will disappear as b gets closer and closer to b=1. Figure 5b illustrates this point. For a=100, we need before the asymptotic theory becomes valid.

5.  Discussion

One of the compelling questions in the study of rhythms in the nervous system is how neurons synchronize. Weak coupling theory and other approximations have led to some general principles about what leads to neural synchronization (Van Vreeswijk, Abbott, & Ermentrout, 1994; Hansel et al., 1995; Pfeuty, Mato, Golomb, & Hansel, 2003; and Wang, 2010, for a recent review). One such principle is that nonnegative phase-resetting curves lead to synchrony with inhibitory synapses and are less likely to synchronize with excitatory synapses. Thus, one of the interesting biophysical questions is what determines the shape of the PRC. Excitatory synaptic coupling between neurons with a bimodal PRC encourages synchrony. Starting with the work in Hansel et al. (1995) and Ermentrout (1996), a common assumption has been that neurons that become rhythmic through a saddle node on a limit cycle (SNIC) will have a nonnegative PRC. This is based on theory that is strictly true only very close to the bifurcation when the stimulus will lead to a delay of the next action potential.

In this letter, we have described a very simple model that undergoes a SNIC bifurcation, and yet, quite close to the bifurcation, the PRC can be almost a perfect sine wave and thus strongly bimodal. If we regard the distance from the bifurcation as proportional to the drive, and thus the frequency of the oscillation, we see that the PRC of this model is very sensitive to a change in frequency (Gutkin, Ermentrout, & Reyes, 2005) and can switch from class I to class II over a very narrow range of parameters. Three parameters play a large role in determining the shape of the PRC; the attraction to the limit cycle a, the distance from the SNIC b, and finally the position of the stable equilibria for b>1, determined by . In this work, we find that for , the PRC will be largely positive, whereas for , the PRC will be largely negative. Consequently, the statement that systems with a SNIC will generally have a positive PRC is not true, since systems with a SNIC can also be bimodal or largely negative, as shown in Figures 3 and 5. However, it is nevertheless possible that real biological oscillators that have a largely positive PRC arise in general from an oscillation by a SNIC bifurcation.

Equations 3.2 are not suitable for detailed modeling of neural systems. However, the trace in Figure 4 is similar to the shape of an action potential. Further, the related model, in which there is uniform velocity on the unit circle (Hoppensteadt & Keener, 1982; Guevara & Glass, 1982), predicts many of the qualitative features observed when cardiac tissue is subjected to periodic stimulation. Similarly, we are optimistic that equations 3.2 might be useful in understanding the dynamics of nerve cells under periodic stimulation or when coupled to other neurons.

Our main conclusions can be summarized as follows. If the location of the incipient SNIC bifurcation is at a position of the limit cycle at which a stimulus will lead to a delay of the next action potential, the PRC will be primarily negative; if the location of the incipient SNIC bifurcation is at a position of the limit cycle at which a stimulus will lead to an advance of the next action potential, the PRC will be primarily positive; and if the location of the incipient SNIC bifurcation is at a position of the limit cycle lying between the region of advance and delay of the next action potential (e.g., for in Figure 3 here), the PRC will be primarily bimodal. However, in the last case, due to the curvature of the isochrons, the PRC may nevertheless be primarily positive in the asymptotic limit (see Figure 5).

Acknowledgments

L. G. thanks the Natural Sciences Engineering and Research Council (Canada) and B. E. thanks the National Science Foundation (USA) for financial support.

References

Achuthan
,
S.
,
Butera
,
R. J.
, &
Canavier
,
C. C.
(
2011
).
Synaptic and intrinsic determinants of the phase resetting curve for weak coupling
.
J. Comput. Neuro.
,
16
(
30
),
373
390
.
Brown
,
E.
,
Moehlis
,
J.
, &
Holmes
,
P.
(
2004
).
On the phase reduction and response dynamics of neural oscillator populations
.
Neural Computation
,
16
(
4
),
673
715
.
Canavier
,
C. C.
(
2006
).
Phase response curve
.
Scholarpedia
,
1
(
12
),
1332
.
Doedel
,
E. J.
,
Kooi
,
B. W.
,
van Voorn
,
G.A.K.
, &
Kuznetsov
,
Y. A.
(
2008
).
Continuation of connecting orbits in 3D-ODEs (I): Point-to-cycle connections
.
International Journal of Bifurcation and Chaos
,
18
(
7
),
1889
1903
.
Doedel
,
E. J.
, &
Oldeman
,
B. E.
(
2012
).
AUTO-07P: Continuation and bifurcation software for ordinary differential Equations
. http://cmvl.cs.concordia.ca/auto/.
Ermentrout
,
G. B.
(
1996
).
Type I membranes, phase resetting curves, and synchrony
.
Neural Computation
,
8
(
5
),
979
1001
.
Ermentrout
,
G. B.
,
Galán
,
R. F.
, &
Urban
,
N. N.
(
2008
).
Reliability, synchrony and noise
.
Trends in Neurosciences
,
31
(
8
),
428
434
.
Ermentrout
,
G. B.
, &
Kopell
,
N.
(
1984
).
Frequency plateaus in a chain of weakly coupled oscillators, I.
SIAM J. Math. Anal.
,
15
,
215
237
.
Ermentrout
,
G. B.
, &
Terman
,
D. H.
(
2010
).
Mathematical foundations of neuroscience
.
New York
:
Springer
.
Glass
,
L.
, &
Mackey
,
M. C.
(
1988
).
From clocks to chaos: The rhythms of life
.
Princeton, NJ
:
Princeton University Press
.
Guckenheimer
,
J.
(
1975
).
Isochrones and phaseless sets
.
Journal of Mathematical Biology
,
1
,
259
273
.
Guevara
,
M. R.
, &
Glass
,
L.
(
1982
).
Phase locking, periodic doubling bifurcations and chaos in a mathematical model of a periodically driven oscillator: A theory for the entrainment of biological oscillators and the generation of cardiac dysrhythmias
.
Journal of Mathematical Biology
,
14
,
1
23
.
Guevara
,
M. R.
,
Glass
,
L.
, &
Shrier
,
A.
(
1981
).
Phase locking, period-doubling bifurcations, and irregular dynamics in periodically stimulated cardiac cells
.
Science
,
214
,
1350
1353
.
Gutkin
,
B. S.
,
Ermentrout
,
G. B.
, &
Reyes
,
A. D.
(
2005
).
Phase-response curves give the responses of neurons to transient inputs
.
Journal of Neurophysiology
,
94
(
2
),
1623
1635
.
Hansel
,
D.
,
Mato
,
G.
, &
Meunier
,
C.
(
1995
).
Synchrony in excitatory neural networks
.
Neural Computation
,
7
(
2
),
307
337
.
Hindmarsh
,
J. L.
, &
Rose
,
R. M.
(
1984
).
A model of neuronal bursting using three coupled first order differential equations
.
Proceedings of the Royal Society of London, Series B, Biological Sciences
,
221
(
1222
),
87
102
.
,
F. C.
, &
Keener
,
J. P.
(
1982
).
Phase locking of biological clocks
.
Journal of Mathematical Biology
,
15
,
339
349
.
Izhikevich
,
E.
(
2007
).
Dynamical systems in neuroscience: The geometry of excitability and bursting
.
Cambridge, MA
:
MIT Press
.
,
T.
,
Butera
,
R.
,
Ermentrout
,
G. B.
, &
Glass
,
L.
(
2012
).
Phase resetting neural oscillators: Topological theory versus the real world
. In
R.B.N.W. Schultheiss & A. Prinz
(Eds.),
Phase response curves in neuroscience: Theory, experiment, and analysis
(pp.
33
51
).
New York
:
Springer-Verlag
.
Marella
,
S.
, &
Ermentrout
,
G. B.
(
2008
).
Class-II neurons display a higher degree of stochastic synchronization than class-I neurons
.
Physical Review E
,
77
(
4
),
041918
.
Osinga
,
H. M.
, &
Moehlis
,
J.
(
2010
).
Continuation-based computation of global isochrons
.
SIAM Journal of Applied Dynamical Systems
,
9
(
4
),
1201
1228
.
Pfeuty
,
B.
,
Mato
,
G.
,
Golomb
,
D.
, &
Hansel
,
D.
(
2003
).
Electrical synapses and synchrony: The role of intrinsic currents
.
Journal of Neuroscience
,
23
(
15
),
6280
6294
.
Rinzel
,
J.
, &
Ermentrout
,
G. B.
(
1998
).
Analysis of neural excitability and oscillations
. In
C. Koch & I. Segev (Eds.)
,
Methods in neuronal modeling: From ions to networks
(2nd ed., pp.
251
292
).
Cambridge, MA
:
MIT Press
.
Strogatz
,
S. H.
(
1994
).
Nonlinear dynamics and chaos
.
:
.
Van Vreeswijk
,
C.
,
Abbott
,
L. F.
, &
Ermentrout
,
G. B.
(
1994
).
When inhibition not excitation synchronizes neural firing
.
Journal of Computational Neuroscience
,
1
,
313
321
.
Wang
,
X.-J.
(
2010
).
Neurophysiological and computational principles of cortical rhythms in cognition
.
Physiological Reviews
,
90
(
3
),
1195
1268
.
Winfree
,
A. T.
(
2000
).
The geometry of biological time
(2nd ed.).
New York
:
Springer
.

Note

1

Sometimes this is called type I rather than class I. Since type 1 and type 0 have also been used to refer to the topology of PRCs (Krogh-Madsen et al., 2012), we use the class to classify the shapes of the PRCs and type to classify the topology.