企业官网建设流程全解析

简介：《极大似然法》压缩包面向系统辨识初学者与工程实践者，围绕极大似然估计及递归极大似然（RML）在动态系统参数估计中的应用展开，帮助读者理解如何利用输入输出数据估计线性时不变系统的传递函数或状态空间模型，并兼顾非高斯噪声与在线更新场景。资源共19个文件，以m脚本、mat数据、jpg结果图和pptx讲义为主，压缩包约865KB；脚本可复现算法流程，数据与图像便于核对参数估计误差、白噪声与有色噪声下的输出表现，讲义补充理论脉络。已有189人学习，适合希望把概率统计基础、极大似然优化与递归在线学习串联起来的读者，通过激励信号设计、参数收敛观察和误差对比，逐步掌握RML的建模思路与实验验证方法。

1. 极大似然法到底在“最大”什么：从一个摸球游戏说起

桌上有个不透明的罐子，里面黑白两色球的数量比例未知。你伸手摸了十次，记下七黑三白。这时有人问：罐子里黑球占多少？多数人会说“大概七成”。这个朴素的回答背后，就是极大似然法的全部直觉——在所有可能的比例里，挑出让“七黑三白”这个已经发生的结果出现概率最大的那一个。放到系统辨识里，罐子换成待估参数的数学模型，摸球结果换成采集到的输入输出数据，“比例”换成模型里的系数。极大似然法要做的就是：给定这批数据，反推哪组参数最像“生成这批数据的真凶”。搞懂它，你就能看清参数估计里一大票方法的来路与边界。

2. 系统辨识为什么需要极大似然：当噪声不再听话

2.1 系统辨识在做什么：给黑匣子配一副数学眼镜

一个动态系统摆在面前，你只能往里面灌输入信号，再从另一头读输出信号，内部结构看不见摸不着，这就是典型的黑匣子。系统辨识干的事，说白了是给这个黑匣子配一副度数合适的数学眼镜：先假定它属于某一类差分方程描述的模型家族，再用实测数据把方程里的系数估出来，最后拿另一批数据检验这副眼镜看得准不准。整个流程有四个绕不开的要素——输入输出数据、模型结构、参数估计准则、模型验证，缺一个都会让结论站不住脚。

模型结构这一步就已经埋下了后面所有分歧。常见的四类结构可以写成下面这样，其中 q⁻¹ 是延迟算子，表示“往前取一个采样点”：

通俗地说，ARX 认为干扰是“干净的白噪声”，ARMAX 认为干扰经过了一层染色滤镜，OE 压根不管干扰怎么来的、只关心输入到输出的传递关系，BJ 则把“过程”和“干扰”当成两条独立的通道分别建模。模型越贴近真实干扰的产生方式，参数估得越准，代价是计算越麻烦。极大似然法正是那种“愿意为噪声多写几行方程”的做法。

2.2 最小二乘的软肋：噪声一进回归量就“带偏”

最小二乘的思路极其讨喜：让残差平方和最小。写成矩阵就是 θ̂=(ΦᵀΦ)⁻¹ΦᵀY，其中 Φ 是由历史输出和历史输入拼成的回归矩阵，Y 是当前输出。它有个漂亮的闭式解，一行代码就能算完，因此在工程里被用得像家常便饭。但它无偏的前提相当苛刻：回归矩阵里的每一列，都必须和方程误差不相关。白噪声满足这个条件，所以 ARX 结构下最小二乘几乎没有争议。

麻烦出在真实系统的噪声往往是有色的。假设真实过程是 ARMAX，干扰 e(k) 经过 C(q) 染色后才进入输出，此时输出 y(k-1) 里就残留着 e(k-1) 的成分，而回归矩阵恰好包含 y(k-1)。回归量和误差相关，最小二乘的解就会系统性地偏离真值，偏差量大致是 (ΦᵀΦ)⁻¹Φᵀv，v 是那部分没被建模的染色噪声。这里最容易踩的坑是：加大数据量并不能救它。偏差属于模型设定错误带来的系统性偏移，采样点从一千加到十万，估出来的错误值只会越来越稳定地停在错误的位置上，方差小了， bias 却纹丝不动。这种“数据越多错得越自信”的现象，是很多人第一次做辨识时最想不通的地方。

2.3 极大似然的解法：把噪声模型一起请进来

极大似然换了个问法：不是“让残差平方和最小”，而是“让这批观测数据出现的概率最大”。把 N 个采样点拼成一个向量 Y，在参数 θ 给定的条件下，Y 有一个联合概率密度 p(Y|θ)，这个函数看成 θ 的函数就是似然函数 L(θ)。求解就是找一个 θ，让 L(θ) 取最大值。由于连乘容易数值下溢，实际都取对数，变成求负对数似然的最小值。

关键的转折在这里：一旦你把噪声模型 C(q) 也当作待估参数的一部分，似然函数就自动“知道”噪声是有色的。对 ARMAX 来说，一步预测误差不再等于方程残差，而是经过白化滤波后的 ε(k)=(A(q)/C(q))y(k)−(B(q)/C(q))u(k)。负对数似然写成：

−log L = N·log σ + (1/2σ²)·Σε(k)² + 常数

优化这个目标，等于在“白化”之后再做一次加权最小二乘：先把有色噪声滤成白的，再让白化后的残差平方和最小。最小二乘在这里只是 C(q)=1、σ 固定时的一个特例。换句话说，极大似然不是最小二乘的对手，而是它的父集——高斯白噪声条件下两者给出同一个解，噪声一有色，父集才显出威力。

import numpy as np from scipy.optimize import minimize np.random.seed(0) N = 2000 u = np.random.randn(N) e = 0.3 * np.random.randn(N) y = np.zeros(N) # 真实系统 ARMAX(2,2,1)：a=[1,-1.5,0.7], b=[0,1.0,0.5], c=[1,0.8] for k in range(2, N): y[k] = 1.5*y[k-1] - 0.7*y[k-2] + 1.0*u[k-1] + 0.5*u[k-2] + e[k] + 0.8*e[k-1] # 方案一：ARX 结构的最小二乘，只估 a1,a2,b1,b2 Phi = np.column_stack([y[1:-1], y[:-2], u[1:-1], u[:-2]]) Y = y[2:] theta_ls = np.linalg.lstsq(Phi, Y, rcond=None)[0] print("LS 估计 a1,a2,b1,b2 =", np.round(theta_ls, 3)) # 真实值应为 [-1.5, 0.7, 1.0, 0.5]，观察 a1 是否被拉偏 # 方案二：极大似然，把 c1 与噪声强度一并估计 def predict_error(p, y, u): a1, a2, b1, b2, c1 = p n = len(y) v = np.zeros(n) # v(k) = A(q)y(k) - B(q)u(k) = C(q)e(k) eps = np.zeros(n) # 白化后的一步预测误差 for k in range(2, n): v[k] = y[k] + a1*y[k-1] + a2*y[k-2] - b1*u[k-1] - b2*u[k-2] eps[k] = v[k] - c1*eps[k-1] # 用 1/C(q) 把 v 白化 return eps[2:] def neg_log_like(p): eps = predict_error(p, y, u) sigma2 = np.mean(eps**2) + 1e-12 return len(eps)*0.5*np.log(sigma2) + 0.5*np.sum(eps**2)/sigma2 x0 = np.array([*theta_ls, 0.1]) # 用 LS 结果作初值，c1 给个小值 res = minimize(neg_log_like, x0, method="Nelder-Mead", options={"maxiter": 4000, "xatol": 1e-6}) print("ML 估计 a1,a2,b1,b2,c1 =", np.round(res.x, 3))

上面这段代码里，predict_error完成的是“先算方程残差 v，再用 1/C(q) 把它白化”这两步，neg_log_like把白化残差换算成负对数似然。初值取最小二乘结果是个实用习惯：似然面在有色噪声下往往比较崎岖，从 LS 解出发能少走很多弯路。参数maxiter设大一点，是因为单纯形法收敛偏慢；xatol控制的是参数变化小于多少就停，别设得太松，否则 c1 常常停在半路上。运行时你会看到 LS 估的 a1 明显偏离 −1.5，而 ML 的解基本贴着真值，这个对比就是这一章最直观的注脚。

2.4 一致性与有效性，以及为此付出的代价

极大似然之所以在统计估计里地位稳固，靠的是两条性质。一致性说的是采样点趋于无穷时，估计值依概率收敛到真值——只要模型结构包含真实系统、噪声分布假设没错、输入信号足够“有料”。有效性说的是在同样一批数据下，它的估计方差渐近达到克拉美-罗下界，也就是所有无偏估计里方差最小的那一档。这两条加起来，解释了为什么工程师愿意为它多写几行滤波器代码、多花几十倍算力。

代价也实打实。第一，你得假设噪声的概率分布，高斯假设好用但未必对；第二，负对数似然通常对参数非线性，没有闭式解，只能靠数值优化，初值不好就可能掉进局部极小；第三，阶次选错时它的表现未必比简单方法好，模型错了再精巧的准则也救不回来。所以工程上的常见做法是先跑一遍最小二乘看残差，残差自相关明显不白，再上极大似然。这套顺序背后是血泪经验：直接上重炮，往往连问题出在结构还是参数上都分不清。

3. 似然函数是怎么搭出来的：从概率密度到能求导的对数

3.1 联合概率拆成连乘：独立性假设帮了大忙

观测数据是一个向量，要写出它的联合概率密度，最省事的假设是各采样点条件独立。动态系统里严格独立不成立，但可以走条件分解这条路：把联合密度写成第一个点的密度乘上后续每个点在“已知历史”下的条件密度，一层层乘下去。对 ARMAX 这类模型，条件密度只依赖有限阶历史，于是似然函数变成 N 个一元密度的连乘。这一步是整套方法的地基——它把“一整批数据的概率”降解成“每个点预测误差的概率”，计算上才可行。

3.2 取对数：连乘变连加，数值也稳了

连乘在数值上极不友好，几百个小于 1 的数乘下来立刻下溢成零。取对数之后连乘变连加，单调性不变，最优解位置也不变。高斯情形展开来特别清爽：

log L = −(N/2)·log(2πσ²) − (1/2σ²)·Σε(k)²

对 σ² 求导令其为零，得到 σ̂²=Σε²/N，代回去就只剩“让白化残差平方和最小”这一项在起作用。这就解释了一个常被问到的问题：为什么极大似然和高斯假设绑得这么紧？因为在高斯下，似然的对数正好把平方和暴露出来，优化目标和最小二乘长得几乎一样，只多了个白化滤波器。噪声不是高斯时，似然里会出现绝对值、对数或其他非线性项，性质就复杂得多。

3.3 写成预测误差形式：似然其实是一个滤波器

把 ARMAX 的似然写成递推形式，会看到它天然对应一个滤波器结构。先由 A(q)y(k)−B(q)u(k) 得到中间量 v(k)，再让 v(k) 减去历史预测误差的加权和，得到白化误差 ε(k)。这个递推式的系数就是 C(q)，所以整个似然评估过程等价于“用 1/C(q) 对残差做一次白化滤波”。理解这一点，后面看 RML 的递推公式就不会觉得是天上掉下来的——它不过是在线地估计这个滤波器的参数，同时估计过程模型的参数。

3.4 参数化与数值细节：别让优化器替你做决定

有几处细节不处理，优化器给出的结果会似是而非。首项系数必须固定，A(q) 和 C(q) 的常数项都锁成 1，否则分子分母同乘一个常数解不唯一，优化会在一条平坦的脊上乱走。稳定性要有约束，估出来的 A(q) 根若跑到单位圆外，模型虽然拟合得好，但物理上讲不通，用它做预测会发散。初值优先用最小二乘或 ARX 的结果，比随机初始化靠谱得多。最后一定要看残差：把 ε(k) 画出来、算它的自相关，如果还有明显拖尾，说明阶次不够或结构选错，此时调优化器的容差是白费力气。

# 残差白化检查：优化完别急着收工 eps = predict_error(res.x, y, u) acf = np.correlate(eps, eps, mode="full")[len(eps)-1:] acf = acf / acf[0] # 归一化自相关 print("前 5 个自相关系数:", np.round(acf[:5], 3)) # 理想情况下除 0 阶外都应接近 0；若 1 阶仍大于 0.1，考虑升阶或换结构

这段检查只有三行，却比任何拟合优度指标都直白。自相关在 1 阶、2 阶仍有明显峰值，说明白化没做干净，噪声模型阶次偏低；全都很小但参数偏差依旧大，那问题多半在输入信号激励不足，跟估计算法无关。

4. RML 的完整含义：递推极大似然法如何边跑边估

4.1 从一批数据到每一个采样点

批处理极大似然要等数据全部采完，一次性优化，适合离线分析。可很多场景等不起：自适应控制要在每个控制周期内更新模型，在线监测要跟着工况漂移走，数据量大到内存装不下时也不可能一次算完。于是有了递推版本——每来一个新采样点，就用它修正一次参数，历史信息压缩在一个协方差矩阵里。递推最小二乘（RLS）是这条路上的老前辈，公式简单、收敛快；递推极大似然（RML）则是把 RLS 的目标函数从“残差平方和”换成“白化残差平方和”，代价是每步还要多算一个滤波。

4.2 三个动作：算预测误差、算灵敏度、更新参数

RML 的每一次迭代拆成三步。第一步，用当前参数算这一步的一步预测误差 ε(k)，它同时是“模型错得多离谱”的度量。第二步，算预测误差对参数的灵敏度 ψ(k)=−∂ε(k)/∂θ，这一项决定了参数该往哪个方向调、调多少；对 ARMAX 而言，ψ 由历史输入、历史输出和历史预测误差拼成，再经 1/C(q) 滤波。第三步，用类似 RLS 的增益更新参数：

θ(k) = θ(k−1) + K(k)·ε(k) K(k) = P(k−1)ψ(k) / (λ + ψ(k)ᵀP(k−1)ψ(k)) P(k) = (P(k−1) − K(k)ψ(k)ᵀP(k−1)) / λ

λ 是遗忘因子，取 1 表示无限记忆，取 0.95 到 0.99 表示更看重近期数据。P 矩阵的物理含义是“参数还有多少不确定性”，初值给大一点表示一无所知，随着数据进来它逐步收缩。这套结构和 RLS 几乎同构，区别只在 ψ 不是原始回归向量，而是白化后的灵敏度向量——这正是“极大似然”四个字落到递推形式上的全部差别。

# RML 递推骨架（ARMAX，参数顺序 a1,a2,b1,b2,c1），仅示意核心循环 theta = np.array([-1.0, 0.3, 0.8, 0.3, 0.1]) # 初值，可来自一段离线 LS P = 1000.0 * np.eye(5) # 大初值 = 一开始什么都不信 lam = 0.99 # 遗忘因子，越接近 1 记忆越长 hist_y, hist_u, hist_eps = [], [], [] for k in range(2, N): # 动作一：用当前参数算一步预测误差 v = y[k] + theta[0]*y[k-1] + theta[1]*y[k-2] \ - theta[2]*u[k-1] - theta[3]*u[k-2] eps = v - theta[4]*(hist_eps[-1] if hist_eps else 0.0) # 动作二：灵敏度向量 ψ = -∂ε/∂θ，由历史量经白化得到 psi = np.array([-y[k-1], -y[k-2], u[k-1], u[k-2], hist_eps[-1] if hist_eps else 0.0]) # 动作三：增益、参数、协方差三连更新 denom = lam + psi @ P @ psi K = (P @ psi) / denom theta = theta + K * eps P = (P - np.outer(K, psi @ P)) / lam hist_eps.append(eps) print("RML 最终参数:", np.round(theta, 3))

循环里denom那行是标量，K是 5 维列向量，np.outer负责把增益和灵敏度卷成矩阵更新。lam调小能让模型更快跟上参数突变，但也让 P 收缩变慢、估计抖动变大；调大则稳但迟钝，这是在线辨识里最典型的取舍。实际部署时通常还会给 P 加一个上限，防止它在长期无激励的工况下膨胀到数值溢出，这类小保护在教科书公式里看不到，却是跑通线上系统的必要条件。

4.3 和批处理 ML 的关系：同一座山的两种爬法

批处理极大似然每一步都用全部数据算一次梯度方向，步子稳、方向准，但一次只能走一步。RML 每来一个点就朝着“当前这个点的负对数似然”走一小步，方向噪声大，靠时间平均来抵消。数学上它属于随机近似那一类方法，收敛要靠步长衰减或者遗忘因子来保证。两者的最优解是同一个：模型结构一致、数据充分激励时，RML 收敛到的参数和批处理 ML 估出来的应当吻合。差别在于过程——批处理给你一个终值，RML 给你一条轨迹，轨迹本身在在线监测里就是有用信息。

4.4 参数怎么定：阶次、遗忘因子、协方差初值

表里最后一行的“激励”经常被忽略，却是在线辨识翻车的头号原因。输入信号如果长时间停留在某个常量附近，某些参数组合对输出的影响无法区分，P 矩阵会在那个方向上持续膨胀，估出来的数看着在变，其实全是噪声。稳妥做法是给输入叠加一点小幅度的探测信号，或者至少在监测到激励不足时冻结更新。

5. 易混概念盘点：极大似然、最小二乘、贝叶斯与 EM 的区别

5.1 易混点一：极大似然与最小二乘是对立的两种方法

初学者常把两者摆在二选一的位置上。之所以会混，是因为很多教材先把最小二乘讲透，再另起一章讲极大似然，看起来像两套平行体系。其实最小二乘是极大似然在高斯白噪声假设下的一个特例：噪声白且方差固定时，负对数似然里唯一和参数有关的项就是残差平方和，最大化似然等价于最小化平方和。区分的关键看噪声：噪声白、结构是 ARX，直接最小二乘，又快又对；噪声有色、结构是 ARMAX 或 BJ，必须上极大似然或其近似，否则偏差会一直挂着。判断该用哪个，看一眼残差自相关比翻书快。

5.2 易混点二：似然和概率是同一个东西

同一个函数 p(Y|θ)，固定 θ 看 Y 叫概率密度，固定 Y 看 θ 叫似然。会混是因为符号完全一样，公式也一模一样。区分的关键是“谁在变”：概率描述的是数据在给定参数下的分布，所有可能结果的概率加起来为一；似然描述的是在已观测数据下，不同参数的可信程度，它对 θ 积分并不等于 1，也不是一个概率分布。这个区别决定了似然函数不能直接当成“参数是某个值的概率”来读，那一步需要贝叶斯公式加上先验才能完成。

5.3 易混点三：极大似然和贝叶斯估计是竞争关系

两者确实常被拿来比较，但差别只有一个词：先验。极大似然只看数据，找让观测最可能出现的参数；贝叶斯把参数也当成随机变量，先给一个先验分布，再用数据更新成后验，取后验最大就是最大后验估计。之所以容易混，是因为数据量一大，先验的影响被冲淡，两者结果趋于一致。区分场景看样本量：数据充足、先验信息弱，用极大似然省事；样本稀少但你对参数范围有可靠经验，先验能防止估计跑到荒谬区域，这时贝叶斯更稳。

5.4 易混点四：EM 算法是另一种估计准则

EM 经常和极大似然并列出现，让人以为它是和最大似然平行的另一套标准。实际上 EM 是求解极大似然（或最大后验）的一种迭代算法，专门对付“有隐变量导致似然不好直接优化”的情形。之所以混，是因为很多教程把 EM 放在极大似然之后讲，标题并列，读起来像兄弟。区分看层次：极大似然是目标，EM 是工具；能直接对似然求导的场景根本用不上 EM，只有数据缺失、状态不可观测时才需要它那套“先猜隐变量、再更新参数”的交替迭代。

5.5 易混点五：递推极大似然就是递推最小二乘的别名

两个缩写只差一个字母，公式结构又极其相似，混起来毫不奇怪。区别藏在更新向量里：RLS 用的是原始回归向量，本质在最小化方程残差平方和；RML 用的是白化后的灵敏度向量，最小化的是预测误差平方和，参数里多了一组噪声模型系数。区分方法很直接：看估计器有没有在估 C(q)。如果代码里只有 a 和 b，那不管名字叫什么，都是最小二乘家族；只有当 c 系数也进入参数向量并参与递推，才算真正走到了极大似然这条线上。

6. 什么时候该用极大似然法：一个仿真例子给出的边界

拿一组 ARMAX 数据做过对比：真实 a1 是 −1.5，噪声经过 C(q)=1+0.8q⁻¹ 染色。用 ARX 结构跑最小二乘，两千个点估出来 a1 落在 −1.38 附近，加到两万个点还是 −1.38，偏差稳稳地待着；换成极大似然把 c1 一起估，a1 回到 −1.49，残差自相关也压到了 0.05 以下。这个对比划出了第一条边界：噪声明显有色、你又需要参数本身有物理意义时，极大似然值得那点算力。

第二条边界在数据质量上。同样这套代码，把输入换成一个几乎不变的常值，两种方法都开始给出乱七八糟的结果，ML 甚至因为优化面太平坦而报出收敛失败。激励不足时，任何估计准则都是巧妇难为无米之炊，这时该做的是改实验设计，而不是换算法。

也有明确不该用的时候。算力紧张的嵌入式场景，如果残差检查显示噪声接近白，实在没必要为了理论上的优雅去跑非线性优化；模型阶次还没定下来、正在做结构筛选的阶段，用最小二乘快速扫一遍阶次曲线效率高得多；数据里存在明显异常点又没做剔除时，极大似然对分布假设敏感，一个离群点就能把优化拖偏，这时候先清洗数据比换准则有用。

从那以后我每次拿到一组辨识数据，都强制自己先做三件事：画残差自相关、看输入频谱有没有覆盖关心的频段、用最小二乘跑一遍当基准线。三件事做完再决定要不要上极大似然，比一上来就调优化器省下太多时间。希望这套顺序也能帮到你，资源里的完整推导和可运行脚本，正好可以拿它对着自己的数据走一遍。

本文还有配套的精品资源，点击获取