SSM 优化（上）：为什么不能直接算，以及 HiPPO 理论
前面几篇文章建立了 SSM 的基础知识和 S4 Block 的架构理解。但我们刻意回避了一个核心问题：矩阵 $\mathbf{A}$ 到底该怎么初始化？
如果 $\mathbf{A}$ 随机初始化，SSM 在长序列上几乎必然失败——要么数值爆炸，要么瞬间遗忘。HiPPO 理论从数学上证明了什么样的 $\mathbf{A}$ 才是"理论上最优的"。本文将完整介绍 HiPPO 的推导过程和核心结论，为后续的 S4 工程化优化做好铺垫。
一、为什么不能随机初始化 $\mathbf{A}$
1. 随机初始化的后果
回顾离散化后的 SSM 递推公式：$h_k = \bar{\mathbf{A}}h_{k-1} + \bar{\mathbf{B}}x_k$，其中 $\bar{\mathbf{A}} = e^{\mathbf{A}\Delta}$。
如果 $\mathbf{A}$ 像普通神经网络权重那样随机初始化（例如标准正态分布），会发生三种情况：
$\mathbf{A}$ 的特征值实部为正： $\bar{\mathbf{A}}$ 的特征值大于 1。递推 $h_k = \bar{\mathbf{A}}h_{k-1} + \ldots$ 中，矩阵连乘 $\bar{\mathbf{A}}^k$ 会迅速趋于无穷大，状态 $h_k$ 数值爆炸。
$\mathbf{A}$ 的特征值实部为负且绝对值很大： $\bar{\mathbf{A}}$ 的特征值远小于 1。矩阵连乘 $\bar{\mathbf{A}}^k$ 迅速趋于 0，模型瞬间"遗忘"所有历史输入。
$\mathbf{A}$ 是一般随机矩阵： 特征值在复平面上近似均匀分布在某个圆内（Circular Law），无法保证全部实部为负，更无法控制衰减速率的分布。在这种情况下，训练几乎不可能稳定进行。
2. 问题的本质
SSM 需要一个特殊的 $\mathbf{A}$ 矩阵，它必须同时满足两个看似矛盾的要求：
- 稳定性： 特征值实部全为负，系统不发散。
- 长程记忆： 至少有一些特征值的绝对值足够小（衰减足够慢），使得模型能够记住序列开头的信息。
传统 RNN 随机初始化的权重无法同时满足这两个要求，这是 RNN 训练困难的根本原因之一。HiPPO 理论的目标就是从数学上推导出一个理论上最优的 $\mathbf{A}$ 矩阵初始化方案。
二、HiPPO 的核心思想：从"死记硬背"到"函数逼近"
1. 两种记忆策略的对比
传统的 RNN 记忆策略类似于"死记硬背"：试图用隐藏状态 $h_t$ 存储尽可能多的历史输入的原始信息。这需要极大的状态维度 $N$，且信息容易在传递中丢失。
HiPPO 的记忆策略完全不同：它把"记忆"重新定义为"用一组基函数的系数来逼近历史信号"。隐藏状态 $h(t)$ 存储的不是原始数据，而是这组基函数的展开系数。
2. 具体设定
设历史输入信号为 $f(\tau)$，$\tau \in [0, t]$。在时刻 $t$，我们用 $N$ 个基函数 ${p_n(\tau, t)}_{n=0}^{N-1}$ 的线性组合来逼近 $f$ 在 $[0, t]$ 上的行为：
$$\hat{f}(\tau; t) = \sum_{n=0}^{N-1} c_n(t) \cdot p_n(\tau, t)$$
其中 $c_n(t)$ 是第 $n$ 个基函数的系数，也就是隐藏状态 $h(t)$ 的第 $n$ 个分量。
目标是最小化加权 $L^2$ 误差：
$$\min_{c(t)} \int_0^t \left(f(\tau) - \hat{f}(\tau; t)\right)^2 d\mu_t(\tau)$$
其中 $\mu_t$ 是定义在 $[0, t]$ 上的测度，决定了对历史不同时刻的关注权重。
HiPPO 实际拟合的对象是什么？（重要的语义澄清）
在 NLP 语言模型中，输入并不是直接的文字，而是先经过 Tokenization → Embedding 映射，得到每个 token 的特征向量 $x_t \in \mathbb{R}^d$。HiPPO 及其衍生的 SSM 层，各自对"沿序列位置变化的数值信号"进行拟合。
具体来说，对每个特征维度 $j$（$j=1,\ldots,d$），都会有一个独立的数值序列 $$f^{(j)}(\tau) = \text{第 }, \tau \text{个 token 的第 }, j \text{个特征维度的值}$$ 该数值序列被视作 HiPPO 框架中的"连续信号"$f(\tau)$，并用勒让德多项式基对其历史进行逼近。状态向量 $h(t)$ 的每个分量因此存储的是该特征维度信号的多项式展开系数。
例如 $d=768$、$N=64$ 的常见配置意味着：模型共有 768 个并行的状态空间"通道"，每个通道用 64 维的状态向量去压缩对应特征维度沿时间轴的历史信息。
澄清：勒让德多项式本身不理解语义。 基函数 $P_n$ 只是对"时间轴上数值的变化模式"进行数学拟合，它不具备对语言含义的理解。语义能力来自于多层 SSM Block 纵向堆叠后的非线性变换与表达学习，而非某一层状态空间本身。
3. 勒让德多项式（Legendre Polynomials）的显式性质
HiPPO 选择的基函数是勒让德多项式 $P_n(x)$，定义在 $[-1, 1]$ 上。为了不依赖"读者应当已经熟悉"，在此列出推导过程中实际用到的勒让德多项式的显式性质。
定义（Rodrigues 公式）： $$P_n(x) = \frac{1}{2^n n!} \frac{d^n}{dx^n}(x^2 - 1)^n, \quad x \in [-1, 1]$$
前几项的显式表达式：
- $P_0(x) = 1$
- $P_1(x) = x$
- $P_2(x) = \frac{1}{2}(3x^2 - 1)$
- $P_3(x) = \frac{1}{2}(5x^3 - 3x)$
正交性： $$\int_{-1}^{1} P_m(x) P_n(x) , dx = \frac{2}{2n+1} \delta_{mn}$$
递推关系： $$(2n+1)x P_n(x) = (n+1)P_{n+1}(x) + n P_{n-1}(x)$$
微分展开式（后续推导的核心工具）： $$P_n’(x) = \sum_{\substack{k=0 \ k \equiv n-1 \pmod{2}}}^{n-1} (2k+1) P_k(x)$$ 即 $P_n’(x)$ 可以表示为所有与 $n$ 奇偶性不同的低阶 $P_k(x)$ 的线性组合。
边界值： $$P_n(1) = 1, \quad P_n(-1) = (-1)^n$$
选择勒让德多项式的三个硬性原因：
- 在有限区间上正交：保证系数可以独立求解，不会出现病态矩阵（Gram 矩阵条件数恒为 1）。
- 数值稳定：相比普通幂函数基 ${1, x, x^2, \ldots}$（在 $[0,1]$ 上高度线性相关），勒让德正交基从一开始就是良态的。
- 有封闭的递推公式和微分公式：后续推导需要将积分 $\int f \cdot P_n’$ 转化为低阶系数 $c_k$ 的线性组合，而这正需要上面列出的微分展开式。
三、HiPPO 矩阵 $\mathbf{A}$ 和 $\mathbf{B}$ 的推导
1. 构造正交基
将历史区间 $[0, t]$ 映射到 $[-1, 1]$，令 $x = \frac{2\tau - t}{t}$。定义缩放后的基函数：
$$p_n(\tau, t) = P_n!\left(\frac{2\tau - t}{t}\right)$$
在测度 $d\mu_t = \frac{d\tau}{t}$ 下，正交性变为：
$$\int_0^t P_m!\left(\frac{2\tau-t}{t}\right) P_n!\left(\frac{2\tau-t}{t}\right) \frac{d\tau}{t} = \frac{1}{2n+1}\delta_{mn}$$
这里读者可能会问：为什么非要把 $[0,t]$ 换到 $[-1,1]$ 上来？主要有两个理由。
- 勒让德多项式 $P_n(x)$ 本来就是定义在 $[-1,1]$ 上的正交多项式，直接套进去就能用。
- 换到 $[-1,1]$ 之后，正交性、递推关系、微分展开式这些已知公式可以直接使用，不用重新推一遍。
所以这个变量替换本质上是把"在 $[0,t]$ 上的缩放勒让德基"，翻译回"标准的 $[-1,1]$ 上的勒让德基"，方便调公式。
2. 最优系数
由正交性，最小二乘最优系数为：
$$c_n(t) = (2n+1) \int_0^t f(\tau) P_n!\left(\frac{2\tau - t}{t}\right) \frac{d\tau}{t}$$
3. 对 $c_n(t)$ 求时间导数
这是推导的核心。$c_n(t)$ 对 $t$ 的依赖来自三处：积分上限 $t$、被积函数中 $P_n$ 的参数 $\frac{2\tau-t}{t}$、以及前面的 $\frac{1}{t}$。使用 Leibniz 积分法则：
$$\dot{c}_n(t) = (2n+1)\left[\frac{f(t)P_n(1)}{t} - \frac{1}{t^2}\int_0^t f(\tau)P_n!\left(\frac{2\tau-t}{t}\right)d\tau + \frac{1}{t}\int_0^t f(\tau)\frac{\partial}{\partial t}P_n!\left(\frac{2\tau-t}{t}\right)d\tau\right]$$
逐项处理：
第一项： 由 $P_n(1) = 1$，得 $\frac{2n+1}{t}f(t)$。
第二项： 由 $c_n(t)$ 的定义，$\int_0^t f(\tau)P_n(\cdots)d\tau = \frac{t}{2n+1}c_n(t)$，所以第二项为 $-\frac{1}{t}c_n(t)$。
第三项： 先计算偏导数 $$\frac{\partial}{\partial t}P_n!\left(\frac{2\tau-t}{t}\right) = P_n’!\left(\frac{2\tau-t}{t}\right)\cdot \left(-\frac{2\tau}{t^2}\right)$$ 代入得第三项 $$\frac{1}{t}\int_0^t f(\tau) \cdot P_n’!\left(\frac{2\tau-t}{t}\right)\left(-\frac{2\tau}{t^2}\right) d\tau = -\frac{2}{t^3}\int_0^t f(\tau), \tau\cdot P_n’!\left(\frac{2\tau-t}{t}\right) d\tau$$ 使用变量替换 $x = \frac{2\tau-t}{t}$（则 $\tau = \frac{t(x+1)}{2}$，$d\tau = \frac{t}{2}dx$，积分区间从 $\tau\in[0,t]$ 变为 $x\in[-1,1]$），该项可转化为关于 $P_n’(x)$ 的积分，再利用上述勒让德微分展开式 $P_n’(x) = \sum_{k\equiv n-1 !!\pmod{2}}^{n-1}(2k+1)P_k(x)$ 将其拆解为低阶 $P_k$ 的线性组合。利用每个 $P_k$ 项对应的系数定义 $$\int_0^t f(\tau) P_k!\left(\frac{2\tau-t}{t}\right)\frac{d\tau}{t} = \frac{c_k(t)}{2k+1}$$ 最终第三项表示为 ${c_k(t)}_{k<n}$ 的线性组合，从而得到闭合的矩阵形式。
注意：上述变量替换与展开正是 HiPPO 原文中"第三项代数化简"的核心步骤，完整展开依赖勒让德微分关系 $P’n = \sum{k\equiv n-1(2)}(2k+1)P_k$。这也是为什么 HiPPO 必须选取具有显式微分展开式的正交基的原因之一。
这一步如果只看公式，会觉得有点像"突然出现了一个展开式"。它的实际用途只有一个：把 $\int f \cdot P_n’$ 这种"导数形式的积分"，变成 $\sum c_k$ 这种"系数的线性组合"，从而让整个导数方程闭合在已知的系数 $c_n(t)$ 上面。否则你就会得到一堆没法用已知变量表示的积分项，方程就不可能写成状态空间的线性形式。
再多说一句，为什么勒让德多项式能做到这一点，而普通幂函数基 $1, x, x^2, \dots$ 做不到：
- 勒让德的递推式和微分展开式，保证了 $P_n’$ 可以写成只含比 $n$ 低阶的 $P_k$ 的线性组合。
- 这一点让第三项的积分结果天然落在已有的系数集 ${c_k}_{k<n}$ 内部，不需要引入新的信号信息。
- 这也是为什么 HiPPO 选择勒让德的一个隐藏理由：它不只是"因为正交"，更因为它的代数结构能让这类推导闭合。
4. 最终结果
将三项合并，得到 HiPPO 微分方程：
$$\dot{c}n(t) = \frac{1}{t}\left[\sum{k=0}^{n} A_{nk} c_k(t) + B_n f(t)\right]$$
写成矩阵形式：
$$\dot{\mathbf{c}}(t) = \frac{1}{t}\left(\mathbf{A}\mathbf{c}(t) + \mathbf{B}f(t)\right)$$
其中矩阵 $\mathbf{A}$ 和向量 $\mathbf{B}$ 的具体元素为：
$$A_{nk} = \begin{cases} 0 & \text{若 } k > n \ -(2n+1)/2 & \text{若 } k = n \ \sqrt{(2n+1)(2k+1)} & \text{若 } k < n \end{cases}$$
$$B_n = \sqrt{2n+1}$$
这个结果对应于 HiPPO-LegS（Scaled Legendre）变体。
四、$\mathbf{A}$ 矩阵的定量性质
1. 矩阵结构
$\mathbf{A}$ 是一个 $N \times N$ 的实数密集下三角矩阵。以 $N = 4$ 为例：
$$\mathbf{A} = \begin{pmatrix} -1/2 & 0 & 0 & 0 \ \sqrt{3} & -3/2 & 0 & 0 \ \sqrt{5} & \sqrt{15} & -5/2 & 0 \ \sqrt{7} & \sqrt{35} & \sqrt{105} & -7/2 \end{pmatrix}$$
2. 特征值
由于 $\mathbf{A}$ 是下三角矩阵，其特征值就是对角线元素：
$$\lambda_n = -\frac{2n+1}{2}, \quad n = 0, 1, \ldots, N-1$$
即 $\lambda_0 = -1/2$，$\lambda_1 = -3/2$，$\lambda_2 = -5/2$，$\ldots$，$\lambda_{N-1} = -(2N-1)/2$。
3. 特征值与记忆的关系
回忆一维线性系统的解析解：若 $\dot{h} = \lambda h$，则 $h(t) = e^{\lambda t} h(0)$。
稳定性： 所有 $\lambda_n < 0$，系统渐近稳定，不会爆炸。这是随机初始化无法保证的。
多尺度记忆： 每个特征值对应的衰减时间常数为 $\tau_n = 1/|\lambda_n|$。
- 第 0 个模式：$|\lambda_0| = 1/2$，衰减时间常数 $\tau_0 = 2$，衰减最慢，负责长期记忆（全局趋势）。
- 第 $N-1$ 个模式：$|\lambda_{N-1}| = (2N-1)/2$，衰减时间常数 $\tau_{N-1} = 2/(2N-1)$，衰减最快，负责短期记忆（局部细节）。
- 特征值均匀间隔分布：$|\lambda_{n+1}| - |\lambda_n| = 1$，模型自动拥有从短到长的全频段记忆能力。
与随机初始化的对比： 随机矩阵的特征值在复平面上近似均匀分布在圆内（Circular Law），无法保证全部实部为负，更无法控制衰减时间尺度的分布。
4. $\mathbf{B}$ 向量的物理意义
$B_n = \sqrt{2n+1}$。它决定了新输入 $f(t)$ 写入第 $n$ 阶系数的强度。
- 高阶（$n$ 大）的 $B_n$ 更大 → 新输入对高阶多项式系数（捕捉局部快速变化）的影响更强。
- 低阶（$n$ 小）的 $B_n$ 更小 → 新输入对低阶多项式系数（捕捉全局慢速趋势）的影响更温和。
这符合直觉：一个新的数据点主要影响局部细节，对全局趋势的影响应该被"稀释"。
五、HiPPO 是初始化方案，不是固定参数
1. 训练时 $\mathbf{A}$ 和 $\mathbf{B}$ 会变化
HiPPO 推导出的 $\mathbf{A}$ 和 $\mathbf{B}$ 是初始值。在训练过程中，它们会通过梯度下降被更新：
- 初始化： 用 HiPPO 公式计算出的精确值。
- 前向传播： 用当前的 $\mathbf{A}$、$\mathbf{B}$ 计算状态和输出。
- 反向传播： 计算损失对 $\mathbf{A}$ 和 $\mathbf{B}$ 的梯度。
- 参数更新： $\mathbf{A} \leftarrow \mathbf{A} - \eta \frac{\partial \mathcal{L}}{\partial \mathbf{A}}$，$\mathbf{B}$ 同理。
训练结束后，$\mathbf{A}$ 和 $\mathbf{B}$ 的值不再严格等于 HiPPO 公式的值。
2. HiPPO 初始化的意义
HiPPO 的意义是提供一个理论最优的起点，而非固定不变的参数。
- 随机初始化： $\mathbf{A}$ 的特征值随机分布，大概率导致训练初期就爆炸或消失，模型学不到长程依赖。
- HiPPO 初始化： $\mathbf{A}$ 的特征值从一开始就精确分布在负实轴上，系统天然稳定且具备多尺度记忆。训练只是在这个已经非常好的起点上做微调。
可以类比为：用 ImageNet 预训练的 ResNet 权重初始化图像分类网络，而不是随机初始化。HiPPO 就是给 $\mathbf{A}$ 做的"理论预训练"。
3. 为什么有了闭式解，仍然需要训练
以上一节的解析解论述为基础，这里回答一个常见困惑：线性 ODE 天然具备矩阵指数的闭式解，既然数学上可以"直接求解"，为什么 SSM 还需要通过梯度下降进行训练？
答案可以拆成三层：
-
闭式解的前提是所有量都已知。 线性常微分方程 $h’(t) = \mathbf{A}h(t) + \mathbf{B}x(t)$ 的精确闭式解确实存在，形式为 $$h(t) = e^{\mathbf{A}t}h(0) + \int_0^t e^{\mathbf{A}(t-\tau)}\mathbf{B}x(\tau),d\tau$$ 但这个表达式里，$\mathbf{A},\mathbf{B},\mathbf{C}$ 的具体数值、初始状态 $h(0)$、以及输入信号 $x(\tau)$ 都必须已知才能代入。
-
在深度学习场景中，这些量最初都不知。
- 输入 $x(\tau)$（例如一段文本的 Embedding 序列）只有在前向运行时才逐步可见，无法事先写出闭式积分。
- 矩阵 $\mathbf{A},\mathbf{B},\mathbf{C}$ 本身是模型参数——如果它们已知且最优，就不需要 SSM 这一整层架构去"学习"了。
- 因此闭式解在这里不是可直接调用的计算公式，而是离散化方案（矩阵指数、ZOH、双线性变换等）的理论依据。
-
SSM 的训练，本质上是在参数空间里搜索最优的 $\mathbf{A,B,C}$。 离散化公式（例如 $h_k = \bar{\mathbf{A}}h_{k-1} + \bar{\mathbf{B}}x_k$）将连续闭式解的结构"嵌入"到网络的前向计算中，然后通过 loss 反传更新 $\mathbf{A,B,C}$。换句话说：
- 不是"用解析解算出答案"，
- 而是"用解析解启发的离散结构去搭建可微分的前向计算，并在此基础上学习参数"。
这也是为什么博客系列在介绍完离散化方法之后，重点转向 HiPPO 初始化和 DPLR 加速：前者回答"参数从什么最优位置出发"，后者回答"考虑到参数未知且需学习，怎样把矩阵运算变得高效以便训练"。两者共同服务于"可训练的 SSM"这一目标，而不是去替代训练本身。
小结：线性 ODE 的解析解在 SSM 语境里是结构与推导的来源，不是训练的替代品。SSM 之所以还要训练，是因为输入未知、参数未知，且任务目标（如下一个 token 的预测损失）无法由通用的 ODE 通解自动满足，需要在具体数据上把参数学出来。
4. 多层间的参数分化
每一层 SSM Block 的 $\mathbf{A}$ 和 $\mathbf{B}$ 都用相同的 HiPPO 公式来初始化。但在训练过程中，每一层的梯度不同（因为每层处理的特征抽象级别不同），所以训练结束后每一层的 $\mathbf{A}$ 和 $\mathbf{B}$ 会自然分化：浅层偏向捕捉局部模式，深层偏向捕捉全局语义。
六、HiPPO 面临的工程挑战
HiPPO 在理论上给出了完美的 $\mathbf{A}$ 矩阵，但当我们要把它塞进深度学习框架进行训练时，会遇到两个严重的工程问题。
1. 时变性问题
HiPPO-LegS 的方程中含有 $\frac{1}{t}$ 因子：
$$\dot{\mathbf{c}}(t) = \frac{1}{t}(\mathbf{A}\mathbf{c}(t) + \mathbf{B}f(t))$$
这里的 $\frac{1}{t}$ 让系统成为时变系统。这意味着：
- 状态转移的速率随着时间 $t$ 的增长而变化（越往后，衰减越慢）。
- 不能直接使用卷积视角——卷积要求时不变系统。
这是因为测度 $d\mu_t = \frac{d\tau}{t}$ 中的归一化因子 $\frac{1}{t}$ 依赖于当前时间 $t$。随着 $t$ 增大，越早的输入被"稀释"得越厉害，这合理：旧信息不应该拥有和新信息一样的权重。但这种时变性给工程实现带来了困难。
2. 计算复杂度问题
HiPPO-LegS 的 $\mathbf{A}$ 是一个 $N \times N$ 的密集下三角矩阵。计算矩阵指数 $e^{\mathbf{A}\Delta}$ 和卷积核 $\bar{K}_i = \mathbf{C}\bar{\mathbf{A}}^i\bar{\mathbf{B}}$ 的复杂度是 $O(N^3)$。当 $N = 256$ 时，$256^3 \approx 1.6 \times 10^7$，在深层网络中每一步都算一遍，计算量不可接受。
七、总结
本文介绍了 HiPPO 理论的完整推导过程及其理论意义：
- HiPPO 将"记忆"重新定义为"用正交多项式逼近历史信号"，从而推导出理论上最优的 $\mathbf{A}$ 矩阵。
- HiPPO-LegS 的 $\mathbf{A}$ 是一个密集下三角矩阵，特征值均匀分布在负实轴上，天然具备稳定性和多尺度记忆能力。
- HiPPO 只是一个初始化方案，训练中参数会通过梯度下降更新。
- 线性 ODE 的闭式解在 SSM 中是推导与结构的来源，不能替代训练——因为输入未知、参数未知，且任务目标需要在数据上学习参数。
- 但 HiPPO 有两个工程问题：时变性（含 $\frac{1}{t}$）和 $O(N^3)$ 的计算复杂度。
下一篇文章将介绍 S4 如何通过 DPLR 结构、双线性变换和 Z 变换，解决这两个工程问题，让 HiPPO 理论真正能在 GPU 上高效运行。从而引出本文第一节遗留的问题：为什么单纯知道 $\mathbf{A}$ 的理论形式还不够，还需要 S4 的一系列工程化改造。