简介这份系统辨识教案面向自动化专业学生与相关课程学习者围绕从输入输出数据中建立系统模型这一核心问题梳理辨识理论、实验设计与编程计算的学习路径。压缩包共1个文件为doc格式文档大小约1.42MB内容以课时讲义形式组织涵盖参考书目、学习要求、考核办法与概论章节。文档从实体与模型的区分讲起介绍直觉模型、物理模型与数学模型等分类并展开Zadeh对系统辨识的定义及数据、模型类、准则三要素随后结合长网造纸过程的质量控制案例完整演示明确辨识目的、收集先验知识、设计试验、数据采集与预处理、结构辨识与参数估计、模型检验与转换等一般步骤。已有150人学习适合希望减少课堂依赖、通过自学掌握辨识方法并独立完成数据分析报告的读者参考。1. 拿到一份 40 学时的系统辨识教案先别急着翻公式如果你手上正躺着一份《系统辨识与建模》的教案40 学时、五本参考书、考核方式是“给一个数据文件编程辨识写报告”那这份资源的价值不在“读完”而在“照着它把一条辨识链路跑通”。我带过几届自动化方向的课程设计最深的体会是系统辨识这门课教材讲得再全学生卡住的地方永远是那几个——M 序列怎么设步宽、采样间隔和扰动周期怎么配、数据预处理做完参数怎么反演回去。这份教案把 Zadeh 的三要素数据、模型类、准则、辨识十一大步、M 序列设计公式、标准化反演公式都写进去了属于少见的“理论够用、步骤能抄”的讲义型资料。它适合两类人一是要独立完成一次辨识实验并交报告的学生二是需要快速把辨识流程落到代码里的初级工程师。下面我按“先立理论骨架、再拆实验设计、最后抠预处理和避坑”的顺序把它拆成能复现的东西。2. 三要素与十一大步把辨识流程拆成可执行清单2.1 为什么先锁死“目的—模型类—准则”这三件事教案开篇就引 Zadeh 1962 年的定义辨识是在输入输出数据基础上从一组给定模型类中确定一个与所观测系统等价的模型。这句话的实操含义是数据、模型类、准则三者必须同时固定否则不同方法之间根本没有可比性。我见过太多报告前面用差分方程拟合后面又拿传递函数去比误差结论自然站不住。选型上有个直接对应关系目的是验证理论模型就要求零极点、结构、时延都准模型类必须和理论模型同构目的是设计数字控制器动态响应和时延要准模型类要便于计算机运算差分方程比连续传递函数更顺手目的是监视过程参数做故障诊断那参数准确比动态特性重要模型类要能直观体现被监视量。教案里造纸过程的例子把这点讲透了——控制定量和水份先调查出车速同步困难要恒定、流量到定量约 60 秒延迟而响应只有 2 秒、浓度到定量约 120 秒延迟响应约 80 秒、蒸汽压力到水份约 45 秒延迟响应约 60 秒。这些先验知识直接决定了后面哪些输入加扰动、采样频率取多少。准则这块教案给了三种误差定义我建议新手先吃透方程误差和输出误差的区别误差类型准则形式对参数是否线性预测能力适用模型参数误差J∑/n是仅算法研究理论验证方程误差J∑e(k)²/Ne(k)y(k)ay(k-1)-bu(k-1)-cu(k-2)线性只能一步预测AYB1U1B2U2输出误差J∑e(k)²/Ne(k)y(k)-ŷ(k)非线性输入已知可多步预测B/A 结构方程误差是参数的线性函数最小二乘能直接解但只能做一步预测输出误差是参数的非线性函数得迭代但能做多步预测。教案还提醒了一个容易被忽略的点模型阶次越高 J 值越小会导致过估计所以不能拿 J 值单独比模型好坏采样间隔增大一倍模型阶次大致减半跨采样步长比较时误差公式要改成 J∑e(k)ᵀWe(k)·T0²/N。这条在写报告做模型对比时是硬约束。2.2 十一大步里真正会返工的是哪几步教案把流程画成明确目的→收集先验知识→设计辨识试验→现场准备→数据采集→数据预处理→选择模型类→结构辨识与参数估计→模型检验→模型转换→应用评价检验和应用评价不通过要回退。实际做下来返工最频繁的是第 3 步试验设计和第 6 步预处理因为这两步的错误在参数估计阶段才暴露而那时你已经以为快结束了。设计辨识试验的核心原则教案写得很直白在安全前提下尽可能激励系统、保持输入输出关系、适当解耦。落到操作上就是先划边界——从边界外连进来的其它信号尽量保持稳定并当作噪声再定哪些输入叠加扰动、哪些保持稳定再看过渡过程时间常数是否有明显差异、噪信比多大、噪声是白色还是有色。造纸例子里保持车速和纸浆浓度不变对流量和蒸汽压力分别施加伪随机序列扰动幅度以不引起断纸为限采样频率设为 2试验 1000 秒采集定量、水份、纸浆流量和蒸汽压力四路信号。这套参数不是拍脑袋是前面先验知识推出来的。现场准备那步教案列了五件事申请试验时段、准备扰动信号发生器并预发检验、测验现场干扰必要时设计模拟滤波器、准备模数转换设备调好零迁和放大、接线。这几件事在课程设计里常被简化成“直接采数”但零迁和放大没调好采回来的数据要么削顶要么分辨率不够后面预处理再努力也救不回来。3. M 序列与采样设计扰动信号参数怎么算才不翻车3.1 M 序列的周期、游程和自相关性质教案花了不小篇幅讲 M 序列因为它是辨识实验里最常用的扰动信号。n 阶 M 序列周期 N2ⁿ-1是一个周期内移位寄存器能表示的最多状态数一个周期里 1 的个数比 0 多 1游程长度从 1 到 n长度为 1 的游程有 2ⁿ⁻² 个长度每加一数量减半长度为 n 的游程只有 1 个。把 0 换成 -1 后自相关函数在 τiN 时为 1在 τiNkk1…N-1时为 -1/N。N 充分大时这个 -1/N 趋近 0近似 δ 函数所以以 M 序列为输入的线性系统互相关函数序列就等于脉冲响应序列——这是用相关分析法辨识的数学基础。多输入单输出系统有个实用技巧同一 M 序列的不同移位序列M(k)、M(kJ)、M(k2J)分别作为各输入扰动只要输出对各输入的过渡过程小于 J就可认为输入之间正交。这解决了多路扰动互相污染的问题但前提是 J 选得够大J 太小正交性就破了这是后面避坑章要展开的点。逆重复 M 序列是 M 序列乘方波得到的周期 2N一个周期里 1 和 -1 个数相等自相关在 τi·2N 时为 1、τi·2N1 时为 -1、其余为 0。它保留了白噪声的优点而幅值可控直流分量问题比普通 M 序列更干净。3.2 步宽 Ts 和级数 n 的两套算法教案给了两套 M 序列设计方法我建议新手用第一套逻辑链清楚# M 序列参数设计第一套算法 # 输入系统过渡过程时间 T秒系数 kfmax k/T常取 3~15 T 60.0 # 系统过渡过程时间由阶跃预实验估计 k 5 # 最高工作频率系数fmax k/T fmax k / T # 系统最高工作频率 # 步骤1由 1/Ts 3*fmax 得 Ts 上限 Ts_max 1.0 / (3 * fmax) Ts Ts_max * 0.8 # 留余量实际取略小于上限 # 步骤2由 1.5*T (2^n - 1)*Ts 反解 n import math n math.ceil(math.log2(1.5 * T / Ts 1)) N 2**n - 1 print(fTs{Ts:.3f}s, n{n}, N{N}, 周期长度{N*Ts:.1f}s)这段代码的逻辑是先由过渡过程时间 T 估出最高工作频率 fmax(3~15)/T再由 1/Ts3·fmax 定出步宽上限最后用 1.5·T(2ⁿ-1)·Ts 保证一个周期覆盖至少 1.5 倍过渡过程反解出最小阶次 n。参数说明k 取 3 偏保守、取 15 偏激进课程设计里取 5 左右比较稳Ts 乘 0.8 是给实际采样留余量因为现场信号发生器的时间精度未必理想。第二套算法是反过来先定 n 再算 Ts(0.1~0.8)T/n适合已经知道要用几阶移位寄存器的场合。多输入时还要加一条约束(2ⁿ-1-1)·Ts1.5·∑(输出对各输入的过渡过程时间)否则移位序列之间的正交性不成立。采样间隔单独设计满足香农定理 1/T02·fmax有扰动信号时 T0Ts且 Ts 必须是 T0 的整数倍。实验时间取 3~5 个扰动信号周期预实验 1~2 个周期。这几条是硬约束T0 和 Ts 不整除后面重抽样和对齐就会出小数索引血泪经验。3.3 白噪声、有色噪声和成形滤波器教案指出白噪声谱密度为常数、频谱宽度无限、不同时刻互不相关用途包括激发系统所有模态、防止数据病态、以输出估计误差是否具有白色性判断辨识方法优劣。物理上白噪声不存在常见的是有色噪声而有色噪声表示定理说任何谱密度为实函数的平稳噪声都可由白噪声经过一个渐近稳定的线性环节成形滤波器得到。这条定理的实操意义是如果你的残差不是白的说明模型没把噪声的成形环节吸收进去要么加阶次要么换模型结构而不是简单说“噪声太大”。4. 数据预处理标准化、滤波、重抽样和参数反演4.1 标准化为什么必须做反演公式怎么推教案把预处理拆成四件事传感器非线性校正、滤波、重抽样、标准化。前三件偏工程标准化是数学上最容易做错又最影响结果的。目的写得很清楚降低舍入误差、防止病态方程、提高计算精度。方法是对输入输出各自零均值、单位方差化辨识完再把参数反演回原尺度。以二阶差分方程为例原始模型是 y(k)a·y(k-1)b·u(k-1)c·u(k-2)d。先减均值消掉常数项 d得到 Y(k)a·Y(k-1)b·U(k-1)c·U(k-2)再除以标准差得到 Yσ(k)a·Yσ(k-1)B·Uσ(k-1)C·Uσ(k-2)其中 Bb·σu/σy、Cc·σu/σy。辨识得到 B、C 后反演bB·σy/σucC·σy/σu。a 不变因为它作用在同尺度的 Y 上。import numpy as np def normalize_io(u, y): 输入输出零均值单位方差化返回标准化数据与反演系数 u_mean, y_mean u.mean(), y.mean() U u - u_mean Y y - y_mean sigma_u, sigma_y U.std(), Y.std() U_sigma U / sigma_u Y_sigma Y / sigma_y return U_sigma, Y_sigma, sigma_u, sigma_y, u_mean, y_mean def denormalize_params(B, C, sigma_u, sigma_y): 把标准化域辨识出的 B、C 反演回原尺度 b、c b B * sigma_y / sigma_u c C * sigma_y / sigma_u return b, c逻辑说明先减均值是为了消掉模型里的常数项 d否则 d 会和参数耦合导致有偏再除标准差是为了让方程系数落在同一量级避免最小二乘的法方程条件数爆炸。参数说明σu、σy 用样本标准差ddof0 还是 1 影响不大但全流程要统一反演时 b、c 都乘 σy/σu因为输入被缩小了 σu 倍、输出被缩小了 σy 倍系数要补回这个比例。注意 a 不需要反演它对应的是 Y(k) 和 Y(k-1)同尺度相消。4.2 滤波、重抽样和非线性校正的先后顺序教案说滤波时对输入输出加相同的滤波器输入输出关系不变——这句话的前提是同一个滤波器、同样的相位响应否则输入输出之间的时移对不上辨识出的时延就是错的。重抽样针对过渡过程较长的子系统重新选采样间隔以降低模型阶次配合前面“采样间隔翻倍、阶次减半”的结论用。传感器非线性校正要按厂家给的静态非线性曲线编程序校正这一步在课程设计里常被跳过但如果你的传感器在量程两端明显非线性不校正的话辨识出的模型会把非线性折进参数里换个工作点就失效。顺序上我一般按非线性校正→滤波→重抽样→标准化。理由是非线性校正是逐点映射必须先做滤波和重抽样都改变数据长度或频带放在标准化之前保证标准化统计量是基于最终数据算的。5. 避坑与排查五个让辨识报告返工的具体问题5.1 现象辨识出的时延总是差一两个采样点原因输入输出加了不同滤波器或者滤波器的群延迟没对齐导致互相关峰值偏移。教案明确要求输入输出加相同滤波器但实操里常有人给输出加低通、给输入不加理由是“输入干净”。解决输入输出走同一个滤波器对象、同一组系数如果必须分别处理先算滤波器的群延迟并在数据上做补偿再用互相关找时延。5.2 现象最小二乘解出的参数发散或大得离谱原因没做标准化输入输出量级差几个数量级法方程条件数过大或者模型里保留了常数项 d 又做了零均值化造成共线。解决先零均值单位方差化再辨识常数项单独估计或直接消掉检查数据里有没有削顶、跳变、缺失异常点对最小二乘的杠杆效应极大。5.3 现象多输入实验里各路扰动互相干扰辨识结果对不上单输入实验原因移位序列的移位量 J 取得太小不满足“输出对各输入的过渡过程小于 J”正交性假设破了。解决按 (2ⁿ-1-1)·Ts1.5·∑(各输入过渡过程时间) 重新核算 n 和 TsJ 至少取最大过渡过程对应的采样点数宁大勿小。5.4 现象模型检验残差不是白的但加大阶次后 J 值还在降原因把“J 值下降”当成“模型变好”忽略了教案说的阶次越高 J 越小会导致过估计残差非白可能是噪声成形环节没建模不是阶次不够。解决用残差白性检验自相关是否落在置信带内而不是只看 J 值跨采样步长比较时用 J·T0²/N 归一化必要时给噪声加成形滤波器结构。5.5 现象报告里模型转换后应用效果和辨识阶段对不上原因参考模型转应用模型时离散化方法或采样率变了但参数没跟着重算或者标准化反演漏了某一路系数。解决模型转换后重新在验证数据上跑一遍预测对比转换前后的输出误差反演公式逐项核对b、c 乘 σy/σua 不动常数项按 (1a)·ȳ-(bc)·ū 回填。6. 用一份数据文件把整条链路验一遍我的收尾习惯教案的考核方式是给一个数据文件、编程辨识、写报告这其实是最好的自检方式。我一般会拿一份已知结构的仿真数据先跑通全链路再换真实数据。具体做法是用二阶差分方程加有色噪声生成 2000 点数据按第 3 章算 M 序列参数、按第 4 章做预处理、用最小二乘估参数、反演回原尺度最后对比真值。import numpy as np # 生成已知结构的仿真数据y(k)a*y(k-1)b*u(k-1)c*u(k-2)e(k) np.random.seed(0) N 2000 a_true, b_true, c_true -0.7, 0.5, 0.3 u np.random.randn(N) e np.random.randn(N) * 0.05 y np.zeros(N) for k in range(2, N): y[k] -a_true*y[k-1] b_true*u[k-1] c_true*u[k-2] e[k] # 构造回归矩阵并最小二乘 Phi np.column_stack([-y[1:-1], u[1:-1], u[:-2]]) Y y[2:] theta, *_ np.linalg.lstsq(Phi, Y, rcondNone) print(估计参数:, theta, 真值:, [a_true, b_true, c_true]) # 残差白性粗检滞后1~5阶自相关应落在 ±2/sqrt(N) 内 res Y - Phi theta for lag in range(1, 6): r np.corrcoef(res[:-lag], res[lag:])[0, 1] print(flag{lag} 自相关{r:.3f}, 阈值{2/np.sqrt(len(res)):.3f})这段代码的价值不在估参数本身而在给你一个“真值已知”的对照系参数估偏了能立刻定位是预处理、回归矩阵构造还是噪声水平的问题残差自相关逐阶打印超过 2/√N 就说明还有结构没吸收。参数说明噪声标准差 0.05 是故意设小先验证算法正确性再逐步加大到噪信比接近真实工况回归矩阵第一列取 -y[k-1] 是为了让 θ 直接对应 a 的正负号和教案公式一致。从那以后我每次交辨识报告前都强制先用仿真数据把“参数估计—反演—残差检验”走一遍确认链路无误再碰实测数据这个习惯帮我省掉了至少三次整篇返工。希望这份拆解能帮到你把教案里的公式真正变成能跑出结果的代码。本文还有配套的精品资源点击获取