第12章：系统辨识

系统辨识是从输入输出数据中建立动态系统数学模型的科学与艺术。在实际控制工程中，精确的系统模型往往难以从第一性原理推导，或者系统参数会随时间变化。系统辨识技术使我们能够从实验数据中提取系统的动态特性，为控制器设计提供可靠的模型基础。本章将介绍主要的辨识方法，从经典的参数估计到现代的子空间方法，并特别关注实际应用中的挑战，如闭环辨识和非线性系统辨识。

12.1 系统辨识概述

12.1.1 辨识的基本问题

系统辨识的核心任务是找到一个数学模型 $M$，使其能够最好地描述观测到的输入输出数据：

$$y(t) = G(q, \theta)u(t) + H(q, \theta)e(t)$$ 其中：

$G(q, \theta)$ 是系统传递函数
$H(q, \theta)$ 是噪声模型
$\theta$ 是待辨识的参数向量
$e(t)$ 是白噪声序列
$q$ 是移位算子

12.1.2 辨识流程

实验设计 → 数据采集 → 模型结构选择 → 参数估计 → 模型验证
    ↑                                                    ↓
    ←────────────── 不满足要求则重新设计 ←─────────────

12.1.3 模型类型选择

线性模型族 - ARX (AutoRegressive with eXogenous input) - ARMAX (ARX with Moving Average) - Output Error (OE) - Box-Jenkins (BJ)
状态空间模型 - 连续时间：$\dot{x} = Ax + Bu$ - 离散时间：$x_{k+1} = Ax_k + Bu_k$
非参数模型 - 频率响应 - 脉冲响应 - 阶跃响应

12.2 参数估计方法

12.2.1 最小二乘法（LS）

考虑ARX模型： $$A(q)y(t) = B(q)u(t) + e(t)$$ 其中： $$A(q) = 1 + a_1q^{-1} + \cdots + a_{n_a}q^{-n_a}$$ $$B(q) = b_1q^{-1} + \cdots + b_{n_b}q^{-n_b}$$ 写成线性回归形式： $$y(t) = \phi^T(t)\theta + e(t)$$ 其中回归向量： $$\phi(t) = [-y(t-1), \ldots, -y(t-n_a), u(t-1), \ldots, u(t-n_b)]^T$$ 参数向量： $$\theta = [a_1, \ldots, a_{n_a}, b_1, \ldots, b_{n_b}]^T$$ 最小二乘估计： $$\hat{\theta}_{LS} = \arg\min_\theta \sum_{t=1}^N [y(t) - \phi^T(t)\theta]^2$$ 解析解： $$\hat{\theta}_{LS} = \left[\sum_{t=1}^N \phi(t)\phi^T(t)\right]^{-1} \sum_{t=1}^N \phi(t)y(t)$$

12.2.2 递推最小二乘（RLS）

实时更新参数估计： $$\hat{\theta}(t) = \hat{\theta}(t-1) + K(t)[y(t) - \phi^T(t)\hat{\theta}(t-1)]$$ 其中增益矩阵： $$K(t) = P(t-1)\phi(t)[1 + \phi^T(t)P(t-1)\phi(t)]^{-1}$$ 协方差矩阵更新： $$P(t) = [I - K(t)\phi^T(t)]P(t-1)$$ 遗忘因子：处理时变系统 $$P(t) = \frac{1}{\lambda}[I - K(t)\phi^T(t)]P(t-1)$$ 其中 $0.95 \leq \lambda < 1$ 是遗忘因子。

12.2.3 极大似然估计（MLE）

假设噪声服从高斯分布 $e(t) \sim \mathcal{N}(0, \sigma^2)$，似然函数： $$L(\theta) = \prod_{t=1}^N \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left(-\frac{[y(t) - \hat{y}(t|\theta)]^2}{2\sigma^2}\right)$$ 对数似然： $$\log L(\theta) = -\frac{N}{2}\log(2\pi\sigma^2) - \frac{1}{2\sigma^2}\sum_{t=1}^N [y(t) - \hat{y}(t|\theta)]^2$$ MLE估计： $$\hat{\theta}_{MLE} = \arg\max_\theta \log L(\theta)$$

12.2.4 预测误差方法（PEM）

一般化的辨识准则： $$\hat{\theta} = \arg\min_\theta V_N(\theta) = \arg\min_\theta \frac{1}{N}\sum_{t=1}^N \ell(\varepsilon(t, \theta))$$ 其中：

$\varepsilon(t, \theta)$ 是预测误差
$\ell(\cdot)$ 是损失函数（如二次损失）

对于不同模型结构，预测误差有不同形式：

ARX: $\varepsilon(t) = A(q)y(t) - B(q)u(t)$
OE: $\varepsilon(t) = y(t) - B(q)/F(q)u(t)$
ARMAX: $\varepsilon(t) = [A(q)y(t) - B(q)u(t)]/C(q)$

12.3 频域辨识

12.3.1 频率响应估计

通过正弦激励获取频率响应：

输入：$u(t) = A\sin(\omega t)$ 稳态输出：$y(t) = |G(j\omega)|A\sin(\omega t + \angle G(j\omega))$

频率响应函数（FRF）： $$\hat{G}(j\omega) = \frac{Y(j\omega)}{U(j\omega)} = \frac{\mathcal{F}\{y(t)\}}{\mathcal{F}\{u(t)\}}$$

12.3.2 谱分析方法

经验传递函数估计（ETFE）： $$\hat{G}_{ETFE}(e^{j\omega}) = \frac{Y_N(e^{j\omega})}{U_N(e^{j\omega})}$$ 其中 $Y_N$ 和 $U_N$ 是DFT变换。

谱分析估计： $$\hat{G}(e^{j\omega}) = \frac{\hat{\Phi}_{yu}(e^{j\omega})}{\hat{\Phi}_u(e^{j\omega})}$$ 其中：

$\hat{\Phi}_{yu}$ 是输入输出互功率谱
$\hat{\Phi}_u$ 是输入自功率谱

12.3.3 频域加权最小二乘

给定频率响应数据 $\{G_k, \omega_k\}_{k=1}^M$，拟合参数模型： $$\hat{\theta} = \arg\min_\theta \sum_{k=1}^M W_k |G_k - G(e^{j\omega_k}, \theta)|^2$$ 其中 $W_k$ 是频率相关的权重。

12.4 子空间辨识方法

12.4.1 基本原理

子空间方法直接从输入输出数据辨识状态空间模型： $$\begin{aligned} x_{k+1} &= Ax_k + Bu_k + w_k \\ y_k &= Cx_k + Du_k + v_k \end{aligned}$$ 无需迭代优化，通过奇异值分解（SVD）获得模型。

12.4.2 数据方程

构造Hankel矩阵： $$Y_{p} = \begin{bmatrix} y_0 & y_1 & \cdots & y_{j-1} \\ y_1 & y_2 & \cdots & y_j \\ \vdots & \vdots & \ddots & \vdots \\ y_{i-1} & y_i & \cdots & y_{i+j-2} \end{bmatrix}$$ 类似构造 $U_p$（过去输入）、$Y_f$（未来输出）、$U_f$（未来输入）。

12.4.3 MOESP算法

构造数据矩阵并进行QR分解： $$\begin{bmatrix} U_p \\ U_f \\ Y_p \\ Y_f \end{bmatrix} = \begin{bmatrix} R_{11} & 0 & 0 & 0 \\ R_{21} & R_{22} & 0 & 0 \\ R_{31} & R_{32} & R_{33} & 0 \\ R_{41} & R_{42} & R_{43} & R_{44} \end{bmatrix} \begin{bmatrix} Q_1^T \\ Q_2^T \\ Q_3^T \\ Q_4^T \end{bmatrix}$$
对 $R_{32}$ 进行SVD： $$R_{32} = U\Sigma V^T = \begin{bmatrix} U_1 & U_2 \end{bmatrix} \begin{bmatrix} \Sigma_1 & 0 \\ 0 & \Sigma_2 \end{bmatrix} \begin{bmatrix} V_1^T \\ V_2^T \end{bmatrix}$$
确定系统阶次 $n$（通过奇异值的跳变）
提取系统矩阵： - 可观性矩阵：$\mathcal{O} = U_1\Sigma_1^{1/2}$ - 从 $\mathcal{O}$ 提取 $C$ 和 $A$ - 通过最小二乘估计 $B$ 和 $D$

12.4.4 N4SID算法

N4SID（Numerical algorithms for Subspace State Space System IDentification）是另一种流行的子空间方法：

斜投影： $$Y_f / U_f \perp W_p$$ 其中 $W_p = \begin{bmatrix} U_p \\ Y_p \end{bmatrix}$
通过SVD估计状态序列
使用最小二乘估计系统矩阵

12.5 非线性系统辨识

12.5.1 Hammerstein-Wiener模型

级联结构：静态非线性 + 线性动态 + 静态非线性

u → [静态非线性 f(·)] → v → [线性系统 G(q)] → w → [静态非线性 h(·)] → y

辨识策略：

分离辨识：先辨识线性部分，再辨识非线性部分
迭代优化：交替优化线性和非线性参数

12.5.2 NARMAX模型

非线性自回归移动平均模型： $$y(t) = f(y(t-1), \ldots, y(t-n_y), u(t-1), \ldots, u(t-n_u), e(t-1), \ldots, e(t-n_e)) + e(t)$$ 常用非线性函数形式：

多项式
神经网络
径向基函数（RBF）
小波网络

12.5.3 神经网络辨识

循环神经网络（RNN）模型： $$\begin{aligned} h_t &= \sigma(W_h h_{t-1} + W_u u_t + b_h) \\ y_t &= W_y h_t + b_y \end{aligned}$$ LSTM模型：处理长期依赖 $$\begin{aligned} f_t &= \sigma(W_f [h_{t-1}, u_t] + b_f) \\ i_t &= \sigma(W_i [h_{t-1}, u_t] + b_i) \\ \tilde{C}_t &= \tanh(W_C [h_{t-1}, u_t] + b_C) \\ C_t &= f_t \odot C_{t-1} + i_t \odot \tilde{C}_t \\ o_t &= \sigma(W_o [h_{t-1}, u_t] + b_o) \\ h_t &= o_t \odot \tanh(C_t) \end{aligned}$$

12.5.4 稀疏辨识（SINDy）

系统动力学的稀疏辨识： $$\dot{x} = \Theta(x, u)\xi$$ 其中：

$\Theta(x, u)$ 是候选函数库（多项式、三角函数等）
$\xi$ 是稀疏系数向量

通过L1正则化促进稀疏性： $$\hat{\xi} = \arg\min_\xi ||\dot{X} - \Theta\xi||_2^2 + \lambda||\xi||_1$$

12.6 闭环辨识问题

12.6.1 闭环辨识的挑战

在许多实际系统中，出于安全或稳定性考虑，必须在闭环条件下进行辨识：

      r(t) →⊕→ [控制器 C] → u(t) → [系统 G] → y(t)
            ↑                                    ↓
            ←────────────────────────────────────

主要问题：

输入输出相关性：$u(t)$ 不再是独立激励
可辨识性降低：某些模态可能被控制器抑制
偏差问题：直接方法会产生有偏估计

12.6.2 直接法

忽略反馈回路，直接使用开环辨识方法： $$\hat{G} = \arg\min_G \sum_{t=1}^N [y(t) - G(q)u(t)]^2$$ 适用条件：

参考信号 $r(t)$ 充分激励
噪声较小
控制器不过于激进

12.6.3 间接法

同时辨识闭环传递函数，然后推导开环模型：

闭环传递函数： $$T(q) = \frac{GC}{1 + GC}, \quad S(q) = \frac{1}{1 + GC}$$ 从测量数据辨识 $T$ 和 $S$，然后计算： $$G = \frac{T}{C \cdot S}$$

12.6.4 联合输入输出法

使用输入输出的联合信息构造无偏估计： $$\hat{\theta} = \arg\min_\theta \left|\left|\begin{bmatrix} Y \\ U \end{bmatrix} - \begin{bmatrix} G(\theta) \\ I \end{bmatrix} \begin{bmatrix} U \\ E \end{bmatrix}\right|\right|_F^2$$

12.6.5 双级法（Two-Stage Method）

第一级：辨识闭环灵敏度函数 $$\hat{S} = \arg\min_S ||Y - S \cdot R||^2$$ 第二级：使用模拟开环数据 $$\tilde{u} = \hat{S} \cdot r, \quad \tilde{y} = \hat{S} \cdot y$$ 然后用开环方法辨识： $$\hat{G} = \arg\min_G ||\tilde{y} - G \cdot \tilde{u}||^2$$

12.6.6 预测误差方法的闭环扩展

考虑闭环下的预测误差： $$\varepsilon(t, \theta) = H^{-1}(\theta)[y(t) - G(\theta)u(t)]$$ 其中 $u(t) = C[r(t) - y(t)]$

优化准则： $$\hat{\theta} = \arg\min_\theta \sum_{t=1}^N \varepsilon^2(t, \theta)$$

12.7 实践考虑

12.7.1 实验设计

激励信号选择：

PRBS（伪随机二进制序列）： - 宽频谱特性 - 实现简单 - 适合线性系统
多正弦信号： $$u(t) = \sum_{k=1}^M A_k \sin(\omega_k t + \phi_k)$$

精确控制频率内容
避免非线性失真

Chirp信号： $$u(t) = A \sin(\omega(t) \cdot t)$$ 其中 $\omega(t)$ 线性增长

采样率选择：

至少10倍于系统带宽
考虑抗混叠滤波器

数据长度：

参数数量的10-20倍（最小）
考虑计算复杂度

12.7.2 预处理

去趋势：消除直流偏移和线性趋势
滤波：去除高频噪声
归一化：改善数值条件
异常值处理：检测和处理离群点

12.7.3 模型验证

残差分析：

白度检验：$R_{\varepsilon\varepsilon}(\tau) \approx \delta(\tau)$
独立性检验：$R_{\varepsilon u}(\tau) \approx 0$

交叉验证：

将数据分为训练集和验证集
计算验证误差： $$V_{val} = \frac{1}{N_{val}} \sum_{t \in \text{val}} [y(t) - \hat{y}(t|\theta_{train})]^2$$ FIT值： $$\text{FIT} = 100 \times \left(1 - \frac{||y - \hat{y}||}{||y - \bar{y}||}\right) \%$$

12.7.4 模型阶次选择

信息准则：

AIC: $\text{AIC} = N \log(\hat{\sigma}^2) + 2p$
BIC: $\text{BIC} = N \log(\hat{\sigma}^2) + p \log(N)$
MDL: $\text{MDL} = N \log(\hat{\sigma}^2) + p \log(N)/2$

其中 $p$ 是参数数量，$\hat{\sigma}^2$ 是残差方差。

12.8 案例研究：锂电池等效电路模型辨识

12.8.1 问题背景

锂离子电池的精确建模对电池管理系统（BMS）至关重要。等效电路模型（ECM）因其简单性和实用性被广泛采用。我们需要从充放电数据中辨识模型参数，用于SOC（State of Charge）估计和健康状态监测。

12.8.2 二阶RC等效电路模型

     R0
I →──/\/\/──┬──────────────────→ V
            │
           ┌┴┐    R1        R2
           │ │OCV  ║         ║
           └┬┘    ┌┴┐       ┌┴┐
            │     │C1│      │C2│
            │     └┬┘       └┬┘
            │      │         │
            └──────┴─────────┘

模型方程： $$V(t) = OCV(SOC) - I(t)R_0 - V_1(t) - V_2(t)$$ 其中RC网络动态： $$\begin{aligned} \dot{V}_1 &= -\frac{1}{R_1C_1}V_1 + \frac{1}{C_1}I \\ \dot{V}_2 &= -\frac{1}{R_2C_2}V_2 + \frac{1}{C_2}I \end{aligned}$$

12.8.3 辨识流程

步骤1：OCV-SOC曲线辨识

通过低倍率充放电实验，获取准静态OCV-SOC关系： $$OCV(SOC) = a_0 + a_1 \cdot SOC + a_2 \cdot SOC^2 + \cdots + a_n \cdot SOC^n$$ 使用最小二乘拟合多项式系数。

步骤2：内阻R0辨识

通过脉冲放电测试，测量瞬时电压跌落： $$R_0 = \frac{\Delta V_{instant}}{\Delta I}$$ 步骤3：RC参数辨识

离散化RC网络方程（采样时间 $T_s$）： $$\begin{aligned} V_1(k+1) &= e^{-T_s/(R_1C_1)} V_1(k) + R_1(1 - e^{-T_s/(R_1C_1)}) I(k) \\ V_2(k+1) &= e^{-T_s/(R_2C_2)} V_2(k) + R_2(1 - e^{-T_s/(R_2C_2)}) I(k) \end{aligned}$$ 定义参数： $$\begin{aligned} \alpha_1 &= e^{-T_s/(R_1C_1)}, \quad \beta_1 = R_1(1 - \alpha_1) \\ \alpha_2 &= e^{-T_s/(R_2C_2)}, \quad \beta_2 = R_2(1 - \alpha_2) \end{aligned}$$ 构造回归模型： $$V(k) = OCV(SOC_k) - R_0 I(k) - V_1(k) - V_2(k)$$ 使用递推最小二乘（RLS）在线辨识 $\alpha_i, \beta_i$。

12.8.4 实验设计

HPPC测试（Hybrid Pulse Power Characterization）：

电流 I
  │     10s放电          40s休息        10s充电
  │  ┌─────────┐                    ┌─────────┐
  │  │         │                    │         │
──┼──┘         └────────────────────┘         └──
  │
  └─────────────────────时间 t

这种脉冲序列能够激发不同时间常数的动态响应。

12.8.5 参数随SOC和温度的变化

实际应用中，参数是SOC和温度的函数： $$\theta = f(SOC, T)$$ 建立查找表或拟合函数： $$R_0(SOC, T) = r_{00} + r_{01} \cdot SOC + r_{10} \cdot T + r_{11} \cdot SOC \cdot T$$

12.8.6 验证结果

使用动态应力测试（DST）或联邦城市驾驶循环（FUDS）验证模型：

电压误差统计：

- 最大误差: < 50mV
- 均方根误差(RMSE): < 20mV  
- 平均绝对误差(MAE): < 15mV

12.8.7 在线参数更新

考虑电池老化，使用自适应算法在线更新参数： $$\begin{aligned} \hat{\theta}(k) &= \hat{\theta}(k-1) + K(k) \varepsilon(k) \\ K(k) &= \frac{P(k-1)\phi(k)}{\lambda + \phi^T(k)P(k-1)\phi(k)} \\ P(k) &= \frac{1}{\lambda}[P(k-1) - K(k)\phi^T(k)P(k-1)] \end{aligned}$$ 其中 $\lambda = 0.98$ 是遗忘因子，适应缓慢的老化过程。

12.9 历史人物：Lennart Ljung (1987)

Lennart Ljung 是系统辨识理论体系化的关键人物。他在1987年出版的经典著作《System Identification: Theory for the User》成为该领域的权威参考。Ljung的主要贡献包括：

预测误差方法的统一框架：将各种辨识方法统一在预测误差最小化的框架下
System Identification Toolbox：开发了MATLAB系统辨识工具箱，使辨识技术广泛应用
渐近理论：建立了参数估计的统计性质理论基础
实用算法：提出了许多实用的数值算法，如Gauss-Newton迭代

他的工作将系统辨识从分散的技术集合发展成系统的学科，极大推动了控制工程的实践应用。

12.10 前沿专题：稀疏辨识与符号回归在控制中的应用

12.10.1 稀疏辨识的动机

传统辨识方法需要预先假定模型结构，而实际系统的真实结构往往未知。稀疏辨识通过数据驱动的方式自动发现系统的控制方程，特别适合：

复杂非线性系统
部分可观测系统
高维系统的降阶建模

12.10.2 SINDy算法扩展

带控制输入的SINDy（SINDyc）： $$\dot{x} = \Theta(x)\Xi_x + \Theta(u)\Xi_u$$ 噪声鲁棒的SINDy：使用积分形式减少噪声影响： $$x(t) - x(0) = \int_0^t \Theta(x(\tau), u(\tau))\xi d\tau$$ 约束SINDy：加入物理约束（如能量守恒）： $$\begin{aligned} \min_\xi &\quad ||\dot{X} - \Theta\xi||_2^2 + \lambda||\xi||_1 \\ \text{s.t.} &\quad C\xi = d \end{aligned}$$

12.10.3 符号回归方法

遗传编程（GP）：通过进化算法搜索函数空间：

初始化：随机生成表达式树群体
评估：计算每个表达式的拟合误差
选择：选择适应度高的个体
交叉变异：生成新的表达式
迭代直到收敛

深度符号回归：结合神经网络和符号回归：

神经网络学习特征表示
符号层提取可解释的方程

12.10.4 在控制中的应用

模型预测控制的符号化：将复杂的MPC策略蒸馏为简单的符号表达式： $$u = f_{symbolic}(x, r)$$ 优势：

计算效率高
可解释性强
易于嵌入式实现

未知动力学的自适应控制：

在线稀疏辨识系统动力学
基于辨识模型设计控制器
持续更新模型和控制器

12.11 本章小结

系统辨识是连接理论与实践的桥梁，本章介绍了从经典到现代的主要辨识方法：

核心概念：

参数估计：最小二乘、极大似然、预测误差方法
频域方法：谱分析、频率响应拟合
子空间方法：直接从数据辨识状态空间模型
非线性辨识：Hammerstein-Wiener、NARMAX、神经网络
闭环辨识：处理反馈存在下的辨识问题

关键公式：

最小二乘：$\hat{\theta} = (X^TX)^{-1}X^Ty$
递推更新：$\hat{\theta}(k) = \hat{\theta}(k-1) + K(k)\varepsilon(k)$
预测误差：$V(\theta) = \frac{1}{N}\sum_{t=1}^N \varepsilon^2(t,\theta)$
信息准则：$\text{AIC} = N\log(\hat{\sigma}^2) + 2p$

实践要点：

实验设计决定辨识质量
模型验证不可或缺
闭环辨识需要特殊处理
在线更新应对时变系统

12.12 练习题

基础题

习题12.1 最小二乘估计的性质给定线性回归模型 $y = X\theta + e$，其中 $e \sim \mathcal{N}(0, \sigma^2I)$。证明最小二乘估计 $\hat{\theta}_{LS}$ 是无偏的，并求其协方差矩阵。

提示：利用 $E[e] = 0$ 和 $\text{Cov}(e) = \sigma^2I$

答案

无偏性证明： $$E[\hat{\theta}_{LS}] = E[(X^TX)^{-1}X^Ty] = (X^TX)^{-1}X^T E[X\theta + e] = \theta$$ 协方差矩阵： $$\text{Cov}(\hat{\theta}_{LS}) = (X^TX)^{-1}X^T \text{Cov}(y) X(X^TX)^{-1} = \sigma^2(X^TX)^{-1}$$

习题12.2 持续激励条件考虑一阶系统 $y(k) = ay(k-1) + bu(k-1) + e(k)$。如果输入信号是常数 $u(k) = c$，能否唯一辨识参数 $a$ 和 $b$？说明原因。

提示：考虑信息矩阵的秩

答案

不能唯一辨识。当输入为常数时，稳态下 $y(k) = y_{ss} = \frac{bc}{1-a}$。只能得到一个约束方程，无法唯一确定两个参数。信息矩阵秩亏，需要输入信号至少包含两个不同频率成分。

习题12.3 模型阶次选择给定一组数据，使用不同阶次的ARX模型拟合，得到：

2阶：$\hat{\sigma}^2 = 1.2$，参数数4个
3阶：$\hat{\sigma}^2 = 1.0$，参数数6个
4阶：$\hat{\sigma}^2 = 0.95$，参数数8个

数据长度 $N = 100$。使用AIC准则选择最佳阶次。

提示：$\text{AIC} = N\log(\hat{\sigma}^2) + 2p$

答案

计算各阶AIC值：

2阶：$\text{AIC}_2 = 100\log(1.2) + 2(4) = 18.23 + 8 = 26.23$
3阶：$\text{AIC}_3 = 100\log(1.0) + 2(6) = 0 + 12 = 12$
4阶：$\text{AIC}_4 = 100\log(0.95) + 2(8) = -5.13 + 16 = 10.87$

选择4阶模型（AIC最小）。

挑战题

习题12.4 闭环可辨识性分析系统 $G(s) = \frac{K}{s+a}$ 在单位负反馈下运行。参考输入 $r(t)$ 是白噪声。分析： (a) 能否同时辨识 $K$ 和 $a$？ (b) 如果控制器改为 $C(s) = K_p + \frac{K_i}{s}$，可辨识性如何变化？

提示：分析闭环传递函数的参数化

答案

(a) 闭环传递函数：$T(s) = \frac{K}{s+a+K}$。只能辨识 $a+K$ 的组合，不能分别辨识 $K$ 和 $a$。

(b) 使用PI控制器后，闭环特征方程变为二阶，包含 $s^2 + (a+KK_p)s + KK_i = 0$。如果 $K_p$ 和 $K_i$ 已知，可以唯一确定 $K$ 和 $a$。

习题12.5 子空间方法的数值实现给定输入输出数据，实现MOESP算法的关键步骤。设系统阶次未知，如何从奇异值确定阶次？给出判断准则。

提示：观察奇异值的"膝点"

答案

阶次选择准则：

奇异值比：$\sigma_{n+1}/\sigma_n < \epsilon$（如0.01）
累积能量：$\sum_{i=1}^n \sigma_i^2 / \sum_{i=1}^{2i} \sigma_i^2 > 0.99$
奇异值间隙：$\max_i(\sigma_i - \sigma_{i+1})$ 对应的 $i$ 为阶次

实践中常结合多个准则，并通过验证数据确认。

习题12.6 电池模型的实时辨识设计一个实时算法，在电动汽车运行中持续更新电池模型参数。要求： (a) 计算复杂度 $O(n^2)$ 以内 (b) 能处理参数突变（如温度变化） (c) 保证数值稳定性

提示：考虑带遗忘因子的RLS或滑动窗口方法

答案

推荐方案：带遗忘因子和协方差重置的RLS

初始化：θ(0), P(0) = αI
循环：

  1. 计算预测误差：ε = y - φ'θ
  2. 更新增益：K = Pφ/(λ + φ'Pφ)
  3. 更新参数：θ = θ + Kε
  4. 更新协方差：P = (P - Kφ'P)/λ
  5. 监测突变：if |ε| > threshold, P = αI
  6. 数值稳定：每N步做P = (P + P')/2

复杂度：$O(n^2)$ per sample

习题12.7 非线性系统的SINDy辨识考虑Lorenz系统的部分观测（只测量 $x$ 和 $y$）。设计稀疏辨识方案恢复完整动力学。讨论库函数选择和正则化参数调节。

提示：利用时延嵌入重构状态空间

答案

方案：

时延嵌入：$\tilde{x} = [x(t), y(t), x(t-τ), y(t-τ), x(t-2τ), y(t-2τ)]$
候选库：多项式到3阶 + 交叉项
分步辨识：先辨识 $\dot{x}, \dot{y}$，再通过数值微分估计 $z$
正则化：使用交叉验证选择 $\lambda$，典型范围 $[0.001, 0.1]$

关键：选择合适的延迟 $τ$ 使嵌入空间展开吸引子。

习题12.8 开放性问题：深度学习与经典辨识的结合讨论如何结合深度学习的表达能力和经典辨识的理论保证。设计一个混合框架，用于辨识具有未知非线性但需要稳定性保证的控制系统。考虑：

模型结构设计
训练策略
稳定性验证
实时更新能力

答案

建议框架：

结构：线性骨架 + 神经网络残差 $$\dot{x} = Ax + Bu + f_{NN}(x,u)$$ 其中 $A$ 保证稳定，$f_{NN}$ 有界
训练： - 第一阶段：经典方法辨识线性部分 - 第二阶段：固定线性，训练NN拟合残差 - 第三阶段：联合微调
稳定性： - Lipschitz约束：谱归一化 - Lyapunov正则化：$\mathcal{L} = \mathcal{L}_{fit} + \lambda \max(0, \dot{V})$
在线更新： - 线性部分：RLS - 非线性部分：小批量SGD - 切换逻辑：基于预测误差

12.13 常见陷阱与错误（Gotchas）

数据相关陷阱

数据质量问题 - 陷阱：忽视传感器噪声、量化误差、采样抖动 - 解决：预滤波、过采样后降采样、时间戳校正
激励不充分 - 陷阱：使用单频正弦或阶跃信号辨识高阶系统 - 解决：PRBS、chirp或多正弦激励，验证持续激励条件
非线性影响 - 陷阱：线性模型辨识时忽视饱和、死区等非线性 - 解决：选择工作点、分段线性化、使用非线性模型

算法相关陷阱

初值敏感性 - 陷阱：非线性优化陷入局部最优 - 解决：多初值尝试、全局优化方法、物理直觉指导
过拟合 - 陷阱：模型阶次过高，训练误差小但泛化差 - 解决：交叉验证、信息准则、正则化
数值问题 - 陷阱：信息矩阵病态导致数值不稳定 - 解决：数据归一化、正则化、SVD代替直接求逆

实施相关陷阱

闭环忽视 - 陷阱：闭环数据用开环方法导致偏差 - 解决：使用专门的闭环辨识方法
时变忽视 - 陷阱：固定参数模型用于时变系统 - 解决：滑动窗口、遗忘因子、自适应算法
延迟处理 - 陷阱：忽略系统延迟导致模型错位 - 解决：同时辨识延迟、互相关分析确定延迟

调试技巧

残差诊断：始终检查残差的白度和独立性
物理合理性：检查辨识参数的物理意义
渐进复杂：从简单模型开始，逐步增加复杂度
仿真验证：用辨识模型仿真，对比实际输出

12.14 最佳实践检查清单

实验设计阶段

[ ] 确定辨识目标（控制设计、仿真、故障诊断）
[ ] 选择合适的工作点和工作范围
[ ] 设计激励信号（频谱覆盖、幅值范围）
[ ] 确定采样率（10倍带宽以上）
[ ] 规划数据长度（参数数的20倍以上）
[ ] 考虑安全约束和操作限制

数据采集阶段

[ ] 传感器校准和零点检查
[ ] 同步采样和时间戳记录
[ ] 数据完整性检查（丢包、异常值）
[ ] 记录实验条件（温度、湿度等）
[ ] 备份原始数据

预处理阶段

[ ] 去除趋势和直流偏移
[ ] 处理异常值和缺失数据
[ ] 必要的滤波（注意相位影响）
[ ] 数据分段（训练集、验证集、测试集）
[ ] 归一化处理

模型选择阶段

[ ] 从简单模型开始（ARX、低阶）
[ ] 考虑物理知识（质量守恒、能量守恒）
[ ] 比较不同模型结构
[ ] 使用信息准则选择阶次
[ ] 考虑模型的可解释性

参数估计阶段

[ ] 选择合适的损失函数
[ ] 设置合理的参数约束
[ ] 多初值尝试（非线性情况）
[ ] 监控优化收敛性
[ ] 记录计算时间和资源消耗

验证阶段

[ ] 残差白度检验
[ ] 输入输出互相关检验
[ ] 交叉验证（新数据）
[ ] 仿真预测对比
[ ] 不确定性量化

应用阶段

[ ] 模型简化（必要时）
[ ] 鲁棒性分析
[ ] 实时性评估
[ ] 文档化（假设、限制、使用条件）
[ ] 版本控制和更新机制