第2章 生物运动力学：行走与奔跑
原课程: Robots That Learn, UC Berkeley CS294-291, Spring 2026 本章基于 Lecture 2A: The Neuromechanics of Locomotion — Walking and Running 核心论文: Pavan Ramdya & Auke Jan Ijspeert, The neuromechanics of animal locomotion: From biology to robotics and back, Science Robotics, Vol. 8, eadg0279 (2023) 补充参考: A.J. Ijspeert et al., Dynamical movement primitives: Learning attractor models for motor behaviors, Neural Computation 25(2):328-373 (2013)
一、本章概要
1.1 核心问题
一只家蝇 (Drosophila melanogaster) 仅凭约 \(10^4\) 个腹神经索神经元，就能在垂直墙壁上行走、倒挂在天花板上奔跑，在遭遇扰动时瞬间恢复稳定姿态。一只猎豹能以超过 100 km/h 的速度奔跑，每一步都涉及数十块肌肉的精确协调激活。这些惊人的运动能力与当前最先进的仿人机器人形成鲜明对比——后者即使在平坦地面上行走，也需要兆瓦级的计算资源来实时求解非线性模型预测控制问题。
核心问题: 动物运动控制的根本原理是什么？这些原理能否用数学语言精确描述，并通过代码实现为机器人控制器？
本章从第一性原理出发，建立从神经力学 (neuromechanics) 到计算模型的完整数学框架。我们将深入到五个核心控制原则，推导中枢模式发生器 (CPG)、动态运动基元 (DMP) 和 Tegotae 反馈的数学公式，并给出可运行的 Python 实现。
1.2 为什么从生物学开始？
传统机器人控制依赖于精确的逆动力学求解——你需要知道每个关节的质量、惯量、摩擦系数，然后求解非线性微分方程来规划轨迹。这是第 3 章的内容。然而，动物显然不这么做。动物的大脑不下发关节力矩指令，而是下发低维指令（”加速”、”左转”、”跳跃”），由脊髓/腹神经索中的低层回路自动协调产生节律运动。这一架构的关键优势是：
- 鲁棒性: 对未建模的扰动有内在的弹性恢复力
- 效率: 低维控制信号 + 被动动力学 = 远低于传统方法所需的计算量
- 泛化性: 同一套 CPG 网络可以通过调节幅度和频率生成行走、小跑、奔跑的不同步态
1.3 五大生物运动控制原则
Ramdya & Ijspeert (2023) 总结了跨越脊椎动物和无脊椎动物的数十年比较研究，提炼出五条高度保守的运动控制组织原则。尽管小鼠脊髓约有 \(10^6\) 个神经元而果蝇腹神经索仅有约 \(10^4\) 个神经元——相差两个数量级——这些原则仍然深具同源性 (deep homology) [Strausfeld & Hirth, Science, 2013]。
| 原则 | 英文术语 | 定义 | 神经轴位置 |
|---|---|---|---|
| 下行驱动 | Descending Drive | 大脑通过少量下行神经元发放低维指令，激活脊髓/腹神经索中的运动回路 | 脑 \(\to\) SC/VNC |
| 中枢模式发生器 | Central Pattern Generators (CPGs) | 无需节律输入即可产生周期性输出的神经元或神经回路 | SC/VNC |
| 感觉反馈 | Sensory Feedback | 机械感受器和本体感受器将身体状态反馈给运动回路，实现反射控制 | 身体 \(\to\) SC/VNC |
| 上行投射 | Ascending Projections | 运动回路将自身活动状态报告给大脑，使大脑获知正在进行的行为 | SC/VNC \(\to\) 脑 |
| 具身化 | Embodiment | 主动神经控制经由被动生物力学（组织粘弹性、骨骼约束、环境物理）过滤后才产生可观测的运动学 | 全身 |
这五个原则不是孤立存在的——它们嵌套在复杂的多层闭环控制中，构成双重控制回路：
脑 (下行驱动/上行感知)
↕ ↕
SC/VNC (CPG + 反射回路)
↕ ↕
本体感觉/机械感受 ← 身体 (具身化)
↕
环境物理交互
1.4 与前后章节的逻辑关系
Ch01 导论 ──→ Ch02 生物运动力学 ──→ Ch03 机器人机构学
(为什么学) (生物如何动) (机器如何动)
本章的 CPG/DMP 框架是后续章节的基石:
- Ch03 (机器人机构学): 动物的可变肌肉刚度 vs. 机器人的确定性运动学链——理解差异从本章开始
- Ch07 (运动控制的发展视角): CPG 和运动基元是发展认知中”运动词汇”的神经基础
- Ch09 (计算神经科学): 前向/逆向模型的生物学实现直接基于 CPG-反射交互
- Ch14 (运动控制案例研究): 四足机器人的 RMA 控制器本质上是用 RL 学出了一个 CPG 调制策略
二、核心概念与数学基础
2.1 中枢模式发生器 (CPG) 的数学理论
2.1.1 定义与生物学动机
定义 2.1 (中枢模式发生器, CPG)：一个 CPG 是一个动力系统
\[\dot{\mathbf{x}} = f(\mathbf{x}, \mathbf{u}; \boldsymbol{\theta}), \quad \mathbf{y} = h(\mathbf{x})\]其中 \(\mathbf{x} \in \mathbb{R}^n\) 是内部状态（可理解为神经元膜电位），\(\mathbf{u} \in \mathbb{R}^m\) 是外部输入（下行指令 + 感觉反馈），\(\boldsymbol{\theta}\) 是参数集，\(\mathbf{y}\) 是输出（可理解为运动指令），满足以下性质：当 \(\mathbf{u} = \mathbf{0}\) (或 \(\mathbf{u}\) 为常数) 时，系统存在一个稳定的极限环——即 \(\exists \mathbf{x}^*(\cdot): \mathbf{x}^*(t + T) = \mathbf{x}^*(t)\) 对所有 \(t\) 成立，其中 \(T > 0\) 是周期，且该极限环是渐近稳定的。
“中央” (central) 意味着节律产生于回路内部，不需要节律输入；”模式发生” (pattern generation) 意味着输出是结构化的时空模式。生物 CPG 的经典实验证据来自去传入实验——切断所有感觉反馈后，动物仍能产生节律运动模式 [Grillner & Wallen, Annu. Rev. Neurosci., 1985]。
2.1.2 Matsuoka 振荡器：神经元层面的 CPG 模型
Matsuoka 振荡器 [Matsuoka, Biol. Cybern., 1985, 1987] 是应用最广泛的 CPG 数学模型，因为它同时捕捉了三个关键生物物理特性：(1) 适应 (adaptation) — 持续输入下频率自适应；(2) 夹带 (entrainment) — 被外部节律信号同步；(3) 互抑制耦合 (reciprocal inhibition) — 成对振荡器自然地反相振荡。
Matsuoka 振荡器由两个相互抑制的”神经元”组成 (flexor, \(\phi_f\) 和 extensor, \(\phi_e\))。每个神经元的动力学由膜电位 \(u_i\)、适应变量 \(v_i\) 和输出 \(y_i\) 描述：
膜电位方程:
\[\tau_u \dot{u}_i = -u_i - w_{ji}y_j - \beta v_i + s_i + \text{feed}_i\]适应方程:
\[\tau_v \dot{v}_i = -v_i + y_i\]输出（发放率）:
\[y_i = \max(0, u_i)\]其中 \(i \in \{f, e\}\) (\(f\) = flexor, \(e\) = extensor, \(j\) 表示另一个神经元):
| 参数 | 含义 | 典型范围 | 生物对应 |
|---|---|---|---|
| \(\tau_u\) | 膜电位时间常数 | 10–50 ms | 神经膜电容充电 |
| \(\tau_v\) | 适应时间常数 | 100–1000 ms | 钙依赖钾通道 |
| \(w_{ji}\) | 互抑制权重 | 2–5 | 抑制性突触强度 |
| \(\beta\) | 适应增益 | 2–5 | 后超极化电流 |
| \(s_i\) | 紧张性驱动输入 | 0–10 | 下行/感觉输入 |
| \(\text{feed}_i\) | 外部反馈信号 | 可变 | 来自 Tegotae 规则 |
振荡存在的必要条件: 在对称参数下 (\(w_{fe} = w_{ef} = w\), \(s_f = s_e = s\))，稳定极限环存在要求:
\[\sqrt{w_{fe}w_{ef}} > 1 + \frac{\tau_u}{\tau_v}\]这来自系统在不动点处的特征值分析。直观上，条件说明：(1) 互抑制强度 \(w\) 必须足够大（超过 \(1 + \tau_u/\tau_v\)），以产生”胜者为王”的竞争动力学；(2) 适应时间常数 \(\tau_v\) 必须显著大于膜电位时间常数 \(\tau_u\)，产生快-慢系统的弛豫振荡。
相平面分析: 考虑对称情况。系统在 \((u_f, u_e)\) 平面上的零倾线为:
\(\dot{u}_f = 0: \; u_f = -w y_e - \beta v_f + s\) \(\dot{u}_e = 0: \; u_e = -w y_f - \beta v_e + s\)
慢变量 \(v_f, v_e\) 使这些零倾线随时间漂移，从而产生振荡。当代表点在相空间中穿越激活和不激活区域时，输出交替出现。
2.1.3 相位振荡器：简化描述与耦合
当 Matsuoka 振荡器的内在动力学不是关注重点，而主要关心多个振荡器之间的同步时，可以将每个 CPG 单元约化为相位振荡器。这种简化将 \(n\) 维非线性系统的极限环约化为 1 维相位变量。
定义 2.2 (相位约化)：设 CPG 系统有一个极限环 \(\mathbf{x}^*(\phi)\), \(\phi \in [0, 2\pi)\), 其中相位沿极限环匀速演化: \(\dot{\phi} = \omega\) (固有频率)。
第 \(k\) 个肢体的相位振荡器动力学为：
\[\dot{\phi}_k = \omega_k + \sum_{j \neq k} \Gamma_{kj}(\phi_j - \phi_k) + \text{sensory}_k\]其中 \(\omega_k\) 是固有频率，\(\Gamma_{kj}(\Delta\phi)\) 是相位耦合函数，\(\text{sensory}_k\) 是感觉反馈项。
Kuramoto 模型形式（最简耦合）:
\[\dot{\phi}_k = \omega_k + \frac{K}{N}\sum_{j=1}^{N} \sin(\phi_j - \phi_k - \Delta_{kj})\]其中 \(K\) 是耦合强度，\(\Delta_{kj}\) 是期望相位差，\(N\) 是振荡器总数。
步态作为相位差的稳定解: 对于四足动物，定义期望相位差模式：
| 步态 | 相位差 \((\phi_{LF}, \phi_{RF}, \phi_{LH}, \phi_{RH})\) | 占空比 | 速度 | 稳定性 |
|---|---|---|---|---|
| 行走 (walk) | \((0, \frac{\pi}{2}, \pi, \frac{3\pi}{2})\) | \(>0.5\) | 低 | 静态稳定 |
| 小跑 (trot) | \((0, \pi, \pi, 0)\) | \(\approx 0.5\) | 中 | 动态稳定 |
| 疾驰 (gallop) | \((0, \approx 0.6\pi, \pi, \approx 1.6\pi)\) | \(<0.5\) | 高 | 动态稳定 |
这些相位差模式是耦合函数 \(\Gamma_{kj}\) 的稳定不动点。增大 \(\omega_k\) 时，相位差模式发生分岔，出现自发的步态转换——在猫的电刺激实验 [Shik et al., Biophysics, 1966] 和四足机器人实验 [Owaki & Ishiguro, Sci. Rep., 2017] 中都得到了验证。
2.1.4 CPG 输出的映射：从相位到关节角
CPG 输出必须映射为关节角度指令。常见方法：
方法 1: 基于相位的分段节律函数
对于第 \(k\) 个肢体的第 \(l\) 个关节:
\[\theta_{kl}^{\text{des}}(\phi_k) = \theta_{kl}^{\text{offset}} + A_{kl} \cdot f_{kl}(\phi_k)\] \[f_{kl}(\phi_k) = \begin{cases} g_{\text{swing}}(\phi_k), & 0 \le \phi_k < 2\pi \cdot d \\ g_{\text{stance}}(\phi_k), & 2\pi \cdot d \le \phi_k < 2\pi \end{cases}\]其中 \(d\) 是占空比 (duty factor)。
方法 2: 连续可微基函数 (用于梯度优化):
\[\theta_{kl}^{\text{des}}(\phi_k) = \sum_{m=1}^{M} w_{klm} \cdot \psi_m(\phi_k)\]其中 \(\psi_m\) 是傅里叶基或 B-spline 基函数，\(w_{klm}\) 是可学习参数。这种表示使 CPG 参数可以被 RL 策略或梯度下降直接优化 [Bellegarda & Ijspeert, IEEE RA-L, 2022]。
2.2 动态运动基元 (DMPs) 的数学推导
2.2.1 动机：从 CPG 到 DMP
CPG 擅长产生节律运动（行走、奔跑、游泳）。然而，动物还需要产生离散运动（伸手、抓取、踢腿）。DMP 正是为了在数学上统一这两种运动类型而提出的 [Ijspeert et al., Neural Computation, 2013]。
DMP 的核心洞察是：将运动的”时间调度” (timing) 与”空间形状” (spatial shape) 解耦。时间调度由正则系统管理，空间形状由转换系统中的可学习强迫项编码——而强迫项本身通过局部加权回归 (Locally Weighted Regression, LWR) 从示教数据中学习。
2.2.2 离散 DMP
正则系统 (一维跨临界分岔动力学):
\[\tau \dot{s} = -\alpha_s s\]其中 \(s\) 是相位变量，\(s(0) = 1\), \(\tau > 0\) 是时间缩放常数，\(\alpha_s > 0\) 是衰减速率。解为 \(s(t) = e^{-\alpha_s t / \tau}\)。随着 \(t \to \infty\), \(s \to 0\), 强迫项在运动结束时消失——转换系统收敛于一个稳定的点吸引子。
转换系统 (二阶线性动力学 + 非线性强迫项):
\[\tau \dot{z} = \alpha_z (\beta_z (g - y) - z) + f(s)\] \[\tau \dot{y} = z\]其中 \(y \in \mathbb{R}\) 是运动变量，\(z\) 是速度尺度变量，\(g\) 是目标/吸引子位置，\(\alpha_z, \beta_z\) 是正增益。通常设 \(\beta_z = \alpha_z / 4\) 以实现临界阻尼（无振荡收敛）。
强迫项的构造:
\[f(s) = \frac{\sum_{i=1}^{N} \Psi_i(s) \cdot w_i}{\sum_{i=1}^{N} \Psi_i(s)} \cdot s \cdot (g - y_0)\]其中 \(\Psi_i(s) = \exp\left(-\frac{1}{2\sigma_i^2}(s - c_i)^2\right)\) 是第 \(i\) 个高斯基函数，\(c_i\) 是基函数中心（均匀分布在 \([0, 1]\)），\(\sigma_i\) 是基函数宽度，\(w_i\) 是可学习权重。因子 \((g - y_0)\) 保证空间缩尺不变性；因子 \(s\) 确保 \(s \to 0\) 时 \(f(s) \to 0\)。
从示教数据学习权重:
给定示教轨迹 \(\{y_{\text{demo}}(t), \dot{y}_{\text{demo}}(t), \ddot{y}_{\text{demo}}(t)\}_{t=0}^{T}\)，目标强迫项通过求解转换系统反问题得到:
\[f_{\text{target}}(s) = \tau^2 \ddot{y}_{\text{demo}} - \alpha_z\left(\beta_z(g - y_{\text{demo}}) - \tau\dot{y}_{\text{demo}}\right)\]通过局部加权回归 (LWR) 学习权重:
\[w_i = \frac{\mathbf{s}^T \mathbf{\Gamma}_i \mathbf{f}_{\text{target}}}{\mathbf{s}^T \mathbf{\Gamma}_i \mathbf{s}}\]其中 \(\mathbf{\Gamma}_i = \text{diag}(\Psi_i(s_1), \ldots, \Psi_i(s_T))\)。
核心性质:
- 时间缩尺: 改变 \(\tau\) 快进或慢放相同的运动形状
- 目标缩尺: 改变 \(g\) 泛化到新目标位置，保持运动形状
- 扰动鲁棒性: 动力学本身是吸引子系统，外部扰动被渐近恢复
2.2.3 节律 DMP
正则系统替换为相位振荡器:
\[\tau \dot{\phi} = 1\]转换系统使用 von Mises 基函数（环形高斯）:
\[\Psi_i(\phi) = \exp\left(h_i(\cos(\phi - c_i) - 1)\right)\]\(h_i\) 控制基函数的集中度。与离散 DMP 不同，节律 DMP 的 \(f(\phi)\) 是周期为 \(2\pi\) 的周期函数，没有 \(s\) 衰减因子。这一统一性意味着同一个代码框架可以处理从”伸手拿杯子”到”行走步态”的所有运动类型。
2.2.4 DMP 的多维推广
对于 \(d\) 维运动，每个维度使用独立的转换系统，但共享同一个正则系统：
\(\tau \dot{s} = -\alpha_s s \quad \text{(共享时间)}\) \(\tau \dot{\mathbf{z}} = \alpha_z(\beta_z(\mathbf{g} - \mathbf{y}) - \mathbf{z}) + \mathbf{f}(s)\) \(\tau \dot{\mathbf{y}} = \mathbf{z}\)
这种结构确保不同自由度之间的时间同步——例如，机械臂的三个关节同时到达目标。
2.3 Tegotae 原理：感觉反馈驱动的去中心化协调
2.3.1 Tegotae 的定义
定义 2.3 (Tegotae) [Owaki et al., J. R. Soc. Interface, 2013]: Tegotae 是一个日语概念，描述”感知到的反应与预期的运动指令之间的匹配程度”。在数学上，第 \(k\) 个肢体的局部 Tegotae 函数 \(T_k\) 定义为感觉信号 \(s_k\) 与运动指令 \(m_k\) 的匹配度量：
\[T_k(t) = s_k(t) \cdot m_k(t)\]更一般地，正的 Tegotae 意味着感觉反馈与运动指令一致（”感知到的反应符合预期”），负的 Tegotae 意味着不一致（”遭遇意外阻力”）。
2.3.2 基于 Tegotae 的相位调制
在四足机器人控制中，Tegotae 被用于调制 CPG 相位振荡器的频率：
\[\dot{\phi}_k = \omega + \frac{\partial T_k}{\partial \phi_k}\]直觉上：如果感觉反馈告诉你”当前的运动是有效的”（正的 Tegotae），则保持速度；如果感觉反馈告诉你”遭遇了阻力”（负的 Tegotae），则减慢或调整。
负载反馈的 Tegotae 设计:
定义感觉信号 \(s_k = N_k\) (第 \(k\) 个肢体的法向地面反力)，运动指令 \(m_k = \cos(\phi_k)\) (期望负载模式——\(\phi_k = 0\) 时腿在支撑相，期望大负载):
\[T_k(\phi_k) = \sigma \cdot N_k \cdot \cos(\phi_k)\]相位动力学:
\[\dot{\phi}_k = \omega + \frac{\partial T_k}{\partial \phi_k} = \omega - \sigma \cdot N_k \cdot \sin(\phi_k)\]当 \(N_k > 0\) (腿在地上) 且 \(\phi_k \approx 0\) (腿在支撑相) 时，\(\sin(\phi_k) \approx 0\)，频率不变。当 \(N_k \gg 0\) (腿负荷过大) 且 \(\phi_k\) 接近支撑末期时，\(-\sin(\phi_k)\) 为负，相位减慢——在移动之前多停留一会儿，等负荷减轻。
2.3.3 去中心化协调的核心洞察
Tegotae 的核心贡献在于：肢间协调不需要直接的振荡器-振荡器耦合。通过使每个振荡器对局部感觉反馈做出反应——而这种感觉反馈又受到其他肢体与环境交互的物理影响——协调自然地涌现出来。这一原则与传统的 CPG 耦合模型形成互补：传统模型依赖显式的 \(\Gamma_{kj}(\phi_j - \phi_k)\) 项，而 Tegotae 依赖物理世界的隐式通信 (“information through the environment”)。
Owaki 等人展示了仅用 Tegotae 驱动的振荡器就可以在四足机器人上涌现出稳定的行走、小跑和疾驰步态——仅需增大全局频率参数 \(\omega\) [Owaki et al., Sci. Rep., 2017]。
2.4 生物力学基础：步态周期与肌肉模型
2.4.1 步态周期的数学描述
生物行走的步态周期 (gait cycle) 从一侧脚跟触地开始，到同侧脚跟再次触地结束。定义相位变量 \(\phi \in [0, 2\pi)\):
\[\text{stance phase}: \phi \in [0, 2\pi \cdot d), \quad \text{swing phase}: \phi \in [2\pi \cdot d, 2\pi)\]其中 \(d\) 是占空比 (duty factor)。对于行走 \(d \approx 0.6\)，对于奔跑 \(d \approx 0.3\)。
2.4.2 Hill 肌肉模型
Hill 型肌肉 [Hill, Proc. R. Soc. B, 1938] 由三个元素组成:
- 收缩单元 (CE): 主动力产生，\(F_{\text{CE}} = a \cdot F_{\text{max}} \cdot f_l(l) \cdot f_v(v)\)，其中 \(a \in [0,1]\) 是激活度
- 串联弹性单元 (SEE): 肌腱弹性
- 并联弹性单元 (PEE): 结缔组织被动刚度
肌肉总力:
\[F = F_{\text{max}} \left[a \cdot f_l(l) \cdot f_v(v) + f_{\text{PEE}}(l)\right]\]这一模型是 NeuroMechFly 中使用的肌肉力计算基础。
三、论文精读：Ramdya & Ijspeert (2023)
3.1 论文概况
Ramdya 和 Ijspeert 2023 年发表于 Science Robotics 的综述论文是本章的基石文献。全文分为三大板块：(1) 用仿真和机器人测试已知的生物学运动控制原则；(2) 将这些原则工程化为机器人控制器；(3) 展望神经科学与机器人学深度融合的未来。
论文的核心论证结构：
\[\begin{aligned} &\text{大前提: } &\text{脊椎动物和无脊椎动物共享一套保守的运动控制架构}\\ &\text{小前提: } &\text{纯工程控制方法在非结构化环境中缺乏鲁棒性和效率}\\ &\hline &\text{结论: } &\text{嵌入生物原则的机器人控制可显著提升性能} \end{aligned}\]3.2 核心方法论框架：四阶段互动
论文 Fig. 2 提出的四阶段框架是本章最重要的视觉信息：
| 阶段 | 描述 | 当前状态 |
|---|---|---|
| (A) 经典方法 | 人类科学家解读生物实验，提取抽象控制原则 | 历史上的主要途径 |
| (B) 当前主流 | 用仿真和机器人测试这些原则的有效性 | 活跃中 |
| (C) 近期前沿 | 利用可学习 ANN（可能受生物连接组约束）发现新原则 | 初步探索 |
| (D) 远期愿景 | AI 科学家在人类定义的框架下自动化科学推理 | 未来方向 |
3.3 案例研究精读
3.3.1 飞蝇粘附与步态选择 (Ramdya et al., Nat. Commun., 2017)
核心发现：在物理仿真中，有腿尖粘附的垂直行走条件下，优化速度得出的最佳步态是昆虫中广泛观察到的静态稳定三脚架步态；而在无粘附的地面行走条件下，优化出的是自然界中不存在的动态稳定”双脚架”步态。
解读：昆虫在进化中保留了三脚架步态，不是因为它在”纯速度”上最优，而是因为昆虫需要在三维空间（包括倒挂）中行走——这需要持续的附着力。快速奔跑的”腾空期”（动态稳定）在攀爬时是不可能的。这是一个优雅的例子，展示了生态位如何通过物理约束塑造了神经控制策略。
代价函数形式化:
\[J(\text{gait}) = \alpha \cdot v_{\text{forward}} - \beta \cdot P_{\text{fall}}(\text{gait})\]在地面上，\(P_{\text{fall}}(\text{tripod}) \approx P_{\text{fall}}(\text{bipod}) \approx 0\)，优化器选择更快的 bipod 步态。在墙面上，\(P_{\text{fall}}(\text{bipod})\) 显著大于 \(P_{\text{fall}}(\text{tripod})\)，优化器选择三脚架步态。
3.3.2 RHex：开环 CPG 步态的工程化 (Saranli et al., Int. J. Robot. Res., 2001)
核心洞察：蟑螂在快速奔跑时，肌电图记录显示运动模式几乎不随地形变化——暗示快速步态本质上是开环控制。RHex 的设计哲学：用极简的正弦信号驱动旋转腿，依靠腿刚度的被动力学性质吸收扰动。
RHex 的控制器:
\[\theta_i(t) = A \cdot \sin(\omega t + \phi_i) + \theta_0\]其中 \((A, \omega, \phi_i)\) 是常数。没有反馈。没有模型。没有在线优化。这种极简控制器使 RHex 能在泥地、沙地和密集植被中稳定行走。RHex 后来被 Boston Dynamics 商业化。
3.3.3 Geyer-Herr 人类步行模型 (Geyer & Herr, IEEE TNSRE, 2010)
核心发现：一个仅由脊髓反射驱动的 2D 肌肉骨骼模型（7 块肌肉，3 旋转自由度/腿，零 CPG），能产生极为逼真的人类步态——关节运动学、地面反作用力、肌肉激活模式都与实验数据高度吻合。
步态是一个极限环——系统在状态空间中自组织出一个稳定的周期性轨道，完全由反射回路和被动动力学维持。扰动使状态点偏离极限环，但反馈动力学逐渐将其收回。
工程影响：Geyer-Herr 模型的反射回路设计被移植到踝假肢 [Eilenberg et al., 2010]、膝-踝假肢 [Thatte & Geyer, 2016] 和外骨骼 [Wu et al., 2017]。
3.3.4 反射 vs. CPG：关于运动起源的论战
论文中一个引人入胜的隐性论战：
反射派证据：
- 线虫爬行：绝大多数神经肌肉方案可不使用 CPG，仅凭牵张感受器反馈产生振荡 [Izquierdo & Beer, Phil. Trans. R. Soc. B, 2018]
- 人类步行：28 个反射参数即可产生极限环步态 [Geyer & Herr, 2010]
CPG 派证据：
- CPG 可独立于外周输入产生周期性活动 [Grillner & Wallen, 1985]
- 七鳃鳗游泳：纯 CPG 和纯感觉反馈均可独立产生，但两者结合对神经损伤的鲁棒性最强 [Thandiackal et al., Sci. Robot., 2021]
合流：Tegotae 揭示了第三种可能——物理交互而非神经耦合驱动肢间协调。
3.3.5 NeuroMechFly：果蝇的数字孪生 (Lobato-Rios et al., Nat. Methods, 2022)
研究流程：成年雌蝇 CT 扫描 → 数字化重建生物力学外壳 → 基于真实运动学测量启用关节自由度 → PyBullet 物理引擎仿真。
关键能力：
- 通过回放运动学数据，推断无法直接测量的关节力矩和地面反力
- 改变前肢几何形状会破坏 grooming 运动学的回放——证明了形态逼真度的重要性
- 用遗传算法优化简单神经控制器实现快速稳定行走
NeuroMechFly 2.0 [Wang-Chen et al., Nat. Methods, 2024] 进一步升级：MuJoCo 引擎、视觉/嗅觉输入、上行动作反馈、FlyGym 发布。
3.3.6 OroBOT：逆向工程灭绝物种的运动 (Nyakatura et al., Nature, 2019)
研究目标：重建二叠纪早期干羊膜动物 Orobates pabsti（同时保存完整骨骼和足迹化石的最古老物种）的行走步态。
四管齐下的方法：
- Micro-CT 扫描化石，校正埋藏学变形，建立 3D 骨骼模型
- 分析四种现存四足爬行动物的 X 射线运动学和地面反力
- 在约束足部位置与化石足迹匹配的条件下，测试参数组合
- OroBOT 物理机器人验证——某些步态导致机器人摔倒甚至折断肋骨
关键发现：Orobates 的步态比此前预期更”先进”——接近凯门鳄式的相对直立姿态，将高级行走能力的演化与羊膜卵的出现解耦。
四、算法与代码实现
以下代码实现了本章的核心数学模型：Matsuoka CPG 振荡器、离散 DMP 和 Tegotae 负载反馈。代码基于纯 NumPy/SciPy，可直接运行。
4.1 Matsuoka CPG 振荡器
"""
matsuoka_cpg.py — Matsuoka 双神经元 CPG 振荡器
参考文献:
Matsuoka, K. (1987) "Mechanisms of frequency and pattern control
in the neural rhythm generators", Biological Cybernetics, 56: 345-353.
Ramdya & Ijspeert (2023), Science Robotics, 8, eadg0279.
振荡条件: sqrt(a12*a21) > 1 + tau_x/tau_y 且 tau_y > 3*tau_x
"""
import numpy as np
from typing import Tuple
class MatsuokaOscillator:
"""
Matsuoka 双神经元振荡器 — 生物 CPG 的最小模型
动力学 (Section 2.1.2):
tau_x * dx_i/dt = -x_i - b*y_i - sum_j(a_ij * max(x_j,0)) + s_i
tau_y * dy_i/dt = -y_i + max(x_i, 0)
Parameters
----------
tau_x : float 膜电位时间常数 (默认 1.0)
tau_y : float 适应时间常数 (默认 12.0, >> tau_x)
a12, a21 : float 相互抑制权重 (默认 2.5)
b : float 适应强度 (默认 2.5)
s1, s2 : float 强直输入 (默认 1.0)
dt : float 仿真步长 (默认 0.01)
"""
def __init__(self, tau_x=1.0, tau_y=12.0, a12=2.5, a21=2.5,
b=2.5, s1=1.0, s2=1.0, dt=0.01):
self.tau_x = tau_x
self.tau_y = tau_y
self.a12 = a12
self.a21 = a21
self.b = b
self.s1 = s1
self.s2 = s2
self.dt = dt
# 验证振荡条件
cond = np.sqrt(a12 * a21)
threshold = 1.0 + tau_x / tau_y
if cond <= threshold:
print(f"警告: sqrt(a12*a21)={cond:.2f} <= 1+tau_x/tau_y={threshold:.2f}, "
f"可能不满足振荡条件")
# 初始化状态: [x1, y1, x2, y2]
self.state = np.array([0.5, 0.0, -0.5, 0.0])
def step(self) -> Tuple[float, float]:
"""执行一步 Euler 积分，返回 (z1, z2) = (屈肌输出, 伸肌输出)"""
x1, y1, x2, y2 = self.state
z1 = max(x1, 0) # ReLU 非线性
z2 = max(x2, 0)
dx1 = (-x1 - self.b * y1 - self.a12 * z2 + self.s1) / self.tau_x
dy1 = (-y1 + z1) / self.tau_y
dx2 = (-x2 - self.b * y2 - self.a21 * z1 + self.s2) / self.tau_x
dy2 = (-y2 + z2) / self.tau_y
self.state += np.array([dx1, dy1, dx2, dy2]) * self.dt
return max(self.state[0], 0), max(self.state[2], 0)
def run(self, duration: float) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
"""运行指定时长 (s)，返回 (t, flex, exten) 数组"""
n_steps = int(duration / self.dt)
t = np.arange(n_steps) * self.dt
flex = np.zeros(n_steps)
exten = np.zeros(n_steps)
for i in range(n_steps):
flex[i], exten[i] = self.step()
return t, flex, exten
def set_tonic_drive(self, s1: float, s2: float):
"""调节强直输入 — 控制振幅 (频率基本不变)"""
self.s1 = max(0.0, s1)
self.s2 = max(0.0, s2)
def add_feedback(self, ds1: float, ds2: float):
"""添加感觉反馈调制"""
self.s1 = max(0.0, self.s1 + ds1)
self.s2 = max(0.0, self.s2 + ds2)
# ======== 演示 ========
if __name__ == "__main__":
cpg = MatsuokaOscillator(tau_x=1.0, tau_y=12.0, a12=2.5, a21=2.5,
b=2.5, s1=1.0, s2=1.0, dt=0.01)
t, flex, exten = cpg.run(duration=5.0)
# 反相关性: 屈肌和伸肌应交替激活
cross_corr = np.correlate(flex - flex.mean(),
exten - exten.mean(), mode='full')
peak_lag = np.argmax(cross_corr) - len(flex) + 1
period_samples = abs(2 * peak_lag)
print(f"估计振荡周期: {period_samples * cpg.dt:.2f} s")
print(f"估计振荡频率: {1.0 / (period_samples * cpg.dt):.2f} Hz")
print(f"输出振幅范围: flex [{flex.min():.2f}, {flex.max():.2f}], "
f"exten [{exten.min():.2f}, {exten.max():.2f}]")
4.2 动态运动基元 (DMP)
"""
dmp.py — 离散动态运动基元 (DMP)
参考文献:
Ijspeert et al. (2013) "Dynamical Movement Primitives",
Neural Computation, 25(2): 328-373.
"""
import numpy as np
from typing import Optional, Tuple
class DiscreteDMP:
"""
离散 DMP: 从示教轨迹学习点吸引子动力学 (Section 2.2.2)
规范系统: tau * ds/dt = -alpha_s * s
转换系统: tau * dz/dt = alpha_z*(beta_z*(g - y) - z) + f(s)
tau * dy/dt = z
强迫项: f(s) = [sum(psi_i(s)*w_i) / sum(psi_i(s))] * s * (g - y0)
Parameters
----------
n_bfs : int 基函数数量 (默认 50)
alpha_s : float 规范系统衰减率 (默认 3.0)
alpha_z : float 转换系统刚度 (默认 25.0)
beta_z : float 转换系统阻尼 (默认 alpha_z/4, 临界阻尼)
dt : float 仿真步长 (默认 0.001)
"""
def __init__(self, n_bfs=50, alpha_s=3.0, alpha_z=25.0,
beta_z=None, dt=0.001):
self.n_bfs = n_bfs
self.alpha_s = alpha_s
self.alpha_z = alpha_z
self.beta_z = beta_z if beta_z is not None else alpha_z / 4.0
self.dt = dt
# 基函数中心 (均匀分布在 [0, 1])
self.c = np.exp(-alpha_s * np.linspace(0, 1, n_bfs))
# 基函数宽度 (自适应)
self.h = np.ones(n_bfs) * (0.5 / (0.5 * (self.c[0] - self.c[-1]) / n_bfs)**2
if n_bfs > 1 else 1.0)
# 确保第一个和最后一个基函数有足够的宽度覆盖边界
self.h[0] = self.h[0] * 2
self.h[-1] = self.h[-1] * 2
self.w = np.zeros(n_bfs)
# 状态
self.y = 0.0
self.z = 0.0
self.s = 1.0
self.y0 = 0.0
self.g = 1.0
self.tau = 1.0
def reset(self, y0: float, g: float, tau: float = 1.0):
"""重置 DMP 状态并设置新目标。"""
self.y = y0
self.z = 0.0
self.s = 1.0
self.y0 = y0
self.g = g
self.tau = tau
def _basis(self, s_val: float) -> np.ndarray:
"""计算所有基函数在 s 处的激活值。"""
return np.exp(-self.h * (s_val - self.c)**2)
def _forcing_term(self, s_val: float) -> float:
"""计算强迫项 f(s)。"""
psi = self._basis(s_val)
psi_sum = np.sum(psi)
if psi_sum < 1e-10:
return 0.0
f = np.dot(psi, self.w) / psi_sum
return f * s_val * (self.g - self.y0)
def step(self) -> Tuple[float, float, float]:
"""推进一个时间步。返回 (y, z, s)。"""
# 规范系统
ds = -self.alpha_s * self.s / self.tau
# 转换系统
f_val = self._forcing_term(self.s)
dz = (self.alpha_z * (self.beta_z * (self.g - self.y) - self.z) + f_val) / self.tau
dy = self.z / self.tau
# Euler 积分
self.s += ds * self.dt
self.z += dz * self.dt
self.y += dy * self.dt
self.s = max(self.s, 0.0)
return self.y, self.z, self.s
def rollout(self, y0: float, g: float, tau: float = 1.0,
max_steps: int = 10000) -> np.ndarray:
"""执行完整 rollout，返回轨迹。"""
self.reset(y0, g, tau)
traj = [self.y]
for _ in range(max_steps):
y, _, s = self.step()
traj.append(y)
if s <= 1e-6:
break
return np.array(traj)
def imitate(self, y_demo: np.ndarray, tau: float = 1.0):
"""从单条示教轨迹学习 DMP 权重 (LWR, Section 2.2.2)。"""
T = len(y_demo)
dt_demo = tau / T
yd = np.asarray(y_demo)
dyd = np.gradient(yd) / dt_demo
ddyd = np.gradient(dyd) / dt_demo
# 规范系统轨迹
t_demo = np.linspace(0, tau, T)
s_demo = np.exp(-self.alpha_s * t_demo / tau)
self.y0 = yd[0]
self.g = yd[-1]
self.tau = tau
# 反算目标强迫项
f_target = (tau**2 * ddyd
- self.alpha_z * (self.beta_z * (self.g - yd) - tau * dyd))
# LWR
for i in range(self.n_bfs):
psi_i = np.exp(-self.h[i] * (s_demo - self.c[i])**2)
num = np.sum(psi_i * f_target)
den = np.sum(psi_i) + 1e-12
self.w[i] = num / den
# ======== 演示 ========
if __name__ == "__main__":
# 生成示教轨迹: 最小急动度点对点运动
T_demo = 1.0
dt_demo = 0.001
t = np.arange(0, T_demo, dt_demo)
s_frac = t / T_demo
y_demo = 10*s_frac**3 - 15*s_frac**4 + 6*s_frac**5 # y: 0 -> 1
dmp = DiscreteDMP(n_bfs=50, dt=0.001)
dmp.imitate(y_demo, tau=T_demo)
# 验证: 原始目标
traj_orig = dmp.rollout(y0=0.0, g=1.0, tau=1.0)
mse = np.mean((traj_orig[:len(y_demo)] - y_demo[:len(traj_orig)])**2)
print(f"重构 MSE: {mse:.6f}")
# 目标泛化
traj_new = dmp.rollout(y0=0.0, g=2.0, tau=1.0)
print(f"泛化到 g=2.0: 终点={traj_new[-1]:.4f} (期望 2.0)")
# 时间缩尺
traj_fast = dmp.rollout(y0=0.0, g=1.0, tau=0.5)
print(f"2x 速度 (tau=0.5): 步数={len(traj_fast)} vs 原始 ({len(traj_orig)})")
4.3 Tegotae 负载反馈 + 四足步态生成
"""
tegotae_quadruped.py — Tegotae 驱动的分散式四足步态
参考文献:
Owaki et al. (2013) J. R. Soc. Interface, 10, 20120669.
Owaki & Ishiguro (2017) Sci. Rep., 7, 277.
Ramdya & Ijspeert (2023) Science Robotics, 8, eadg0279.
"""
import numpy as np
from typing import Tuple, Dict
class TegotaeQuadruped:
"""
基于 Tegotae 的去中心化四足相位振荡器 (Section 2.3)
相位动力学 (Eq. 2.24):
dphi_k/dt = omega - sigma1 * N_k * cos(phi_k)
+ sigma2 * avg(N_neighbors) * cos(phi_k)
Parameters
----------
omega : float 固有角频率 (rad/s), 默认 2*pi*1.2 (1.2 Hz)
sigma1 : float 局部负载反馈增益, 默认 2.0
sigma2 : float 肢间负载反馈增益, 默认 1.0
coupling : float 辅助振荡器耦合强度, 默认 1.0
dt : float 仿真步长 (s), 默认 0.001
"""
def __init__(self, omega=2.0*np.pi*1.2, sigma1=2.0, sigma2=1.0,
coupling=1.0, dt=0.001):
self.omega = omega
self.sigma1 = sigma1
self.sigma2 = sigma2
self.coupling = coupling
self.dt = dt
# 四足相位 (LF=0, RF=1, LH=2, RH=3)
self.phi = np.array([0.0, np.pi, np.pi, 0.0]) # 小跑初始
# 腿名映射
self.names = ['LF', 'RF', 'LH', 'RH']
# 邻接关系 (每条腿连接其他三条)
self.neighbors = {
0: [1, 2, 3], # LF
1: [0, 2, 3], # RF
2: [0, 1, 3], # LH
3: [0, 1, 2], # RH
}
def step(self, grf: np.ndarray) -> np.ndarray:
"""
执行一个控制步。
Parameters
----------
grf : ndarray shape=(4,) — 每条腿的地面反力 (N)
Returns
-------
phi : ndarray shape=(4,) — 更新后的相位
"""
dphi = np.full(4, self.omega)
# 振荡器间弱耦合 (辅助同步)
for k in range(4):
for j in range(4):
if k != j:
# 小跑步态期望相位差
if (k, j) in [(0, 1), (1, 0), (2, 3), (3, 2),
(0, 2), (2, 0), (1, 3), (3, 1)]:
delta = np.pi
dphi[k] += self.coupling * np.sin(
self.phi[j] - self.phi[k] - delta)
else:
dphi[k] += self.coupling * np.sin(
self.phi[j] - self.phi[k])
# Tegotae 负载反馈 (Eq. 2.24)
for k in range(4):
# 局部反馈: f_k1 = -sigma1 * N_k * cos(phi_k)
dphi[k] -= self.sigma1 * grf[k] * np.cos(self.phi[k])
# 肢间反馈: f_k2 = +sigma2 * mean(N_neighbors) * cos(phi_k)
nbr_grf = np.mean([grf[j] for j in self.neighbors[k]])
dphi[k] += self.sigma2 * nbr_grf * np.cos(self.phi[k])
# Euler 积分
self.phi += dphi * self.dt
self.phi = self.phi % (2.0 * np.pi)
return self.phi.copy()
def get_joint_angle(self, leg_idx: int, amplitude: float = 0.6,
offset: float = 0.1, duty_factor: float = 0.6
) -> Tuple[float, str]:
"""将相位映射为关节角度 (分段函数, Section 2.1.4)。"""
phi_leg = self.phi[leg_idx]
phi_swing = 2.0 * np.pi * duty_factor
if phi_leg < phi_swing:
frac = phi_leg / max(phi_swing, 1e-10)
angle = offset + amplitude * (1.0 - 2.0 * frac)
return angle, 'stance'
else:
frac = (phi_leg - phi_swing) / max(2.0*np.pi - phi_swing, 1e-10)
angle = offset + amplitude * (-1.0 + 2.0 * frac)
return angle, 'swing'
def analyze_gait(self, phi_history: np.ndarray) -> str:
"""分析稳态步态类型 (Section 2.1.3 步态表)。"""
steady = phi_history[-500:] if len(phi_history) > 500 else phi_history
mean_phi = np.mean(steady, axis=0) % (2*np.pi)
sorted_indices = np.argsort(mean_phi)
sorted_phi = mean_phi[sorted_indices]
diffs = np.diff(sorted_phi)
diffs = np.append(diffs, 2*np.pi - sorted_phi[-1] + sorted_phi[0])
mean_diff = np.mean(diffs)
if mean_diff > 2.5: # 一对近 0，一对近 pi → trot
return "trot"
elif mean_diff > 2.0: # 近似均匀分布 → walk
return "walk"
else:
return "gallop/bound"
# ======== 演示 ========
if __name__ == "__main__":
# 在不同频率下测试步态
for omega in [np.pi, 2*np.pi, 4*np.pi]:
robot = TegotaeQuadruped(omega=omega)
n_steps = 3000
phi_hist = np.zeros((n_steps, 4))
for t in range(n_steps):
# 简化的地面反力模型 (支撑相有力, 摆动相无力)
grf = np.zeros(4)
for k in range(4):
if robot.phi[k] < np.pi: # 前半周期为支撑相
grf[k] = max(2.0 * np.cos(robot.phi[k]) + 0.3, 0.0)
phi_hist[t] = robot.step(grf)
gait = robot.analyze_gait(phi_hist)
freq_hz = omega / (2*np.pi)
print(f"omega={omega:.1f} rad/s ({freq_hz:.1f} Hz): 步态 → {gait}")
4.4 整合演示
"""
integrate_demo.py — 整合 CPG + Tegotae + DMP 的完整演示
运行:
python integrate_demo.py
"""
import numpy as np
from matsuoka_cpg import MatsuokaOscillator
from dmp import DiscreteDMP
from tegotae_quadruped import TegotaeQuadruped
def demo_cpg_comparison():
"""比较不同强直输入下的 Matsuoka 振荡器行为。"""
print("=" * 60)
print("演示 1: Matsuoka CPG — 强直输入与振幅/频率的关系")
print("=" * 60)
for s in [0.8, 1.0, 1.5, 2.0]:
cpg = MatsuokaOscillator(s1=s, s2=s, dt=0.01)
t, flex, _ = cpg.run(duration=5.0)
# 分析稳定后的频率和振幅
steady = flex[-1000:] # 最后 10s
zero_crossings = np.where(np.diff(np.signbit(steady - np.mean(steady))))[0]
if len(zero_crossings) >= 2:
period = np.mean(np.diff(zero_crossings[::2])) * cpg.dt
freq = 1.0 / period
amp = (steady.max() - steady.min()) / 2
print(f" s={s:.1f}: 频率={freq:.2f} Hz, 振幅={amp:.3f}")
def demo_dmp_generalization():
"""演示 DMP 的目标泛化和时间缩尺。"""
print("\n" + "=" * 60)
print("演示 2: DMP — 目标泛化与时间缩尺")
print("=" * 60)
# 示教: 最小急动度轨迹
T = 1.0
dt = 0.001
t = np.arange(0, T, dt)
s = t / T
y_demo = 10*s**3 - 15*s**4 + 6*s**5 # 0 -> 1
dmp = DiscreteDMP(n_bfs=50, dt=dt)
dmp.imitate(y_demo, tau=T)
# 测试 1: 原始目标
traj1 = dmp.rollout(y0=0.0, g=1.0)
print(f" 原始目标 g=1.0: 终点={traj1[-1]:.4f}")
# 测试 2: 目标泛化
for g in [0.5, 1.5, 2.0]:
traj = dmp.rollout(y0=0.0, g=g)
print(f" 泛化目标 g={g:.1f}: 终点={traj[-1]:.4f}")
# 测试 3: 时间缩尺
for tau in [0.5, 2.0]:
traj = dmp.rollout(y0=0.0, g=1.0, tau=tau)
print(f" 时间缩尺 tau={tau:.1f}: 步数={len(traj)}")
def demo_tegotae_gait_transition():
"""演示 Tegotae 频率扫描导致的步态转换。"""
print("\n" + "=" * 60)
print("演示 3: Tegotae — 频率驱动步态转换")
print("=" * 60)
for omega in [np.pi, 2*np.pi, 4*np.pi, 6*np.pi]:
robot = TegotaeQuadruped(omega=omega, dt=0.001)
n_steps = 4000
phi_hist = np.zeros((n_steps, 4))
for t in range(n_steps):
grf = np.zeros(4)
for k in range(4):
if robot.phi[k] < np.pi:
grf[k] = max(2.0 * np.cos(robot.phi[k]) + 0.3, 0.0)
phi_hist[t] = robot.step(grf)
gait = robot.analyze_gait(phi_hist)
freq = omega / (2*np.pi)
print(f" omega={omega:.1f} ({freq:.1f} Hz): 步态 → {gait}")
if __name__ == '__main__':
demo_cpg_comparison()
demo_dmp_generalization()
demo_tegotae_gait_transition()
五、核心 Takeaway 与延伸阅读
5.1 五个核心洞察
1. “傀儡牵线”架构是动物运动的统一范式。 大脑通过少量下行通路控制脊髓/腹神经索中的 CPG 网络——就像傀儡师通过几根线控制傀儡的复杂动作一样 [Loeb, J. Physiol., 2001]。这种架构将高维的运动规划问题约化为低维的”动机信号”（方向、速度、步态选择）。工程上，这意味着机器人控制策略不需要直接输出关节力矩——它只需要学习如何在 CPG 参数空间中导航。Bellegarda & Ijspeert (2022) 的 CPG-RL 工作正是利用了这一洞察。
2. 被动动力学是最便宜的控制器。 RHex 和 Geyer-Herr 模型的成功都指向一个反直觉的结论：正确设计的身体力学本身就能产生稳定行为。RHex 的旋转腿 + 合适的刚度 = 无需反馈即可穿越复杂地形。人类的跟腱弹性 = 无代谢代价的行走能量回收。在设计机器人时，”让物理替你计算”比”编写显式控制器来对抗物理”更高效、更鲁棒。
3. 感觉反馈和 CPG 的关系不是”选择”，而是”集成”。 Thandiackal et al. (2021) 的最强实验证据：七鳃鳗游泳可以通过纯 CPG 或纯感觉反馈实现，但混合控制器对神经元损伤的鲁棒性远高于任一单独机制。这一”并联冗余”是生物系统鲁棒性的关键来源——机器人控制系统应效仿这一点。
4. Tegotae 原则揭示了物理交互作为通信通道。 Owaki 等人的工作表明：当每个 CPG 振荡器仅依赖局部感觉反馈时，全局协调步态仍然可以通过身体-环境物理交互涌现。这从根本上挑战了传统的”强耦合 CPG 网络”观点，为去中心化的分布式机器人控制开辟了新路径。
5. NeuroMechFly 和 OroBOT 代表了具身神经科学的新范式。 这些工作不仅是”用仿真验证假说”(A/B 阶段)，更是”用仿真生成假说”(C/D 阶段)。NeuroMechFly 可以推断无法测量的关节力矩；OroBOT 逆向工程了灭绝动物的步态。未来的 Fig. 2D 愿景——AI 科学家自动设计实验、更新模型、生成可测试的神经科学假说——虽然遥远，但技术上已初具雏形。
5.2 开放问题与研究前沿
| 问题 | 描述 |
|---|---|
| 从连接组到控制策略 | 果蝇腹神经索连接组提供了约 \(10^4\) 个神经元的完整连接图谱 [Phelps et al., Cell, 2021]。如何从”接线图”推导出可运行的控制策略？ |
| 可微 CPG | 如何将 Matsuoka 振荡器的 \(\max(0, u)\) 替换为可微替代（如 Softplus），以便端到端梯度优化？ |
| Sim-to-Real 中的动力学差距 | NeuroMechFly 的仿真动力学在多大程度上忠实反映真实果蝇的物理？如何量化”仿真保真度”？ |
| 多模态运动统一框架 | 如何设计一个 CPG 网络，使其能通过改变耦合参数在行走、奔跑、攀爬、游泳之间连续切换？[Ijspeert et al., Science, 2007] |
| 形态与控制协同设计 | 如果身体形态和神经控制是协同进化的，在机器人设计中是否应该同时优化形态学和控制器？ |
| 敏捷性基准统一 | 如何为动物和机器人定义可量化的”敏捷性”指标，在同一尺度上比较？[Eckert & Ijspeert, IEEE TRO, 2019] |
5.3 推荐延伸阅读
| 优先级 | 文献 | 核心内容 |
|---|---|---|
| ★★★★★ | Ramdya & Ijspeert (2023) Science Robotics | 本章核心论文，124 篇参考文献的全面综述 |
| ★★★★★ | Ijspeert et al. (2013) Neural Computation | DMP 的权威数学推导 |
| ★★★★☆ | Geyer & Herr (2010) IEEE TNSRE | 基于反射的人类步行极限环模型 |
| ★★★★☆ | Owaki et al. (2013) J. R. Soc. Interface | Tegotae 概念的原始提出 |
| ★★★★☆ | Lobato-Rios et al. (2022) Nature Methods | NeuroMechFly — 果蝇神经力学数字孪生 |
| ★★★☆☆ | Grillner & Wallen (1985) Annu. Rev. Neurosci. | CPG 的经典综述，七鳃鳗游泳网络 |
| ★★★☆☆ | Dickinson et al. (2000) Science | “How animals move” — 整合视角经典 |
| ★★★☆☆ | Nyakatura et al. (2019) Nature | OroBOT — 逆向工程灭绝物种步态 |
| ★★★☆☆ | Bellegarda & Ijspeert (2022) IEEE RA-L | CPG-RL: 用 RL 学习高层 CPG 调制策略 |
| ★★★☆☆ | Saranli et al. (2001) Int. J. Robot. Res. | RHex 机器人 — 开环 CPG 控制的工程化 |
六、练习与思考
6.1 数学推导
练习 2.1 (Matsuoka 振荡器) 证明：对于对称的 Matsuoka 振荡器 (\(s_f = s_e = s\), \(w_{fe} = w_{ef} = w\)), 当 \(s > 0\) 且 \(\tau_v > 3\tau_u\) 时，系统存在稳定极限环。提示：分析 \((u_f, u_e)\) 相平面中的零倾线和慢变量 \(v_f, v_e\) 的漂移方向。
练习 2.2 (Kuramoto 同步) 考虑两个耦合相位振荡器: \(\dot{\phi}_1 = \omega_1 + K \sin(\phi_2 - \phi_1)\), \(\dot{\phi}_2 = \omega_2 + K \sin(\phi_1 - \phi_2)\)。推导相位差 \(\Delta = \phi_1 - \phi_2\) 的演化方程。证明当 \(K > \vert \omega_1 - \omega_2\vert /2\) 时，两个振荡器达到频率锁定。这个不等式有什么物理含义？
练习 2.3 (DMP 时间缩尺) 证明对于离散 DMP，如果将以时间缩尺 \(\tau\) 学习的轨迹与以 \(\tau'\) 回放的轨迹比较，空间形状被保持，时间演化被均匀缩尺。提示：对转换系统进行变量替换。
练习 2.4 (Tegotae 设计) 考虑 Tegotae 函数 \(T_k(\phi_k) = \sigma \cdot N_k \cdot \cos(\phi_k - \phi_0)\)。解释 \(\phi_0\) 的物理含义。如果希望肢体在支撑相末端被负载时加速（而非减速），应如何修改 Tegotae 函数？
6.2 编程练习
练习 2.5 扩展 MatsuokaOscillator 类，实现六足 (hexapod) 的三脚架步态。需要 6 个振荡器，耦合定义期望相位差使 \((LF, LM, LH, RF, RM, RH)\) 的稳态相位为 \((0, \frac{4\pi}{3}, \frac{2\pi}{3}, \pi, \frac{\pi}{3}, \frac{5\pi}{3})\)。
练习 2.6 为 DiscreteDMP 实现多维支持。对于二维示教轨迹 \((y_1(t), y_2(t))\)，使用共享正则系统但独立转换系统的架构。测试：学习一个圆弧轨迹并泛化到不同起始点。
练习 2.7 修改 TegotaeQuadruped 的 analyze_gait 方法，使其能精确区分四种步态：walk, trot, pace, bound。
练习 2.8 在 demo_tegotae_gait_transition 中添加对 \(\text{sigma}_1\) 和 \(\text{sigma}_2\) 的扫描，观察反馈增益对比步态的敏感度。画出步态-参数相图。
6.3 思考题
思考 2.1 (CPG vs. 传统轨迹规划) 在什么场景下选择 CPG 驱动而非传统轨迹规划？在什么场景下选择传统轨迹规划？列出至少三个场景并给出理由。讨论：随着深度 RL 的兴起 (Ch11, Ch14), CPG 是否正在变得过时，还是获得了新的使用方式？
思考 2.2 (形态与控制协同设计) 你正在设计一个搜救机器人。它需要穿越碎石堆、爬楼梯、挤过窄缝。从 RHex 的设计中你获得了什么启示？如果在设计阶段可以同时选择身体形态和控制架构，你会如何利用”具身化”来简化控制？
思考 2.3 (从果蝇到人形机器人) NeuroMechFly 项目建立了与真实果蝇一一对应的数字孪生。如果我们对 Atlas/ANYmal 做同样的事：在 MuJoCo 中建立高精度模型，然后用模仿学习将人类行走数据迁移到仿真中。这种方法相比于从零用 RL 训练步态策略 (如 Ch14 的 RMA) 的优势和劣势分别是什么？
思考 2.4 (动物的最优性) CPG-反射-具身化架构在动物中如此普遍，但从控制理论的角度来看，它是最优的吗？或者说，CPG 只是进化在经济性约束下的一个”足够好”的解决方案？讨论：如果让你用现代最优控制理论（LQR, MPC）从头设计一个四足机器人控制器，你会得到类似 CPG + 反射的结构吗？
下一章预告: 第 3 章《机器人机构学：运动学与动力学》将从生物力学转入机器人机构学。你将学习 SO(3)/SE(3) 李群-李代数对应、指数坐标的 Product-of-Exponentials 公式、Denavit-Hartenberg 参数、以及拉格朗日与牛顿-欧拉两种动力学推导方法——这是将 CPG 输出映射到真实关节力矩所必需的工程基础。