第九章 相位重置与耦合振子

读书笔记

---

9.1 振子中的相位重置

9.1.1 背景与生物学动机

生物振子广泛存在于自然界中，理解外部扰动对其后续振荡的影响具有重要的理论和实际意义。A.T. Winfree 在20世纪60年代关于昼夜节律的开创性研究中，首次在生物学背景下提出了这一基本问题。他通过对果蝇（Drosophila melanogaster）蛹期周期性羽化现象的实验研究，揭示了外部刺激对生物时钟的深远影响。

心脏中的起搏器是一个重要的振子，Jalife 和 Antzelevitch（1979）研究了心脏组织中的起搏器活动。神经元细胞在一定条件下表现出规律性的周期性放电，理解外部刺激对这种振荡的影响具有重要的神经生物学意义。Best（1979）将 Hodgkin-Huxley 模型（神经动作电位传播的标准模型）置于周期性脉冲下进行扰动实验，展示了本章将讨论的一些重要现象。

9.1.2 果蝇羽化实验

在果蝇的蛹期发育过程中，会发生变态，最终羽化为成虫果蝇。如果将正在变态的蛹放置在完全黑暗的环境中，果蝇仍然以大约每24小时一次的频率羽化，持续约6至8小时。如果在黑暗环境中对蛹施加短暂的光脉冲，羽化的时序或相位会发生偏移。相移的大小取决于光脉冲的时序 \(T\) 和剂量 \(D\)。Winfree（1975）的实验表明，存在一个临界剂量 \(D^*\)，如果在特定的时序 \(T^*\) 施加，会导致周期性羽化行为完全消失。

9.1.3 摆的类比

考虑一个以周期 \(\omega\) 摆动的摆。如果我们给摆锤一个脉冲，可以明确地扰乱规律的周期性摆动。脉冲之后，摆最终再次表现出简谐运动，但摆锤不再在 \(t = n\omega\) 时刻到达位置 \(S\)，而是在 \(t = t_s + n\omega\) 时刻到达，其中 \(t_s\) 是某个常数——这称为相位重置。如果我们恰好在摆摆至垂直位置时给予适当大小的刺激，可以使摆完全停止振荡，这就是刺激-相位-响应空间中的奇异点。

9.1.4 极限环振子的数学描述

假设振子由向量状态变量 \(\mathbf{u}\) 描述，满足微分方程系统：

\[\frac{d\mathbf{u}}{dt} = \mathbf{f}(\mathbf{u}, \lambda) \qquad(9.1)\]

其中 \(\mathbf{f}\) 是非线性速率函数，\(\lambda\) 表示振子的参数。为简单起见，假设（9.1）描述一个只涉及2个变量 \(x\) 和 \(y\) 的极限环振子。极限环轨迹 \(\gamma\) 是相平面中的一条简单闭合曲线。通过适当的变量变换，可以将极限环转化为圆，其中振子的状态本质上由角度 \(\theta\)（相位）描述，原点在圆上的某个任意位置。极限环以速度 \(v = d\theta/dt\) 被遍历，一周遍历中 \(\theta\) 增加 \(2\pi\)。

一个简单的极限环系统例子是：

\[\frac{dr}{dt} = R(r), \quad \frac{d\theta}{dt} = \Phi(r) \qquad(9.2)\]

其中：

\[R(r) \begin{cases} >0 & 0 < r < r_0 \\ =0 & r = r_0 \\ <0 & r > r_0 \end{cases}, \quad R(r_0) = 0, \quad \Phi(r_0) = 1 \qquad(9.3)\]

这些条件意味着系统（9.2）具有唯一的吸引极限环 \(r = r_0\)，\(d\theta/dt = 1\)。

我们取最简单的非平凡极限环振子系统：

\[\frac{dr}{dt} = r(1 - r), \quad \frac{d\theta}{dt} = 1 \qquad(9.4)\]

此时相位 \(\theta(t) = \theta_0 + t\)（模 \(2\pi\)）。我们用这个系统来讨论两类基本的相位重置：1型（Type 1）和0型（Type 0）。

---

9.2 相位重置曲线

9.2.1 1型相位重置曲线

首先只扰动相位，使控制方程变为：

\[r = 1, \quad \frac{d\theta}{dt} = 1 + v(\theta, I) \qquad(9.5)\]

其中 \(v(\theta, I)\) 表示施加在角速度上的速度变化（即刺激），\(I\) 是表示脉冲强度的参数。为便于说明，取一个简单的但非平凡的 \(v\)：

\[\frac{d\theta}{dt} = 1 + I\cos 2\theta \qquad(9.6)\]

其中 \(I\) 可正可负。如果刺激在 \(t = 0\) 施加并持续时间 \(T\)，积分（9.6）给出新相位 \(\phi\) 与旧相位 \(\theta\) 的关系：

\[\int_\theta^\phi (1 + I\cos 2s)^{-1} ds = \int_0^T dt = T \qquad(9.7)\]

积分结果为：

\[|I| < 1: \quad \tan\phi = A\tan[TB + \tan^{-1}(A^{-1}\tan\theta)]\]

\[I = 1: \quad \tan\phi = 2T + \tan\theta\]

\[I = -1: \quad \tan\phi = \frac{\tan\theta}{1 - 2T\tan\theta}\]

其中：

\[A = \sqrt{\frac{|1+I|}{|1-I|}}, \quad B = |1 - I^2|^{1/2} \qquad(9.9)\]

这些公式明确给出了新相位 \(\phi\) 作为旧相位 \(\theta\) 以及刺激强度 \(I\) 和持续时间 \(T\) 的函数。

对于像（9.5）中的刺激，一个重要性质是 \(d\phi/d\theta > 0\)（对所有 \(I, T, \theta\)）。这意味着较晚的旧相位导致较晚的新相位。如果我们绘制新相位 \(\phi\) 关于旧相位 \(\theta\) 的曲线（即相位重置曲线），得到的是典型的递增曲线。

1型相位重置曲线的特征是 \(d\phi/d\theta > 0\)（对所有 \(0 < \theta \leq 2\pi\)），一个周期内的平均斜率为1——这就是名称的由来。

9.2.2 0型相位重置曲线

考虑相同的极限环（9.4），但现在施加一个将解移出极限环 \(r = 1\) 的刺激 \(I\)。具体来说，取 \(I\) 为平行于 \(y\) 轴的脉冲。从图中分析可得：

\[\rho\cos\phi = \cos\theta, \quad \rho\sin\phi + I = \sin\theta \qquad(9.10)\]

消去 \(\rho\) 得到：

\[\tan\phi = \tan\theta - \frac{I}{\cos\theta} \qquad(9.11)\]

对 \(\theta\) 求导：

\[\frac{d\phi}{d\theta} = \frac{1 - I\sin\theta}{\cos^2\theta} \begin{cases} >0 & \text{对所有 } 0 < \theta < 2\pi, \text{ 如果 } |I| < 1 \\ <0 & \text{对于 } \sin\theta > \frac{1}{I} \end{cases} \qquad(9.12)\]

当 \(|I| > 1\) 时，对于某些 \(\theta\) 范围有 \(d\phi/d\theta < 0\)。随着 \(\theta\) 取遍完整周期 \(2\pi\)，新相位 \(\phi\) 只取其中一部分值，这种曲线称为0型（Type 0）相位重置曲线。其平均斜率为0，这就是名称的来源。

0型相位重置曲线不能仅由相位刺激获得。当 \(I\) 从0增加穿过 \(I = 1\) 时，发生从1型到0型的分岔。

---

9.3 黑洞

9.3.1 相位重置的分岔

从上节分析可见，当刺激 \(I\) 从0增加时，在 \(I = 1\) 处存在相位重置类型的分岔。为清楚理解物理意义，考虑将 \(\phi = \phi(\theta, I)\) 曲面投影到 \((\theta, I)\) 平面上。对于各种 \(\phi\) 值，构造曲线：

\[I = \sin\theta - \cos\theta\tan\phi \qquad(9.13)\]

结果显示，存在两个奇异点 \(S_1\) 和 \(S_2\)，每个点都汇聚了相位 \(0 \leq \phi \leq 2\pi\) 的完整集合的常数相位曲线——在一个情况下（即 \(S_1\)）曲线以逆时针方式排列，在另一个情况下（即 \(S_2\)）以顺时针方式排列。

9.3.2 黑洞的定义与性质

假设有这样的振子，在给定相位 \(\theta\) 施加刺激 \(I\)。只要 \(|I| < 1\)，给定旧相位 \(\theta\)，可以简单地读出新相位 \(\phi\)，且结果是唯一的。对于 \(|I| > 1\)，给定旧相位 \(\theta\)，新相位同样唯一确定。但对于 \(|I| > 1\)，给定刺激 \(I\)，对于两个不同的旧相位 \(\theta\) 可以得到相同的新相位 \(\phi\)。

现在考虑在相位 \(\theta = \pi/2\) 施加特定刺激 \(I = 1\) 的情况。结果是图中的奇异点 \(S_1\)，它没有单一特定的相位 \(\phi\) 与之关联，而是整个范围 \(0 \leq \phi \leq 2\pi\)。换句话说，这种特定刺激在特定相位下的效果是不确定的。奇异点 \(S_1\) 和 \(S_2\) 是刺激-相位空间中的"黑洞"，是在该点刺激结果未知的点。

从实际角度看，生物振子受到这种刺激的结果是不可预测的。数学上，如果精确的刺激 \(I = 1\) 恰好在精确的相位 \(\theta = \pi/2\) 施加，就没有结果性的新相位 \(\phi\)。这就是单摆情况中，当摆恰好摆过垂直位置时给予正确脉冲所能发生的情况。

生物振荡器的一个关键特征是，在其旧相位-刺激空间中存在与黑洞相对应的脉冲和相位。生物振荡器最重要的应用之一是，对于这样的振子，存在一种刺激，如果在特定相位施加，可以完全消灭振荡。黑洞存在的连续性论证是：如果刺激连续增加时从1型重置转变为0型重置，那么在相位的转变值和刺激处就存在一个黑洞。

---

9.4 真实生物振子中的黑洞

9.4.1 Hodgkin-Huxley 模型中的黑洞

Best（1979）对 Hodgkin 和 Huxley（1952）的空间夹钳神经元放电模型进行了数值研究。当振子受到电压变化（这些是刺激）扰动时，他发现了预期的1型和0型相位重置曲线。在 Hodgkin-Huxley 模型中，黑洞或零空间的存在表现为从1型到0型重置曲线的转变。由于数值模拟中固有的近似，不可能确定单个奇异点，而是在奇点周围存在一个区域，即黑洞或零空间，在适当的旧相位范围内经过适当扰动后，新相位是不确定的。

9.4.2 心脏组织中的实验验证

Jalife 和 Antzelevitch（1979）在狗、猫和小牛的心脏组织上进行了类似的研究，将基本振荡施加电刺激。他们从实验中获得相位重置曲线，展示了1型和0型重置曲线。从重置曲线可以预期存在转变值或值，用于刺激持续时间，从而消灭振荡，即在内源心脏振子的零空间或黑洞中。他们确实发现了这一点。当刺激时间施加在周期的130毫秒处时，第三个轨迹中振荡完全被抑制。在狗心脏组织中，通过刺激消灭振荡也得到了证明。

9.4.3 与心脏病理的关联

Winfree 多年来一直在研究心脏骤停的可能原因及其与起搏器振子拓扑的关联。虽然心脏收缩涉及起搏器，但偏离正常的情况通常涉及循环收缩波的出现，而不是起搏机制的中断。在心室颤动的情况下，可能对黑洞或奇异点的深入理解有助于阐明这个问题。

---

9.5 耦合振子：动机与模型系统

9.5.1 耦合振子的重要性

生物振子和周期性过程在生态学、流行病学、发育生物学等领域普遍存在。在许多情况下，振子以某种方式耦合以获得所需的输出。耦合极限环振子已被广泛研究，它们共同表现出的现象范围比任何单个振子都要大得多。这里我们主要关注同步过程及其崩溃条件。

9.5.2 Belousov-Zhabotinskii 反应实验

Marek 和 Stuchl（1975）研究了将两个具有不同参数（即不同周期振荡）的 Belousov-Zhabotinskii 反应系统耦合的效果。他们通过在分离的搅拌槽反应器中进行每个反应，并通过有孔壁材料交换将它们耦合。他们观察到，如果内源振子具有几乎相同的频率，则相位差趋于常数——这称为相位锁定。然而，如果内源频率差异太大，相位锁定不持续，耦合系统在缓慢变化的相位差长间隔之间穿插着非常短的快速波动。

9.5.3 耦合振子模型

假设极限环振子是相同的，每个单独由以下方程控制：

\[\frac{dx_i}{dt} = F(x_i, y_i), \quad \frac{dy_i}{dt} = G(x_i, y_i), \quad i = 1, 2 \qquad(9.14)\]

假设（9.14）的解表现出稳定极限环行为，周期为 \(T\)：

\[x_i = X(t + \psi_i), \quad y_i = Y(t + \psi_i), \quad i = 1, 2 \qquad(9.15)\]

为便于表述弱耦合，我们考虑非维模型系统：

\[\frac{dx_1}{dt} = F(x_1, y_1) + \varepsilon\{k(x_2 - x_1) + \lambda f(x_1, y_1)\}\]

\[\frac{dy_1}{dt} = G(x_1, y_1) + \varepsilon\{k(y_2 - y_1) + \lambda g(x_1, y_1)\}\]

\[\frac{dx_2}{dt} = F(x_2, y_2) + \varepsilon k(x_1 - x_2)\]

\[\frac{dy_2}{dt} = G(x_2, y_2) + \varepsilon k(y_1 - y_2) \qquad(9.16)\]

其中 \(0 < \varepsilon \ll 1\)，\(k > 0\) 是耦合常数。当 \(\varepsilon = 0\) 时振子未耦合，这些方程简化为（9.14）。

---

9.6 振荡相位锁定：萤火虫同步

9.6.1 萤火虫同步现象

耦合生物振子的一个显著而众所周知的例子是萤火虫（Pteroptyx malaccae）大规模群体的周期性闪光同步。雄性萤火虫闪光以吸引雌性。Hanson（1978）的早期研究表明，单个萤火虫可以改变其光发射振子的相位，只要其周期在萤火虫自然周期（约0.9秒）的邻域内，就可以同步或 entrain 到闪光。如果人工光刺激的周期离自然周期太远，则无法实现 entrainment。Pteroptyx malaccae 似乎最为擅长，能够改变其频率近15%。

9.6.2 萤火虫同步的数学模型

令萤火虫振子在时刻 \(t\) 的相位为 \(\theta(t)\)，自然频率为 \(\omega\)。即，在没有外部刺激的情况下：

\[\frac{d\theta}{dt} = \omega \qquad(9.17)\]

萤火虫试图将其频率同步到外部刺激的频率。萤火虫振荡器因此相位锁定到刺激，但如果 \(\delta = 0\)，\(\phi = 0\) 是（9.22）的稳定解，则振子闪光没有相位差，因此是一致的。

同步条件由外源频率 \(\omega_e\) 满足：

\[\omega - I \leq \omega_e \leq \omega + I \qquad(9.23)\]

这就是 entrainment 的刺激频率范围。当 \(-1 \leq \delta \leq 1\) 时，相位锁定的萤火虫振荡器的维周期 \(T\) 由（9.22）得到，作为 \(\phi\) 变化 \(2\pi\) 的时间：

\[T = \frac{1}{I}\int_0^{2\pi} \frac{d\phi}{\delta - \sin\phi} = \frac{2\pi}{I(\delta^2 - 1)^{1/2}} = \frac{2\pi}{\sqrt{(\omega_e - \omega)^2 - I^2}} \qquad(9.25)\]

当 \(\delta \to \pm 1\) 时，周期变得无限大——即没有 entrainment。

---

9.7-9.10 奇异摄动分析：初步变换、变换系统、双时间展开与相位偏移方程分析

9.7.1 相位-幅度变量变换

方程（9.16）通常难以分析。由于在许多感兴趣的情况下耦合是弱的，我们利用 \(0 < \varepsilon \ll 1\) 并使用奇异摄动理论。

每个振子有其自己的极限环解，可以用闭合轨迹 \(\gamma\) 表示。我们可以用相 \(\theta\)（绕 \(\gamma\) 完整一圈时 \(\theta\) 从0到 \(T\)）和垂直于 \(\gamma\) 的扰动幅度 \(A\) 来表征。适当的变量变换为：

\[x_i = X(\theta_i) + \varepsilon A_i Y'(\theta_i), \quad y_i = Y(\theta_i) - \varepsilon A_i X'(\theta_i), \quad i = 1, 2 \qquad(9.37)\]

这里 \(X(\theta_i)\) 和 \(Y(\theta_i)\) 是 \(\theta_i\) 的 \(T\) 周期函数，\(\varepsilon A_i\) 强调轨道扰动很小。

9.8 变换后的系统

使用变量变换（9.37）在耦合系统（9.16）上，经过泰勒级数展开，得到以 \(\theta_i\) 和 \(A_i\) 表示的方程：

\[\frac{dA_1}{dt} = \Phi(\theta_1)A_1 + kU(\theta_1, \theta_2) + \lambda\phi(\theta_1) + \varepsilon\Psi_1(A, \theta) + O(\varepsilon^2) \qquad(9.45)\]

\[\frac{d\theta_1}{dt} = 1 + \varepsilon\{[\Omega(\theta_1) - \Gamma(\theta_1)]A_1 + \lambda\gamma(\theta_1) + kr(\theta_1) + kV(\theta_1, \theta_2)\} + O(\varepsilon^2)\]

类似的方程适用于 \(A_2\) 和 \(\theta_2\)。

所有函数 \(\Phi, \phi, \Gamma, \gamma, \Omega, \Psi_1, \Psi_2, r\) 都是 \(\theta_1\) 和 \(\theta_2\) 的 \(T\) 周期函数，且 \(\Phi\) 满足：

\[\int_0^T \Phi(\sigma) d\sigma < 0 \qquad(9.47)\]

这来自原始未耦合振子的极限环解是稳定的这一事实。

9.9 双时间展开

\(O(\varepsilon)\) 项在长时间 \(O(1/\varepsilon)\) 后才会有影响。这提示我们寻找 \(A_i\) 和 \(\theta_i\) 的渐近解，形式为：

\[A_i \sim {}_0A_i + \varepsilon {}_1A_i, \quad \theta_i \sim {}_0\theta_i + \varepsilon {}_1\theta_i \qquad(9.50)\]

其中 \(A\) 和 \(\theta\) 是快时间 \(t\) 和慢时间 \(\tau = \varepsilon t\) 的函数。所有时间导数：

\[\frac{d}{dt} = \frac{\partial}{\partial t} + \varepsilon\frac{\partial}{\partial\tau} \qquad(9.51)\]

经过代入和按 \(\varepsilon\) 的幂次整理，得到方程层次。经过复杂的分析，最终得到相位偏移 \(\chi = \psi_2 - \psi_1\) 的方程：

\[\frac{d\chi}{d\tau} = P(\chi) - \lambda\beta \qquad(9.69)\]

其中 \(P(\chi) = H(-\chi) - H(\chi)\)。

9.10 相位偏移方程的分析与应用

方程（9.69）中的函数 \(H(\chi)\) 和 \(H(-\chi)\) 是 \(T\) 周期函数，因此 \(P(\chi)\) 也是 \(T\) 周期函数。如果 \(\chi = 0\)，则 \(P(0) = 0\)。

相位锁定是当相位差恒定时的状态，即 \(\chi = \chi_0\) 为常数。相位锁定解由（9.69）令 \(d\chi/d\tau = 0\) 给出，即 \(P(\chi) - \lambda\beta = 0\) 的解。线性稳定性由在 \(\chi_0\) 处线性化（9.69）给出：

\[\frac{d(\chi - \chi_0)}{d\tau} \approx P'(\chi_0)(\chi - \chi_0)\]

因此如果 \(P'(\chi_0) < 0\)，\(\chi_0\) 稳定；如果 \(P'(\chi_0) > 0\)，\(\chi_0\) 不稳定。

当 \(P(\chi_m) < \lambda\beta < P(\chi_M)\) 时，有两个稳态解 \(\chi_1\) 和 \(\chi_2\)，分别是稳定和不稳定的。

当 \(\lambda\beta\) 变得小于临界值 \((\lambda\beta)_c = P(\chi_m)\) 时，稳态解消失，解演变为"节律分裂"现象：解在 \(\chi_m + nT\) 附近停留很长时间（\(\tau = O(1/\delta^2)\)），然后在短时间（\(\tau = O(1)\)）内快速变化到下一个值。这正是 Marek 和 Stuchl（1975）在 Belousov-Zhabotinskii 反应实验中观察到的现象。

---

公式汇总表

方程编号	公式名称	公式内容
(9.1)	基本振子方程	\(\frac{d\mathbf{u}}{dt} = \mathbf{f}(\mathbf{u}, \lambda)\)
(9.2)-(9.3)	极限环系统	\(\frac{dr}{dt} = R(r), \quad \frac{d\theta}{dt} = \Phi(r)\)
(9.4)	简单极限环	\(\frac{dr}{dt} = r(1-r), \quad \frac{d\theta}{dt} = 1\)
(9.5)	相位扰动方程	\(r=1, \quad \frac{d\theta}{dt} = 1 + v(\theta,I)\)
(9.6)	余弦相位速度	\(\frac{d\theta}{dt} = 1 + I\cos 2\theta\)
(9.11)	0型重置关系	\(\tan\phi = \tan\theta - \frac{I}{\cos\theta}\)
(9.13)	投影曲线方程	\(I = \sin\theta - \cos\theta\tan\phi\)
(9.14)	未耦合振子	\(\frac{dx_i}{dt} = F(x_i,y_i), \quad \frac{dy_i}{dt} = G(x_i,y_i)\)
(9.15)	周期解表示	\(x_i = X(t+\psi_i), \quad y_i = Y(t+\psi_i)\)
(9.16)	耦合振子模型	含 \(\varepsilon\) 耦合项的四方程系统
(9.17)-(9.18)	萤火虫振荡器	\(\frac{d\theta}{dt}=\omega, \quad \frac{d\theta_e}{dt}=\omega_e\)
(9.19)	同步方程	\(\frac{d\theta}{dt} = \omega + I\sin(\theta_e-\theta)\)
(9.20)	相位差方程	\(\frac{d\phi}{dt} = \omega_e - \omega - I\sin\phi\)
(9.22)	无量纲相位差	\(\frac{d\phi}{d\tau} = \delta - \sin\phi\)
(9.23)	同步条件	\(\omega - I \leq \omega_e \leq \omega + I\)
(9.25)	相锁周期	\(T = \frac{2\pi}{\sqrt{(\omega_e-\omega)^2-I^2}}\)
(9.27)	相位表示	\(x_i = X(\theta_i), \quad y_i = Y(\theta_i)\)
(9.37)	变量变换	\(x_i = X+\varepsilon A_i Y', \quad y_i = Y-\varepsilon A_i X'\)
(9.45)-(9.46)	变换后系统	\(A_i\) 和 \(\theta_i\) 的 \(O(\varepsilon)\) 方程
(9.47)	稳定性条件	\(\int_0^T \Phi(\sigma)d\sigma < 0\)
(9.50)	双时间展开	\(A_i \sim {}_0A_i + \varepsilon{}_1A_i, \quad \theta_i \sim {}_0\theta_i + \varepsilon{}_1\theta_i\)
(9.51)	时间导数分解	\(\frac{d}{dt} = \frac{\partial}{\partial t} + \varepsilon\frac{\partial}{\partial\tau}\)
(9.69)	相位偏移方程	\(\frac{d\chi}{d\tau} = P(\chi) - \lambda\beta\)
(9.71)	相位锁定条件	\(P(\chi) - \lambda\beta = 0\)
(9.72)	稳定性判据	\(P'(\chi_0) < 0\) 稳定，\(P'(\chi_0) > 0\) 不稳定
(9.73)	双稳态条件	\(P(\chi_m) < \lambda\beta < P(\chi_M)\)

---

关键概念总结

1. 相位重置：外部刺激改变振子的相位但不改变其周期。
2. 1型 vs 0型重置：1型中 \(d\phi/d\theta > 0\)，0型中平均斜率为0。
3. 黑洞：刺激-相位空间中的奇异点，在该点刺激结果不确定或振荡被消灭。
4. 相位锁定：两个耦合振子的相位差趋于常数。
5. 节律分裂：当参数变化时，从相位锁定突然转变为准周期行为。
6. 双时间展开：处理弱耦合振子的奇异摄动技术。

---

本笔记基于 J.D. Murray《Mathematical Biology I》第九章内容整理。