回顾 Art. 14 算子全景 的概念地图：
我们沿着时序算子子线来到最后一站。前两篇处理的是离散状态 ——有限个类别之间的跳转，转移矩阵的每一行是概率分布。现在，我们从表格的上方移到下方 ：状态不再是 { s 1 , … , s N } \{s_1, \ldots, s_N\} { s 1 , … , s N } 连续向量 x ( t ) ∈ R n \mathbf{x}(t) \in \mathbb{R}^n x ( t ) ∈ R n
核心变化：
离散 → 连续 ：差分方程 x k + 1 = A x k \mathbf{x}_{k+1} = A\mathbf{x}_k x k + 1 = A x k x ˙ = A x \dot{\mathbf{x}} = A\mathbf{x} x ˙ = A x 矩阵幂 → 矩阵指数 ：A n A^n A n e A t e^{At} e A t 模 < 1 → 实部 < 0 ：离散稳定性判据变成连续稳定性判据Forward 算法 → Kalman 滤波 ：离散求和变成连续积分，但”预测 + 修正”的核心结构不变
为什么这在 ML 中重要？因为状态空间模型 （State-Space Model, SSM）——Mamba（Art. 26 SSM/Mamba ）的数学基础——正是连续时间线性系统的离散化。本文建立的公式将在 Part 3 中直接复用。
从离散到连续：A n A^n A n e A t e^{At} e A t
离散系统回顾
在 Art. 15 马尔可夫链 中，系统的演化规则是：
x k + 1 = A x k \mathbf{x}_{k+1} = A\mathbf{x}_k x k + 1 = A x k
从初始状态 x 0 \mathbf{x}_0 x 0 k k k
x k = A k x 0 \mathbf{x}_k = A^k \mathbf{x}_0 x k = A k x 0
对角化 A = Q Λ Q − 1 A = Q\Lambda Q^{-1} A = Q Λ Q − 1 Art. 2 特征分解 ）：
x k = Q Λ k Q − 1 x 0 = Q diag ( λ 1 k , λ 2 k , … , λ n k ) Q − 1 x 0 \mathbf{x}_k = Q \Lambda^k Q^{-1} \mathbf{x}_0 = Q\, \text{diag}(\lambda_1^k, \lambda_2^k, \ldots, \lambda_n^k)\, Q^{-1} \mathbf{x}_0 x k = Q Λ k Q − 1 x 0 = Q diag ( λ 1 k , λ 2 k , … , λ n k ) Q − 1 x 0
每个特征方向独立演化，第 i i i λ i k \lambda_i^k λ i k 离散系统是”一步一步走”——每乘一次 A A A
连续系统：微分方程形式
现在想象步长无限缩小。离散系统中 x k + 1 − x k = ( A − I ) x k \mathbf{x}_{k+1} - \mathbf{x}_k = (A - I)\mathbf{x}_k x k + 1 − x k = ( A − I ) x k
x ˙ ( t ) = d x d t = A x ( t ) \dot{\mathbf{x}}(t) = \frac{d\mathbf{x}}{dt} = A\mathbf{x}(t) x ˙ ( t ) = d t d x = A x ( t )
这是一阶线性常微分方程组（linear ODE system）。A ∈ R n × n A \in \mathbb{R}^{n \times n} A ∈ R n × n 系统矩阵 （system matrix），它编码了状态各分量之间的耦合关系。x ˙ \dot{\mathbf{x}} x ˙ 连续系统是”平滑流动”——状态在每个瞬间都在改变，而非等到下一步。
类比 ：离散系统像翻阅照片——每张照片是一个时间步的快照。连续系统像播放视频——状态在时间轴上连续变化。
离散步进 vs 连续流动 离散：x_k = λ^k · x₀ k = 0, 1, 2, ... 像翻阅照片——每张是一步快照 连续：x(t) = e^{at} · x(0) t 连续 像播放视频——状态平滑变化 Aᵏ → eᴬᵗ | |λ| < 1 → Re(λ) < 0 | 单位圆 → 左半平面 解：矩阵指数 e A t e^{At} e A t
标量情形大家很熟悉：x ˙ = a x \dot{x} = ax x ˙ = a x x ( t ) = e a t x ( 0 ) x(t) = e^{at} x(0) x ( t ) = e a t x ( 0 )
x ( t ) = e A t x ( 0 ) \boxed{\mathbf{x}(t) = e^{At}\, \mathbf{x}(0)} x ( t ) = e A t x ( 0 )
其中 e A t e^{At} e A t 矩阵指数 （matrix exponential），定义为幂级数：
e A t = ∑ k = 0 ∞ ( A t ) k k ! = I + A t + ( A t ) 2 2 ! + ( A t ) 3 3 ! + ⋯ e^{At} = \sum_{k=0}^{\infty} \frac{(At)^k}{k!} = I + At + \frac{(At)^2}{2!} + \frac{(At)^3}{3!} + \cdots e A t = ∑ k = 0 ∞ k ! ( A t ) k = I + A t + 2 ! ( A t ) 2 + 3 ! ( A t ) 3 + ⋯
逐项理解这个级数：
k = 0 k = 0 k = 0 I I I k = 1 k = 1 k = 1 A t At A t t t t x ( t ) ≈ ( I + A t ) x ( 0 ) \mathbf{x}(t) \approx (I + At)\mathbf{x}(0) x ( t ) ≈ ( I + A t ) x ( 0 ) k = 2 k = 2 k = 2 A 2 t 2 2 ! \frac{A^2 t^2}{2!} 2 ! A 2 t 2 后续项：越来越高阶的修正
矩阵指数 e^{At} = 幂级数求和 e^{At} = I + At + (At)²/2! + (At)³/3! + ... I 不动 + At 线性近似 + A²t²/2! 二阶修正 + A³t³/3! 三阶修正 + ... 高阶项递减 阶乘 k! 压制增长——级数对任意 A 和 t 都绝对收敛（标量 e^a 的矩阵推广）
为什么级数收敛？ 对于任何有限维矩阵 A A A t t t ∥ A ∥ k t k k ! \frac{\|A\|^k t^k}{k!} k ! ∥ A ∥ k t k k ! k! k ! e a t = ∑ ( a t ) k k ! e^{at} = \sum \frac{(at)^k}{k!} e a t = ∑ k ! ( a t ) k a a a
验证它确实是解 ：对 x ( t ) = e A t x ( 0 ) \mathbf{x}(t) = e^{At}\mathbf{x}(0) x ( t ) = e A t x ( 0 ) t t t
d d t e A t = d d t ( I + A t + A 2 t 2 2 ! + A 3 t 3 3 ! + ⋯ ) = A + A 2 t + A 3 t 2 2 ! + ⋯ = A ⋅ e A t \frac{d}{dt} e^{At} = \frac{d}{dt}\left(I + At + \frac{A^2 t^2}{2!} + \frac{A^3 t^3}{3!} + \cdots\right) = A + A^2 t + \frac{A^3 t^2}{2!} + \cdots = A \cdot e^{At} d t d e A t = d t d ( I + A t + 2 ! A 2 t 2 + 3 ! A 3 t 3 + ⋯ ) = A + A 2 t + 2 ! A 3 t 2 + ⋯ = A ⋅ e A t
所以 x ˙ ( t ) = A e A t x ( 0 ) = A x ( t ) \dot{\mathbf{x}}(t) = A e^{At}\mathbf{x}(0) = A\mathbf{x}(t) x ˙ ( t ) = A e A t x ( 0 ) = A x ( t ) ■ \blacksquare ■
Inline Box：矩阵幂级数与标量 Taylor 展开的正式对应
标量 矩阵 x ˙ = a x \dot{x} = ax x ˙ = a x x ˙ = A x \dot{\mathbf{x}} = A\mathbf{x} x ˙ = A x 解：x ( t ) = e a t x ( 0 ) x(t) = e^{at} x(0) x ( t ) = e a t x ( 0 ) 解：x ( t ) = e A t x ( 0 ) \mathbf{x}(t) = e^{At}\mathbf{x}(0) x ( t ) = e A t x ( 0 ) e a t = ∑ k = 0 ∞ ( a t ) k k ! e^{at} = \sum_{k=0}^{\infty} \frac{(at)^k}{k!} e a t = ∑ k = 0 ∞ k ! ( a t ) k e A t = ∑ k = 0 ∞ ( A t ) k k ! e^{At} = \sum_{k=0}^{\infty} \frac{(At)^k}{k!} e A t = ∑ k = 0 ∞ k ! ( A t ) k 稳定：a < 0 a < 0 a < 0 稳定：所有特征值 Re ( λ i ) < 0 \text{Re}(\lambda_i) < 0 Re ( λ i ) < 0 振荡：a a a 振荡：特征值是纯虚数
从标量到矩阵，公式形式完全相同——“标量 a a a A A A 。这不是类比，而是严格的数学推广。唯一的注意事项是矩阵乘法不满足交换律，所以 e ( A + B ) t ≠ e A t e B t e^{(A+B)t} \neq e^{At}e^{Bt} e ( A + B ) t = e A t e B t A B = B A AB = BA A B = B A
通过对角化计算矩阵指数
直接用幂级数求和需要无限项，实际中不可行。但如果 A A A A = Q Λ Q − 1 A = Q\Lambda Q^{-1} A = Q Λ Q − 1 Art. 2 特征分解 中已经预告过）：
A k = Q Λ k Q − 1 A^k = Q\Lambda^k Q^{-1} A k = Q Λ k Q − 1
代入幂级数：
e A t = ∑ k = 0 ∞ ( Q Λ Q − 1 ) k t k k ! = Q ( ∑ k = 0 ∞ Λ k t k k ! ) Q − 1 = Q e Λ t Q − 1 e^{At} = \sum_{k=0}^{\infty} \frac{(Q\Lambda Q^{-1})^k t^k}{k!} = Q \left(\sum_{k=0}^{\infty} \frac{\Lambda^k t^k}{k!}\right) Q^{-1} = Q\, e^{\Lambda t}\, Q^{-1} e A t = ∑ k = 0 ∞ k ! ( Q Λ Q − 1 ) k t k = Q ( ∑ k = 0 ∞ k ! Λ k t k ) Q − 1 = Q e Λ t Q − 1
而对角矩阵的指数就是逐元素取指数：
e Λ t = diag ( e λ 1 t , e λ 2 t , … , e λ n t ) e^{\Lambda t} = \text{diag}(e^{\lambda_1 t}, e^{\lambda_2 t}, \ldots, e^{\lambda_n t}) e Λ t = diag ( e λ 1 t , e λ 2 t , … , e λ n t )
所以：
e A t = Q diag ( e λ 1 t , e λ 2 t , … , e λ n t ) Q − 1 \boxed{e^{At} = Q\, \text{diag}(e^{\lambda_1 t}, e^{\lambda_2 t}, \ldots, e^{\lambda_n t})\, Q^{-1}} e A t = Q diag ( e λ 1 t , e λ 2 t , … , e λ n t ) Q − 1
这与离散情形 A k = Q diag ( λ 1 k , … , λ n k ) Q − 1 A^k = Q\,\text{diag}(\lambda_1^k, \ldots, \lambda_n^k)\, Q^{-1} A k = Q diag ( λ 1 k , … , λ n k ) Q − 1 λ i k \lambda_i^k λ i k e λ i t e^{\lambda_i t} e λ i t
数值例子
取 A = [ − 1 0.5 0 − 2 ] A = \begin{bmatrix}-1 & 0.5 \\ 0 & -2\end{bmatrix} A = [ − 1 0 0.5 − 2 ]
第一步：特征分解。 A A A λ 1 = − 1 \lambda_1 = -1 λ 1 = − 1 λ 2 = − 2 \lambda_2 = -2 λ 2 = − 2
对 λ 1 = − 1 \lambda_1 = -1 λ 1 = − 1 ( A + I ) q 1 = [ 0 0.5 0 − 1 ] q 1 = 0 (A + I)\mathbf{q}_1 = \begin{bmatrix}0 & 0.5 \\ 0 & -1\end{bmatrix}\mathbf{q}_1 = \mathbf{0} ( A + I ) q 1 = [ 0 0 0.5 − 1 ] q 1 = 0 q 1 = [ 1 0 ] \mathbf{q}_1 = \begin{bmatrix}1 \\ 0\end{bmatrix} q 1 = [ 1 0 ]
对 λ 2 = − 2 \lambda_2 = -2 λ 2 = − 2 ( A + 2 I ) q 2 = [ 1 0.5 0 0 ] q 2 = 0 (A + 2I)\mathbf{q}_2 = \begin{bmatrix}1 & 0.5 \\ 0 & 0\end{bmatrix}\mathbf{q}_2 = \mathbf{0} ( A + 2 I ) q 2 = [ 1 0 0.5 0 ] q 2 = 0 q 2 = [ − 0.5 1 ] \mathbf{q}_2 = \begin{bmatrix}-0.5 \\ 1\end{bmatrix} q 2 = [ − 0.5 1 ]
Q = [ 1 − 0.5 0 1 ] , Q − 1 = [ 1 0.5 0 1 ] Q = \begin{bmatrix}1 & -0.5 \\ 0 & 1\end{bmatrix}, \quad Q^{-1} = \begin{bmatrix}1 & 0.5 \\ 0 & 1\end{bmatrix} Q = [ 1 0 − 0.5 1 ] , Q − 1 = [ 1 0 0.5 1 ]
第二步：矩阵指数。
e A t = [ 1 − 0.5 0 1 ] [ e − t 0 0 e − 2 t ] [ 1 0.5 0 1 ] = [ e − t 0.5 ( e − t − e − 2 t ) 0 e − 2 t ] e^{At} = \begin{bmatrix}1 & -0.5 \\ 0 & 1\end{bmatrix} \begin{bmatrix}e^{-t} & 0 \\ 0 & e^{-2t}\end{bmatrix} \begin{bmatrix}1 & 0.5 \\ 0 & 1\end{bmatrix} = \begin{bmatrix}e^{-t} & 0.5(e^{-t} - e^{-2t}) \\ 0 & e^{-2t}\end{bmatrix} e A t = [ 1 0 − 0.5 1 ] [ e − t 0 0 e − 2 t ] [ 1 0 0.5 1 ] = [ e − t 0 0.5 ( e − t − e − 2 t ) e − 2 t ]
第三步：验证。 取 x ( 0 ) = [ 1 1 ] \mathbf{x}(0) = \begin{bmatrix}1 \\ 1\end{bmatrix} x ( 0 ) = [ 1 1 ]
x ( t ) = e A t x ( 0 ) = [ e − t + 0.5 ( e − t − e − 2 t ) e − 2 t ] = [ 1.5 e − t − 0.5 e − 2 t e − 2 t ] \mathbf{x}(t) = e^{At}\mathbf{x}(0) = \begin{bmatrix}e^{-t} + 0.5(e^{-t} - e^{-2t}) \\ e^{-2t}\end{bmatrix} = \begin{bmatrix}1.5e^{-t} - 0.5e^{-2t} \\ e^{-2t}\end{bmatrix} x ( t ) = e A t x ( 0 ) = [ e − t + 0.5 ( e − t − e − 2 t ) e − 2 t ] = [ 1.5 e − t − 0.5 e − 2 t e − 2 t ]
在 t = 1 t = 1 t = 1 x ( 1 ) = [ 1.5 × 0.368 − 0.5 × 0.135 0.135 ] = [ 0.484 0.135 ] \mathbf{x}(1) = \begin{bmatrix}1.5 \times 0.368 - 0.5 \times 0.135 \\ 0.135\end{bmatrix} = \begin{bmatrix}0.484 \\ 0.135\end{bmatrix} x ( 1 ) = [ 1.5 × 0.368 − 0.5 × 0.135 0.135 ] = [ 0.484 0.135 ]
两个分量都在衰减（因为 λ 1 = − 1 < 0 , λ 2 = − 2 < 0 \lambda_1 = -1 < 0, \lambda_2 = -2 < 0 λ 1 = − 1 < 0 , λ 2 = − 2 < 0 x 2 x_2 x 2 e − 2 t e^{-2t} e − 2 t e − t e^{-t} e − t
稳定性：特征值的实部决定一切
连续系统的稳定性判据
对角化后，每个特征方向的分量独立演化：第 i i i e λ i t e^{\lambda_i t} e λ i t
如果 λ i \lambda_i λ i
λ i < 0 \lambda_i < 0 λ i < 0 e λ i t → 0 e^{\lambda_i t} \to 0 e λ i t → 0 λ i = 0 \lambda_i = 0 λ i = 0 e λ i t = 1 e^{\lambda_i t} = 1 e λ i t = 1 λ i > 0 \lambda_i > 0 λ i > 0 e λ i t → ∞ e^{\lambda_i t} \to \infty e λ i t → ∞
如果 λ i \lambda_i λ i λ i = α + i β \lambda_i = \alpha + i\beta λ i = α + i β
e λ i t = e ( α + i β ) t = e α t ⋅ e i β t = e α t ( cos β t + i sin β t ) e^{\lambda_i t} = e^{(\alpha + i\beta)t} = e^{\alpha t} \cdot e^{i\beta t} = e^{\alpha t}(\cos\beta t + i\sin\beta t) e λ i t = e ( α + i β ) t = e α t ⋅ e i β t = e α t ( cos β t + i sin β t )
α < 0 \alpha < 0 α < 0 α = 0 \alpha = 0 α = 0 α > 0 \alpha > 0 α > 0
稳定性完全由实部 α = Re ( λ i ) \alpha = \text{Re}(\lambda_i) α = Re ( λ i ) ，虚部 β = Im ( λ i ) \beta = \text{Im}(\lambda_i) β = Im ( λ i )
连续系统稳定 ⟺ 所有特征值满足 Re ( λ i ) < 0 \boxed{\text{连续系统稳定} \iff \text{所有特征值满足 } \text{Re}(\lambda_i) < 0} 连续系统稳定 ⟺ 所有特征值满足 Re ( λ i ) < 0
与离散系统的对比
离散系统 x k + 1 = A x k \mathbf{x}_{k+1} = A\mathbf{x}_k x k + 1 = A x k 连续系统 x ˙ = A x \dot{\mathbf{x}} = A\mathbf{x} x ˙ = A x 解 x k = A k x 0 \mathbf{x}_k = A^k \mathbf{x}_0 x k = A k x 0 x ( t ) = e A t x ( 0 ) \mathbf{x}(t) = e^{At}\mathbf{x}(0) x ( t ) = e A t x ( 0 ) 特征方向的行为 λ i k \lambda_i^k λ i k e λ i t e^{\lambda_i t} e λ i t 稳定条件 $ \lambda_i 稳定区域 复平面上的单位圆内部 复平面上的左半平面 振荡条件 λ i \lambda_i λ i λ i \lambda_i λ i
关键联系 ：离散的 ∣ λ ∣ < 1 |\lambda| < 1 ∣ λ ∣ < 1 Re ( λ ) < 0 \text{Re}(\lambda) < 0 Re ( λ ) < 0 z = e s Δ t z = e^{s\Delta t} z = e s Δ t s s s z z z Re ( s ) < 0 \text{Re}(s) < 0 Re ( s ) < 0 ∣ e s Δ t ∣ = e Re ( s ) Δ t < 1 |e^{s\Delta t}| = e^{\text{Re}(s)\Delta t} < 1 ∣ e s Δ t ∣ = e Re ( s ) Δ t < 1
下图并排展示了两种稳定区域——离散系统的单位圆和连续系统的左半平面，以及连接它们的离散化映射 z = e s Δ t z = e^{s\Delta t} z = e s Δ t
稳定区域对比：离散系统 vs 连续系统
离散：|λ| < 1（单位圆内） ↔ 连续：Re(λ) < 0（左半平面） 通过 z = eˢΔᵗ 联系
离散系统 xₖ₊₁ = Axₖ Re Im 稳定 |λ| < 1 不稳定 |λ| = 1 z = eˢΔᵗ 离散化映射 连续系统 ẋ = Ax Re Im 稳定 Re(λ) < 0 不稳定 Re(λ) = 0 Re(s) < 0 ⟹ |eˢΔᵗ| = e^(Re(s)Δt) < 1 — 连续稳定 ⟹ 离散稳定 复特征值 → 振荡行为
当系统矩阵 A A A λ = α ± i β \lambda = \alpha \pm i\beta λ = α ± i β 振荡 （oscillation）叠加指数包络 （exponential envelope）：
振荡频率 ：ω = β \omega = \beta ω = β T = 2 π / β T = 2\pi / \beta T = 2 π / β 包络变化 ：e α t e^{\alpha t} e α t α < 0 \alpha < 0 α < 0 α > 0 \alpha > 0 α > 0
一个典型的物理例子：阻尼弹簧系统 m x ¨ + c x ˙ + k x = 0 m\ddot{x} + c\dot{x} + kx = 0 m x ¨ + c x ˙ + k x = 0
[ x ˙ v ˙ ] = [ 0 1 − k / m − c / m ] [ x v ] \begin{bmatrix}\dot{x}\\\dot{v}\end{bmatrix} = \begin{bmatrix}0 & 1\\-k/m & -c/m\end{bmatrix}\begin{bmatrix}x\\v\end{bmatrix} [ x ˙ v ˙ ] = [ 0 − k / m 1 − c / m ] [ x v ]
当阻尼 c c c
可视化：线性系统的相平面
下面的交互组件展示了六种典型 2D 线性系统的状态演化。每条彩色轨迹从不同的初始条件出发，在系统矩阵 A A A x ( t ) = e A t x ( 0 ) \mathbf{x}(t) = e^{At}\mathbf{x}(0) x ( t ) = e A t x ( 0 )
2D 线性系统状态演化
ẋ = Ax 的相平面轨迹
稳定节点 不稳定节点 稳定螺旋 不稳定螺旋 中心（纯振荡） 鞍点
A = [-1, 0.5; 0, -2] 两个负实特征值：所有轨迹指数衰减到原点
-3 -3 -2 -2 -1 -1 1 1 2 2 3 3 x₁ x₂ t = 3.00 特征值位置
稳定 不稳定 Re Im λ1=-1 λ2=-2 连续系统稳定性判据：
播放 重置 t = 3.00
探索建议 ：
稳定节点 ：观察所有轨迹如何从不同方向收敛到原点。特征值都在左半平面，两个衰减速率不同。不稳定节点 ：轨迹向外发散。特征值都在右半平面。稳定螺旋 ：轨迹螺旋式收敛——这是复特征值（负实部）的标志性行为。不稳定螺旋 ：轨迹螺旋式发散。复特征值（正实部）。中心（纯振荡） ：轨迹在闭合曲线上运动——纯虚特征值，既不收敛也不发散。鞍点 ：一正一负特征值——沿一个方向收敛，另一个方向发散。
离散化：从连续到离散的桥梁
为什么需要离散化？因为数字计算机只能处理离散的时间步。连续系统 x ˙ = A x + B u \dot{\mathbf{x}} = A\mathbf{x} + B\mathbf{u} x ˙ = A x + B u u \mathbf{u} u x k + 1 = A ˉ x k + B ˉ u k \mathbf{x}_{k+1} = \bar{A}\mathbf{x}_k + \bar{B}\mathbf{u}_k x k + 1 = A ˉ x k + B ˉ u k Mamba 的状态矩阵就是通过离散化从连续参数得到的 。
零阶保持（Zero-Order Hold, ZOH）
ZOH 假设输入 u ( t ) \mathbf{u}(t) u ( t ) [ k Δ t , ( k + 1 ) Δ t ) [k\Delta t, (k+1)\Delta t) [ k Δ t , ( k + 1 ) Δ t ) u ( t ) = u k \mathbf{u}(t) = \mathbf{u}_k u ( t ) = u k
从连续方程 x ˙ = A x + B u \dot{\mathbf{x}} = A\mathbf{x} + B\mathbf{u} x ˙ = A x + B u [ k Δ t , ( k + 1 ) Δ t ) [k\Delta t, (k+1)\Delta t) [ k Δ t , ( k + 1 ) Δ t )
x ( ( k + 1 ) Δ t ) = e A Δ t x ( k Δ t ) + ( ∫ 0 Δ t e A τ d τ ) B u k \mathbf{x}((k+1)\Delta t) = e^{A\Delta t}\mathbf{x}(k\Delta t) + \left(\int_0^{\Delta t} e^{A\tau}\, d\tau\right) B\, \mathbf{u}_k x (( k + 1 ) Δ t ) = e A Δ t x ( k Δ t ) + ( ∫ 0 Δ t e A τ d τ ) B u k
定义离散化后的矩阵：
A ˉ = e A Δ t , B ˉ = ( ∫ 0 Δ t e A τ d τ ) B = A − 1 ( e A Δ t − I ) B \bar{A} = e^{A\Delta t}, \qquad \bar{B} = \left(\int_0^{\Delta t} e^{A\tau}\, d\tau\right) B = A^{-1}(e^{A\Delta t} - I) B A ˉ = e A Δ t , B ˉ = ( ∫ 0 Δ t e A τ d τ ) B = A − 1 ( e A Δ t − I ) B
（当 A A A
则离散系统为：
x k + 1 = A ˉ x k + B ˉ u k \mathbf{x}_{k+1} = \bar{A}\mathbf{x}_k + \bar{B}\mathbf{u}_k x k + 1 = A ˉ x k + B ˉ u k
ZOH 离散化的关键性质：
精确性 ：在采样点上，离散系统给出的值与连续系统完全一致（不是近似）稳定性保持 ：如果连续系统稳定（Re ( λ i ) < 0 \text{Re}(\lambda_i) < 0 Re ( λ i ) < 0 ∣ e λ i Δ t ∣ < 1 |e^{\lambda_i \Delta t}| < 1 ∣ e λ i Δ t ∣ < 1
双线性变换用一种不同的方式将连续系统转换为离散系统。核心思想是用梯形法则近似积分：
s ≈ 2 Δ t ⋅ z − 1 z + 1 s \approx \frac{2}{\Delta t} \cdot \frac{z - 1}{z + 1} s ≈ Δ t 2 ⋅ z + 1 z − 1
这给出：
A ˉ = ( I − Δ t 2 A ) − 1 ( I + Δ t 2 A ) \bar{A} = \left(I - \frac{\Delta t}{2}A\right)^{-1}\left(I + \frac{\Delta t}{2}A\right) A ˉ = ( I − 2 Δ t A ) − 1 ( I + 2 Δ t A )
B ˉ = ( I − Δ t 2 A ) − 1 Δ t ⋅ B \bar{B} = \left(I - \frac{\Delta t}{2}A\right)^{-1} \Delta t \cdot B B ˉ = ( I − 2 Δ t A ) − 1 Δ t ⋅ B
双线性变换的关键性质：
保持频率响应 ：将连续系统的频率特性保真地映射到离散域绝对稳定性保持 ：连续系统的左半平面精确映射到离散系统的单位圆内部——这意味着稳定的连续系统离散化后一定 稳定（无论 Δ t \Delta t Δ t 这正是 S4 (Structured State Spaces) 论文中使用的离散化方法 （Gu, Goel & Ré, 2022）
离散化方法对比
方法 离散 A ˉ \bar{A} A ˉ 精确性 稳定性保持 用于 ZOH e A Δ t e^{A\Delta t} e A Δ t 采样点精确 是 控制工程（经典） 双线性变换 ( I − Δ t 2 A ) − 1 ( I + Δ t 2 A ) (I - \frac{\Delta t}{2}A)^{-1}(I + \frac{\Delta t}{2}A) ( I − 2 Δ t A ) − 1 ( I + 2 Δ t A ) 频率保真 无条件 S4 / SSM 欧拉前向 I + A Δ t I + A\Delta t I + A Δ t O ( Δ t ) O(\Delta t) O ( Δ t ) 不保证 （Δ t \Delta t Δ t 教学/简单仿真
离散化方法对比：从连续 A 到离散 Ā 连续 ẋ = Ax + Bu → 离散 x_{k+1} = Āx_k + B̄u_k 欧拉前向 Ā = I + AΔt 精度: O(Δt) 近似 稳定性: 不保证 教学 ZOH Ā = e^{AΔt} 精度: 采样点精确 稳定性: 保持 控制 / Mamba 双线性变换 Ā = (I-½ΔtA)⁻¹(I+½ΔtA) 精度: 频率保真 稳定性: 无条件保持 S4 / SSM 精度和稳定性保证递增 SSM/Mamba 通过离散化获得分辨率不变性——训练和推理可用不同 Δt
为什么 SSM/Mamba 关注离散化？ 因为 SSM 的核心参数（A , B , C , D A, B, C, D A , B , C , D Δ t \Delta t Δ t Δ t \Delta t Δ t
Kalman 滤波：连续世界的”Forward 算法”
问题设置
回顾 Art. 16 HMM 的核心问题：状态隐藏了，只能通过噪声观测来推断。 Kalman 滤波处理的是同一个问题，但在连续状态空间 中：
HMM（Art. 16） Kalman 滤波（本文） 隐状态 离散 X t ∈ { s 1 , … , s N } X_t \in \{s_1, \ldots, s_N\} X t ∈ { s 1 , … , s N } 连续 x t ∈ R n \mathbf{x}_t \in \mathbb{R}^n x t ∈ R n 状态转移 转移矩阵 A A A 状态方程 x k + 1 = A x k + w k \mathbf{x}_{k+1} = A\mathbf{x}_k + \mathbf{w}_k x k + 1 = A x k + w k 观测 发射矩阵 B B B 观测方程 y k = C x k + v k \mathbf{y}_k = C\mathbf{x}_k + \mathbf{v}_k y k = C x k + v k 推断方法 Forward 算法（求和） Kalman 滤波（高斯积分） 核心结构 传播 + 修正 预测 + 更新
State-Space 形式 ( A , B , C , D ) (A, B, C, D) ( A , B , C , D )
一个完整的线性 state-space 模型由四个矩阵描述：
x ˙ ( t ) = A x ( t ) + B u ( t ) \dot{\mathbf{x}}(t) = A\mathbf{x}(t) + B\mathbf{u}(t) x ˙ ( t ) = A x ( t ) + B u ( t ) y ( t ) = C x ( t ) + D u ( t ) \mathbf{y}(t) = C\mathbf{x}(t) + D\mathbf{u}(t) y ( t ) = C x ( t ) + D u ( t )
逐项解释：
x ( t ) ∈ R n \mathbf{x}(t) \in \mathbb{R}^n x ( t ) ∈ R n 状态向量 ——系统的内部状态，维度 n n n u ( t ) ∈ R m \mathbf{u}(t) \in \mathbb{R}^m u ( t ) ∈ R m 输入向量 ——外部驱动y ( t ) ∈ R p \mathbf{y}(t) \in \mathbb{R}^p y ( t ) ∈ R p 输出向量 ——可观测的量A ∈ R n × n A \in \mathbb{R}^{n \times n} A ∈ R n × n 状态矩阵 ——描述状态的自演化（没有输入时系统如何自行演化）B ∈ R n × m B \in \mathbb{R}^{n \times m} B ∈ R n × m 输入矩阵 ——描述输入如何影响状态（本文首次引入）C ∈ R p × n C \in \mathbb{R}^{p \times n} C ∈ R p × n 输出矩阵 ——描述状态如何被观测D ∈ R p × m D \in \mathbb{R}^{p \times m} D ∈ R p × m 直通矩阵 ——描述输入对输出的直接影响（很多系统中 D = 0 D = 0 D = 0
第一个方程是状态方程 ——描述系统的内部动态。第二个是观测方程 ——描述我们能看到什么。
与 HMM 的对应 ：HMM 的转移矩阵 A HMM A_{\text{HMM}} A HMM A A A B HMM B_{\text{HMM}} B HMM C C C B B B B B B
离散 Kalman 滤波
实际中，Kalman 滤波通常在离散时间上实现。离散化后的 state-space 模型为：
x k + 1 = A x k + B u k + w k , w k ∼ N ( 0 , Q ) \mathbf{x}_{k+1} = A\mathbf{x}_k + B\mathbf{u}_k + \mathbf{w}_k, \quad \mathbf{w}_k \sim \mathcal{N}(\mathbf{0}, Q) x k + 1 = A x k + B u k + w k , w k ∼ N ( 0 , Q ) y k = C x k + v k , v k ∼ N ( 0 , R ) \mathbf{y}_k = C\mathbf{x}_k + \mathbf{v}_k, \quad \mathbf{v}_k \sim \mathcal{N}(\mathbf{0}, R) y k = C x k + v k , v k ∼ N ( 0 , R )
其中 w k \mathbf{w}_k w k v k \mathbf{v}_k v k Q Q Q R R R
Kalman 滤波在每个时间步执行两步操作——与 HMM 的 Forward 算法结构完全对应：
步骤 1：预测（Predict） —— 对应 HMM Forward 的”传播”
x ^ k ∣ k − 1 = A x ^ k − 1 ∣ k − 1 + B u k − 1 \hat{\mathbf{x}}_{k|k-1} = A\hat{\mathbf{x}}_{k-1|k-1} + B\mathbf{u}_{k-1} x ^ k ∣ k − 1 = A x ^ k − 1∣ k − 1 + B u k − 1 P k ∣ k − 1 = A P k − 1 ∣ k − 1 A T + Q P_{k|k-1} = A P_{k-1|k-1} A^T + Q P k ∣ k − 1 = A P k − 1∣ k − 1 A T + Q
第一个方程用状态方程将上一步的最优估计 x ^ k − 1 ∣ k − 1 \hat{\mathbf{x}}_{k-1|k-1} x ^ k − 1∣ k − 1 先验估计 x ^ k ∣ k − 1 \hat{\mathbf{x}}_{k|k-1} x ^ k ∣ k − 1 估计协方差矩阵 P k − 1 ∣ k − 1 P_{k-1|k-1} P k − 1∣ k − 1 A A A Q Q Q
步骤 2：更新（Update） —— 对应 HMM Forward 的”修正”
K k = P k ∣ k − 1 C T ( C P k ∣ k − 1 C T + R ) − 1 K_k = P_{k|k-1} C^T (C P_{k|k-1} C^T + R)^{-1} K k = P k ∣ k − 1 C T ( C P k ∣ k − 1 C T + R ) − 1 x ^ k ∣ k = x ^ k ∣ k − 1 + K k ( y k − C x ^ k ∣ k − 1 ) \hat{\mathbf{x}}_{k|k} = \hat{\mathbf{x}}_{k|k-1} + K_k(\mathbf{y}_k - C\hat{\mathbf{x}}_{k|k-1}) x ^ k ∣ k = x ^ k ∣ k − 1 + K k ( y k − C x ^ k ∣ k − 1 ) P k ∣ k = ( I − K k C ) P k ∣ k − 1 P_{k|k} = (I - K_kC) P_{k|k-1} P k ∣ k = ( I − K k C ) P k ∣ k − 1
逐项解释：
K k K_k K k Kalman 增益 （Kalman gain）——核心权重矩阵y k − C x ^ k ∣ k − 1 \mathbf{y}_k - C\hat{\mathbf{x}}_{k|k-1} y k − C x ^ k ∣ k − 1 新息 （innovation）——实际观测与预测观测的差异x ^ k ∣ k \hat{\mathbf{x}}_{k|k} x ^ k ∣ k 后验估计 ——用新息修正先验估计P k ∣ k P_{k|k} P k ∣ k 后验协方差 ——更新后的不确定性（因为获得了观测信息，不确定性降低）
Kalman 增益的直觉
K k K_k K k
K k = P k ∣ k − 1 C T ⏟ 状态预测的不确定性投影到观测空间 ⋅ ( C P k ∣ k − 1 C T + R ) − 1 ⏟ 总不确定性（预测 + 观测）的逆 K_k = \underbrace{P_{k|k-1} C^T}_{\text{状态预测的不确定性投影到观测空间}} \cdot \underbrace{(C P_{k|k-1} C^T + R)^{-1}}_{\text{总不确定性（预测 + 观测）的逆}} K k = 状态预测的不确定性投影到观测空间 P k ∣ k − 1 C T ⋅ 总不确定性（预测 + 观测）的逆 ( C P k ∣ k − 1 C T + R ) − 1
直觉上：
如果观测噪声 R R R （传感器很准）：K k → C − 1 K_k \to C^{-1} K k → C − 1
如果观测噪声 R R R （传感器很差）：K k → 0 K_k \to 0 K k → 0
如果预测不确定性 P k ∣ k − 1 P_{k|k-1} P k ∣ k − 1 （模型不确定）：K k K_k K k
如果预测不确定性 P k ∣ k − 1 P_{k|k-1} P k ∣ k − 1 （模型很确定）：K k K_k K k
Kalman 增益 K：模型预测 vs 传感器观测的信任平衡 完全信模型 完全信观测 K → 0 K → C⁻¹ K R 小（传感器准） K R ≈ P（平衡） K R 大（传感器差） K_k = P_k|k-1 C^T (C P_k|k-1 C^T + R)^{-1} ← 预测不确定性 P 大则 K 大，观测噪声 R 大则 K 小
“预测 + 更新”与 HMM “传播 + 修正”的统一结构 ：
步骤 HMM Forward Kalman 滤波 传播/预测 乘转移矩阵：α ′ = α k − 1 A \alpha' = \alpha_{k-1} A α ′ = α k − 1 A 状态方程推进：x ^ k ∣ k − 1 = A x ^ k − 1 ∣ k − 1 \hat{x}_{k \mid k-1} = A \hat{x}_{k-1 \mid k-1} x ^ k ∣ k − 1 = A x ^ k − 1 ∣ k − 1 修正/更新 乘发射概率（逐元素） 用观测新息修正：加 K k K_k K k 数学本质 离散概率的乘法 高斯分布的贝叶斯更新
结构完全一致：先用系统模型推进一步 ，再用观测数据修正估计 。HMM 用离散概率的乘法实现修正，Kalman 用高斯分布的贝叶斯公式实现修正。这就是 Art. 16 HMM 末尾提到的”传播 + 修正范式”在连续状态空间的自然延伸。
下图展示了 Kalman 滤波的”预测-更新”循环结构，以及它与 HMM Forward 算法的对应关系。
Kalman 滤波的"预测-更新"循环
与 HMM Forward 的"传播-修正"结构完全对应
x̂ₖ₋₁|ₖ₋₁ Pₖ₋₁|ₖ₋₁ 上一步后验 ① 预测 (Predict) x̂ₖ|ₖ₋₁ = A·x̂ₖ₋₁|ₖ₋₁ Pₖ|ₖ₋₁ = A·P·Aᵀ + Q 不确定性增大 ② 更新 (Update) Kₖ = P·Cᵀ·(C·P·Cᵀ+R)⁻¹ x̂ₖ|ₖ = x̂ₖ|ₖ₋₁ + Kₖ·(yₖ − C·x̂) 观测修正，不确定性减小 yₖ 观测数据 后验成为下一步的先验 → 循环 HMM 对应：① = 乘转移矩阵 A（传播） ② = 乘发射概率 b(oₜ)（修正）
Inline Box：Cholesky 分解在协方差计算中的角色
Kalman 增益的计算涉及矩阵求逆 ( C P k ∣ k − 1 C T + R ) − 1 (C P_{k|k-1} C^T + R)^{-1} ( C P k ∣ k − 1 C T + R ) − 1 Cholesky 分解 。
因为 S k = C P k ∣ k − 1 C T + R S_k = C P_{k|k-1} C^T + R S k = C P k ∣ k − 1 C T + R 对称正定矩阵 （协方差矩阵之和），它可以唯一分解为：
S k = L k L k T S_k = L_k L_k^T S k = L k L k T
其中 L k L_k L k S k K k T = ( C P k ∣ k − 1 ) T S_k K_k^T = (C P_{k|k-1})^T S k K k T = ( C P k ∣ k − 1 ) T
Cholesky 分解的优势：
计算量是一般 LU 分解的一半（1 3 n 3 \frac{1}{3}n^3 3 1 n 3 2 3 n 3 \frac{2}{3}n^3 3 2 n 3
保证数值稳定性（利用了正定性）
保持协方差矩阵的正定性（避免舍入误差导致”负方差”）
这种技巧在高维 Kalman 滤波（如 GPS 定位、SLAM）中至关重要。“square-root Kalman filter”就是将整个算法改写为 Cholesky 因子的传播，而非协方差矩阵本身的传播。
Kalman 滤波的最优性
Kalman (1960) 证明了一个深刻的结果：在线性系统 + 高斯噪声的假设下，Kalman 滤波是最小均方误差（MMSE）最优估计器。 即：
x ^ k ∣ k = E [ x k ∣ y 0 , y 1 , … , y k ] \hat{\mathbf{x}}_{k|k} = \mathbb{E}[\mathbf{x}_k \mid \mathbf{y}_0, \mathbf{y}_1, \ldots, \mathbf{y}_k] x ^ k ∣ k = E [ x k ∣ y 0 , y 1 , … , y k ]
这不是某种启发式——它是给定所有过去观测后对状态的条件期望 ，在均方误差意义下不可能被任何其他估计器超越。这个最优性源于高斯分布在线性变换下的封闭性：高斯的先验乘以高斯的似然，后验仍然是高斯——所以只需要追踪均值和协方差就足够了。
与 SSM/Mamba 的 State-Space 公式对照
Kalman 滤波的 state-space 公式与 SSM/Mamba（Art. 26 SSM/Mamba ）使用的公式是同一套公式 ，只是语境和目标不同：
Kalman 滤波 SSM / Mamba 连续形式 x ˙ = A x + B u \dot{\mathbf{x}} = A\mathbf{x} + B\mathbf{u} x ˙ = A x + B u y = C x \mathbf{y} = C\mathbf{x} y = C x x ˙ = A x + B u \dot{\mathbf{x}} = A\mathbf{x} + B\mathbf{u} x ˙ = A x + B u y = C x + D u \mathbf{y} = C\mathbf{x} + D\mathbf{u} y = C x + D u 离散形式 x k + 1 = A ˉ x k + B ˉ u k \mathbf{x}_{k+1} = \bar{A}\mathbf{x}_k + \bar{B}\mathbf{u}_k x k + 1 = A ˉ x k + B ˉ u k x k + 1 = A ˉ x k + B ˉ u k \mathbf{x}_{k+1} = \bar{A}\mathbf{x}_k + \bar{B}\mathbf{u}_k x k + 1 = A ˉ x k + B ˉ u k A A A 由物理系统决定 可学习参数 （约束为对角矩阵）B , C B, C B , C 由物理系统决定 可学习参数 离散化 ZOH（经典） 双线性变换（S4）/ ZOH（Mamba） 目标 从噪声观测中推断 隐状态 从输入序列预测 输出序列 核心关注点 最优估计（最小方差） 高效序列建模（线性时间复杂度）
关键差异 ：Kalman 滤波中 ( A , B , C ) (A, B, C) ( A , B , C ) ( A , B , C ) (A, B, C) ( A , B , C ) 学习 的参数——这正是 Part 3 “汇——学习算子”的核心转变。
关键连接 ：SSM/Mamba 将 A A A 对角矩阵 （A = diag ( a 1 , … , a n ) A = \text{diag}(a_1, \ldots, a_n) A = diag ( a 1 , … , a n ) e A t = diag ( e a 1 t , … , e a n t ) e^{At} = \text{diag}(e^{a_1 t}, \ldots, e^{a_n t}) e A t = diag ( e a 1 t , … , e a n t ) A ˉ \bar{A} A ˉ Art. 2 特征分解 中”对角化简化一切”这个洞察的终极应用。
总结与展望
本文从离散步进走向连续流动，完成了 Part 2 时序算子子线的最后一站：
连续线性系统 x ˙ = A x \dot{\mathbf{x}} = A\mathbf{x} x ˙ = A x x ( t ) = e A t x ( 0 ) \mathbf{x}(t) = e^{At}\mathbf{x}(0) x ( t ) = e A t x ( 0 ) 矩阵指数的计算 依赖对角化：e A t = Q diag ( e λ 1 t , … , e λ n t ) Q − 1 e^{At} = Q\,\text{diag}(e^{\lambda_1 t}, \ldots, e^{\lambda_n t})\, Q^{-1} e A t = Q diag ( e λ 1 t , … , e λ n t ) Q − 1 A k = Q Λ k Q − 1 A^k = Q\Lambda^k Q^{-1} A k = Q Λ k Q − 1 稳定性判据 从离散的 ∣ λ i ∣ < 1 |\lambda_i| < 1 ∣ λ i ∣ < 1 Re ( λ i ) < 0 \text{Re}(\lambda_i) < 0 Re ( λ i ) < 0 复特征值 带来振荡行为——实部控制包络增长/衰减，虚部控制振荡频率离散化 （ZOH、双线性变换）将连续系统转换为可在计算机上实现的离散系统——这是 SSM/Mamba 的关键步骤Kalman 滤波 是 HMM Forward 算法在连续状态空间的自然延伸——“预测 + 更新”对应”传播 + 修正”，是线性高斯系统的 MMSE 最优估计State-Space 公式 ( A , B , C , D ) (A, B, C, D) ( A , B , C , D )
统一线索：特征分解贯穿始终。 A k = Q Λ k Q − 1 A^k = Q\Lambda^k Q^{-1} A k = Q Λ k Q − 1 e A t = Q e Λ t Q − 1 e^{At} = Qe^{\Lambda t}Q^{-1} e A t = Q e Λ t Q − 1 A ˉ = e A Δ t \bar{A} = e^{A\Delta t} A ˉ = e A Δ t Art. 2 特征分解 → 本文 → Art. 26 SSM/Mamba 。
至此，Part 2 的时序算子子线画上句号。从最简单的马尔可夫链（离散、可观测）到 HMM（离散、隐藏）到本文的连续系统和 Kalman 滤波——我们沿着概念地图走完了”离散 → 连续”和”可观测 → 隐藏”的两个维度。
下一篇 ，我们进入 Part 2 的第二条子线——图/空间算子。PageRank 将马尔可夫链的思想应用于互联网的链接图：把整个万维网看作一个巨大的马尔可夫链，网页是状态，超链接定义转移概率，稳态分布就是每个网页的”重要性”。它是时序算子子线和图算子子线的交汇点 ——同时身兼两个身份，连接我们已经走过的路和即将踏入的新领域。