Phase Model

• Eugene M. Izhikevich, Editor-in-Chief of Scholarpedia, the peer-reviewed open-access encyclopedia

• Dr. Bard Ermentrout, Dept of Mathematics, Univ Pittsburgh, Pittsburgh PA

Coupled oscillators interact via mutual adjustment of their amplitudes and phases. When coupling is weak, amplitudes are relatively constant and the interactions could be described by phase models.

Figure 1: Phase of oscillation (denoted by \(\vartheta\) in the rest of the article) of the FitzHugh-Nagumo model with I=0.5. The phase here is measured in units of time. The zero-phase point \(x_0\) is chosen to correspond to the peak of the potential (the peak of spike).

Phase of oscillation

Many physical, chemical, and biological systems can produce rhythmic oscillations (Winfree 2001), which can be represented mathematically by a nonlinear dynamical system \[ x' = f(x) \] having a periodic orbit \(\gamma\ .\) Let \(x_0\) be an arbitrary point on \(\gamma\ ,\) then any other point on the periodic orbit can be characterized by the time, \(\vartheta\ ,\) since the last passing of \(x_0\ ;\) see Figure 1. The variable \(\vartheta\) is called phase of oscillation, and it is bounded by the period of oscillation \(T\ .\) The phase is often normalized by \(T\) or \(T/2\pi\ ,\) so that it is bounded by \(1\) or \(2\pi\ ,\) respectively.

The phase of oscillation can also be defined outside \(\gamma\) using the notion of isochrons. The change of variables \(x(t) = \gamma(\vartheta(t))\) transforms the nonlinear system in a neighborhood of \(\gamma\) into an equivalent but simpler phase model \[ \vartheta' = 1. \] Such a change of variables removes the amplitude but saves the phase of oscillation (Ermentrout 1986). It is often convenient to assume that the phase \(\vartheta\) is defined on the unit circle \(S^1\ .\)

The phase of oscillation could also be defined for chaotic systems (Pikovsky et al. 2001) using the observation that many chaotic attractors, as in Rossler oscillator, look like smeared limit cycles.

Figure 2: Examples of function \(Q=(Q_1, Q_2)\) for Andronov-Hopf oscillator \(z'=(1+i)z-z|z|^2\) and van der Pol oscillator \(x'=x-x^3-y\ ,\) \(y'=x\ ;\) modified from Izhikevich (2007).

Weak forcing

The same change of variables transforms a weakly forced oscillator \[ x' = f(x) + \varepsilon s(t) \] into the phase model of the form \[ \vartheta' = 1 + \varepsilon Q(\vartheta) \cdot s(t), \] where, the term \(\varepsilon s(t)\) denotes a weak time-dependent (and possibly \(x\)-dependent) input, e.g., from other oscillators in a network; the dot, "\(\cdot\)", denotes the scalar (dot) product of two vectors; The function \(Q(\vartheta)\ ,\) illustrated in the Figure 2, is called linear response function, sensitivity function, or infinitesimal PRC. This impulse response function, kernel, or Green's function can then be convolved with the actual input received by each oscillator in order to compute the total phase resetting of the oscillator received over one cycle of the network oscillation. \(Q(\vartheta)\) satisfies three equivalent conditions (Izhikevich 2007):

• Winfree\[Q(\vartheta)\] is a normalized phase response curve (PRC) to infinitesimal pulsed perturbations. That is, one measures PRC of the oscillator \(x'=f(x)\) by perturbing each component of state vector \(x\) with brief pulses of small amplitude with area of the pulse \(A\ ,\) and then takes \(Q=\)PRC\(/A\) in the limit \(A\rightarrow 0 \ .\)
• Kuramoto\[Q(\vartheta) = \]grad\(\Theta(x)\ ,\) where \(\Theta(x)\) is the isochron function defined in a neighborhood of the periodic orbit \(\gamma\ .\) That is, one starts from every point \(x\) in a neighborhood of \(\gamma\) and determines its asymptotic phase, \(\Theta(x)\ ,\) relative to the phase of the solution starting with \(x_0\ .\)
• Malkin\[Q(\vartheta)\] is the solution to the adjoint problem \( dQ/d\vartheta=-\{Df(\gamma(\vartheta))\}^\top Q\;, \) with the normalization condition \(Q(\vartheta) \cdot f(\gamma(\vartheta))=1\) for any \(\vartheta\ .\) That is, one determines the Jacobian matrix \(Df\) along the periodic orbit and then solves, usually numerically, the adjoint problem.

Malkin's condition, though least intuitive, is the most useful in applications.

Examples of reduction

The infinitesimal PRC function \(Q(\vartheta)\) can be found analytically in a few simple cases.

Phase oscillators

A nonlinear phase oscillator \(\dot{x} = f(x)\) with periodic phase variable \(x \in [0, 1]\) and \(f>0\) has \(Q(\vartheta) = 1/f(\gamma(\vartheta))\ .\) Indeed, the function can be found from Malkin's normalization condition \(Q(\vartheta) \cdot f(\gamma(\vartheta))=1\ .\)

SNIC oscillators

A system near saddle-node on invariant circle (SNIC) bifurcation has \(Q(\vartheta)\) proportional to \(1-\cos^2 \vartheta\ .\)

Andronov-Hopf oscillators

A system near supercritical Andronov-Hopf bifurcation has \(Q(\vartheta)\) proportional to \(\sin (\vartheta-\psi)\ ,\) where \(\psi\) is a constant phase shift.

Other interesting cases

Izhikevich (2000) derived the phase model for weakly coupled relaxation oscillators. Brown et al. (2004) consider other interesting cases, including homoclinic oscillators. Coupled bursters are considered by Izhikevich (2007). Pulse coupled oscillators provide many other analytically solvable examples.

Weakly coupled oscillators

Let us treat \(s(t)\) in \(x' = f(x) + \varepsilon s(t)\) as the input from the network, and consider weakly coupled oscillators \[ x_i' = f_i(x_i) + \varepsilon \sum_{j=1}^n g_{ij}(x_i, x_j). \] The corresponding phase model \[ \vartheta_i' = 1 + \varepsilon \, Q_i(\vartheta_i) \cdot \sum_{j=1}^n g_{ij}(\gamma_i(\vartheta_i), \gamma_j(\vartheta_j)), \] has the form \[ \vartheta_i' = 1 + \varepsilon \sum_{j=1}^n h_{ij}(\vartheta_i, \vartheta_j), \] where \(h_{ij} = Q_i g_{ij}\) describes the influence of phase of the j-th oscillator on the i-th oscillator, and each \(\vartheta_i\) has its own period \(T_i\ .\)

The implicit assumption of weak coupling is that the relative position (phase) of the oscillators changes slowly with respect to their motion around the limit cycle (absolute phase). This implies slow convergence to a steady state phase locking. If the coupling is not sufficiently weak but is pulsatile in nature, the methods for pulse-coupled oscillators can be utilized, otherwise there are no general methods.

Phase model

Introducing phase deviation variables \(\vartheta_i = t + \varphi_i\ ,\) one can transform the system above into the form \[ \varphi_i' = \varepsilon \sum_{j=1}^n h_{ij}(t + \varphi_i, \ t + \varphi_j). \] This system can be averaged to the phase model \[ \varphi_i' = \varepsilon \sum_{j=1}^n H_{ij}(\varphi_j-\varphi_i)\;, \] where each function \[ H_{ij}(\chi) = \lim_{T\rightarrow \infty} \frac{1}{T} \int_0^T h_{ij}(t, \ t +\chi)\, dt \] describes the interaction between oscillators. This function is just a constant unless the oscillators have nearly resonant periods, i.e., the ratio \(T_i/T_j\) is \(\varepsilon\)-close to a low-order rational number \(p/q\) (\(p+q\) is small). Since the dynamics of two coupled non-resonant oscillators is described by an uncoupled phase model (\(H=\)const), such oscillators cannot phase lock. That is, the phase of one of them cannot change the phase of the other one even on the long time scale of order \(1/\varepsilon\ .\)

Computational neuroscience provides an important application of phase models. In this case, the state variables \(x_i\) and \(x_j\) describe activities of the post-synaptic (forced) and pre-synaptic (forcing) periodically firing neurons, and the function \(g_{ij}\) describes the time course of synaptic input. The phase variables \(\varphi_i\) and \(\varphi_j\) describe the timings of firings of the neurons, and the function \(H_{ij}\) describes normalized phase resetting curve (Netoff et al 2005).

Analysis

Two coupled oscillators

Consider two mutually coupled oscillators with nearly identical periods \[ \varphi'_1 = 1 + \varepsilon \omega_1 + \varepsilon H_{12}(\varphi_2-\varphi_1) \] \[ \varphi'_2 = 1 + \varepsilon \omega_2 + \varepsilon H_{21}(\varphi_1-\varphi_2)\;, \] where \(\omega_i = H_{ii}(0)\) are small frequency deviations. Let \(\chi = \varphi_2-\varphi_1\) denote the phase difference between the oscillators, then \[ \chi' = \varepsilon \omega + \varepsilon H(\chi), \]

Figure 3: Examples of the connection function H; modified from Izhikevich (2007). The units on the y-axis are "phase/time".

where \[ \omega = \omega_2 - \omega_1 \] and \( H(\chi) = H_{21}(-\chi)-H_{12}(\chi), \) is the frequency mismatch and the anti-symmetric part of the coupling, respectively, illustrated in the Figure 3, dashed curves. A stable equilibrium of this system corresponds to a stable limit cycle of the phase model.

All equilibria of this system are solutions to \(H(\chi) = -\omega\ ,\) provided that \(\omega\) is small enough. Geometrically, the equilibria are intersections of the horizontal line \(-\omega\) with the graph of \(H\ .\) They are stable if the slope of the graph is negative at the intersection, i.e., \(H'(\chi)<0\ .\) If oscillators are identical, then \(H(\chi)\) is an odd function (i.e., \(H(-\chi)=-H(\chi)\)), and \(\chi=0\) and \(\chi=\pi\) are always equilibria, possibly unstable, corresponding to the in-phase and anti-phase synchronized solutions. The in-phase synchronization of coupled oscillators in the figure is stable because the slope of \(H\) (dashed curves) is negative at \(\chi=0\ .\) The max and min values of the function \(H\) determine the tolerance of the network to the frequency mismatch \(\omega\ ,\) since there are no equilibria outside this range.

Chains of oscillators

The behavior of chains of phase models is considerably more complex than that of pairs, even for nearest neighbor coupling. The reason for this is that when coupling is local, oscillators at the ends get different inputs from those in the middle so that phase locking may not even exist. However, in a large class of models, chains can be analyzed either by direct calculation or by letting the size of the chains tend to infinity. In the former case, Cohen et al. (1982) examined a linear chain of nearest neighbor oscillators with a frequency gradient: \[ \theta_i' = \omega_i + \sin (\theta_{i+1}-\theta_i) + \sin(\theta_{i-1}-\theta_i). \] As long as the differences in the frequencies are small enough, there will be a phase-locked solution. Interestingly, if the length of the chain is \( N \) and the frequency gradient is linear with slope \(b\) then \(b = O(1/N)\) as \(N\to\infty\ .\) That is, nearest neighbor chains can support very small gradients when the coupling is sinusoidal (and, in fact, any odd periodic function). However, if the coupling function contains any even components (that is, replace \( sin\theta \) with \( \sin (\theta+\beta) \ ,\) then frequency gradients as that are \( O(1) \) can be supported in nearest neighbor chains of coupled phase oscillators. Kopell & Ermentrout (1986,1990) derived a set of continuum equations from which general phase-locked solutions could be found.

Networks of oscillators that are not arranged in a ring have "boundaries" that can lead to patterns of phase difference that look like waves. If the coupling is isotropic, the waves take the form of one-dimensional target waves, either originating at the center and propagating symmetrically to the edges, or starting at the edges and propagating to the center. For example, an array of nearest neighbor-coupled oscillators of the form \[ \theta_i' = \sin(\theta_{i+1}-\theta_i - \alpha) + \sin(\theta_{i-1}-\theta_i -\alpha) \quad i=0,\ldots,N \] with no "i-1" term for oscillator 0 and no "i+1" term for oscillator N, will generate a wave with phase which are roughly of the form \( \theta_i = \alpha |N/2-i|. \) With anisotropic coupling, the waves are almost straight lines, \( \theta_i = \pm \alpha i. \) These results were proven in a series of papers by Ermentrout and Kopell.

Linear arrays of oscillators

Now consider a network of \(n>2\) weakly all to all coupled oscillators. To determine the existence and stability of synchronized states in the network, we need to study equilibria of the corresponding phase model \[ \varphi_i' = \varepsilon \omega_i + \varepsilon \sum_{j\neq i}^n H_{ij}(\varphi_j-\varphi_i). \] Existence of one equilibrium of the phase model above implies the existence of the entire circular family of equilibria, since translation of all \(\varphi_i\) by a constant phase shift does not change the phase differences \(\varphi_i-\varphi_j\) and hence the form of the phase model. This family corresponds to a periodic orbit, on which all oscillators have equal frequencies and constant phase shifts, i.e., they are synchronized, possibly out-of-phase.

Vector \(\phi=(\phi_1,\ldots,\phi_n)\) is an equilibrium when \[ 0 = \omega_i + \sum_{j\neq1}^n H_{ij}(\phi_j-\phi_i) \] for all \(i\ .\) It is stable when all eigenvalues of the linearization matrix (Jacobian) at \(\phi\) have negative real parts, except one zero eigenvalue corresponding to the eigenvector along the circular family of equilibria (\(\phi\) plus a phase shift is a solution too since the phase shifts \(\phi_j-\phi_i\) are not affected).

In general, determining the stability of equilibria is a difficult problem. Ermentrout (1992) found a simple sufficient condition. If

• \(a_{ij} = H_{ij}'(\phi_j-\phi_i) \geq 0\ ,\) and
• the directed graph defined by the matrix \(a = (a_{ij})\) is connected, (i.e., each oscillator is influenced, possibly indirectly, by every other oscillator),

then the equilibrium \(\phi\) is neutrally stable, and the corresponding limit cycle \(x(t+\phi)\) of the phase model is asymptotically stable.

Another sufficient condition was found by Hoppensteadt and Izhikevich (1997). If the phase model satisfies

• \(\omega_1=\cdots = \omega_n = \omega\) (identical frequencies)
• \(H_{ij}(-\chi) = - H_{ji}(\chi)\) (pair-wise odd coupling)

for all \(i\) and \(j\ ,\) then the network dynamics converge to a limit cycle. On the cycle, all oscillators have equal frequencies \(1+\varepsilon\omega\) and constant phase deviations. The proof follows from the observation that the phase model is a gradient system in a rotating coordinate system.

2D Arrays of oscillators

Two-dimensional arrays provide a much richer class of dynamics. Consider the two-dimensional analogues of the one-dimensional chain of nearest-neighbor coupled oscillators: \[ \theta_{i,j}' = H_N(\theta_{i+1,j}-\theta_{i,j}) + H_S(\theta_{i-1,j}-\theta_{i,j}) + H_E(\theta_{i,j+1}-\theta_{i,j}) + H_W(\theta_{i,j-1}-\theta_{i,j}) \] It can be shown that the patterns of phase-locked activity are of the form \[ \theta_{i,j} = \Psi_i + \Phi_j \] where \( \Psi_i,\Phi_j \) are solutions to the one-dimensional chains using \( H_{N,S} \) or \( H_{E,W} \) respectively.

Beyond these cross product patterns, it is possible to find non-trivial spatial patterns which are not a consequence of boundary effects. For example, consider the simple sinusoidal array \[ \theta_{i,j}' = \omega + \sum_{k,l} \sin(\theta_{k,l}-\theta_{i,j}) \] where \( (k,l) \) are the 2,3, or 4 nearest neighbors of \( (i,j) \) in a \( 2N \times 2N \) array of oscillators. Paullet and Ermentrout proved that there is a stable rotating wave solution to this in addition to the stable synchronous solution. For example, when N=2, the pattern of phases is: \[ \begin{array}{cccc} 0 & x & \pi/2-x & \pi/2 \\ -x & 0 & \pi/2 & \pi/2+x \\ 3\pi/2+x & 3\pi/2 & \pi & \pi-x \\ 3\pi/2 & 3\pi/2-x & \pi+x & \pi \end{array} \] where \( \cos 2x = 2 \sin x. \) This generalizes to any odd interaction function. For non-odd interactions, the rotating waves have a twist and look more like classic spiral waves.

References

• Brown E., Moehlis J., and Holmes P. (2004) On the phase reduction and response dynamics of neural oscillator populations. Neural Computation, 16:673-715.
• Cohen, A.H., Holmes, P.J. and Rand, R.H., (1982) The nature of the coupling between segmental oscillators of the lamprey spinal generator for locomotion: A mathematical model, J. Mathematical Biology 13:345-369.
• Ermentrout, G. B. (1986) Losing amplitude and saving phase. Lecture Notes in Biomath., 66, Springer, Berlin-New York.
• Ermentrout G. B. (1992) Stable periodic solutions to discrete and continuum arrays of weakly coupled nonlinear oscillators. SIAM Journal on Applied Mathematics 52:1665-1687.
• Glass L. and Mackey M.C. (1988) From Clocks to Chaos. Princeton University Press.
• Hoppensteadt F.C. and Izhikevich E.M. (1997) Weakly Connected Neural Networks. Springer-Verlag, NY
• Izhikevich E.M. (2007) Dynamical Systems in Neuroscience: The Geometry of Excitability and Bursting. The MIT Press.
• Izhikevich E.M. (2000) Phase Equations For Relaxation Oscillators. SIAM Journal on Applied Mathematics, 60:1789-1805
• Kopell, N. and Ermentrout G.B. (1986) Symmetry and phaselocking in chains of weakly coupled oscillators. Comm. Pure Appl. Math. 39: 623-660.
• Kopell N. and Ermentrout G.B. (1990) Phase transitions and other phenomena in chains of coupled oscillators. SIAM J. Appl. Math. 50:1014-1052
• Kuramoto Y. (1984) Chemical Oscillations, Waves, and Turbulence. Springer-Verlag, New York.
• Netoff T. I., Acker C. D., Bettencourt J.C. and White J.A. (2005) Beyond two-cell networks: Experimental measurement of neuronal responses to multiple synaptic inputs. J. Comp. Neurosci. 18:287-295
• Pikovsky A., Rosenblum M., Kurths J. (2001) Synchronization: A Universal Concept in Nonlinear Science. CUP, Cambridge.
• Paullet J. E. and Ermentrout G. B. (1994) Stable rotating waves in two-dimensional discrete active media. SIAM J. Appl. Math. 54: 1720-1744
• Ren L. and Ermentrout G. B. (1998) Monotonicity of phaselocked solutions in chains and arrays of nearest-neighbor coupled oscillators. SIAM J. Math. Anal. 29:208-234
• Winfree A. (2001) The Geometry of Biological Time. Springer-Verlag, New York, second edition.

Internal references

• Yuri A. Kuznetsov (2006) Andronov-Hopf bifurcation. Scholarpedia, 1(10):1858.

• John W. Milnor (2006) Attractor. Scholarpedia, 1(11):1815.

• John Guckenheimer (2007) Bifurcation. Scholarpedia, 2(6):1517.

• Jonathan E. Rubin (2007) Burst synchronization. Scholarpedia, 2(10):1666.

• James Meiss (2007) Dynamical systems. Scholarpedia, 2(2):1629.

• Eugene M. Izhikevich (2007) Equilibrium. Scholarpedia, 2(10):2014.

• Bard Ermentrout (2008) Ermentrout-Kopell canonical model. Scholarpedia, 3(3):1398.

• Eugene M. Izhikevich and Richard FitzHugh (2006) FitzHugh-Nagumo model. Scholarpedia, 1(9):1349.

• Kresimir Josic, Eric T. Shea-Brown, Jeff Moehlis (2006) Isochron. Scholarpedia, 1(8):1361.

• Rodolfo Llinas (2008) Neuron. Scholarpedia, 3(8):1490.

• Jeff Moehlis, Kresimir Josic, Eric T. Shea-Brown (2006) Periodic orbit. Scholarpedia, 1(7):1358.

• Carmen C. Canavier (2006) Phase response curve. Scholarpedia, 1(12):1332.

• Carmen C. Canavier and Srisairam Achuthan (2007) Pulse coupled oscillators. Scholarpedia, 2(4):1331.

• Philip Holmes and Eric T. Shea-Brown (2006) Stability. Scholarpedia, 1(10):1838.

• Emmanuil E. Shnol (2007) Stability of equilibria. Scholarpedia, 2(3):2770.

• Arkady Pikovsky and Michael Rosenblum (2007) Synchronization. Scholarpedia, 2(12):1459.

• Frank Hoppensteadt (2006) Voltage-controlled oscillations in neurons. Scholarpedia, 1(11):1599.

See Also

Ermentrout-Kopell canonical model, Isochron, Kuramoto model, Periodic orbit, Phase response curve, Pulse coupled oscillators, Relaxation oscillator, Synchronization, Voltage-controlled oscillations in neurons

Sponsored by: Eugene M. Izhikevich, Editor-in-Chief of Scholarpedia, the peer-reviewed open-access encyclopedia

Reviewed by: Dr. Kresimir Josic, University of Houston, Houston, Texas

Reviewed by: Dr. Sorinel Oprisan, College of Charleston

Reviewed by: Carmen C. Canavier, Louisiana State University Health Sciences Ctr, New Orleans, Louisiana

Accepted on: 2008-09-30 00:33:45 GMT