第 2 章介绍的参数估计数学模型描述了观测值与参数之间的关系，是静态估计模型。静态估计模型不考虑研究对象自身的运动规律，仅用观测值来估计参数。而事实上，我们的研究对象是一个动态系统，有自身的运动特性和规律。如果用数学的语言来描述这些运动规律，就可以在估计时利用这些信息。本章首先介绍描述运动规律的微分方程并将其线性化，得到随机线性连续系统的数学模型。然后，本章将介绍状态转移矩阵和随机连续线性系统的解，并将随机连续线性系统离散化为特定时间点上的状态。最后，本章给出随机线性系统的可控性和可测性的概念和判定条件。
连续线性系统的数学模型
微分方程是描述某对象连续变化的数学语言，也被称为连续动态系统的状态方程。本节先以两个例子来说明如何利用微分方程对运动系统进行描述，得到连续动态系统状态方程的一般表达式，然后对连续动态系统状态方程进行线性化，最后得到线性连续系统的数学模型。
连续系统的数学模型
例 3.1如图 3.1 所示，设一个质量为 \(m\) 的长方体被弹簧吊着，从平衡状况开始，在空气中作垂直方向的振动，用 \(y(t)\) 表示质点在时刻 \(t\) 的位置，质点在运动中受到的力有；弹簧的恢复力 \(F_K\)、阻尼器阻力 \(F_V(t)\) 和外力 \(F(t)\)。弹簧的恢复力 \(F_K\) 与位移成正比，方向与位移相反，大小为 \(-K_y(t)\)（\(K>0\) 为常数）；阻尼器阻力 \(F_V(t)\) 与速度成正比，方向相反，大小为 \(-f\dot{y}(t)\)（\(f>0\) 为常数）。试建立弹簧阻尼系统的运动方程。
解：由牛顿第二定律得到弹簧的运动方程为 \[\frac{-K}{m}y(t)-\frac{f}{m}\dot{y}(t)+\frac{1}{m}F(t)=\ddot{y}(t) \tag{3.1.1}\] 此外，小球还可能受到其他未知的作用力，考虑这些未知或者不确定的作用力，式 (3.1.1) 表示为 \[\frac{-K}{m}y(t)-\frac{f}{m}\dot{y}(t)+\frac{1}{m}F(t)+e(t)=\ddot{y}(t) \tag{3.1.2}\] 随机变量 \(e(t)\) 用来补偿这些未知或者不确定的作用力。二阶微分方程 (3.1.2) 描述了弹簧阻尼系统的运动规律。这样的二阶线性微分方程可以转换为一阶线性微分方程。设状态变量 \[\bm{X}(t)=\left[\begin{array}{ll}X_1(t) & X_2(t)\end{array}\right]^{\mathrm{T}}=\left[\begin{array}{ll}y(t) & \dot{y}(t)\end{array}\right]^{\mathrm{T}} \tag{3.1.3}\]
状态变量的本质是系统的“记忆”：对弹簧阻尼系统，只需记录位置 \(y(t)\) 和速度 \(\dot{y}(t)\) 两个量，就足以依据运动方程推算此后任意时刻的整个运动过程——再多就是冗余，再少就信息不全。这和第 2 章的静态参数估计有本质区别：那里 \(\bm{X}\) 是一组不变的数字，估计一次即可；这里 \(\bm{X}(t)\) 本身随 \(t\) 演化，估计要逐时刻反复进行，这正是动态估计与静态估计的分水岭。把二阶微分方程改写为 \(\bm{X}(t)=\left[\begin{array}{ll}y(t) & \dot{y}(t)\end{array}\right]^{\mathrm{T}}\) 的形式，等于把系统的“记忆”显式地装进了一个向量里，也为下面将最高阶导数“藏”进一阶方程组做了准备。
式 (3.1.2) 可以表示为 \[\begin{aligned} \dot{X}_1(t)&=X_2(t)\\ \dot{X}_2(t)&=-\frac{K}{m}X_1(t)-\frac{f}{m}X_2(t)+\frac{1}{m}F(t)+e(t) \end{aligned} \tag{3.1.4}\] 将上式表达为向量形式 \[\frac{\mathrm{d}}{\mathrm{d}t}\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix} =\begin{bmatrix}0 & 1\\[4pt] -\dfrac{K}{m} & -\dfrac{f}{m}\end{bmatrix} \begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix} +\begin{bmatrix}0\\[4pt] \dfrac{1}{m}\end{bmatrix}F(t)+\begin{bmatrix}0\\ 1\end{bmatrix}e(t) \tag{3.1.5}\] 这样就将一个动态系统的二阶微分方程转化为了一组一阶微分方程组，它表达了外力 \(F(t)\)（即输入）与状态变量 \(\bm{X}(t)\)，以及状态变量之间的关系。
补出“二阶微分方程化为一阶微分方程组”的关键跳步。设 \(X_1(t)=y(t)\)、\(X_2(t)=\dot{y}(t)\) 后，第一式 \(\dot{X}_1(t)=X_2(t)\) 只是恒等式（速度是位置的导数），真正的动力学内容全部压进第二式 \(\dot{X}_2(t)=-\dfrac{K}{m}X_1(t)-\dfrac{f}{m}X_2(t)+\dfrac{1}{m}F(t)+e(t)\)。换言之，换元并不改变物理内容，只是把“谁是谁的导数”这层关系显式写成方程组，使最高阶导数只出现在一个方程的左端。式 (3.1.8) 中 \(\bm{A}(t)\) 的“友矩阵”形式正是这种换元的一般化：主对角线上方一排 1 就是一连串 \(\dot{X}_i=X_{i+1}\) 的恒等式，最后一行才承载真正的微分方程。理解这一点，后面的状态转移矩阵（3.2 节）就是“解”这个换元后的方程组。
例 3.2设有 \(n\) 阶线性微分方程 \[\frac{\mathrm{d}^ny(t)}{\mathrm{d}t^n}=-a_0(t)y-a_1(t)\frac{\mathrm{d}y(t)}{\mathrm{d}t}-\cdots-a_{n-1}(t)\frac{\mathrm{d}^{n-1}y}{\mathrm{d}t^{n-1}}+e(t) \tag{3.1.6}\] 其中 \(e(t)\) 为补偿系统不确定因素的随机部分。设置状态变量，将其转化为一组一阶微分方程。
解：设 \[\begin{aligned} \bm{X}(t)&=\left[\begin{array}{llll}y(t) & \dot{y}(t) & \cdots & y^{(n-1)}(t)\end{array}\right]^{\mathrm{T}}\\ &=\left[\begin{array}{llll}X_1(t) & X_2(t) & \cdots & X_n(t)\end{array}\right]^{\mathrm{T}} \end{aligned} \tag{3.1.7}\] 那么 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\\ \vdots\\ \vdots\\ \dot{X}_n(t)\end{bmatrix} =\begin{bmatrix} 0 & 1 & 0 & \cdots & 0\\ 0 & 0 & 1 & \cdots & 0\\ \vdots & \vdots & \vdots & & \vdots\\ 0 & 0 & 0 & \cdots & 1\\ -a_0(t) & -a_1(t) & -a_2(t) & \cdots & -a_{n-1}(t) \end{bmatrix} \begin{bmatrix}X_1(t)\\ X_2(t)\\ \vdots\\ X_{n-1}(t)\\ X_n(t)\end{bmatrix} +\begin{bmatrix}0\\ 0\\ \vdots\\ 0\\ 1\end{bmatrix}e(t) \tag{3.1.8}\]
综合例 3.1 和例 3.2，如果设控制输入向量为 \[\underset{p\times 1}{\bm{u}(t)}=\left[\begin{array}{llll}u_1(t) & u_2(t) & \cdots & u_p(t)\end{array}\right]^{\mathrm{T}} \tag{3.1.9}\] 微分方程中不确定的随机部分为 \[\underset{q\times 1}{\bm{e}(t)}=\left[\begin{array}{llll}e_1(t) & e_2(t) & \cdots & e_q(t)\end{array}\right]^{\mathrm{T}} \tag{3.1.10}\] 并且设 \(\bm{X}(t)\) 前的系数矩阵为 \(\bm{A}(t)\)，控制输入 \(\bm{u}(t)\) 前的系数矩阵为 \(\bm{B}(t)\)，\(\bm{e}(t)\) 前的系数为 \(\bm{C}(t)\)，那么微分方程的一般表达形式为 \[\underset{n\times 1}{\dot{\bm{X}}(t)}=\underset{n\times n}{\bm{A}(t)}\ \underset{n\times 1}{\bm{X}(t)} +\underset{n\times p}{\bm{B}(t)}\ \underset{p\times 1}{\bm{u}(t)}+\underset{n\times q}{\bm{C}(t)}\ \underset{q\times 1}{\bm{e}(t)} \tag{3.1.11}\] 随机变量部分 \(\bm{e}(t)\) 也称为系统噪声或者状态噪声。假设 \(\bm{e}(t)\) 为白噪声过程，即 \[E\left[\,\bm{e}(t)\,\right]=0 \tag{3.1.12}\] \[\mathrm{Cov}\left[\,\bm{e}(t),\ \bm{e}(\tau)\,\right]=\bm{D}_e(t)\delta(t-\tau) \tag{3.1.13}\] 其中 \(\bm{D}_e(t)\) 为 \(q\times q\) 维对称非负定矩阵，是系统噪声 \(\bm{e}(t)\) 的均方值；\(\delta(t-\tau)\) 为狄拉克-\(\delta\) 函数。式 (3.1.11) (3.1.13) 描述了系统的运动规律，它不仅给出了外部输入 \(\bm{u}(t)\) 与状态 \(\bm{X}(t)\) 的关系，也表明状态之间的关系以及不确定因素 \(\bm{e}(t)\) 对状态的影响。
白噪声假定的两处易错点。其一，式 (3.1.13) 中的 \(\delta(t-\tau)\) 是 Dirac-\(\delta\) 函数，它严格表达“不同时刻的噪声互不相关”——只要 \(t\neq\tau\)，协方差为零；\(\bm{D}_e(t)\) 刻画同一时刻噪声的强度，即均方值。其二，只有当 \(\bm{e}(t)\) 是白噪声时，后续方差传播（如 3.3 节 \(\bm{D}_w(k-1)\) 的积分式）中交叉项才全部消失；若 \(\bm{e}(t)\) 为有色噪声，\(\bm{D}_e(t)\delta(t-\tau)\) 的简洁形式不再成立，交叉协方差无法消去。还要注意 \(\bm{D}_e(t)\) 只要求对称非负定（允许某些方向没有噪声激励），这与观测噪声方差阵要求正定（滤波中需取逆）形成对比，详见 3.3 节的随机模型。
本节的状态方程一般形式 (3.1.11)，与《广义测量平差》§4-1“连续线性系统的状态方程和观测方程”中 (4-1-3)、(4-1-4) 两式完全同构：控制输入项 \(\bm{B}(t)\bm{u}(t)\) 对应该书 (4-1-3) 的 \(\bm{C}(t)\bm{U}(t)\)，过程噪声项 \(\bm{C}(t)\bm{e}(t)\) 对应该书 (4-1-3) 的 \(\bm{F}(t)\bm{\varOmega}(t)\)。该书用“卫星轨道角位置偏差”作例，本节用弹簧阻尼与卫星轨道作例，两者都说明同一件事：状态方程与观测方程合起来构成动态系统的函数模型。《广义测量平差》§4-1 同时给出随机模型 (4-1-28) (4-1-32)，与本节的 (3.1.12)、(3.1.13) 对应，可互相参照。
如果希望了解系统从某一初始时刻 \(t_0\) 之后任意时刻的状态，就需要知道系统在 \(t_0\) 时刻的状态值 \(\bm{X}(t_0)=\bm{X}_0\)，即状态的初始条件。式 (3.1.11)、式 (3.1.12) 和式 (3.1.13) 加上初始值 \(\bm{X}_0\) 称为随机连续系统的状态方程。对于一个 \(n\) 阶微分方程的动态系统至少需要 \(n\) 个状态变量来描述，这些状态变量应选取最少但能够描述系统必需的变量。当给出某时间点上的状态变量值，就可以确定其他任何时间点上的状态值。
如果式 (3.1.11) 中的 \(\bm{A}(t)\)、\(\bm{B}(t)\) 和 \(\bm{C}(t)\) 不随时间变化，就退化为常数矩阵 \(\bm{A}\)、\(\bm{B}\) 和 \(\bm{C}\)，式 (3.1.11) 成为 \[\dot{\bm{X}}(t)=\bm{A}\bm{X}(t)+\bm{B}\bm{u}(t)+\bm{C}\bm{e}(t) \tag{3.1.14}\] 上式被称为连续时不变系统（定常系统）的状态方程。如果式 (3.1.11) 的控制输入为零，即动态系统无输入，或者不考虑系统的输入，那么无输入的状态方程为 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{C}(t)\bm{e}(t) \tag{3.1.15}\] 它表明系统本身在无外力的作用下自由运动。在没有控制输入，并且 \(\bm{A}(t)\) 和 \(\bm{C}(t)\) 不随时间变化的情况下，状态方程为 \[\dot{\bm{X}}(t)=\bm{A}\bm{X}(t)+\bm{C}\bm{e}(t) \tag{3.1.16}\] 上式即为无控制输入的时不变系统。
式 (3.1.11) (3.1.16) 展示了状态方程的三级“退化”：一般时变 \(\bm{A}(t)\) 降为时不变 \(\bm{A}\)，有控制输入 \(\bm{u}(t)\) 降为无输入。从“处处不同”到“处处相同”、从“可驾驶”到“自由滑行”，每降一级计算就简化一档。其中无输入时不变形式 \(\dot{\bm{X}}(t)=\bm{A}\bm{X}(t)+\bm{C}\bm{e}(t)\) 是最常用的一档：它刻画系统在无外力作用下靠自身惯性运动的方式，3.2 节的状态转移矩阵正是为求解这一档方程而引入的。实际估计问题（如卫星导航、惯性导航）中，常把未模型化的力统一打包进 \(\bm{e}(t)\)，而不必显式建模 \(\bm{u}(t)\)——这正是为什么很多应用直接采用无输入模型。
以上给出的都是线性的微分方程，但在现实应用中，对系统运动规律的描述更多的是非线性的微分方程。
例 3.3卫星的轨迹可以用 \(r(t)\) 和 \(\theta(t)\) 两个极坐标变量表示。其中 \(r(t)\) 是卫星到地心的距离，\(\theta(t)\) 是卫星和地心的连线相对于参考坐标轴的角度。假定卫星具有在轨道径向和切向的推力控制 \(u_r(t)\) 和 \(u_l(t)\)，试建立卫星轨迹控制系统的模型。
解：根据力学规律，卫星的运动方程可以写为 \[\begin{aligned} \ddot{r}(t)&=r(t)\left[\,\dot{\theta}(t)\,\right]^2-\frac{G}{r^2(t)}+u_r(t)+e_r(t)\\[4pt] \ddot{\theta}(t)&=-\frac{2}{r(t)}\dot{\theta}(t)\dot{r}(t)+\frac{1}{r(t)}u_l(t)+e_\theta(t) \end{aligned} \tag{3.1.17}\]
原书式 (3.1.17) 第一式引力项排印为 \(-G/r^2(t)\)，而式 (3.1.19)、例 3.4 的 \(\bm{A}(t)\) 矩阵中相应项分别为 \(-GM/X_1^2(t)\)、\(2GM/[\,X_1^{*}(t)\,]^3\)（\(G\) 为万有引力常数，\(M\) 为地球质量），前后不一致（引力加速度应为 \(GM/r^2\)）。此处均照原样排印。
其中，\(G\) 为万有引力常数，\(M\) 为地球的质量；随机变量 \(e_r(t)\) 和 \(e_\theta(t)\) 用来补偿卫星在轨道径向和切向方向上未知或者难以描述的作用力。选择状态变量为 \[\bm{X}(t)=\begin{bmatrix}X_1(t)\\ X_2(t)\\ X_3(t)\\ X_4(t)\end{bmatrix} =\begin{bmatrix}r(t)\\ \dot{r}(t)\\ \theta(t)\\ \dot{\theta}(t)\end{bmatrix} \tag{3.1.18}\] 可将式 (3.1.17) 表达成向量的形式 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\\ \dot{X}_3(t)\\ \dot{X}_4(t)\end{bmatrix} =\begin{bmatrix} X_2(t)\\[6pt] X_1(t)X_4^2(t)-\dfrac{GM}{X_1^2(t)}+u_r(t)\\[10pt] X_4(t)\\[6pt] -\dfrac{2}{X_1(t)}X_4(t)X_2(t)+\dfrac{1}{X_1(t)}u_l(t) \end{bmatrix} +\begin{bmatrix}0 & 0\\ 1 & 0\\ 0 & 0\\ 0 & 1\end{bmatrix} \begin{bmatrix}e_r(t)\\ e_\theta(t)\end{bmatrix} \tag{3.1.19}\] 显然式 (3.1.19) 是非线性的微分方程式。将非线性的微分方程表达为更一般形式 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\\ \vdots\\ \dot{X}_n(t)\end{bmatrix} =\begin{bmatrix} g_1\left(X_1(t),\ \cdots X_n(t),\ u_1(t),\ \cdots u_p(t),\ e_1(t),\ \cdots e_q(t)\right)\\ g_2\left(X_1(t),\ \cdots X_n(t),\ u_1(t),\ \cdots u_p(t),\ e_1(t),\ \cdots e_q(t)\right)\\ \vdots\\ g_n\left(X_1(t),\ \cdots X_n(t),\ u_1(t),\ \cdots u_p(t),\ e_1(t),\ \cdots e_q(t)\right) \end{bmatrix} \tag{3.1.20}\] 或表示为 \[\underset{n\times 1}{\dot{\bm{X}}(t)}=\underset{n\times 1}{\bm{g}}\left[\,\underset{n\times 1}{\bm{X}(t)},\ \underset{p\times 1}{\bm{u}(t)},\ \underset{q\times 1}{\bm{e}(t)}\,\right] \tag{3.1.21}\] 其中 \[\underset{n\times 1}{\bm{g}\left[\,\cdot\,\right]}=\begin{bmatrix}g_1(\cdot)\\ g_2(\cdot)\\ \vdots\\ g_n(\cdot)\end{bmatrix} \tag{3.1.22}\] 和 \[\begin{aligned} E\left[\,\bm{e}(t)\,\right]&=0\\ \mathrm{Cov}\left[\,\bm{e}(t),\ \bm{e}(\tau)\,\right]&=\bm{D}_e(t)\delta(t-\tau) \end{aligned} \tag{3.1.23}\] 如果还对系统进行观测，观测方程为 \[\underset{\ell\times 1}{\bm{Z}(t)}=\underset{\ell\times 1}{\bm{F}\left[\,\bm{X}(t)\,\right]}+\underset{\ell\times 1}{\bm{\Delta}(t)} \tag{3.1.24}\] 其中 \(\bm{\Delta}(t)\) 为白噪声过程，与过程噪声 \(\bm{e}(t)\) 无关： \[\begin{aligned} E\left[\,\bm{\Delta}(t)\,\right]&=0\\ \mathrm{Cov}\left[\,\bm{\Delta}(t),\ \bm{\Delta}(\tau)\,\right]&=\bm{D}_{\Delta}(t)\delta(t-\tau) \end{aligned} \tag{3.1.25}\] 和 \[\mathrm{Cov}\left[\,\bm{\Delta}(t),\ \bm{e}(\tau)\,\right]=\bm{0} \tag{3.1.26}\] 其中，\(\bm{D}_{\Delta}(t)\) 为 \(\bm{\Delta}(t)\) 的均方值。以上的状态方程和观测方程一起构成了随机连续系统的数学模型。对于这样的非线性微分方程和观测方程，首先需要将其线性化，得到线性的数学模型。
连续系统数学模型的线性化
在许多情况下，动态系统的数学模型都是状态变量的非线性函数，有随机干扰的状态方程为 \[\dot{\bm{X}}(t)=\bm{g}\left[\,\bm{X}(t),\ \bm{u}(t),\ \bm{e}(t)\,\right] \tag{3.1.27}\] 在第 2 章中我们讨论了用泰勒级数法展开观测方程，并舍去高阶项来得到线性的观测方程。这里采用同样的方法将随机线性系统模型线性化。
假设 \(\bm{X}(t)\) 的近似值（参考值）为 \(\bm{X}^{*}(t)\)，系统噪声 \(\bm{e}(t)\) 的近似值为 \(\bm{e}^{*}(t)\)，用泰勒级数将式 (3.1.27) 在近似值处展开，并忽略高阶项（二阶和二阶以上项）： \[\dot{\bm{X}}(t)=\bm{g}\left[\begin{array}{lll}\bm{X}^{*}(t) & \bm{u}(t) & \bm{e}^{*}(t)\end{array}\right] +\left[\frac{\partial\bm{g}}{\partial\bm{X}(t)}\right]^{*}\left[\,\bm{X}(t)-\bm{X}^{*}(t)\,\right] +\left[\frac{\partial\bm{g}}{\partial\bm{e}(t)}\right]^{*}\left[\,\bm{e}(t)-\bm{e}^{*}(t)\,\right] \tag{3.1.28}\] 这里的 \(\left[\,\cdot\,\right]^{*}\) 表示将近似值代入得到的雅可比矩阵。由于 \(\bm{e}(t)\) 的期望为零，所以这里取得近似值为 \(\bm{e}^{*}(t)=0\)，得到 \[\dot{\bm{X}}(t)=\bm{g}\left[\begin{array}{lll}\bm{X}^{*}(t) & \bm{u}(t) & 0\end{array}\right] +\left[\frac{\partial\bm{g}}{\partial\bm{X}(t)}\right]^{*}\left[\,\bm{X}(t)-\bm{X}^{*}(t)\,\right] +\left[\frac{\partial\bm{g}}{\partial\bm{e}(t)}\right]^{*}\bm{e}(t) \tag{3.1.29}\] 令 \[\begin{aligned} \bm{A}(t)&=\left[\frac{\partial\bm{g}}{\partial\bm{X}(t)}\right]^{*}\\ \bm{C}(t)&=\left[\frac{\partial\bm{g}}{\partial\bm{e}(t)}\right]^{*} \end{aligned} \tag{3.1.30}\] 得到 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{g}\left[\begin{array}{lll}\bm{X}^{*}(t) & \bm{u}(t) & 0\end{array}\right]-\bm{A}(t)\bm{X}^{*}(t)+\bm{C}(t)\bm{e}(t) \tag{3.1.31}\] 上式的 \(\bm{g}\left[\begin{array}{lll}\bm{X}^{*}(t) & \bm{u}(t) & 0\end{array}\right]-\bm{A}(t)\bm{X}^{*}(t)\) 是关于 \(\bm{X}^{*}(t)\) 和 \(\bm{u}(t)\) 的函数，也是已知的常数项，为简单起见，在后面的推导中用 \(\bm{G}(t)\) 来代替所有的与状态和系统噪声无关的项 \[\bm{G}(t)\rightarrow\bm{g}\left[\begin{array}{lll}\bm{X}^{*}(t) & \bm{u}(t) & 0\end{array}\right]-\bm{A}(t)\bm{X}^{*}(t) \tag{3.1.32}\] 式 (3.1.31) 为 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{G}(t)+\bm{C}(t)\bm{e}(t) \tag{3.1.33}\] 这样就得到了线性的微分方程。
补出线性化的两个关键跳步。第一，泰勒展开在近似值 \(\bm{X}^{*}(t)\) 处进行并舍去二阶及以上项，得到式 (3.1.28)；\(\left[\,\cdot\,\right]^{*}\) 表示雅可比矩阵在 \(\bm{X}^{*}(t)\) 处取值，所以 \(\bm{A}(t)\)、\(\bm{C}(t)\) 是“当时当地的切平面斜率”，随参考轨迹 \(\bm{X}^{*}(t)\) 的变化而变化。第二，把只含 \(\bm{X}^{*}(t)\)、\(\bm{u}(t)\) 的已知项 \(\bm{g}\left[\bm{X}^{*}(t),\ \bm{u}(t),\ 0\right]-\bm{A}(t)\bm{X}^{*}(t)\) 合并为 \(\bm{G}(t)\)（式 (3.1.32)），这一确定性驱动项在 3.3 节离散化时被吸收为 \(\bm{\Omega}(k-1)\)。注意展开点取 \(\bm{e}^{*}(t)=0\) 是因为 \(E[\bm{e}(t)]=\bm{0}\)——零均值白噪声的期望就是它的“最可能取值”。
如果 \(\bm{G}(t)\) 是 \(\bm{u}(t)\) 的线性函数，那么可以将式 (3.1.33) 表示为 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t) \tag{3.1.34}\]
同样，对于非线性的观测方程 (3.1.24) 来说，以 \(\bm{X}^{*}(t)\) 作为近似值将观测方程展开，舍去二阶和二阶以上的高阶项，得到 \[\bm{Z}(t)=\bm{F}\left[\,\bm{X}^{*}(t)\,\right] +\left.\frac{\partial\bm{F}\left[\,\bm{X}(t)\,\right]}{\partial\bm{X}(t)}\right|_{\bm{X}(t)=\bm{X}^{*}(t)} \left[\,\bm{X}(t)-\bm{X}^{*}(t)\,\right]+\bm{\Delta}(t) \tag{3.1.35}\] 设 \[\bm{H}(t)=\left.\frac{\partial\bm{F}\left[\,\bm{X}(t)\,\right]}{\partial\bm{X}(t)}\right|_{\bm{X}(t)=\bm{X}^{*}(t)} \tag{3.1.36}\] 并将已知常量部分都移到方程的左边 \[\bm{Z}(t)-\bm{F}\left[\,\bm{X}^{*}(t)\,\right]+\bm{H}(t)\bm{X}^{*}(t)=\bm{H}(t)\bm{X}(t)+\bm{\Delta}(t) \tag{3.1.37}\] 令 \[\bm{z}(t)=\bm{Z}(t)-\bm{F}\left[\,\bm{X}^{*}(t)\,\right]+\bm{H}(t)\bm{X}^{*}(t) \tag{3.1.38}\] 最后，有观测方程 \[\bm{z}(t)=\bm{H}(t)\bm{X}(t)+\bm{\Delta}(t) \tag{3.1.39}\]
线性化模型是“局部”模型：雅可比矩阵在 \(\bm{X}^{*}(t)\) 处取值，只有当真实状态始终靠近参考轨迹 \(\bm{X}^{*}(t)\) 时，舍去高阶项才成立；一旦偏差积累过大（如卫星轨道偏差逐步增大），一阶近似的误差将不可忽略。这正是后续滤波中“线性化+再线性化”迭代的动机——每一轮用新的估计当展开点重算 \(\bm{A}(t)\)、\(\bm{C}(t)\)、\(\bm{H}(t)\)。还要注意式 (3.1.34) 中 \(\bm{B}(t)\) 的出现是有条件的：\(\bm{G}(t)\) 必须是 \(\bm{u}(t)\) 的线性函数；若控制输入以非线性方式进入模型，则不能直接写成 \(\bm{B}(t)\bm{u}(t)\) 的形式，只能保留 \(\bm{G}(t)\) 的一般形式。
本节“泰勒级数展开舍高阶项”的线性化手法，与《广义测量平差》§1-9“广义测量平差原理”中观测方程线性化的处理一脉相承，也是第 2 章 §2-1 参数估计线性化的直接推广——区别仅在于这里 \(\bm{X}(t)\) 随时间变化，近似值 \(\bm{X}^{*}(t)\) 也必须随时间更新。线性化后的观测方程 (3.1.39) \(\bm{z}(t)=\bm{H}(t)\bm{X}(t)+\bm{\Delta}(t)\) 与第 2 章线性模型 (2.1.8) 形式一致，这是第 4 章 Kalman 滤波能直接套用静态估计结论的根基：一旦把“时变参数”换成“时变状态”，静态最小二乘的整套数学就搬进了动态估计。
算例分析
例 3.4例 3.3 中的微分方程为 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\\ \dot{X}_3(t)\\ \dot{X}_4(t)\end{bmatrix} =\begin{bmatrix} X_2(t)\\[6pt] X_1(t)X_4^2(t)-\dfrac{GM}{X_1^2(t)}+u_r(t)\\[10pt] X_4(t)\\[6pt] -\dfrac{2}{X_1(t)}X_4(t)X_2(t)+\dfrac{1}{X_1(t)}u_l(t) \end{bmatrix} +\begin{bmatrix}0 & 0\\ 1 & 0\\ 0 & 0\\ 0 & 1\end{bmatrix} \begin{bmatrix}e_r(t)\\ e_\theta(t)\end{bmatrix}\] 已知状态的近似值为 \(\bm{X}^{*}(t)=\left[\begin{array}{llll}X_1^{*}(t) & X_2^{*}(t) & X_3^{*}(t) & X_4^{*}(t)\end{array}\right]\)，将此微分方程线性化。
解：设原微分方程为 \[\dot{\bm{X}}(t)=\bm{g}\left[\,\bm{X}(t),\ \bm{u}(t)\,\right]+\bm{C}\bm{e}(t)\] 将上式在 \(\bm{X}^{*}\left[\,t\,\right]\) 处用泰勒级数展开，略去高阶项后为 \[\dot{\bm{X}}(t)=\bm{g}\left[\,\bm{X}^{*}(t),\ \bm{u}(t)\,\right]+\bm{A}(t)\left[\,\bm{X}(t)-\bm{X}^{*}(t)\,\right]+\bm{C}\cdot\bm{e}(t)\] 其中 \[\bm{g}\left[\,\bm{X}^{*}(t),\ \bm{u}(t)\,\right] =\begin{bmatrix} X_2^{*}(t)\\[6pt] X_1^{*}(t)\left[\,X_4^{*}(t)\,\right]^2-\dfrac{G}{\left(X_1^{*}\right)^2}+u_r(t)\\[10pt] X_4^{*}(t)\\[6pt] -\dfrac{2}{X_1^{*}(t)}X_4^{*}(t)X_2^{*}(t)+\dfrac{1}{X_1^{*}(t)}u_l(t) \end{bmatrix}\] \[\bm{A}(t)=\begin{bmatrix} \dfrac{\partial g_1}{\partial X_1} & \dfrac{\partial g_1}{\partial X_2} & \dfrac{\partial g_1}{\partial X_3} & \dfrac{\partial g_1}{\partial X_4}\\[10pt] \vdots & \vdots & \vdots & \vdots\\[4pt] \dfrac{\partial g_4}{\partial X_1} & \dfrac{\partial g_4}{\partial X_2} & \dfrac{\partial g_4}{\partial X_3} & \dfrac{\partial g_4}{\partial X_4} \end{bmatrix}_{\bm{X}(t)=\bm{X}^{*}(t)}\] 从而 \[\bm{A}(t)=\begin{bmatrix} 0 & 1 & 0 & 0\\[6pt] X_4^{*2}(t)+\dfrac{2GM}{\left[\,X_1^{*}(t)\,\right]^3} & 0 & 0 & 2X_1^{*}(t)X_4^{*}(t)\\[10pt] 0 & 0 & 0 & 1\\[6pt] \dfrac{2X_2^{*}(t)X_4^{*}(t)-u_l(t)}{\left[\,X_1^{*}(t)\,\right]^2} & -\dfrac{2X_4^{*}(t)}{X_1^{*}(t)} & 0 & -\dfrac{2X_2^{*}(t)}{X_1^{*}(t)} \end{bmatrix}\] \[\bm{C}=\begin{bmatrix}0 & 0\\ 1 & 0\\ 0 & 0\\ 0 & 1\end{bmatrix}\]
状态转移矩阵和连续线性系统的解
求解微分方程的方法有解析解法和数值解法。解析解法有多种，如初等积分法、级数展开法、拉氏变换法、对角化方法和利用 Caylay-Hamilton 定律计算（待定系数法）等方法。但实际上很多微分方程的解是不能用初等函数来表示的，有时即便是形式上非常简单的微分方程也不能用解析方法得到，所以在现实应用中，通常需要用数值解法得到微分方程的解。用数值解法得到的微分方程的解是通过递推的方法得到的在离散点上的状态近似数值，其解算近似程度与步长和阶数有关。
本节首先介绍状态转移矩阵及其性质，在此基础上，推导给出有控制输入和噪声影响下的连续线性系统的解，为下一步连续系统离散化做准备。
状态转移矩阵
设某系统的微分方程为 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t) \tag{3.2.1}\] 初始值为 \(\bm{X}(t_0)=\bm{X}_0\)。式 (3.2.1) 是一阶线性非齐次微分方程。当不考虑控制输入和噪声影响时 \[\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t)=0 \tag{3.2.2}\] 式 (3.2.1) 的就成为齐次微分方程 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t) \tag{3.2.3}\] 设上式的解为 \[\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{X}(t_0) \tag{3.2.4}\] 只要得到了 \(\bm{\Phi}(t,\ t_0)\)，那么就可以通过常数变异法来求式 (3.2.1) 的解。
式 (3.2.4) 也表明微分方程的解实质上是 \(\bm{\Phi}(t,\ t_0)\) 将状态从初始时刻 \(t_0\) 到时刻 \(t\) 的转移，转移过程和转移的结果完全由 \(\bm{\Phi}(t,\ t_0)\) 和初始状态 \(\bm{X}(t_0)\) 所决定，所以称 \(\bm{\Phi}(t,\ t_0)\) 为状态转移矩阵，\(\bm{\Phi}(t,\ t_0)\) 也是求解微分方程 (3.2.1) 的关键。
状态转移矩阵 \(\bm{\Phi}(t,\ t_0)\) 就是“时间推进器”：它把 \(t_0\) 时刻的状态向量 \(\bm{X}(t_0)\) 搬运到 \(t\) 时刻，搬运规则完全由系统矩阵 \(\bm{A}(t)\) 决定。把状态方程比作“轨道”，\(\bm{X}(t_0)\) 是出发时刻的坐标，\(\bm{\Phi}(t,\ t_0)\) 则是沿轨道的“乘车券”——给定初始坐标，凭券即可换算出任意后续时刻的坐标。本节后面给出的三条性质（分段转移、转移可逆、时间不变）分别对应“换乘”、“回程”、“原地不动”，都是这张“乘车券”应具备的基本规则；第 4 章 Kalman 滤波的预测步 \(\hat{\bm{X}}(k)=\bm{\Phi}_{k,\ k-1}\hat{\bm{X}}(k-1)\) 正是把它当作“一步推进器”来使用。
将式 (3.2.4) 代入式 (3.2.3) 得到 \[\dot{\bm{\Phi}}(t,\ t_0)\bm{X}(t_0)=\bm{A}(t)\bm{\Phi}(t,\ t_0)\bm{X}(t_0) \tag{3.2.5}\] 由于有初值问题的微分方程的解是唯一的，所以 \[\dot{\bm{\Phi}}(t,\ t_0)=\bm{A}(t)\bm{\Phi}(t,\ t_0) \tag{3.2.6}\] 可见，满足式 (3.2.6) 的 \(\bm{\Phi}(t,\ t_0)\) 就是齐次方程 (3.2.3) 的解。
1. 差分法求状态转移矩阵
差分法是最简单的一种解微分方程的数值方法。通过这样简单的差分方法，我们可以看到如何用数值方法得到状态转移矩阵。
对于微分方程 \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t) \tag{3.2.7}\] 如果 \(\Delta t=t-t_0\) 较小，那么 \[\dot{\bm{X}}(t)=\frac{1}{\Delta t}\left[\,\bm{X}(t)-\bm{X}(t_0)\,\right] \tag{3.2.8}\] 将式 (3.2.8) 代入式 (3.2.7) \[\frac{1}{\Delta t}\left[\,\bm{X}(t)-\bm{X}(t_0)\,\right]=\bm{A}(t_0)\bm{X}(t_0) \tag{3.2.9}\] 得到 \[\bm{X}(t)=\left(\bm{I}+\bm{A}(t_0)\Delta t\right)\bm{X}(t_0) \tag{3.2.10}\] 设 \[\bm{\Phi}(t,\ t_0)=\bm{I}+\bm{A}(t_0)\Delta t \tag{3.2.11}\] 式 (3.2.10) 成为 \[\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{X}(t_0) \tag{3.2.12}\] 如果动态系统为时不变系统，那么式 (3.2.11) 为 \[\bm{\Phi}(t-t_0)=\bm{I}+\bm{A}\Delta t \tag{3.2.13}\]
差分法 (3.2.8) 是一阶 Euler 前向格式：用 \(\dfrac{1}{\Delta t}\left[\,\bm{X}(t)-\bm{X}(t_0)\,\right]\) 近似 \(\dot{\bm{X}}(t)\)，并把 \(\bm{A}(t_0)\bm{X}(t_0)\) 固定在区间左端，得到的 \(\bm{\Phi}(t,\ t_0)\approx\bm{I}+\bm{A}(t_0)\Delta t\) 只是精确转移矩阵 \(e^{\bm{A}\Delta t}\) 的一阶截断。它只有在 \(\Delta t\) 足够小、\(\bm{A}\Delta t\) 的高阶项可忽略时才成立；\(\Delta t\) 越大误差越大，且局部误差是 \(\Delta t^2\) 量级。若要更高精度，应改用 Euler 改进法或四阶 Runge-Kutta 法，代价是计算量增加——这就是数值求解中“精度与计算量”权衡的第一课，也是正文强调“步长越小、阶数越高、近似程度越高但计算量也越大”的原因。
差分方法是最简单的数值解法，它舍去了 \(\Delta t\) 的高阶项。除了差分法，还有很多其他数值解算微分方程的方法，其数值计算结果与原模型的近似程度与步长 \(\Delta t\) 和阶数有关，步长越小，阶数越高，近似程度越高，但计算量也越大。在实际应用中，需要综合考虑计算负荷和近似程度，当近似造成的误差可以忽略，就不需要无谓地增加阶数和减小步长了。常用的数值解法有 Euler 改进方法和四阶 Runge-Kutta 方法。此外，在实际的数值计算过程中，是通过 \[\dot{\bm{\Phi}}(t,\ t_0)=\bm{A}(t)\bm{\Phi}(t,\ t_0)\] 直接求解 \(\bm{\Phi}(t_k,\ t_{k-1})\) 来实现 \(\bm{X}(t_{k-1})\) 与 \(\bm{X}(t_k)\) 的递推关系，\(\bm{\Phi}(t_k,\ t_{k-1})\) 也将用于后面的滤波计算。
2. 时不变系统状态转移矩阵的解析解
当系统为时不变系统时，式 (3.2.3) 成为 \[\dot{\bm{X}}(t)=\bm{A}\bm{X}(t) \tag{3.2.14}\] 它的解为 \[\bm{X}(t)=e^{\int\!A\mathrm{d}t}\times C=e^{At}\times C \tag{3.2.15}\] 根据初值 \(\bm{X}(t_0)=\bm{X}_0\) 可以得到 \[C=e^{-At_0}\bm{X}(t_0) \tag{3.2.16}\] 因此式 (3.2.14) 的解为 \[\bm{X}(t)=e^{A\times(t-t_0)}\bm{X}(t_0) \tag{3.2.17}\] 这时状态转移矩阵为 \[\bm{\Phi}(t,\ t_0)=e^{\int_{t_0}^{t}A\mathrm{d}\tau} \tag{3.2.18}\] 也为 \[\bm{\Phi}(t-t_0)=e^{A\times(t-t_0)} \tag{3.2.19}\] 上式表明时不变系统的状态转移只与 \(\Delta t=t-t_0\) 有关。
将 \(\bm{\Phi}(t-t_0)\) 用级数展开 \[\begin{aligned} \bm{\Phi}(t-t_0)&=\bm{I}+\bm{A}\times(t-t_0)+\frac{1}{2!}\bm{A}^2\times(t-t_0)^2+\cdots+\frac{1}{k!}\bm{A}^k\times(t-t_0)^k+\cdots\\ &=\sum_{k=0}^{\infty}\frac{1}{k!}\bm{A}^k\times(t-t_0)^k \end{aligned} \tag{3.2.20}\] 当 \(t_0=0\) 时，有 \[\bm{\Phi}(t)=e^{At}=\bm{I}+\bm{A}t+\frac{1}{2!}\bm{A}^2t^2+\cdots+\frac{1}{k!}\bm{A}^kt^k+\cdots=\sum_{k=0}^{\infty}\frac{1}{k!}\bm{A}^kt^k \tag{3.2.21}\]
补出时不变系统解析解的两处推导细节。第一，式 (3.2.15) 的 \(e^{\int\!A\mathrm{d}t}\) 是矩阵指数记号：对常矩阵 \(\bm{A}\)，\(\int\bm{A}\,\mathrm{d}t=\bm{A}t\)，\(e^{\bm{A}t}\) 定义为级数 (3.2.20)（逐项收敛），它保持标量指数最重要的性质 \(\dfrac{\mathrm{d}}{\mathrm{d}t}e^{\bm{A}t}=\bm{A}e^{\bm{A}t}\)，因而代入齐次方程 (3.2.14) 恒等成立。第二，用初值定常数时用到 \(\bm{A}\) 与 \(e^{\bm{A}t}\) 可交换：\(\bm{X}(t)=e^{\bm{A}t}C=e^{\bm{A}t}e^{-At_0}\bm{X}(t_0)=e^{\bm{A}(t-t_0)}\bm{X}(t_0)\)——把两个指数合并成 \(e^{\bm{A}(t-t_0)}\) 正是式 (3.2.17) 的结果，也说明了时不变系统的转移只依赖时间差 \(\Delta t\)。
对矩阵指数要警惕三处“标量直觉”陷阱。其一，\(e^{\bm{A}t}\) 是对矩阵逐项作幂级数求和，不能对每个元素单独取指数——例 3.5 的 \(\bm{A}=\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\) 满足 \(\bm{A}^2=\bm{0}\)，级数只剩 \(\bm{I}+\bm{A}t\) 两项，若逐元素取指数就会得到错误结果。其二，一般有 \(e^{(\bm{A}+\bm{B})t}\neq e^{\bm{A}t}e^{\bm{B}t}\)，只有当 \(\bm{A}\bm{B}=\bm{B}\bm{A}\)（可交换）时等式才成立。其三，式 (3.2.36) 的行列式公式实为 Jacobi 公式 \(|\,\bm{\Phi}(t,\ t_0)\,|=\exp\left(\int_{t_0}^{t}\mathrm{trace}\left[\,\bm{A}(\tau)\,\right]\mathrm{d}\tau\right)\)（原书漏印指数函数，见编者注），它说明转移矩阵恒非奇异，从而性质 (2) 的“可逆”对任意 \(\bm{A}(t)\) 总成立。
3. 状态转移矩阵的性质
状态转移矩阵有如下性质：
(1) 分段转移： \[\bm{\Phi}(t_2,\ t_0)=\bm{\Phi}(t_2,\ t_1)\bm{\Phi}(t_1,\ t_0) \tag{3.2.22}\]
(2) 转移可逆： \[\bm{\Phi}(t,\ t_0)=\bm{\Phi}^{-1}(t_0,\ t) \tag{3.2.23}\]
(3) 时间不变： \[\bm{\Phi}(t_0,\ t_0)=\bm{I} \tag{3.2.24}\]
现给出以上性质的证明。
根据式 (3.2.12)，有 \[\bm{X}(t_1)=\bm{\Phi}(t_1,\ t_0)\bm{X}(t_0) \tag{3.2.25}\] \[\bm{X}(t_2)=\bm{\Phi}(t_2,\ t_1)\bm{X}(t_1) \tag{3.2.26}\] 将式 (3.2.25) 代入式 (3.2.26) \[\bm{X}(t_2)=\bm{\Phi}(t_2,\ t_1)\bm{\Phi}(t_1,\ t_0)\bm{X}(t_0) \tag{3.2.27}\] 又由于 \[\bm{X}(t_2)=\bm{\Phi}(t_2,\ t_0)\bm{X}(t_0) \tag{3.2.28}\] 根据微分方程解的唯一性 \[\bm{\Phi}(t_2,\ t_0)=\bm{\Phi}(t_2,\ t_1)\bm{\Phi}(t_1,\ t_0) \tag{3.2.29}\] 这说明了状态可以从 \(t_0\) 转移到 \(t_2\)，也可以分步进行：先从 \(t_0\) 转移到 \(t_1\)，再从 \(t_1\) 转移到 \(t_2\)。图 3.2 表示了以上的转移过程。
将式 (3.2.4) 的两边同时左乘 \(\bm{\Phi}^{-1}(t,\ t_0)\)，得到 \[\bm{X}(t_0)=\bm{\Phi}^{-1}(t,\ t_0)\bm{X}(t) \tag{3.2.30}\] 因为有 \[\bm{X}(t_0)=\bm{\Phi}(t_0,\ t)\bm{X}(t) \tag{3.2.31}\] 所以 \[\bm{\Phi}(t,\ t_0)=\bm{\Phi}^{-1}(t_0,\ t) \tag{3.2.32}\] 上式也说明了 \(\bm{\Phi}(t,\ t_0)\) 必须是可逆矩阵。
如果时间 \(t=t_0\)，那么式式 (3.2.4) 为 \[\bm{X}(t_0)=\bm{\Phi}(t_0,\ t_0)\bm{X}(t_0) \tag{3.2.33}\] 因此， \[\bm{\Phi}(t_0,\ t_0)=\bm{I} \tag{3.2.34}\]
性质 (1) (3) 可以概括成一句话：转移可以分段、可以回退、原地不动就是恒等变换。性质 (1) 分段转移 \(\bm{\Phi}(t_2,\ t_0)=\bm{\Phi}(t_2,\ t_1)\bm{\Phi}(t_1,\ t_0)\) 说的是“先到 \(t_1\) 再到 \(t_2\)，等价于直达 \(t_2\)”，滤波预测中常把它拆成多步递推；性质 (2) 转移可逆 \(\bm{\Phi}(t,\ t_0)=\bm{\Phi}^{-1}(t_0,\ t)\) 说的是“\(t_0\to t\) 的乘车券的逆，就是 \(t\to t_0\) 的回程票”，后续 Kalman 平滑需要把滤波值“往回搬”，正是靠它实现；性质 (3) \(\bm{\Phi}(t_0,\ t_0)=\bm{I}\) 则是“不花时间等于没动”。对时不变系统，\(\bm{\Phi}(t-t_0)=e^{\bm{A}(t-t_0)}\) 使这些性质归结为指数运算的常规法则，例 3.5 的 \(\begin{bmatrix}1 & t\\ 0 & 1\end{bmatrix}\) 就是 \(e^{\bm{A}t}\) 的直观形态。
状态转移矩阵 \(\bm{\Phi}(t,\ t_0)\) 除了有以上的性质，它的行列式还满足 \[\frac{\mathrm{d}}{\mathrm{d}t}\left|\,\bm{\Phi}(t,\ t_0)\,\right|=\mathrm{trace}\left[\,\bm{A}(t)\,\right]\,\left|\,\bm{\Phi}(t,\ t_0)\,\right| \tag{3.2.35}\] 和 \[\left|\,\bm{\Phi}(t,\ t_0)\,\right|=\int_{t_0}^{t}\mathrm{trace}\left[\,\bm{A}(\tau)\,\right]\mathrm{d}\tau \tag{3.2.36}\]
原书式 (3.2.36) 右端漏印了指数函数，按 Jacobi 公式应为 \(|\,\bm{\Phi}(t,\ t_0)\,|=\exp\left(\int_{t_0}^{t}\mathrm{trace}\left[\,\bm{A}(\tau)\,\right]\mathrm{d}\tau\right)\)。此处照原样排印。
时不变系统的状态转移矩阵 \(\bm{\Phi}(t-t_0)\) 除了具有 \(\bm{\Phi}(t,\ t_0)\) 的性质外，还具有以下性质 \[\left[\,\bm{\Phi}(\Delta t)\,\right]^k=\bm{\Phi}(k\cdot\Delta t) \tag{3.2.37}\]
状态转移矩阵的定义与性质，在《广义测量平差》§4-1“状态方程的解”中完整对应：定义见该书 (4-1-20)，性质 (1) 分段转移对应 (4-1-26)，性质 (2) 转移可逆对应 (4-1-27)，初值条件 (3.2.24) 对应 (4-1-24)。该书 4-6 节预测、4-7 节平滑正是反复运用这两条性质；本节例 3.5 求得的 \(\bm{\Phi}(t)=\begin{bmatrix}1 & t\\ 0 & 1\end{bmatrix}\) 与《广义测量平差》例 4-1-1 的状态转移矩阵完全相同，两书此处的推导过程也可互相参照。此外，时不变系统 \(\left[\,\bm{\Phi}(\Delta t)\,\right]^k=\bm{\Phi}(k\Delta t)\)（(3.2.37)）为 3.3 节把一步转移阵递推为 \(k\) 步提供了依据。
连续线性系统的解
连续线性系统的微分方程为： \[\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t) \tag{3.2.38}\] 其中 \(\bm{e}(t)\) 为白噪声，且 \[E\left[\,\bm{e}(t)\,\right]=0 \tag{3.2.39}\] \[\mathrm{Cov}\left[\,\bm{e}(t),\ \bm{e}(\tau)\,\right]=\bm{D}_e(t)\delta(t-\tau) \tag{3.2.40}\] 如式 (3.2.38) 的非齐次微分方程，可以用“常数变易法”求解。设状态方程 (3.2.38) 的解为 \[\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{\xi}(t) \tag{3.2.41}\] 当 \(t=t_0\) 时，有 \[\bm{X}(t_0)=\bm{\Phi}(t_0,\ t_0)\bm{\xi}(t_0) \tag{3.2.42}\] 由于 \(\bm{\Phi}(t_0,\ t_0)=\bm{I}\)，得到 \[\bm{X}(t_0)=\bm{\xi}(t_0) \tag{3.2.43}\] 对式 (3.2.41) 求微分 \[\dot{\bm{X}}(t)=\bm{\Phi}(t,\ t_0)\dot{\bm{\xi}}(t)+\dot{\bm{\Phi}}(t,\ t_0)\bm{\xi}(t) \tag{3.2.44}\] 将式 (3.2.41) 和 (3.2.44) 代入式 (3.2.38)，并考虑式 (3.2.6) 得到 \[\bm{\Phi}(t,\ t_0)\dot{\bm{\xi}}(t)+\bm{A}(t)\bm{\Phi}(t,\ t_0)\bm{\xi}(t) =\bm{A}(t)\bm{\Phi}(t,\ t_0)\bm{\xi}(t)+\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t) \tag{3.2.45}\] 消去共同项并左乘以 \(\bm{\Phi}^{-1}(t,\ t_0)\) 后 \[\dot{\bm{\xi}}(t)=\bm{\Phi}^{-1}(t,\ t_0)\bm{B}(t)\bm{u}(t)+\bm{\Phi}^{-1}(t,\ t_0)\bm{C}(t)\bm{e}(t) \tag{3.2.46}\] 由于 \(\bm{\Phi}^{-1}(t,\ t_0)=\bm{\Phi}(t_0,\ t)\)，所以有 \[\dot{\bm{\xi}}(t)=\bm{\Phi}(t_0,\ t)\bm{B}(t)\bm{u}(t)+\bm{\Phi}(t_0,\ t)\bm{C}(t)\bm{e}(t) \tag{3.2.47}\] 将上式积分并根据初始条件，得到 \[\bm{\xi}(t)-\bm{\xi}(t_0)=\int_{t_0}^{t}\bm{\Phi}(t_0,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau +\int_{t_0}^{t}\bm{\Phi}(t_0,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \tag{3.2.48}\] 由于 \(\bm{\xi}(t_0)=\bm{X}(t_0)\)， \[\bm{\xi}(t)=\bm{X}(t_0)+\int_{t_0}^{t}\bm{\Phi}(t_0,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau +\int_{t_0}^{t}\bm{\Phi}(t_0,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \tag{3.2.49}\] 将式 (3.2.49) 代入式 (3.2.41) 得到微分方程的解 \[\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{X}(t_0)+\int_{t_0}^{t}\bm{\Phi}(t,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau +\int_{t_0}^{t}\bm{\Phi}(t,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \tag{3.2.50}\] 式 (3.2.50) 即为连续线性系统的解。
补出“常数变易法”的完整链条：(3.2.41) 设 \(\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{\xi}(t)\) 是齐次解 \(\bm{\Phi}(t,\ t_0)\bm{X}(t_0)\) 的推广——把常向量 \(\bm{X}(t_0)\) 换成待定向量函数 \(\bm{\xi}(t)\)；由于 \(\bm{\Phi}(t_0,\ t_0)=\bm{I}\)，所以 \(\bm{\xi}(t_0)=\bm{X}(t_0)\)，“变易”只发生在 \(\bm{\xi}\) 上。求导代入后，\(\dot{\bm{\Phi}}\bm{\xi}=\bm{A}\bm{\Phi}\bm{\xi}\) 两项相消（这正是选取转移矩阵作“基础解”的好处），只剩 \(\bm{\Phi}\dot{\bm{\xi}}=\bm{B}\bm{u}+\bm{C}\bm{e}\)，解得 \(\dot{\bm{\xi}}=\bm{\Phi}(t_0,\ t)\left[\,\bm{B}\bm{u}+\bm{C}\bm{e}\,\right]\)（用到 \(\bm{\Phi}^{-1}(t,\ t_0)=\bm{\Phi}(t_0,\ t)\)），积分并代回即得 (3.2.50)。该式结构为“齐次解 + 输入与噪声的卷积积分”，与第 4 章离散滤波中“预测 + 修正”的结构一脉相承。
连续线性系统的解 (3.2.50)，与《广义测量平差》§4-1“状态方程的解”中 (4-1-21) 完全一致：右端均为“转移矩阵作用初始状态 + 对控制与噪声的卷积积分”。两书记号略有差异（本书噪声项写 \(\bm{C}(t)\bm{e}(t)\)，该书写作 \(\bm{F}(t)\bm{\varOmega}(t)\)），结构完全相同。这一结果也是 3.3 节离散化的出发点：令 \(t_0=t_{k-1}\)、\(t=t_k\) 即得一步递推 (3.3.4)。另外，例 3.6 中被积函数出现 \(e^{-\beta(t-t_0)}\)，其来源一阶高斯-马尔可夫过程的解可见第 1 章 §1.10；《广义测量平差》§4-1 的随机模型 (4-1-28) (4-1-32) 则给出解中噪声项的统计特性。
算例分析
例 3.5已知系统的微分方程为 \(\begin{bmatrix}\dot{x}_1(t)\\ \dot{x}_2(t)\end{bmatrix}=\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\begin{bmatrix}x_1(t)\\ x_2(t)\end{bmatrix}\)，初始时刻 \(t_0=0\) 的状态为 \(\bm{x}(0)=\begin{bmatrix}x_1(0)\\ x_2(0)\end{bmatrix}\)，试求微分方程的解。
解：此系统为时不变系统 \[\bm{A}=\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\] 并且 \(t_0=0\)。可设微分方程的解为 \[\bm{x}(t)=\bm{\Phi}(t)\bm{x}(0) \tag{3.2.51}\] 状态转移矩阵 \(\bm{\Phi}(t)\) 为 \[\bm{\Phi}(t)=e^{At}=\bm{I}+\bm{A}t+\frac{1}{2!}\bm{A}^2t^2+\cdots+\frac{1}{k!}\bm{A}^kt^k+\cdots =\sum_{k=0}^{\infty}\frac{1}{k!}\bm{A}^kt^k\] 由于 \[\bm{A}^2=\bm{A}^3=\cdots=\bm{A}^k=\begin{bmatrix}0 & 0\\ 0 & 0\end{bmatrix}\] 所以，状态转移矩阵为 \[\bm{\Phi}(t)=\bm{I}+\bm{A}t=\begin{bmatrix}1 & 0\\ 0 & 1\end{bmatrix}+\begin{bmatrix}0 & t\\ 0 & 0\end{bmatrix} =\begin{bmatrix}1 & t\\ 0 & 1\end{bmatrix}\] 微分方程的解为 \[\begin{bmatrix}x_1(t)\\ x_2(t)\end{bmatrix}=\begin{bmatrix}1 & t\\ 0 & 1\end{bmatrix}\times\begin{bmatrix}x_1(0)\\ x_2(0)\end{bmatrix} \tag{3.2.52}\]
在此问题中，如果视 \(x_1(t)\) 为位移，\(x_2(t)\) 为速度，从式 (3.2.52) 来看，此系统显然处于匀速运动状态。如果已知初始时刻的位移和速度，由式 (3.2.52) 可以求得任何时刻的状态 \(\bm{x}(t)\)。更一般的，如果 \[\underset{n\times n}{\bm{A}}=\begin{bmatrix} 0 & 1 & 0 & & 0\\ 0 & 0 & 1 & & 0\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ 0 & 0 & 0 & \cdots & 1\\ 0 & 0 & 0 & \cdots & 0 \end{bmatrix} \tag{3.2.53}\] 矩阵 \(\bm{A}\) 的 \(k\)（\(k\leqslant n\)）次幂为 \[\underset{n\times n}{\bm{A}^k}=\begin{bmatrix} \underset{(n-k)\times k}{\bm{0}} & \underset{(n-k)\times(n-k)}{\bm{I}}\\[6pt] \underset{k\times k}{\bm{0}} & \underset{k\times(n-k)}{\bm{0}} \end{bmatrix} \tag{3.2.54}\] 当 \(k=n\) 时， \[\underset{n\times n}{\bm{A}^n}=\bm{0} \tag{3.2.55}\] 代入式 (3.30) 得到状态转移矩阵
原书此处“代入式 (3.30)”疑为“式 (3.2.20)”之排印笔误（指状态转移矩阵的级数展开式），此处照原样排印。
\[\bm{\Phi}(t-t_0)=\begin{bmatrix} 1 & t-t_0 & \dfrac{1}{2}(t-t_0)^2 & \dfrac{1}{1\cdot 2\cdot 3}(t-t_0)^3 & \cdots & \dfrac{1}{(n-1)!}(t-t_0)^{n-1}\\[10pt] 0 & 1 & t-t_0 & \dfrac{1}{2}(t-t_0)^2 & \cdots & \dfrac{1}{(n-2)!}(t-t_0)^{n-2}\\[10pt] 0 & 0 & 1 & t-t_0 & \cdots & \dfrac{1}{(n-3)!}(t-t_0)^{n-3}\\[10pt] 0 & 0 & 0 & 1 & \cdots & \dfrac{1}{(n-4)!}(t-t_0)^{n-4}\\[8pt] \vdots & \vdots & \vdots & \vdots & & \vdots\\ 0 & 0 & 0 & 0 & 0 & 1 \end{bmatrix} \tag{3.2.56}\] 当 \(\bm{A}\) 为 \(2\times 2\) 维矩阵时，状态转移矩阵就是其中的前两行和前两列；当 \(\bm{A}\) 为 \(3\times 3\) 维矩阵时，状态为位置、速度和加速度，转移矩阵就是其中的前三行和前三列，依次类推。
例 3.6已知 \(\bm{X}(t_0)\)，试给出连续时间系统 \(\dot{\bm{X}}(t)=\begin{bmatrix}0 & 1\\ 0 & -\beta\end{bmatrix}\bm{X}(t)\) 的解。
解：状态变量为 \(\bm{X}(t)=\left[\begin{array}{ll}X_1(t) & X_2(t)\end{array}\right]^{\mathrm{T}}\)，只要求得 \(\bm{\Phi}(t,\ t_0)\)，就得到了此微分方程的解。由系统的微分方程可知 \[\begin{aligned} \dot{X}_1(t)&=X_2(t) & \text{①}\\ \dot{X}_2(t)&=-\beta X_2(t) & \text{②} \end{aligned}\] 由②和初始值，可得 \[\begin{aligned} X_2(t)&=e^{\int_{t_0}^{t}(-\beta)\mathrm{d}t}X_2(t_0)\\ &=e^{-\beta(t-t_0)}X_2(t_0) & \text{③} \end{aligned}\] 将③代入①得 \[\dot{X}_1(t)=e^{-\beta(t-t_0)}X_2(t_0)\] 积分得到 \[X_1(t)=-\frac{1}{\beta}e^{-\beta(t-t_0)}X_2(t_0)+C\] 当 \(t=t_0\) 时 \[X_1(t_0)=-\frac{1}{\beta}X_2(t_0)+C\] 解得到 \[C=X_1(t_0)+\frac{1}{\beta}X_2(t_0)\] 所以 \[X_1(t)=X_1(t_0)+\frac{1}{\beta}\left(1-e^{-\beta(t-t_0)}\right)X_2(t_0)\] 综上所述 \[\bm{\Phi}(t,\ t_0)=\begin{bmatrix}1 & \dfrac{1}{\beta}\left(1-e^{-\beta(t-t_0)}\right)\\[8pt] 0 & e^{-\beta(t-t_0)}\end{bmatrix}\] 最后，微分方程的解为 \[\bm{X}(t)=\begin{bmatrix}1 & \dfrac{1}{\beta}\left(1-e^{-\beta(t-t_0)}\right)\\[8pt] 0 & e^{-\beta(t-t_0)}\end{bmatrix}\bm{X}(t_0)\]
解微分方程的方法并不唯一，除了上述的方法解此问题中的微分方程，还可以用拉普拉斯变换和约旦规范形等其他方法解得。
离散线性系统的数学模型
式 (3.2.1) 是对系统在时间连续状态下的描述，是连续系统的状态方程。为了研究系统的运动规律，需要对系统进行观测，利用观测值对状态进行估计。对动态系统的观测是某些离散的时间点上的观测，所以在对动态系统的估计前需要将状态方程离散化，得到离散的数学模型，从而便于计算机处理数据。
离散线性系统的函数模型
根据上节推导得出的连续线性系统的解，有 \[\bm{X}(t)=\bm{\Phi}(t,\ t_0)\bm{X}(t_0)+\int_{t_0}^{t}\bm{\Phi}(t,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau +\int_{t_0}^{t}\bm{\Phi}(t,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \tag{3.3.1}\] 其中 \(\bm{e}(t)\) 为白噪声过程，且 \[E\left[\,\bm{e}(t)\,\right]=0 \tag{3.3.2}\] \[\mathrm{Cov}\left[\,\bm{e}(t),\ \bm{e}(\tau)\,\right]=\bm{D}_e(t)\delta(t-\tau) \tag{3.3.3}\] \(\bm{D}_e(t)\) 为 \(q\times q\) 维对称非负定矩阵，是 \(\bm{e}(t)\) 的均方值。现根据需要在时间区间 \(\left[\,t_0,\ t\,\right]\) 中划分不同的时刻，有 \(t_0<t_1<t_2<\cdots<t_{k-1}<t_k\cdots\)。设上式中 \(t_0=t_{k-1}\) 和 \(t=t_k\)，得到 \[\bm{X}(t_k)=\bm{\Phi}(t_k,\ t_{k-1})\bm{X}(t_{k-1}) +\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau +\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \tag{3.3.4}\] 设 \[\begin{aligned} \bm{\Psi}_{k,\ k-1}\bm{u}(t_k)&=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{B}(\tau)\bm{u}(\tau)\,\mathrm{d}\tau\\ \bm{w}(k-1)&=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau \end{aligned} \tag{3.3.5}\] 为了简单起见，这里将原来的变量表示为 \[\begin{aligned} \bm{X}(k)&\rightarrow\bm{X}(t_k)\\ \bm{\Phi}_{k,\ k-1}&\rightarrow\bm{\Phi}(t_k,\ t_{k-1})\\ \bm{u}(k)&=\bm{u}(t_k) \end{aligned} \tag{3.3.6}\] 式 (3.3.4) 成为 \[\underset{n\times 1}{\bm{X}(k)}=\underset{n\times n}{\bm{\Phi}_{k,\ k-1}}\ \underset{n\times 1}{\bm{X}(k-1)} +\underset{n\times p}{\bm{\Psi}_{k,\ k-1}}\ \underset{p\times 1}{\bm{u}(k)}+\underset{n\times 1}{\bm{w}(k-1)} \tag{3.3.7}\] 上式中 \(\bm{\Psi}_{k,\ k-1}\bm{u}(k)\) 为状态方程中确定性的控制输入部分，\(\bm{w}(k-1)\) 为离散化后的状态方程中的随机干扰部分。
离散化相当于在时间轴上“每隔 \(\Delta t\) 拍一张快照”。连续方程 (3.3.1) 描述每时每刻的状态，而观测只在离散时刻 \(t_0<t_1<\cdots<t_k\) 取得，估计也只需这些时刻的状态。令 \(t_0=t_{k-1}\)、\(t=t_k\) 后，(3.3.4) 变成三步结构：状态转移 \(\bm{\Phi}_{k,\ k-1}\bm{X}(k-1)\)（把上个采样点的状态搬过来）、输入驱动 \(\bm{\Psi}_{k,\ k-1}\bm{u}(k)\)（本步内确定性输入的作用）、噪声累积 \(\bm{w}(k-1)\)（本步内白噪声 \(\bm{e}\) 经转移矩阵加权后的总效果）。这三项分别对应连续解 (3.2.50) 中的三项，只是“连续卷积”被压缩成了“步内积分”。
根据式 (3.1.33)，更一般的状态方程表示为 \[\underset{n\times 1}{\bm{X}(k)}=\underset{n\times n}{\bm{\Phi}_{k,\ k-1}}\ \underset{n\times 1}{\bm{X}(k-1)} +\underset{n\times 1}{\bm{\Omega}(k-1)}+\underset{n\times 1}{\bm{w}(k-1)}\] 其中 \[\bm{\Omega}(k-1)=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{G}(\tau)\,\mathrm{d}\tau\] 同样，对观测方程进行离散化。根据式 (3.1.39)，离散化后的线性观测方程为 \[\bm{z}(k)=\bm{H}(k)\bm{X}(k)+\bm{\Delta}(k) \tag{3.3.8}\] 其中 \[\bm{z}(k)=\bm{Z}(k)-\bm{F}\left[\,\bm{X}^{*}(k)\,\right]+\bm{H}_k\bm{X}^{*}(k) \tag{3.3.9}\]
离散线性系统的函数模型 (3.3.7) 与 (3.3.8)，与《广义测量平差》§4-2“离散线性系统的数学模型”同主题对应：该书对连续系统 (4-2-1) 离散化得到离散状态方程 (4-2-8)、(4-2-9)，其中 \(\bm{\varPhi}_{k+1,\ k}\)、\(\bm{\varGamma}_{k+1,\ k}\) 即本书的 \(\bm{\Phi}_{k,\ k-1}\)、\(\bm{\Psi}_{k,\ k-1}\)，可逐式对照。连续到离散的“一步推进 + 噪声积分”思想亦可回溯到《广义测量平差》§4-1 (4-1-21)、(4-1-22) 的连续解。本节 (3.3.7) 与 (3.3.8) 就是第 4 章 Kalman 滤波函数模型的直接基础。
离散线性系统的随机模型
离散化后的状态方程 (3.3.7) 中的系统噪声为 \(\bm{w}(k-1)\)，\(\bm{w}(k-1)\) 为白噪声序列，其期望为 \[E\left[\,\bm{w}(k-1)\,\right]=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)E\left[\,\bm{e}(\tau)\,\right]\mathrm{d}\tau=0 \tag{3.3.10}\] 根据式 (3.3.5) 和方差的定义有 \[\begin{aligned} \bm{D}_w(k-1)=&E\left[\,\bm{w}(k-1)\bm{w}^{\mathrm{T}}(k-1)\,\right]\\ =&E\left\{\left[\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\,\mathrm{d}\tau\right] \left[\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau')\bm{C}(\tau')\bm{e}(\tau')\,\mathrm{d}\tau'\right]^{\mathrm{T}}\right\}\\ =&E\left[\int_{t_{k-1}}^{t_k}\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{e}(\tau)\bm{e}^{\mathrm{T}}(\tau') \bm{C}^{\mathrm{T}}(\tau')\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau')\,\mathrm{d}\tau\mathrm{d}\tau'\right]\\ =&\int_{t_{k-1}}^{t_k}\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\delta(\tau-\tau') \bm{C}^{\mathrm{T}}(\tau')\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau')\,\mathrm{d}\tau'\mathrm{d}\tau\\ =&\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau) \left(\int_{t_{k-1}}^{t_k}\delta(\tau-\tau')\bm{C}^{\mathrm{T}}(\tau')\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau')\,\mathrm{d}\tau'\right)\mathrm{d}\tau \end{aligned} \tag{3.3.11}\] 对 \(\tau'\) 积分得到 \[\bm{D}_w(k-1)=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau)\,\mathrm{d}\tau \tag{3.3.12}\] 式 (3.3.12) 即为离散化后的系统噪声方差矩阵 \(\bm{D}_w(k-1)\) 的严密计算公式。
补出 \(\bm{D}_w(k-1)\) 从双重积分 (3.3.11) 化简到单重积分 (3.3.12) 的关键跳步：由 \(E[\bm{e}(\tau)\bm{e}^{\mathrm{T}}(\tau')]=\bm{D}_e(\tau)\delta(\tau-\tau')\)，对内层积分用 Dirac-\(\delta\) 函数的“筛选性质”\(\int f(\tau')\delta(\tau-\tau')\,\mathrm{d}\tau'=f(\tau)\)，即 \(\int_{t_{k-1}}^{t_k}\delta(\tau-\tau')\,\bm{C}^{\mathrm{T}}(\tau')\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau')\,\mathrm{d}\tau'=\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau)\)，从而把二重积分压成沿 \(\tau\) 的单重积分。注意外层的 \(\bm{\Phi}(t_k,\ \tau)\)、\(\bm{C}(\tau)\) 保留 \(\tau\) 变量，物理含义正是“噪声在 \(\tau\) 时刻产生，再经转移矩阵 \(\bm{\Phi}(t_k,\ \tau)\) 传播到 \(t_k\) 时刻”的累积。
\(\bm{D}_w(k-1)\) 也可以用近似方法计算得到。当 \(\Delta t=t_k-t_{k-1}\rightarrow 0\)，将上式中的 \(\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\) 用 \(\bm{C}(t_{k-1})\bm{D}_e(t_{k-1})\bm{C}^{\mathrm{T}}(t_{k-1})\) 代替，在极短的时间里 \(\bm{\Phi}(t_k,\ \tau)=\bm{I}\)，那么 \[\bm{D}_w(k-1)=\bm{C}(t_{k-1})\bm{D}_e(t_{k-1})\bm{C}^{\mathrm{T}}(t_{k-1})\Delta t \tag{3.3.13}\]
近似公式 (3.3.13) 的成立条件非常明确：\(\Delta t\to 0\)，且把 \(\bm{\Phi}(t_k,\ \tau)\) 取为单位阵 \(\bm{I}\)（忽略步内转移）。它与严密公式 (3.3.12) 的差别在于 \(\Delta t\) 的高阶项——例 3.7 严格积分给出 \(q^2\begin{bmatrix}\dfrac{\Delta t^3}{3} & \dfrac{\Delta t^2}{2}\\[4pt] \dfrac{\Delta t^2}{2} & \Delta t\end{bmatrix}\)，而近似公式只留下对角项 \(\bm{D}_w\approx q^2\begin{bmatrix}0 & 0\\ 0 & \Delta t\end{bmatrix}\)，协方差元素 \(\Delta t^2/2\)、\(\Delta t^3/3\) 被全部丢弃。GNSS/INS 等工程应用中 \(\Delta t\) 通常足够小，近似可接受；但精度要求高（如相位级定位）时必须使用严密积分式。用近似公式前先检查“\(\bm{A}\Delta t\) 是否远小于 1”是一条实用的经验准则。
由于 \(\bm{w}(k)\) 为白噪声序列，所以 \[\mathrm{Cov}\left[\,\bm{w}(k),\ \bm{w}(j)\,\right]=\bm{D}_w(k)\delta(k-j) \tag{3.3.14}\] 对于离散系统来说，\(\delta(k-j)\) 是克罗尼克 \(\delta\) 函数（Kronecker delta function） \[\delta(k-j)=\begin{cases}0, & k\neq j\\ 1, & k=j\end{cases} \tag{3.3.15}\]
如果已知 \(\bm{X}(k-1)\) 和方差 \(\bm{D}_X(k-1)\)，根据状态方程 (3.3.7)，容易得到 \[\bm{D}_X(k)=\bm{\Phi}_{k,\ k-1}\bm{D}_X(k-1)\bm{\Phi}_{k-1}^{\mathrm{T}}+\bm{D}_w(k-1) \tag{3.3.16}\]
方差递推 (3.3.16) \(\bm{D}_X(k)=\bm{\Phi}_{k,\ k-1}\bm{D}_X(k-1)\bm{\Phi}_{k-1}^{\mathrm{T}}+\bm{D}_w(k-1)\) 是第 2 章协方差传播律的动态版本：第一项 \(\bm{\Phi}\bm{D}_X\bm{\Phi}^{\mathrm{T}}\) 是“已知方差经线性变换传播”（与静态平差 \(\bm{D}_L=\bm{B}\bm{D}_X\bm{B}^{\mathrm{T}}\) 同构），第二项 \(\bm{D}_w(k-1)\) 是“本步新注入的噪声”。合起来的意思很直白：状态的不确定度 = 上一步不确定度随运动传播 + 本步新噪声贡献。这个公式将直接成为第 4 章 Kalman 滤波预测步方差更新的原型；而 (3.3.19) 的多步展开则说明“不确定度不会自动消失，只会随传播而增长”——这正是滤波需要不断用观测“压住”发散的根本原因。
\(\bm{X}(k)\) 可由 \(\bm{X}(0)\) 递推得到，即 \[\begin{cases} \bm{X}(1)=\bm{\Phi}_{1,\ 0}\bm{X}(0)+\bm{\Psi}_{1,\ 0}\bm{u}(0)+\bm{w}(0)\\ \bm{X}(2)=\bm{\Phi}_{2,\ 1}\bm{X}(1)+\bm{\Psi}_{2,\ 1}\bm{u}(1)+\bm{w}(1)\\ \quad\vdots\\ \bm{X}(k-1)=\bm{\Phi}_{k-1,\ k-2}\bm{X}(k-2)+\bm{\Psi}_{k-1,\ k-2}\bm{u}(k-2)+\bm{w}(k-2)\\ \bm{X}(k)=\bm{\Phi}_{k,\ k-1}\bm{X}(k-1)+\bm{\Psi}_{k,\ k-1}\bm{u}(k-1)+\bm{w}(k-1) \end{cases} \tag{3.3.17}\] 由于 \(\bm{\Phi}_{k,\ j}=\bm{\Phi}_{k,\ i}\bm{\Phi}_{i,\ j}\)，因此有 \[\bm{X}(k)=\bm{\Phi}_{k,\ 0}\bm{X}(0)+\sum_{i=1}^{k}\bm{\Phi}_{k,\ i}\bm{\Psi}_{i,\ i-1}\bm{u}(i-1) +\sum_{i=1}^{k}\bm{\Phi}_{k,\ i}\bm{w}(i-1) \tag{3.3.18}\] 已知 \(\bm{X}(0)\) 的方差 \(\bm{D}_X(0)\)，那么 \[\bm{D}_X(k)=\bm{\Phi}_{k,\ 0}\bm{D}_X(0)\bm{\Phi}_{k,\ 0}^{\mathrm{T}} +\sum_{i=1}^{k}\bm{\Phi}_{k,\ i}\bm{D}_w(i-1)\bm{\Phi}_{k,\ i}^{\mathrm{T}} \tag{3.3.19}\]
补出 (3.3.17) (3.3.19) 多步递推的归纳跳步。把式 (3.3.17) 逐行代回，并利用分段转移性质 \(\bm{\Phi}_{k,\ i}=\bm{\Phi}_{k,\ j}\bm{\Phi}_{j,\ i}\) 合并同类项，即得 \(\bm{X}(k)=\bm{\Phi}_{k,\ 0}\bm{X}(0)+\sum_{i=1}^{k}\bm{\Phi}_{k,\ i}\bm{\Psi}_{i,\ i-1}\bm{u}(i-1)+\sum_{i=1}^{k}\bm{\Phi}_{k,\ i}\bm{w}(i-1)\)：初始状态乘全程转移 \(\bm{\Phi}_{k,\ 0}\)，第 \(i-1\) 步的输入与噪声经 \(\bm{\Phi}_{k,\ i}\) 传到终点。方差展开 (3.3.19) 是对上式的各项取协方差——因为 \(\bm{X}(0)\) 与各 \(\bm{w}(i)\) 两两不相关（见 (3.3.22)），交叉项全部为零，只剩下两项“平方”。式 (3.3.23) \(\bm{D}_X(k,\ j)=\bm{\Phi}_{k,\ j}\bm{D}_X(j)\) 则是协方差的传播：\(j\) 时刻的状态被搬到 \(k\) 时刻后，它与 \(j\) 时刻自身的协方差随之缩放一个转移矩阵。
更一般的，若已知 \(\bm{X}(j)\) 和 \(\bm{D}_X(j)\)，可以递推得到 \(\bm{X}(k)\ (k>j)\) \[\bm{X}(k)=\bm{\Phi}_{k,\ j}\bm{X}(j)+\sum_{i=j+1}^{k}\bm{\Phi}_{k,\ i}\bm{\Psi}_{i,\ i-1}\bm{u}(i-1) +\sum_{i=j+1}^{k}\bm{\Phi}_{k,\ i}\bm{w}(i-1) \tag{3.3.20}\] 和方差 \(\bm{D}_X(k)\) \[\bm{D}_X(k)=\bm{\Phi}_{k,\ j}\bm{D}_X(j)\bm{\Phi}_{k,\ j}^{\mathrm{T}} +\sum_{i=j+1}^{k}\bm{\Phi}_{k,\ i}\bm{D}_w(i-1)\bm{\Phi}_{k,\ i}^{\mathrm{T}} \tag{3.3.21}\] 从上面的递推关系看出，\(\bm{X}(k)\) 与 \(\bm{w}(j)\)（\(k>j\)）均相关；反之，只要系统噪声 \(\bm{w}(\cdot)\) 的时间超前或等于状态 \(\bm{X}(\cdot)\) 的时间，那么系统噪声与状态无关，因此有 \[\mathrm{Cov}\left[\,\bm{X}(j),\ \bm{w}(k)\,\right]=0\quad (k\geq j) \tag{3.3.22}\] 通过式 (3.3.20)，可以求得 \(\bm{X}(k)\ (k>j)\) 与 \(\bm{X}(j)\) 的协方差为 \[\mathrm{Cov}\left(\bm{X}(k),\ \bm{X}(j)\right)=\bm{D}_X(k,\ j)=\bm{\Phi}_{k,\ j}\bm{D}_X(j) \tag{3.3.23}\]
与离散的状态方程相对应，离散状态下的观测方程为 \[\bm{z}(k)=\bm{H}_k\bm{X}(k)+\bm{\Delta}(k) \tag{3.3.24}\] \(\bm{\Delta}(k)\) 为连续型观测噪声 \(\bm{\Delta}(t)\) 的离散序列。若已知白噪声过程 \(\bm{\Delta}(t)\) 的均方值 \(\bm{D}_{\Delta}(t)\)，在微小的时间段 \(\Delta t\)，将 \(\bm{D}_{\Delta}(t)\) 在 \(\Delta t\) 上的均值作为在 \(t_k\) 时刻的噪声序列 \(\bm{\Delta}(k)\) 的方差 \[\bm{D}_{\Delta}(k)=\bm{D}_{\Delta}(t)/\Delta t \tag{3.3.25}\]
观测噪声的离散化 (3.3.25) \(\bm{D}_{\Delta}(k)=\bm{D}_{\Delta}(t)/\Delta t\) 初看反直觉：采样越密（\(\Delta t\) 越小），离散噪声方差反而越大。其依据是面积等价约定——连续白噪声的均方值 \(\bm{D}_{\Delta}(t)\) 与横轴围成的面积 \(\bm{D}_{\Delta}(k)\Delta t\) 在离散化前后保持不变，即“功率密度不变”。这里 \(\Delta t\) 应与滤波步长一致：若滤波采用别的步长，\(\bm{D}_{\Delta}(k)\) 必须按同一规则换算，否则滤波增益会随采样频率异常漂移。还要区分两类噪声：过程噪声 \(\bm{w}(k-1)\) 的方差按 (3.3.12) 积分得到（与转移矩阵耦合、含 \(\Delta t\) 的高阶项），观测噪声 \(\bm{\Delta}(k)\) 的方差按 (3.3.25) 简单缩放——来源不同，不可混用公式。
也就是说，离散观测噪声序列的方差与横轴所围的面积 \(\bm{D}_{\Delta}(k)\Delta t\) 等于连续噪声的均方值 \(\bm{D}_{\Delta}(t)\)。观测噪声序列仍然为白噪声 \[E\left[\,\bm{\Delta}(k)\,\right]=0\ ,\quad \mathrm{Cov}\left[\,\bm{\Delta}(k),\ \bm{\Delta}(j)\,\right]=\bm{D}_{\Delta}(k)\delta(k-j) \tag{3.3.26}\] 观测值是系统的输出，与系统状态和系统噪声均无关，所以观测噪声 \(\bm{\Delta}(k)\) 与任何时刻的系统噪声 \(\bm{w}(j)\) 和状态 \(\bm{X}(j)\) 互不相关 \[\mathrm{Cov}\left[\,\bm{w}(j),\ \bm{\Delta}(k)\,\right]=0 \tag{3.3.27}\] \[\mathrm{Cov}\left[\,\bm{X}(j),\ \bm{\Delta}(k)\,\right]=0 \tag{3.3.28}\] 此外，初始值 \(\bm{X}(0)\) 一般根据经验或者观测值确定，也是随机量，在估计中假设 \[E\left[\,\bm{X}(0)\,\right]=\hat{\bm{X}}(0)\ ,\quad \mathrm{Var}\left[\,\bm{X}(0)\,\right]=\bm{D}_X(0) \tag{3.3.29}\] \(\bm{D}_X(0)\) 为 \(n\times n\) 维正定矩阵。在估计中还假设初始状态 \(\bm{X}(0)\) 与系统噪声 \(\bm{w}(k)\) 和观测噪声 \(\bm{\Delta}(k)\) 均不相关 \[\mathrm{Cov}\left[\,\bm{X}(0),\ \bm{w}(k)\,\right]=0\ ,\quad \mathrm{Cov}\left[\,\bm{X}(0),\ \bm{\Delta}(k)\,\right]=0 \tag{3.3.30}\]
离散随机模型 (3.3.10) (3.3.30) 与《广义测量平差》§4-2 的随机模型对应，其中白噪声的定义与“噪声与初始状态不相关”的假定可回溯到第 1 章 §1.10“随机过程”与《广义测量平差》§4-1 的随机模型 (4-1-28) (4-1-32)。式 (3.3.22) 的“未来噪声与当前状态不相关”在第 4 章 Kalman 滤波推导 \(\mathrm{cov}[\bm{X}(j),\ \bm{w}(k)]\)、\(\mathrm{cov}[\bm{X}(j),\ \bm{\Delta}(k)]\) 等交叉协方差时会被反复使用，务必记牢。此外，例 3.8 的一阶高斯-马尔可夫过程、例 3.9 的随机常数/随机游走/随机斜坡过程，其连续定义都在第 1 章 §1.10；例 3.10 的 GNSS 钟差模型则是它们组合应用的典型，可对照学习。
算例分析
本节将通过例题给出常用的离散化的状态方程。
例 3.7已知系统的状态方程为 \[\begin{bmatrix}\dot{x}_1(t)\\ \dot{x}_2(t)\end{bmatrix} =\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\begin{bmatrix}x_1(t)\\ x_2(t)\end{bmatrix} +\begin{bmatrix}0\\ 1\end{bmatrix}e(t)\] 系统噪声 \(e(t)\) 为白噪声过程，\(q^2\) 为 \(e(t)\) 的均方值，求离散化的状态方程和它的随机模型。
解：在此微分方程中 \[\bm{A}(t)=\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\ ;\quad \bm{C}(t)=\begin{bmatrix}0\\ 1\end{bmatrix}\] 由例 3.5 知道 \(t=t_k\) 时刻的状态转移矩阵为 \[\bm{\Phi}(t_k,\ \tau)=\begin{bmatrix}1 & t_k-\tau\\ 0 & 1\end{bmatrix}\] 设 \(\Delta t=t_k-t_{k-1}\)，那么 \(t_{k-1}\) 到 \(t_k\) 的状态方程为 \[\bm{X}(k)=\begin{bmatrix}1 & \Delta t\\ 0 & 1\end{bmatrix}\bm{X}(k-1)+\bm{w}(k-1)\] 其中，噪声 \(\bm{w}(k-1)\) 为 \[\bm{w}(k-1)=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)e(\tau)\,\mathrm{d}\tau\] \(\bm{w}(k-1)\) 的方差为 \[\bm{D}_w(k-1)=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau)\,\mathrm{d}\tau\] 对上式积分得到 \[\bm{D}_w(k-1)=q^2\begin{bmatrix}\dfrac{(\Delta t)^3}{3} & \dfrac{(\Delta t)^2}{2}\\[8pt] \dfrac{(\Delta t)^2}{2} & \Delta t\end{bmatrix}\]
| 连续动态系统 | 离散动态系统（\(\Delta t=t_k-t_{k-1}\)） | ||||
|---|---|---|---|---|---|
| \(\dot{\bm{X}}(t)=\bm{A}(t)\bm{X}(t)+\bm{B}(t)\bm{u}(t)+\bm{C}(t)\bm{e}(t)\) | \(\bm{X}(k)=\bm{\Phi}_{k,\ k-1}\bm{X}(k-1)+\bm{\Psi}_{k,\ k-1}\bm{u}(k-1)+\bm{w}(k-1)\) | ||||
| \(\bm{A}(t)\) | \(\bm{C}(t)\) | \(\bm{D}_e(t)\) | \(\bm{\Phi}_{k,\ k-1}\) | \(\bm{D}_w(k-1)\) | |
| \(\left[\,0\,\right]\) | \(\left[\,1\,\right]\) | \(q^2\) | \(\left[\,1\,\right]\) | \(\left[\,q\Delta t\,\right]\) | |
| \(\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\) | \(\begin{bmatrix}0\\ 1\end{bmatrix}\) | \(q^2\) | \(\begin{bmatrix}1 & \Delta t\\ 0 & 1\end{bmatrix}\) | \(q^2\begin{bmatrix}\dfrac{(\Delta t)^3}{3} & \dfrac{(\Delta t)^2}{2}\\[8pt] \dfrac{(\Delta t)^2}{2} & \Delta t\end{bmatrix}\) | |
| \(\begin{bmatrix}0 & 1 & 0\\ 0 & 0 & 1\\ 0 & 0 & 0\end{bmatrix}\) | \(\begin{bmatrix}0\\ 0\\ 1\end{bmatrix}\) | \(q^2\) | \(\begin{bmatrix}1 & \Delta t & \dfrac{1}{2}\Delta t^2\\[6pt] 0 & 1 & \Delta t\\ 0 & 0 & 1\end{bmatrix}\) | \(q^2\begin{bmatrix}\dfrac{(\Delta t)^5}{20} & \dfrac{(\Delta t)^4}{8} & \dfrac{(\Delta t)^3}{6}\\[10pt] \dfrac{(\Delta t)^4}{8} & \dfrac{(\Delta t)^3}{3} & \dfrac{(\Delta t)^2}{2}\\[10pt] \dfrac{(\Delta t)^3}{6} & \dfrac{(\Delta t)^2}{2} & \Delta t\end{bmatrix}\) | |
例 3.7 中的状态方程描述了运动系统在一维空间里处于匀速运动的状态。表 3.1 给出了在一维空间里，运动系统在静止、匀速和匀加速情况下的运动状态，以及离散化后的状态方程和系统噪声方差矩阵。读者可以对表 3.1 进行扩展，给出以上运动状态在三维空间里的状态方程和系统噪声方差矩阵。
例 3.8有连续时间系统 \(\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\end{bmatrix}=\begin{bmatrix}0 & 1\\ 0 & -\beta\end{bmatrix}\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix}+\begin{bmatrix}0\\ 1\end{bmatrix}e(t)\)，\(e(t)\) 为白噪声过程，并且 \(\mathrm{Cov}\left[\,e(t),\ e(\tau)\,\right]=\sigma^2\delta(t-\tau)\)。求采样周期为 \(T\) 的离散化状态方程和状态方程的随机模型。
解：由 1.10.8 节可知，此问题中的 \(X_2(t)\) 为一阶高斯-马尔可夫过程。由例 3.6 可以得到此问题的状态转移矩阵为 \[\bm{\Phi}(t-\tau)=\begin{bmatrix}1 & \dfrac{1}{\beta}\left(1-e^{-\beta(t-\tau)}\right)\\[8pt] 0 & e^{-\beta(t-\tau)}\end{bmatrix}\] 一个周期 \(T=t_k-t_{k-1}\) 的转移矩阵为 \[\bm{\Phi}=\begin{bmatrix}1 & \dfrac{1}{\beta}\left(1-e^{-\beta T}\right)\\[8pt] 0 & e^{-\beta T}\end{bmatrix}\] 离散化后的状态方程为 \[\bm{X}(k)=\begin{bmatrix}1 & \dfrac{1}{\beta}\left(1-e^{-\beta T}\right)\\[8pt] 0 & e^{-\beta T}\end{bmatrix}\bm{X}(k-1)+\bm{w}(k-1)\] 噪声 \(\bm{w}(k-1)\) 为零均值白噪声 \[\bm{w}(k-1)=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)e(\tau)\,\mathrm{d}\tau\] 其中， \[\bm{C}(t)=\begin{bmatrix}0\\ 1\end{bmatrix}\] 根据式 (3.3.12) 求得 \(\bm{w}(k-1)\) 的方差矩阵为 \[\begin{aligned} \bm{D}_w(k)&=\int_{t_k}^{t_{k+1}}\bm{\Phi}(t_{k+1},\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_{k+1},\ \tau)\,\mathrm{d}\tau\\ &=\sigma^2\begin{bmatrix} \dfrac{1}{\beta^2}\left(T-2+\dfrac{1}{2\beta}+2e^{-\beta T}-\dfrac{1}{2\beta}e^{-2\beta T}\right) & \dfrac{1}{2\beta^2}\left(1-2e^{-\beta T}+e^{-2\beta T}\right)\\[10pt] \dfrac{1}{2\beta^2}\left(1-2e^{-\beta T}+e^{-2\beta T}\right) & \dfrac{1}{2\beta}\left(1-e^{-2\beta T}\right) \end{bmatrix} \end{aligned}\]
原书例 3.8 的 \(\bm{D}_w(k)\) 矩阵 (1, 1) 元排印为 \(\dfrac{1}{\beta^2}\left(T-2+\dfrac{1}{2\beta}+2e^{-\beta T}-\dfrac{1}{2\beta}e^{-2\beta T}\right)\)，其中“\(-2\)”与“\(2e^{-\beta T}\)”两项疑漏印分母 \(\beta\)（应为 \(-\dfrac{2}{\beta}\) 与 \(\dfrac{2}{\beta}e^{-\beta T}\)）。此处照原样排印。
例 3.9已知初值 \(x(t_0)\)，给出以下有色噪声过程离散化后的状态方程
(1) \(x(t)\) 为随机常数过程，表示为 \(\dot{x}(t)=0\)；
(2) \(x(t)\) 为一阶高斯-马尔可夫过程，随机过程表示为 \(\dot{x}(t)+\beta x(t)=e(t)\)；
(3) \(x(t)\) 为随机游走过程，随机过程表示为 \(\dot{x}(t)=e(t)\)；
(4) \(\bm{x}(t)\) 为随机斜坡过程，表示为 \(\dot{x}_1(t)=x_2(t)\)，\(\dot{x}_2(t)=0\)。
以上状态方程中 \(e(t)\) 均为白噪声过程，其均方值为 \(\sigma^2\)。
解：已知连续型的微分方程和状态的初始值初值 \(x(t_0)\)，解微分方程，可得到形如式 (3.3.7) 的离散型状态方程，并根据式 (3.3.12) 积分得到系统噪声的方差矩阵。
(1) 解微分方程：\(\dot{x}(t)=0\) \[x(t)=x(t_0)\] 设 \(t=t_k\) 和 \(t_0=t_{k-1}\)，离散化后的状态方程为 \[x(k)=x(k-1)\] 由于此系统无噪声输入，所以 \[\bm{D}_w(k-1)=0\]
(2) 解微分方程：\(\dot{x}(t)+\beta x(t)=e(t)\)，得到 \[x(t)=x(t_0)e^{-\beta(t-t_0)}+\int_{t_0}^{t}e^{-\beta(t-\tau)}e(\tau)\,\mathrm{d}\tau\] 设 \[w(t_0)=\int_{t_0}^{t}e^{-\beta(t-\tau)}e(\tau)\,\mathrm{d}\tau\] 其中 \(w(t_0)\) 为白噪声序列。设 \(t=t_k\) 和 \(t_0=t_{k-1}\)，离散化后的状态方程为 \[\begin{aligned} x(k)&=e^{-\beta(t_k-t_{k-1})}x(k-1)+w(k-1)\\ w(k-1)&=\int_{t_{k-1}}^{t_k}e^{-\beta(t_k-\tau)}e(\tau)\,\mathrm{d}\tau \end{aligned}\] \(w(k-1)\) 的期望为 \[E\left[\,w(k-1)\,\right]=\int_{t_{k-1}}^{t_k}e^{-\beta(t_k-\tau)}E\left[\,e(\tau)\,\right]\mathrm{d}\tau=0\] \(w(k-1)\) 的方差为 \[\begin{aligned} D_w(k-1)&=\sigma_w^2(k-1)\\ &=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau)\,\mathrm{d}\tau\\ &=\int_{t_{k-1}}^{t_k}e^{-\beta(t_k-\tau)}\sigma^2e^{-\beta(t_k-\tau)}\,\mathrm{d}\tau\\ &=\frac{\sigma^2}{2\beta}\left(1-e^{-2\beta(t_k-t_{k-1})}\right) \end{aligned}\]
(3) 解微分方程：\(\dot{x}(t)=e(t)\)，得到 \[x(k)=x(k-1)+w(k-1)\] 其中， \[w(k-1)=\int_{t_{k-1}}^{t_k}e(\tau)\,\mathrm{d}\tau\] \(w(k-1)\) 的方差为 \[D_w(k-1)=\sigma_w^2(k-1)=\int_{t_{k-1}}^{t_k}\sigma^2\,\mathrm{d}\tau=\sigma^2(t_k-t_{k-1})\]
(4) 解随机斜坡过程的微分方程，得到离散化后的状态方程为 \[\begin{bmatrix}x_1(k)\\ x_2(k)\end{bmatrix} =\begin{bmatrix}1 & t_k-t_{k-1}\\ 0 & 1\end{bmatrix}\begin{bmatrix}x_1(k-1)\\ x_2(k-1)\end{bmatrix}\] 由于无噪声输入，所以 \[\bm{D}_w(k-1)=0\]
例 3.10GNSS 接收机时钟与导航系统时间频率不同步导致的时钟偏差 \(t_b\) 是影响 GNSS 测距精度的主要因素，为了消除时钟偏差对定位结果的影响，将接收机钟差模型化并进行估计。接收机时钟频率漂移 \(t_d\)，即振荡器的时间震荡频率与标称频率之间的差异，是引起时钟偏差 \(t_b\) 的根本因素。一般将时钟频率漂移 \(t_d\) 描述为随机游走过程；时钟偏差 \(t_b\) 可认为是由钟频率漂移 \(t_d\) 和白噪声叠加而成，所以将 GNSS 接收机钟差 \(t_b\) 的运动规律描述为 \[\begin{bmatrix}\dot{t}_b(t)\\ \dot{t}_d(t)\end{bmatrix} =\begin{bmatrix}0 & 1\\ 0 & 0\end{bmatrix}\begin{bmatrix}t_b(t)\\ t_d(t)\end{bmatrix} +\begin{bmatrix}1 & 0\\ 0 & 1\end{bmatrix}\begin{bmatrix}e_b(t)\\ e_d(t)\end{bmatrix}\] 其中 \(e_b(t)\) 和 \(e_d(t)\) 为零均值白噪声，谱密度矩阵为 \(\bm{D}_e=\begin{bmatrix}\sigma_{e_b}^2 & 0\\ 0 & \sigma_{e_d}^2\end{bmatrix}\)。将上述的钟差状态方程离散化，并给出系统噪声的方差矩阵。
解：根据表 3.1，可知此系统的状态转移矩阵为 \[\bm{\Phi}(t_k,\ \tau)=\begin{bmatrix}1 & t_k-\tau\\ 0 & 1\end{bmatrix}\] 设 \(\bm{T}(t)=\left[\begin{array}{ll}t_b(t) & t_d(t)\end{array}\right]^{\mathrm{T}}\)，\(\bm{e}(\tau)=\left[\begin{array}{ll}e_b(\tau) & e_d(\tau)\end{array}\right]^{\mathrm{T}}\)，离散化后的状态方程为 \[\bm{T}(k)=\bm{\Phi}_{k,\ k-1}\bm{T}(k-1)+\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k-\tau)\cdot\bm{C}\bm{e}(\tau)\,\mathrm{d}\tau\] 其中， \[\bm{C}=\begin{bmatrix}1 & 0\\ 0 & 1\end{bmatrix}\] 令 \[\bm{w}_{k,\ k-1}=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k-\tau)\cdot\bm{C}\bm{e}(\tau)\,\mathrm{d}\tau\] \[\Delta t=t_k-t_{k-1}\] 那么，状态方程为 \[\bm{T}(k)=\bm{\Phi}(k,\ k-1)\bm{T}(k-1)+\bm{w}_{k,\ k-1}\] 系统噪声的方差为 \[\begin{aligned} \bm{D}_w(k-1)&=\int_{t_{k-1}}^{t_k}\bm{\Phi}(t_k,\ \tau)\bm{C}(\tau)\bm{D}_e(\tau)\bm{C}^{\mathrm{T}}(\tau)\bm{\Phi}^{\mathrm{T}}(t_k,\ \tau)\,\mathrm{d}\tau\\ &=\int_{t_{k-1}}^{t_k}\begin{bmatrix}1 & t_k-\tau\\ 0 & 1\end{bmatrix}\begin{bmatrix}1 & 0\\ 0 & 1\end{bmatrix} \begin{bmatrix}\sigma_{e_b}^2 & 0\\ 0 & \sigma_{e_d}^2\end{bmatrix}\begin{bmatrix}1 & 0\\ 0 & 1\end{bmatrix}^{\mathrm{T}} \begin{bmatrix}1 & t_k-\tau\\ 0 & 1\end{bmatrix}^{\mathrm{T}}\mathrm{d}\tau\\ &=\begin{bmatrix} \sigma_{e_b}^2\Delta t+\sigma_{e_d}^2\dfrac{\Delta t^3}{3} & \sigma_{e_d}^2\dfrac{\Delta t^2}{2}\\[10pt] \sigma_{e_d}^2\dfrac{\Delta t^2}{2} & \sigma_{e_d}^2\Delta t \end{bmatrix} \end{aligned}\]
在以上模型中，谱密度矩阵为 \(\bm{D}_e\) 的取值是关键。参考文献 [10] 通过 Allen 方差分析了接收机钟差的主要噪声成分，并给出了不同类型接收机钟的 \(\sigma_{e_b}^2\) 和 \(\sigma_{e_d}^2\) 的取值： \[\begin{aligned} \sigma_{e_b}^2&\approx\frac{h_0}{2}c^2\\ \sigma_{e_d}^2&\approx 2\pi^2h_{-2}c^2 \end{aligned}\] 其中 \(h_0\) 和 \(h_{-2}\) 为 Allen 方差系数，取值见表 3.2；\(c\) 为光速。
| 记时标准 | \(h_0\) | \(h_{-1}\) | \(h_{-2}\) |
|---|---|---|---|
| TCXO（低质量） | \(2\times 10^{-19}\) | \(7\times 10^{-21}\) | \(2\times 10^{-20}\) |
| TCXO（低质量） | \(2\times 10^{-19}\) | \(1\times 10^{-22}\) | \(3\times 10^{-24}\) |
| OCXO | \(2\times 10^{-25}\) | \(7\times 10^{-25}\) | \(6\times 10^{-25}\) |
| Rubidium | \(2\times 10^{-22}\) | \(4.5\times 10^{-26}\) | \(1\times 10^{-30}\) |
| Cesium | \(2\times 10^{-22}\) | \(5\times 10^{-27}\) | \(1.5\times 10^{-33}\) |
表 3.2 中两行均为“TCXO（低质量）”，按原书排印转录；第二行应为“TCXO（高质量）”（原书排印笔误）。
线性动态系统的可控性和可测性
动态系统数学模型中的状态方程描述了控制输入量及初始状态对系统内部状态的影响，表明了系统内部结构特性，但不是所有的状态方程中的状态变量都受输入量的控制，如状态方程 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\end{bmatrix} =\begin{bmatrix}1 & 0\\ 2 & 3\end{bmatrix}\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix}+\begin{bmatrix}0\\ 2\end{bmatrix}u(t) \tag{3.4.1}\] 输入 \(u(t)\) 对状态的作用如图 3.3 所示。该系统的输入控制 \(u(t)\) 只能影响 \(X_2(t)\) 而不能影响 \(X_1(t)\)，输入控制 \(u(t)\) 无法改变 \(X_1(t)\) 的运动，显然 \(X_1(t)\) 是不可控的，即此系统是不完全可控的。如果改变系统的控制方法，即改变控制输入系数矩阵为 \(\left[\begin{array}{ll}2 & 0\end{array}\right]^{\mathrm{T}}\) 后，那么 \(X_1(t)\) 又是可控的，可控的 \(X_1(t)\) 又可以改变 \(X_2(t)\)，这样就达到了输入 \(u(t)\) 控制 \(X_1(t)\) 和 \(X_2(t)\) 的目的。
如果对该系统进行的测量（输出），观测方程为 \[Z(t)=\left[\begin{array}{ll}1 & 0\end{array}\right]\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix} \tag{3.4.2}\] 此时观测输出 \(Z(t)\) 值有 \(X_1(t)\) 的信息，但从 \(Z(t)\) 的输出中，无法获得 \(X_2(t)\) 的信息，所以这个系统中 \(X_2(t)\) 是不可估计的，这样的系统是不完全可测的。
以上提出的问题就是现代控制论中的可控性和可测性问题。可控性和可测性是卡尔曼于 20 世纪 60 年代首先提出的，是研究线性系统控制问题中重要的概念，在建立动态系统模型时首先需要判断系统中的状态变量是否可控和可测。动态系统的可控性和可测性与噪声无关，所以在下面的推导中均不考虑系统噪声和观测噪声。
可控性问的是“输入 \(u(t)\) 能不能驾驭全部状态”：若某个状态分量无论怎样施加输入都无法改变（如式 (3.4.1) 中的 \(X_1(t)\)），系统就不完全可控。可测性问的是“输出 \(Z(t)\) 能不能看穿全部状态”：若某状态分量在输出中完全不显现（如式 (3.4.2) 中的 \(X_2(t)\)），系统就不完全可测。直观上，可控性考察“控制通道”的覆盖能力（\(\bm{B}\) 矩阵把输入接到了哪些状态上），可测性考察“观测通道”的覆盖能力（\(\bm{H}\) 矩阵让输出看到了哪些状态）。卡尔曼提出这两个概念，正是为了回答滤波开始前必须确认的前提：状态既“能管得着”又“能看得见”，估计才有意义。
可控性
在时间区间 \(\left[\,t_0,\ t_1\,\right]\) 内，如果控制系统的输入 \(u(t)\)（\(t\in\left[\,t_0,\ t_1\,\right]\)，\(t_1>t_0\)）可以将初始值 \(\bm{X}(t_0)\) 转移到预定的 \(\bm{X}(t_1)\)（\(\bm{X}(t_1)\) 为任意的预定状态），在这样的系统中，每一个状态变量都是可以控制的，称这样的系统是完全可控的。如果希望系统能够获得任意的状态值，那么这个系统一定要具备可控性。下面讨论在时不变系统中，连续状态和离散状态下判断系统完全可控的条件。
1. 连续线性时不变系统的可控性条件
连续线性时不变系统的微分方程为 \[\dot{\bm{X}}(t)=\bm{A}\bm{X}(t)+\bm{B}u(t) \tag{3.4.3}\] 根据式 (3.2.50)，它的解为 \[\begin{aligned} \bm{X}(t)&=\bm{\Phi}(t-t_0)\bm{X}(t_0)+\int_{t_0}^{t}\bm{\Phi}(t-\tau)\bm{B}u(\tau)\,\mathrm{d}\tau\\ &=\bm{\Phi}(t)\bm{\Phi}(-t_0)\bm{X}(t_0)+\bm{\Phi}(t)\int_{t_0}^{t}\bm{\Phi}(-\tau)\bm{B}u(\tau)\,\mathrm{d}\tau \end{aligned} \tag{3.4.4}\] 两边同时乘以 \(\bm{\Phi}(-t)\) \[\bm{\Phi}(-t)\bm{X}(t)=\bm{\Phi}(-t_0)\bm{X}(t_0)+\int_{t_0}^{t}\bm{\Phi}(-\tau)\bm{B}u(\tau)\,\mathrm{d}\tau \tag{3.4.5}\] 当 \(t=t_1\) 时， \[\bm{\Phi}(-t_1)\bm{X}(t_1)-\bm{\Phi}(-t_0)\bm{X}(t_0)=\int_{t_0}^{t_1}\bm{\Phi}(-\tau)\bm{B}u(\tau)\,\mathrm{d}\tau \tag{3.4.6}\] 将 \(\bm{\Phi}(-\tau)\) 展开为级数 \[\bm{\Phi}(-\tau)=e^{-\bm{A}\tau}=\sum_{k=0}^{n-1}\bm{A}^k\alpha_k(\tau) \tag{3.4.7}\] 并设 \[\bm{T}=\bm{\Phi}(-t_1)\bm{X}(t_1)-\bm{\Phi}(-t_0)\bm{X}(t_0) \tag{3.4.8}\] 将式 (3.4.7) 和式 (3.4.8) 代入式 (3.4.6)，得到 \[\bm{T}=\bm{B}\int_{t_0}^{t_1}\alpha_0(\tau)u(\tau)\,\mathrm{d}\tau +\bm{A}\bm{B}\int_{t_0}^{t_1}\alpha_1(\tau)u(\tau)\,\mathrm{d}\tau+\cdots +\bm{A}^{n-1}\bm{B}\int_{t_0}^{t_1}\alpha_{n-1}(\tau)u(\tau)\,\mathrm{d}\tau \tag{3.4.9}\] 设 \[\int_{t_0}^{t_1}\alpha_k(\tau)u(\tau)\,\mathrm{d}\tau=\beta_k \tag{3.4.10}\] 那么，式 (3.4.9) 为 \[\left[\,\underset{n\times p}{\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}^2\bm{B}}\ \vdots\ \cdots\ \vdots\ \underset{n\times p}{\bm{A}^{n-1}\bm{B}}\,\right] \begin{bmatrix}\underset{p\times 1}{\bm{\beta}_0}\\ \underset{p\times 1}{\bm{\beta}_2}\\ \vdots\\ \underset{p\times 1}{\bm{\beta}_{n-1}}\end{bmatrix}=\underset{n\times 1}{\bm{T}} \tag{3.4.11}\] 上式中的 \(\bm{\beta}_k\) 为控制输入部分。如果对于任意上的 \(\bm{T}\)，方程都有解，意味着有控制输入 \(u(t)\) 可以将始值 \(\bm{X}(t_0)\) 转移到预定的 \(\bm{X}(t_1)\)。
设 \[\underset{n\times np}{\bm{C}}=\left[\,\underset{n\times p}{\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}^2\bm{B}}\ \vdots\ \cdots\ \vdots\ \underset{n\times p}{\bm{A}^{n-1}\bm{B}}\,\right] \tag{3.4.12}\]
非齐次线性方程有解的条件是 \(\mathrm{Rank}(\bm{C})=\mathrm{Rank}(\bm{C}\vdots\bm{T})\)。由于 \(\mathrm{Rank}(\bm{C})\leq n\)， \(\mathrm{Rank}(\bm{C}\vdots\bm{T})\leq n\)。若 \(\mathrm{Rank}(\bm{C})<n\)，那么 \(\mathrm{Rank}(\bm{C})\) 不一定等于 \(\mathrm{Rank}(\bm{C}\vdots\bm{T})\)，所以只有当 \(\mathrm{Rank}(\bm{C})=n\)， 就一定有 \(\mathrm{Rank}(\bm{C}\vdots\bm{T})=n\)，这时方程一定有解，方程的解即为输入控制 \(\beta_k\ (k=0,\ \cdots,\ n-1)\)，所以线性连续时不变动态系统可控的条件为
\[\mathrm{rank}\left[\,\underset{n\times p}{\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}\bm{B}}\ \vdots\ \underset{n\times p}{\bm{A}^2\bm{B}}\ \vdots\ \cdots\ \vdots\ \underset{n\times p}{\bm{A}^{n-1}\bm{B}}\,\right]=n \tag{3.4.13}\] 线性连续时不变动态系统只要满足式 (3.4.13)，就有控制输入可以将状态由初始值 \(\bm{X}(t_0)\) 转移到 \(\bm{X}(t_1)\)，系统具有可控性，\(\bm{C}\) 也被称为可控性矩阵。
补出“可控性秩判据”推导中的两个关键跳步。第一，式 (3.4.7) 把 \(e^{-\bm{A}\tau}=\bm{\Phi}(-\tau)\) 展开成有限项 \(\sum_{k=0}^{n-1}\bm{A}^k\alpha_k(\tau)\)——这靠的是Cayley-Hamilton 定理：\(n\times n\) 矩阵 \(\bm{A}\) 的幂自 \(n\) 次起可被 \(0\sim n-1\) 次幂线性表示，故矩阵指数只需保留前 \(n\) 项。第二，引入 \(\beta_k=\int_{t_0}^{t_1}\alpha_k(\tau)u(\tau)\,\mathrm{d}\tau\) 后，“存在控制 \(u(t)\) 使任意 \(\bm{T}\) 可解”被转化为非齐次方程 \(\bm{C}\bm{\beta}=\bm{T}\) 对任意 \(\bm{T}\) 有解，即 \(\mathrm{Rank}(\bm{C})=\mathrm{Rank}(\bm{C}\vdots\bm{T})=n\)，于是得到 (3.4.13)。离散情形 (3.4.16) 的推导无需 Cayley-Hamilton，直接由递推展开即可得到完全平行的秩条件。
2. 离散线性时不变系统的可控性条件
设离散的状态方程为 \[\bm{X}(k)=\bm{\Phi}\bm{X}(k-1)+\bm{\Psi}u(k-1) \tag{3.4.14}\] 得到递推方程 \[\begin{cases} \bm{X}(1)=\bm{\Phi}\bm{X}(0)+\bm{\Psi}u(0)\\ \bm{X}(2)=\bm{\Phi}\bm{X}(1)+\bm{\Psi}u(1)=\bm{\Phi}^2\bm{X}(0)+\bm{\Phi}\bm{\Psi}u(0)+\bm{\Psi}u(1)\\ \qquad\vdots\\ \bm{X}(k)=\bm{\Phi}^k\bm{X}(0)+\bm{\Phi}^{k-1}\bm{\Psi}u(0)+\bm{\Phi}^{k-2}\bm{\Psi}u(1)+\cdots+\bm{\Phi}^0\bm{\Psi}u(k-1) \end{cases} \tag{3.4.15}\] 对于任意的 \(\bm{X}(k)\) 和 \(\bm{X}(0)\) 都有 \[\bm{X}(k)-\bm{\Phi}^k\bm{X}(0) =\left[\,\bm{\Phi}^0\bm{\Psi}\ \vdots\ \bm{\Phi}\bm{\Psi}\ \vdots\ \bm{\Phi}^2\bm{\Psi}\ \vdots\ \cdots\ \vdots\ \bm{\Phi}^{k-2}\bm{\Psi}\ \vdots\ \bm{\Phi}^{k-1}\bm{\Psi}\,\right] \begin{bmatrix}u(k-1)\\ u(k-2)\\ u(k-3)\\ \vdots\\ u(1)\\ u(0)\end{bmatrix} \tag{3.4.16}\] 令 \[\underset{n\times np}{\bm{C}}=\left[\,\bm{\Phi}^0\bm{\Psi}\ \vdots\ \bm{\Phi}\bm{\Psi}\ \vdots\ \bm{\Phi}^2\bm{\Psi}\ \vdots\ \cdots\ \vdots\ \bm{\Phi}^{k-1}\bm{\Psi}\,\right] \tag{3.4.17}\] 与式 (3.4.11) 有解的条件一样，式 (3.4.16) 有解的条件为 \[\mathrm{rank}(\bm{C})=n \tag{3.4.18}\] 满足式 (3.4.18) 的线性离散时不变动态系统一定具有可控性，\(\bm{C}\) 也被称为离散系统的可控性矩阵。在判定可控性时，当式 (3.4.17) 中的 \(\bm{\Phi}^{k-1}\bm{\Psi}\) 扩展到 \(k=n\) 时，\(\bm{C}\) 的秩仍然不能为 \(n\)，表明此系统不可控。
用秩判据判可控性有三个易错点。其一，\(\bm{C}=\left[\bm{B}\ \bm{A}\bm{B}\ \cdots\ \bm{A}^{n-1}\bm{B}\right]\) 是 \(n\times np\) 阵，列数远多于行数，判据只问行满秩 \(\mathrm{rank}(\bm{C})=n\)，而不是方阵可逆。其二，连续系统用 \((\bm{A},\ \bm{B})\)，离散系统用 \((\bm{\Phi},\ \bm{\Psi})\)，两者不可混用——3.3 节表 3.1 已给出常用 \(\bm{\Phi}\)、\(\bm{\Psi}\) 的对应形式。其三，\(n\) 是状态维数，只需算到 \(\bm{A}^{n-1}\bm{B}\)；若到这一步秩仍小于 \(n\)，再追加 \(\bm{A}^{n}\bm{B}\) 也不会改变秩（仍是 Cayley-Hamilton 定理），系统必不可控。例 3.11 中 \(\mathrm{rank}(\bm{C})=1<n=2\) 正是典型的不可控演示。
可控性的概念与连续、离散两类判据，与《广义测量平差》§4-8“线性确定系统的能观性和能控性”直接对应：该书同样先给连续系统可控性条件、再给离散系统可控性条件，判定矩阵的构造与本书 (3.4.13)、(3.4.18) 完全一致。卡尔曼滤波稳定性理论（《广义测量平差》§4-9）以能观能控为前提，可作本章概念在该书中的延伸阅读。
可测性
可测性是指从系统的输出中能否能获取状态的信息。在时间区间 \(\left[\,t_0,\ t_1\,\right]\)（\(t_1>t_0\)）内，根据 \(t_0\) 到 \(t_1\) 的观测值 \(\bm{Z}(t)\)（\(t\in\left[\,t_0,\ t_1\,\right]\)）可以唯一地确定系统在初始时刻的状态 \(\bm{X}(t_0)\)，则称系统是完全可测的。如果想通过观测值获得系统的状态，那么这个系统一定要具有可测性。下面分别讨论连续线性时不变系统和离散线性时不变系统完全可测的判断条件。
可测性的定义有一个值得玩味的细节：它要求由观测唯一确定初始状态 \(\bm{X}(t_0)\)。为什么是 \(t_0\) 而不是别的时刻？因为 \(t_0\to t_1\) 之间的运动完全由状态方程决定，一旦 \(\bm{X}(t_0)\) 唯一确定，整段轨迹 \(\bm{X}(t)\) 就全部确定。可控性与可测性是一对对偶概念：可控性研究“输入\(\to\)状态”，可测性研究“输出\(\to\)状态”，两者的判据恰好互为转置——可控性矩阵 \(\bm{C}=\left[\bm{B}\ \bm{A}\bm{B}\ \cdots\ \bm{A}^{n-1}\bm{B}\right]\)，可测性矩阵 \(\bm{M}=\left[\bm{H}^{\mathrm{T}}\ \bm{A}^{\mathrm{T}}\bm{H}^{\mathrm{T}}\ \cdots\ (\bm{A}^{\mathrm{T}})^{n-1}\bm{H}^{\mathrm{T}}\right]^{\mathrm{T}}\)，判据都是满秩。抓住这层对偶关系，两个判据只需记一个，另一个转置即得。
1. 连续线性时不变系统的可测性
不考虑干扰的连续线性时不变系统的数学模型为 \[\begin{cases} \dot{\bm{X}}(t)=\bm{A}\bm{X}(t)+\bm{B}u(t)\\ \bm{Z}(t)=\bm{H}\bm{X}(t) \end{cases} \tag{3.4.19}\] 状态方程的解为 \[\bm{X}(t)=e^{\bm{A}(t-t_0)}\bm{X}(t_0)+\int_{t_0}^{t}e^{\bm{A}(t-\tau)}\bm{B}u(\tau)\,\mathrm{d}\tau \tag{3.4.20}\] 将其代入观测方程 \[\begin{aligned} \bm{Z}(t)&=\bm{H}\bm{X}(t)\\ &=\bm{H}e^{\bm{A}(t-t_0)}\bm{X}(t_0)+\bm{H}\int_{t_0}^{t}e^{\bm{A}(t-\tau)}\bm{B}u(\tau)\,\mathrm{d}\tau \end{aligned} \tag{3.4.21}\] 设 \[\overline{\bm{Z}}(t)=\bm{Z}(t)-\bm{H}\int_{t_0}^{t}e^{\bm{A}(t-\tau)}\bm{B}u(\tau)\,\mathrm{d}\tau \tag{3.4.22}\] 得到关于 \(\bm{X}(t_0)\) 的线性方程组 \[\bm{H}e^{\bm{A}(t-t_0)}\bm{X}(t_0)=\overline{\bm{Z}}(t) \tag{3.4.23}\] 若上式有解，表明 \(\bm{X}(t_0)\) 可以由 \(\overline{\bm{Z}}(t)\) 估计得到。由于 \(\overline{\bm{Z}}(t)\) 可取任意值，所以这等价于研究 \(u(t)=0\) 时由 \(\bm{Z}(t)\) 来估计 \(\bm{X}(t_0)\) 的问题，即可测性与系统输入无关，那么 \[\bm{H}e^{\bm{A}(t-t_0)}\bm{X}(t_0)=\bm{Z}(t) \tag{3.4.24}\] 由于 \[e^{\bm{A}(t-t_0)}=\sum_{k=0}^{n-1}\alpha_k(t)\bm{A}^k \tag{3.4.25}\] 代入式 (3.4.24)，得到 \[\left[\,\alpha_0(t)\bm{I}_{\ell\times\ell}\quad \alpha_1(t)\bm{I}\quad \cdots\quad \alpha_{n-1}(t)\bm{I}\,\right] \begin{bmatrix}\bm{H}\\ \bm{H}\bm{A}^1\\ \vdots\\ \bm{H}\bm{A}^{n-1}\end{bmatrix}_{n\times 1}\bm{X}(t_0)=\bm{Z}(t) \tag{3.4.26}\] 设 \[\bm{M}=\begin{bmatrix}\underset{\ell\times n}{\bm{H}}\\ \underset{\ell\times n\ n\times n}{\bm{H}\bm{A}^1}\\ \vdots\\ \underset{\ell\times n\ n\times n}{\bm{H}\bm{A}^{n-1}}\end{bmatrix} \tag{3.4.27}\] 由于 \(\bm{I}\) 为单位矩阵，并且从 \(e^{\bm{A}\times(t-t_0)}\) 的展开项知道 \(\alpha_k(t)\ (k=0,\ \cdots n-1)\) 互不相等，所以 \(\left[\,\alpha_0(t)\underset{\ell\times\ell}{\bm{I}}\quad \alpha_1(t)\underset{\ell\times\ell}{\bm{I}}\quad \cdots\quad \alpha_{n-1}(t)\underset{\ell\times\ell}{\bm{I}}\,\right]\) 并不改变 \(\bm{M}\) 各列的线性相关性，这样方程 (3.4.26) 有解的条件为 \[\mathrm{rank}(\bm{M})=n \tag{3.4.28}\] \(\bm{M}\) 也被称为可测性矩阵。
补出可测性推导的两个要点。第一，式 (3.4.22) 先把控制输入对观测的贡献 \(\bm{H}\int_{t_0}^{t}e^{\bm{A}(t-\tau)}\bm{B}u(\tau)\,\mathrm{d}\tau\) 从 \(\bm{Z}(t)\) 中扣除，得到 \(\overline{\bm{Z}}(t)=\bm{H}e^{\bm{A}(t-t_0)}\bm{X}(t_0)\)——这正是“可测性与输入无关”、可令 \(u=0\)（(3.4.24)）的依据。第二，对 \(e^{\bm{A}(t-t_0)}=\sum_{k=0}^{n-1}\alpha_k(t)\bm{A}^k\)（同样由 Cayley-Hamilton 定理）作分离变量：\(\overline{\bm{Z}}(t)=\left[\alpha_0(t)\bm{I}\ \ \alpha_1(t)\bm{I}\ \ \cdots\ \ \alpha_{n-1}(t)\bm{I}\right]\bm{M}\bm{X}(t_0)\)，其中 \(\bm{M}=\left[\bm{H}^{\mathrm{T}}\ \ \bm{A}^{\mathrm{T}}\bm{H}^{\mathrm{T}}\ \ \cdots\right]^{\mathrm{T}}\)。由于 \(\alpha_k(t)\) 互不相等且乘 \(\bm{I}\) 不改变列的线性相关性，\(\overline{\bm{Z}}(t)\) 能否携带 \(\bm{X}(t_0)\) 的全部信息就归结为 \(\mathrm{rank}(\bm{M})=n\)。离散情形 (3.4.31) 的推导更直接：把 \(n\) 个时刻的观测方程摞成一个线性方程组，未知量就是 \(\bm{X}(0)\)。
2. 离散线性时不变系统的可测性
前面已经说明了可测性与系统的输入无关，即有无输入并不改变的可测性，所以这里直接分析无输入的离散线性时不变系统的可测条件。离散线性时不变系统为 \[\begin{cases} \bm{X}(k)=\bm{\Phi}\bm{X}(k-1)\\ \bm{Z}(k)=\bm{H}\bm{X}(k) \end{cases} \tag{3.4.29}\] 可以递推得到 \[\begin{cases} \bm{Z}(0)=\bm{H}\bm{X}(0)\\ \bm{Z}(1)=\bm{H}\bm{X}(1)=\bm{H}\bm{\Phi}\bm{X}(0)\\ \bm{Z}(2)=\bm{H}\bm{X}(2)=\bm{H}\bm{\Phi}^2\bm{X}(0)\\ \qquad\vdots\\ \bm{Z}(n-1)=\bm{H}\bm{X}(n-1)=\bm{H}\bm{\Phi}^{n-1}\bm{X}(0) \end{cases} \tag{3.4.30}\] 将上式表示为 \[\begin{bmatrix}\underset{\ell\times 1}{\bm{Z}(0)}\\ \underset{\ell\times 1}{\bm{Z}(1)}\\ \vdots\\ \underset{\ell\times 1}{\bm{Z}(n-1)}\end{bmatrix} =\begin{bmatrix}\underset{\ell\times n}{\bm{H}}\\ \underset{\ell\times n}{\bm{H}\bm{\Phi}}\\ \vdots\\ \underset{\ell\times n}{\bm{H}\bm{\Phi}^{n-1}}\end{bmatrix}\underset{n\times 1}{\bm{X}(0)} \tag{3.4.31}\] 设 \[\bm{M}=\begin{bmatrix}\underset{\ell\times n}{\bm{H}}\\ \underset{\ell\times n}{\bm{H}\bm{\Phi}}\\ \vdots\\ \underset{\ell\times n}{\bm{H}\bm{\Phi}^{n-1}}\end{bmatrix} \tag{3.4.32}\]
显然，当线性方程组 (3.4.31) 系数矩阵满秩时，可以唯一地解得 \(\bm{X}(t_0)\)，所以当 \(\mathrm{rank}(\bm{M})=n\)，根据 \(t_0\) 到 \(t_1\) 的观测值 \(\bm{Z}(t)\)（\(t\in\left[\,t_0,\ t_1\,\right]\)），可以唯一地确定系统在初始时刻的状态 \(\bm{X}(t_0)\)，这样的系统是完全可观测的。
可测性判据的三点提醒。其一，可测性与系统输入无关：由 (3.4.22) 可见控制输入对输出的贡献被完全扣除，所以只要 \((\bm{A},\ \bm{H})\) 满足 \(\mathrm{rank}(\bm{M})=n\) 就完全可测——例 3.11 中改变 \(\bm{B}\) 后可控性改变而可测性不变，正是这个道理。其二，\(\bm{M}\) 是 \(n\ell\times n\) 阵，判据是列满秩 \(\mathrm{rank}(\bm{M})=n\)（与可控性的行满秩互为转置对偶）。其三，秩条件只是“信息充分性”：\(\mathrm{rank}(\bm{M})=n\) 保证 \(\bm{X}(t_0)\) 可唯一确定，但不保证数值上是否病态——实际估计精度还依赖观测噪声水平，那属于第 4 章 Kalman 滤波精度分析的范畴，不要与这里的定性判据混为一谈。
可测性（能观性）的判据与可控性判据构成对偶，在《广义测量平差》§4-8“线性确定系统的能观性和能控性”中有系统论述，含随机系统的能观能控则见该书 §4-9。两书对连续系统均以 \((\bm{A},\ \bm{H})\)、对离散系统均以 \((\bm{\Phi},\ \bm{H})\) 构造判定矩阵。“矩阵满秩 \(\Leftrightarrow\) 方程组有唯一解”这一思路与第 2 章 (2.1.9) 要求设计矩阵列满秩 \(\mathrm{rank}(\bm{H})=n\) 完全同源；例 3.11、例 3.12 的秩计算是与判据配套的直观演练，可与《广义测量平差》§4-8 的例题互相参照。
算例分析
例 3.11某系统的状态方程为 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\end{bmatrix} =\begin{bmatrix}1 & 0\\ 2 & 3\end{bmatrix}\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix}+\begin{bmatrix}0\\ 2\end{bmatrix}u(t)\] 观测方程为 \[Z(t)=\left[\begin{array}{ll}1 & 0\end{array}\right]\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix}\] 判断此系统是否可控和可测。若改变 \(\bm{B}\) 矩阵为 \(\left[\begin{array}{ll}2 & 0\end{array}\right]^{\mathrm{T}}\)，判断其是否可控和可测。
解：在此问题中 \(n=2\)，根据式 (3.4.13) 的判断可控性条件 \[\bm{C}=\left[\,\bm{B}\ \vdots\ \bm{A}\bm{B}\,\right]=\begin{bmatrix}0 & 0\\ 2 & 6\end{bmatrix}\] \[\mathrm{rank}(\bm{C})=1\] 所以系统不可控。
判断其可测性 \[\bm{M}=\begin{bmatrix}\bm{H}\\ \bm{H}\bm{A}\end{bmatrix}=\begin{bmatrix}1 & 0\\ 1 & 0\end{bmatrix}\] 所以此系统不可测。
若改变 \(\bm{B}\) 矩阵为 \(\left[\begin{array}{ll}2 & 0\end{array}\right]^{\mathrm{T}}\)，这时状态方程为 \[\begin{bmatrix}\dot{X}_1(t)\\ \dot{X}_2(t)\end{bmatrix} =\begin{bmatrix}1 & 0\\ 2 & 3\end{bmatrix}\begin{bmatrix}X_1(t)\\ X_2(t)\end{bmatrix}+\begin{bmatrix}2\\ 0\end{bmatrix}u(t)\] 由于 \[\bm{C}=\left[\,\bm{B}\ \vdots\ \bm{A}\bm{B}\,\right]=\begin{bmatrix}2 & 2\\ 0 & 4\end{bmatrix}\] 所以系统可控。
系统的可测性与输入无关，虽然改变了系统的输入，但这并不影响其可测性，所以系统仍然不可测。
例 3.12某离散线性时不变系统为 \[\bm{X}(k)=\begin{bmatrix}1 & -1\\ 1 & 1\end{bmatrix}\bm{X}(k-1)+\begin{bmatrix}2\\ 1\end{bmatrix}u(k-1)\] \[\bm{Z}(k)=\begin{bmatrix}1 & 0\\ -1 & 1\end{bmatrix}\bm{X}(k)\] 判断这个系统的可控性和可测性。
解：此问题中 \(n=2\)，根据式 (3.4.17)，可控性矩阵为 \[\bm{C}=\left[\,\bm{\Psi}\ \ \bm{\Phi}\bm{\Psi}\,\right]=\begin{bmatrix}2 & 1\\ 1 & 3\end{bmatrix}\] \[\mathrm{rank}(\bm{C})=2\] 所以系统是可控的。
根据式 (3.4.32)，可测性矩阵为 \[\bm{M}=\begin{bmatrix}\bm{H}\\ \bm{H}\bm{\Phi}\end{bmatrix}\] 由于 \(\mathrm{rank}(\bm{H})=2\)，所以一定有 \(\mathrm{rank}(\bm{M})=2\)，系统是可测的。