简介面向系统辨识学习者与工程实践者的一套极大似然法专题资料围绕MLE与递归极大似然RML的参数估计展开可用于理解线性时不变系统传递函数或状态空间模型的建模流程也适合作为控制、信号处理课程的实验参考。压缩包共19个文件、约865KB以6个.m算法脚本和9张.jpg结果图为主另含2个.mat数据文件、1份.pptx理论讲义及1个.tmw工程文件分别对应算法实现、误差与输出曲线验证、数据复现和原理讲解。已有189人学习下载。内容覆盖RML递推步骤、参数估计误差对比、白噪声与有色噪声场景及PRBS激励设计读者可据此复现实验、检查估计精度并梳理从理论推导到代码实现的完整思路。1. 极大似然法到底在估什么从「哪种解释最像真的」说起你下班回家发现窗台上的花盆倒了土撒了一地。可能是猫蹿上去碰的可能是窗户没关严被风刮的也可能是你早上出门时胳膊肘蹭到的。你没亲眼看见但脑子里会本能地挑一个「最能解释现场」的原因。极大似然法干的就是这件事只不过它把直觉换成了数学面对一组观测数据在所有可能的参数取值里挑出那个让「这组数据出现的可能性最大」的参数。它不保证找到真相它给出的是在给定模型下最讲得通的那个解释。搞懂它你就明白为什么系统辨识里那么多算法都绕不开它也明白它的天花板卡在哪。2. 为什么需要极大似然法从最小二乘的盲区说起2.1 最小二乘解决了什么又漏掉了什么最小二乘大概是所有人接触参数估计时学到的第一个方法假设模型是 y φθ e那就让残差平方和 Σe² 最小解出来一个闭式公式一行代码就能算完。它好用、稳定、算得快工程上大量场合确实拿它就够了。但它的好用是有代价的——它默认残差是独立同分布、均值零、方差恒定的通俗地说它假设「误差是一堆互不相关、大小差不多的随机抖动」。在这个假设成立时最小二乘非常漂亮。问题在于真实系统很少这么乖巧。传感器采样之间往往有相关性上一拍的扰动会影响下一拍误差不是白噪声而是「有色噪声」。举个例子某实验室的温控台在加热时会有热惯性测量误差和上一时刻的误差正相关这时最小二乘仍然能给出一个解但这个解是有偏的——它系统性地偏向某个方向你再怎么加数据也纠不回来。这就是最小二乘的盲区它只用到了「误差有多小」没用上「误差长什么样」。极大似然法补的正是这一块它把误差的分布形状也写进模型里一起考虑。2.2 似然函数把「像不像」变成一个可以比较的数要理解极大似然先得接受一个有点反直觉的设定数据是已经发生了的参数才是未知的。所以似然函数问的不是「给定参数数据出现的概率是多少」而是反过来——「把数据固定住让参数变化哪种参数设定下这组数据最像是被生成出来的」。数学上这两件事的表达式长得一样但解读完全不同。前者是概率参数固定、数据随机后者是似然数据固定、参数被视为变量。打个比方你手里有一把尺子量了同一个零件五次得到五个略有差异的读数。你不知道这把尺子准不准、零件的真实长度是多少。似然函数做的就是假定真实长度是 L、尺子的误差标准差是 σ那么「量出这五个读数」这件事发生的可能性有多大把 L 和 σ 的所有组合都代进去算一遍画出一条曲线曲线最高点对应的那组 (L, σ)就是极大似然估计。它把一个模糊的「哪个更合理」变成了一个可以排序、可以求极值的数值问题。2.3 噪声分布决定了估计器的天花板这里有个容易被忽略的关键点似然函数的具体形状完全由你假设的噪声分布决定。你假设噪声是高斯的似然函数里就会出现 exp(-e²/2σ²) 这一项取完对数再求导最后推出来的估计式和最小二乘几乎一模一样——这也是为什么很多教材说「高斯噪声下最小二乘就是极大似然」。但如果你假设噪声服从拉普拉斯分布推出来的就不是平方和最小而是绝对值之和最小对离群点更宽容假设噪声服从 t 分布鲁棒性会更强。换句话说极大似然法不是一个单一算法它是一套框架你先声明「我认为误差是这种形状的」框架再帮你推出对应的估计公式。这也是它在系统辨识里地位特殊的原因——真实系统的噪声往往不是白的你想把噪声的统计特性揉进估计过程极大似然是少数能名正言顺做到这件事的路子。代价是你的假设错了它给你的信心也会跟着错而且错得比最小二乘更隐蔽。3. 极大似然法的核心零件似然、对数似然与求解3.1 从联合概率到连乘似然函数怎么搭起来假设你采集了 N 个观测点 y₁ 到 y_N并且认为它们彼此独立、都服从同一个由参数 θ 决定的分布。那么「同时观测到这 N 个值」的联合概率就是每个点概率的乘积L(θ) p(y₁|θ) × p(y₂|θ) × … × p(y_N|θ)。这个连乘就是似然函数最朴素的形态。它有个很直观的性质任何一个点的概率太小整个乘积就会被拉得很低所以极大似然天然地要求模型对每一个样本都不能太离谱。但连乘在工程上很不好用。N 一大每个概率都小于 1几十个数乘下来直接掉到 10 的负几十次方计算机用浮点数表示会下溢成零。更麻烦的是求导一长串乘积求导要反复用乘积法则项数爆炸。所以实际计算从来不用原版似然函数而是走下一步。3.2 为什么要取对数把连乘变成连加对似然函数取自然对数连乘立刻变成连加ℓ(θ) Σ ln p(y_i|θ)。对数函数是单调递增的所以让 ℓ(θ) 最大的 θ 和让 L(θ) 最大的 θ 是同一个argmax 不变。这一步看着像小技巧实际是极大似然能落地的前提——连加求导简单数值上也不会下溢而且每一项的贡献可以单独看哪个样本在「拖后腿」一目了然。用高斯噪声举例单个样本的对数似然是 -0.5·ln(2πσ²) - e_i²/(2σ²)其中 e_i 是残差。把所有样本加起来σ² 固定时最大化对数似然等价于最小化 Σe_i²——最小二乘就这么被「推导」出来了。这个推导过程值得亲手写一遍它让你看清最小二乘不是凭空冒出来的规则而是高斯假设下的一个特例。3.3 求解求导、置零与迭代有了对数似然理论上求极值就是求导置零。但只有极少数情况能解出闭式表达式比如上面那个高斯例子。一旦模型里参数和噪声耦合在一起比如系统的传递函数参数和噪声模型的参数互相纠缠导数方程就是非线性的解不出来。这时候只能上迭代给一个初始猜测算梯度沿着让似然增大的方向挪一小步再算再挪直到挪不动为止。下面的代码示意了一个最简单的一阶系统在高斯噪声假设下的对数似然计算和网格搜索目的是让你看清「似然是一个关于参数的曲面估计就是找山顶」这件事import numpy as np # 模拟一阶系统: y(k) a * y(k-1) b * u(k-1) e(k) # 真实参数 a0.8, b0.5, 噪声 e ~ N(0, 0.1^2) rng np.random.default_rng(7) N 300 a_true, b_true, sigma_true 0.8, 0.5, 0.1 u rng.normal(0, 1, N) # 输入信号 e rng.normal(0, sigma_true, N) # 高斯白噪声 y np.zeros(N) for k in range(1, N): y[k] a_true * y[k-1] b_true * u[k-1] e[k] def log_likelihood(a, b, sigma, y, u): 计算给定参数下的对数似然残差按高斯分布处理 res [] for k in range(1, len(y)): y_hat a * y[k-1] b * u[k-1] # 一步预测 res.append(y[k] - y_hat) # 预测误差 res np.array(res) # 高斯对数似然: -N/2*ln(2πσ²) - Σe²/(2σ²) return -0.5 * len(res) * np.log(2 * np.pi * sigma**2) \ - np.sum(res**2) / (2 * sigma**2) # 在参数网格上搜索找似然最大的那一组 best (-np.inf, None) for a in np.linspace(0.5, 1.1, 61): for b in np.linspace(0.2, 0.8, 61): ll log_likelihood(a, b, sigma_true, y, u) if ll best[0]: best (ll, (a, b)) print(最大对数似然:, best[0]) print(估计参数 a, b:, best[1])这段代码里几个参数值得说清楚a和b是待估的系统参数sigma是噪声标准差实践中它往往也要一起估这里为了简化先固定住y_hat是一步预测值y[k] - y_hat就是预测误差它是似然函数的原料。代码用网格搜索而不是梯度法是因为网格搜索能直接把似然曲面画出来让你看到山顶在哪、周围是不是平坦的。真实场景里参数维度一高网格就爆炸了必须换成梯度上升、牛顿法或者后面要讲的递推形式。提示网格搜索的价值不在于实用而在于建立直觉。第一次跑的时候不妨把似然值打印成矩阵你会发现「山顶」附近往往非常平缓参数稍微动一动似然几乎不变——这正是后面要讲的「可辨识性」问题的源头。4. 从离线到在线RML 递推极大似然法在系统辨识里做什么4.1 系统辨识要解决的是一个「黑匣子」问题系统辨识的典型场景是这样的你手头有一个设备输入输出都能测但内部结构不清楚或者清楚但参数不知道。你想建立一个数学模型让它在这个设备上算出来的输出和真实测到的输出尽量对得上。这就是辨识。它和纯统计估计的区别在于这里的模型往往带有动态结构——这一拍的输出不只取决于这一拍的输入还取决于前几拍的输入和输出误差也会一拍拍地传下去。这个「动态」让事情复杂了不少。静态回归里样本是独立的动态系统里样本前后咬合噪声一旦有色就会顺着反馈路径影响后面的所有预测。极大似然法在这种场景下的优势就体现出来了它可以把噪声模型和系统模型放在一起作为一个整体去最大化似然。代价是整个似然函数变得又长又非线性计算量直线上升。4.2 RML 的核心思路一边估参数一边估噪声RML也就是递推极大似然法解决的是「数据一拍拍进来、没法等全部采完再算」的场合。它的思路可以概括成一句话用一个增广模型把系统参数和噪声参数捆在一起每来一个新样本就更新一次估计同时更新一个近似的信息矩阵来决定每一步该迈多大。通俗地说它像一个不断修正笔记的记账员。每来一笔新账他不会把整本账重算一遍而是根据这笔账和自己的预测差多少按一个增益去微调之前的结论。差得多就多改一点差得少就少改一点随着账目越来越多增益逐渐变小估计趋于稳定。这个增益的计算本质上是对似然函数在当前点做二阶近似用近似的曲率来决定步长这也是 RML 比普通递推最小二乘更「聪明」的地方——它知道自己在参数空间的哪个位置该走多快。RML 通常包含三块一是预测误差的计算用当前参数对下一拍输出做预测二是参数更新用预测误差乘以增益去修正参数向量三是噪声模型参数的同步更新把残差里还残留的相关性吸收掉。三块循环执行参数和噪声模型互相喂饭逐步收敛。4.3 一个最小可运行的示意下面这段代码把 RML 的骨架搭出来用一阶系统加一阶有色噪声做例子。重点是看结构不是看它能不能收敛到多准——真实工程里还要处理数值稳定性、协方差矩阵的正定性、遗忘因子等等那些属于踩坑经验代码里写不下import numpy as np # 增广模型: y(k) a*y(k-1) b*u(k-1) e(k) # e(k) c*e(k-1) v(k), v 为白噪声 # 待估参数向量 theta [a, b, c] rng np.random.default_rng(11) N 500 a_t, b_t, c_t 0.75, 0.6, 0.4 u rng.normal(0, 1, N) v rng.normal(0, 0.2, N) e np.zeros(N) y np.zeros(N) for k in range(1, N): e[k] c_t * e[k-1] v[k] y[k] a_t * y[k-1] b_t * u[k-1] e[k] theta np.zeros(3) # 初始估计 P np.eye(3) * 100 # 初始信息矩阵大一点表示「一开始什么都不信」 for k in range(2, N): # 构造回归向量包含上一拍输出、上一拍输入、上一拍残差 phi np.array([y[k-1], u[k-1], theta[2] * (y[k-1] - theta[0]*y[k-2] - theta[1]*u[k-2])]) y_hat theta phi # 一步预测 eps y[k] - y_hat # 预测误差 # 增益类似递推最小二乘但回归向量随参数更新而变化 K P phi / (1 phi P phi) theta theta K * eps # 参数修正 P (np.eye(3) - np.outer(K, phi)) P # 信息矩阵更新 print(估计参数:, np.round(theta, 3)) print(真实参数:, [a_t, b_t, c_t])代码里有几个点必须解释。phi是回归向量第三个分量用当前估计的噪声参数theta[2]去重构过去的残差这一项是 RML 区别于普通递推最小二乘的关键——普通方法根本没有这一项。P是信息矩阵的逆的近似它的更新公式决定了算法对旧数据的遗忘速度P初始值给大表示不信任初始猜测。K的分母里那个1 phi P phi是归一化项防止新息过大时一步跳飞。跑这段代码你会发现参数估计在头几十拍波动很大之后逐渐稳定但未必精确落在真实值上。这不是代码写错了而是 RML 本身就是渐近有效的——它保证的是样本足够多时接近理论下界不保证小样本下每次都准。实际使用时前几百拍的数据通常直接丢掉等算法「热」起来再开始记录。注意代码里theta[2]出现在phi的构造中意味着回归向量依赖当前参数估计整个更新是耦合的。如果噪声很强或者初值离真值太远这个耦合可能导致发散。工程上常见的做法是先跑一段递推最小二乘拿到粗略初值再切到 RML。5. 常见误区与易混概念极大似然法不是「万能拟合」5.1 易混点一极大似然不等于最大概率很多人第一次听到「似然」就自动翻译成「概率」觉得极大似然就是找概率最大的那个参数。这个翻译在直觉上没错但在数学上会把你带沟里。概率是参数固定、数据随机时某个数据结果出现的可能性所有可能结果的概率加起来等于 1似然是数据固定、参数变化时不同的参数对「这组已发生的数据」的解释力。似然函数对参数积分通常不等于 1它根本不是一个概率分布。怎么区分记住一句话概率描述「未来会看到什么」似然描述「已经看到的东西更像是由谁产生的」。前者可以拿来预测后者只能拿来比较。当你看到某篇材料说「这个参数的概率最大」时它要么是在贝叶斯框架下加了先验要么就是表述不严谨把似然和概率混着用了。5.2 易混点二极大似然不等于最小二乘因为高斯噪声下两者等价很多人就把它们当成一回事这是个隐蔽的坑。等价只在「噪声独立同分布且服从高斯」这个前提下成立。噪声一旦有色、异方差或者有离群点两条路就分道扬镳了最小二乘仍然在最小化平方和极大似然则会根据你声明的噪声分布调整权重离群点的影响被自动压下去或者放大。区分方法很简单看它有没有显式地写出噪声的概率分布。如果只写了残差平方和那就是最小二乘如果写了似然函数或者对数似然并且分布形式可调那才是极大似然。在系统辨识里这个区别直接决定了估计器能不能处理有色噪声——不能处理的遇到相关误差就会有偏。5.3 易混点三似然大不等于模型好似然值越高只说明这组参数让当前数据看起来更合理不说明模型本身抓住了系统的本质。你完全可以把模型搞得非常复杂参数多到能把每个样本的残差都压到接近零似然自然高得离谱但一到新数据上就原形毕露。这就是过拟合在系统辨识里尤其常见因为动态模型加几个极点零点拟合能力就能翻着番往上涨。区分的关键是看「样本外表现」和「参数可辨识性」。似然高但参数估计的置信区间宽得能跑马多半是模型结构冗余了。实操中常用的做法是看不同模型阶次下的准则值或者留一段数据不参与估计专门用来验证。似然是训练集上的成绩不是模型能力的判决书。5.4 易混点四RML 不等于离线极大似然的简单递推有人以为 RML 就是把离线极大似然的公式每一步重新算一遍只是省了存储。不是这样。离线极大似然是在整批数据上一次性优化目标函数是全局的RML 每一步只用当前样本的信息做局部修正它追求的是跟踪而不是一次性最优。这带来两个后果一是 RML 对初值和数据顺序更敏感二是它天然适合时变系统因为旧数据的影响会被逐渐遗忘。区分看目标如果你的系统参数恒定、数据能一次性采齐直接做离线极大似然更稳如果参数可能缓慢漂移、或者必须在线出结果才轮到 RML 上场。把两者当成同一个东西就容易在需要跟踪的时候误用了批处理或者在数据量很小的时候硬上递推结果两头不讨好。6. 什么时候该用极大似然法什么时候别碰一份边界判断6.1 适合用的三种场景第一种噪声明显有色、相关性不能忽略的系统。比如机械臂关节的摩擦扰动、化工过程的慢漂移误差前后咬合最小二乘估计出来的参数会有系统性偏差这时候把噪声模型一起写进似然函数是正解。第二种需要对估计精度有个理论预期。极大似然估计在正则条件下是渐近有效的意思是样本足够多时它的方差能逼近克拉美-罗下界。你要做精度评估、要给参数置信区间这个性质很有用。第三种模型结构可信、样本量充足。似然函数是个曲面样本少的时候曲面崎岖容易卡在局部极值。样本上千、模型阶次合理时曲面才会变得光滑迭代才靠谱。6.2 不建议用的三种情况反过来下面这些场合要克制。噪声分布完全未知、而且明显偏离你假设的形状时似然函数给的是虚假的信心估出来的参数可能比最小二乘还偏。样本量只有几十个点时迭代的方差比参数本身还大不如老老实实用简单的矩估计。实时性要求极高、算力又紧张时RML 每拍要更新信息矩阵、重构回归向量开销不小降阶或者换简单递推更实际。判断维度倾向极大似然倾向简单方法噪声特性有色、相关、需建模近似白噪声样本量数百以上几十以内模型结构有物理依据、阶次确定结构不明、靠试计算资源离线或算力充裕在线、资源受限精度需求需要置信区间只要趋势对我用这套判断吃过亏。早些年做一个模拟项目X数据只有不到一百个点噪声来源也没搞清楚我上来就套了带噪声模型的递推估计跑出来参数漂亮得很似然值也高结果换一批数据重跑参数飘得完全对不上。后来才明白那点样本量根本撑不起噪声模型那几个额外参数似然高只是过拟合的假象。从那以后我每次上手新数据都强制先画残差的自相关图、先跑一遍最小二乘看残差白不白白就直接收工不白再考虑要不要上极大似然这条路。这套习惯救过我不少次希望帮到你。本文还有配套的精品资源点击获取