极大似然法详解：从解析解到神经网络损失函数的完整落地

简介针对系统辨识中的参数估计问题提供极大似然法MLE与递归极大似然RML算法的完整实现资料适合控制理论、自动化、信号处理等领域的学习者与工程人员参考。MLE通过最大化观测数据出现的概率来估计模型参数RML则在数据流中实时递推更新无需重算全部历史数据适用于线性时不变系统的传递函数或状态空间模型辨识。包内共19个文件以MATLAB脚本为主涵盖RML测试、m序列生成、模型递推估计等核心代码并配有参数估计误差、有色噪声与白噪声场景下的输出误差对比图.jpg以及理论讲解PPT、仿真数据.mat和Simulink工程配置文件整体约865KB。已有187人学习下载适合需要在线辨识或处理非高斯噪声的实践者对照实现。结合代码、仿真结果与图示可直观掌握RML递推更新流程、噪声条件下误差分析以及算法收敛性验证便于二次开发、实验复现与教学演示。1. 极大似然法一个压缩包里最值得先吃透的估计方法早几年我处理一批传感器采集的温度数据想判断噪声到底服从什么分布、参数大概是多少。参考资料东一个西一个最后把推导过程、模拟代码和踩坑记录打成一个压缩包命名就叫「极大似然法.zip」。这个标题其实是个很形象的入口极大似然法不是单个知识点而是一整套从数据反推参数的方法论。它回答的核心问题是假设数据来自某个带未知参数的分布如何选出一组参数让「这批数据被观测到」的概率最大。适合做数据分析、算法落地和概率建模的人。下面我按自己重走一遍的顺序把它拆成从手推解析解、数值优化到接入神经网络的完整落地路径。2. 从零推导极大似然拿伯努利和高斯分布手工求一遍2.1 为什么选伯努利和高斯作为起点MLE 的完整流程包含四步写出似然函数、取对数、对参数求导、令导数为零解方程。选伯努利和高斯作为起点是因为它们一个代表离散分布、一个代表连续分布恰好覆盖了你日常遇到的大多数拟合场景而且都能推出闭式解。伯努利分布的参数 p 只有一个求导过程干净能让你把「对数似然为什么比似然更顺手」一次看懂。高斯分布有两个参数而且 μ 和 σ² 的估计之间存在耦合——先估 μ 再估 σ²这个先后关系在多元场景里还会反复出现。更重要的是这两组解析解不是只活在纸面上它们直接对应np.mean和np.var你在任何一本统计教材里看到的「用样本均值估计期望」本质上就是一次极大似然估计。我一般会把这一步当作「体检」如果连这两个分布的推导都顺畅走完说明对似然函数的结构已经建立了手感如果卡住了后面做数值优化和自定义损失时会更容易翻车。2.2 手写对数似然与解析解的完整流程拿伯努利来说观测到一组独立样本 x₁, x₂, …, xₙ每个样本取值 0 或 1。似然函数是每个样本概率的连乘L(p) Π p^{xᵢ} (1-p)^{1-xᵢ}连乘在 n 稍大一点就会下溢成 0所以取对数把它变成连加ℓ(p) (Σxᵢ) ln p (n - Σxᵢ) ln(1-p)对 p 求导并令其为零解出来就是样本均值p̂ (1/n) Σ xᵢ高斯分布同理。设 xᵢ ~ N(μ, σ²)对数似然写作ℓ(μ, σ²) -n/2 ln(2π) - n/2 ln σ² - (1/(2σ²)) Σ (xᵢ - μ)²先对 μ 求导得到 μ̂ x̄再把这个结果代回去对 σ² 求导得到 σ̂² (1/n) Σ (xᵢ - x̄)²。注意这里除的是 n 而不是 n-1所以极大似然给出的方差估计是有偏的样本量小的时候偏差会更明显。教科书里那个 n-1 的无偏版本是在 MLE 基础上做了一次校正。下面用代码验证一下解析解和数值结果是不是对得上顺便看看样本量对估计偏差的影响import numpy as np from scipy import stats # 固定一个真实分布生成模拟数据 np.random.seed(42) mu_true, sigma_true 2.0, 1.5 n 50 data np.random.normal(mu_true, sigma_true, sizen) # 解析解直接用样本均值和样本方差 mu_mle np.mean(data) sigma2_mle np.var(data) # 注意这是除以 n 的 MLE 版本 sigma2_unbiased np.var(data, ddof1) # 除以 n-1 的无偏版本 # 数值对照scipy 的 fit 本质上也是在做极大似然 mu_fit, sigma_fit stats.norm.fit(data) print(fmu_mle {mu_mle:.4f}, mu_fit {mu_fit:.4f}) print(fsigma2_mle {sigma2_mle:.4f}, sigma2_unbiased {sigma2_unbiased:.4f})这段代码里np.var(data)默认除 n正好对应 MLE 的解析解ddof1则对应无偏版本。stats.norm.fit内部用的就是极大似然流程所以mu_fit会和mu_mle完全一致。当你拿到一个不熟悉的分布时用fit的结果和手推解析解对照是验证自己推导有没有写错的最快办法。这里有一个容易被忽略的细节当 n 小到 10 以下σ² 的 MLE 会明显偏离真实值这不是数值 bug而是估计量本身的偏差。所以真实项目中如果样本量小我通常会在 MLE 基础上转用无偏校正或者直接走贝叶斯加先验——MLE 不是万能答案它只是最自然的那一份参考答案。2.3 解析解拿不到时数值优化该怎么接上大多数实际问题没有解析解。比如拟合一个混合高斯、一个威布尔分布或者自定义的物理模型似然函数对参数的导数写不出来或者解方程解不出闭式形式。这时候常规做法是转成数值优化最小化负对数似然NLL因为最小化和最大化是等价的而大多数优化器只做最小化。具体来说就是把参数塞进一个向量 θ定义目标函数为NLL(θ) -Σ ln p(xᵢ | θ)然后交给梯度下降或牛顿法去迭代。这里有一个隐藏的工程决策参数必须做约束处理。比如 σ² 必须大于 0概率 p 必须落在 (0,1) 之间。直接让优化器去搜原始参数很容易在迭代过程中把 σ² 推成负数导致目标函数出现 NaN。我惯用的手法是做变量替换用 log σ 代替 σ用 logit(p) 代替 p。这样参数空间从「带约束的区间」映射到「整个实数轴」优化器可以放心大胆地走回代时再做指数或 sigmoid 变换。这一步虽然简单但能省掉后面 90% 的「优化器报错」类问题。下一章我会直接拿 SciPy 的优化器把这套流程跑通并讲清楚几个必调参数的含义。3. 用现成工具跑通极大似然SciPy 实现与三个必调参数3.1 用 scipy.optimize.minimize 拟合自定义分布假设你已经确认数据来自一个带位置和尺度参数的分布但分布形式是自定义的没法直接用现成fit方法。这时候最通用的做法是写一个负对数似然函数交给minimize去优化。下面这段代码是完整可跑的骨架import numpy as np from scipy.optimize import minimize # 自定义分布把正态分布的均值换成线性模型预测值 # 这里用一个最简版本单变量正态但走通用优化接口 def nll_gauss(theta, x): 负对数似然函数。 theta [mu, log_sigma]用 log_sigma 代替 sigma 保证优化过程中方差恒为正不需要额外约束。 mu, log_sigma theta sigma np.exp(log_sigma) n len(x) # 正态分布 NLL常数项可省略但保留更利于调试 return 0.5 * n * np.log(2 * np.pi * sigma**2) np.sum((x - mu)**2) / (2 * sigma**2) np.random.seed(7) data np.random.normal(loc3.0, scale2.0, size200) # 初始值样本均值 log(样本标准差)通常能更快收敛 theta0 np.array([np.mean(data), np.log(np.std(data))]) # methodNelder-Mead 不依赖梯度适合先跑通 # 正式用推荐 L-BFGS-B见后文参数说明 res minimize(nll_gauss, theta0, args(data,), methodNelder-Mead, options{maxiter: 5000, xatol: 1e-8}) mu_hat, log_sigma_hat res.x print(fmu {mu_hat:.4f}, sigma {np.exp(log_sigma_hat):.4f}) print(fconverged {res.success}, nfev {res.nfev})代码里最关键的一行是log_sigma np.exp(log_sigma)。优化器看到的变量是 log σ范围是整个实数轴怎么迭代都不会产生非法方差。这也解释了为什么在 SciPy 里做 MLE 时我几乎从不给bounds传(0, None)这类约束——变量替换之后根本不需要。options里的两个参数值得注意maxiter控制最大迭代次数默认值在某些复杂模型上偏小容易提前退出还报「成功」xatol是变量变化阈值越小说明收敛判据越严。真正常见的坑是模型很复杂迭代次数不够优化器提前停successTrue但结果明显不对。所以每次跑完我都会打印nfev去看函数评估次数如果刚好等于maxiter基本可以断定是撞到上限了。3.2 参数初始化与边界约束怎么影响结果minimize的默认行为是「从你给的初始点出发沿着下降方向走」所以初始值选得不好结果会直接跑偏。尤其是似然面存在多个局部极小值的时候不同初始点会收敛到完全不同的参数。我一般用两个策略。第一个是矩估计初始化用样本均值初始化位置参数用样本标准差初始化尺度参数。这个初始点通常已经接近真实值优化器只需要少量迭代就能收敛。第二个是多起点策略随机撒几组初始值比较各自的最终 NLL取最小者。下面这张表是我常用的配置参考参数常用取值作用methodL-BFGS-B拟牛顿法适合参数不多、目标函数光滑的场景methodNelder-Mead不依赖梯度目标函数有小噪声时更稳options.maxiter5000~20000防止复杂模型迭代不足提前终止options.xatol1e-8收敛阈值越小越严格bounds尽量不用优先用变量替换代替边界约束这里提一个血泪经验很多人刚接触minimize时喜欢把bounds写成[(None, None), (0, None)]强制 σ 为正。这个思路本身没问题但 σ 一旦接近 0目标函数会变成巨大正数优化器很容易在边界附近来回试探收敛非常慢。换成 log σ 之后同样的模型用 L-BFGS-B 通常几十次迭代就收住了。3.3 收敛判据与数值梯度检查L-BFGS-B 这类梯度法依赖目标函数的梯度但你自己写的 NLL 往往没有解析梯度。SciPy 默认用有限差分去近似这一步本身就有误差。如果你怀疑结果不对先做一次数值梯度检查用中心差分公式计算梯度和优化器内部的梯度做对比。def numerical_gradient(f, theta, args, eps1e-6): 中心差分计算梯度用于验证目标函数实现是否正确 grad np.zeros_like(theta) for i in range(len(theta)): theta_plus theta.copy() theta_minus theta.copy() theta_plus[i] eps theta_minus[i] - eps grad[i] (f(theta_plus, *args) - f(theta_minus, *args)) / (2 * eps) return grad # 在初始点附近比较数值梯度与理论梯度如果有的话 theta_check np.array([np.mean(data), np.log(np.std(data))]) g_num numerical_gradient(nll_gauss, theta_check, args(data,)) print(fnumerical gradient {g_num}) # 也可以用 res.jac 看最终点的梯度L-BFGS-B 会返回该值 print(ffinal gradient {res.jac})如果最终点的梯度各分量都接近 0说明确实停在了驻点。如果梯度的量级还在几千甚至几万说明目标函数本身可能有数值问题最常见的就是数据量纲太大导致目标函数爆炸。此时先做数据标准化再重新拟合往往比调优化器参数更有效。还有一个容易忽略的点检查目标函数是否写错最直接的方法是扫一遍似然面。在最优解附近让 μ 或 log σ 各偏离一点把 NLL 的变化打印出来看是不是在最优解处取最小。这一步虽然土但比对着公式空想有效得多。4. 把极大似然写进神经网络负对数似然作为损失函数4.1 分类任务里 NLL 与交叉熵为什么是同一件事深度学习里你天天用的交叉熵损失本质上就是极大似然。分类模型的输出经过 softmax 之后每个类别的概率可以看作一个 categorical 分布。对这批样本的联合概率取负对数就得到了交叉熵NLL -Σ ln p(yᵢ | xᵢ)展开之后它等于「真实标签 one-hot 向量」与「预测概率对数」的点积负值。PyTorch 里的nn.CrossEntropyLoss内部是先做log_softmax再取负相当于把 softmax 和 NLL 合在一步完成数值稳定性比先 softmax 再 log 更好。这也是为什么框架官方文档总提醒你不要在CrossEntropyLoss之前手动套 softmax否则梯度会翻车。理解这层等价关系最大的实际收益是排查损失异常时有依据。比如训练中 loss 突然变 NaN先检查 logits 是否出现极大值或极小值再查看标签是否有越界。知道了交叉熵的极大似然本质你就不需要背这些规则而是能推导出来。4.2 回归任务里 MSE 的极大似然解释回归任务里常用的均方误差MSE同样可以视为高斯噪声假设下的负对数似然。假设模型的预测值 μ(x) 是真实均值观测值 y 服从 y ~ N(μ(x), σ²)那么单个样本的 NLL 是NLL 0.5 ln(2πσ²) (y - μ(x))² / (2σ²)当 σ² 是常数时忽略与参数无关的项最小化 NLL 就等价于最小化 MSE。这也是为什么 MSE 在大多数回归任务里够用它隐含了「误差同方差、服从正态」的假设。但真实场景里这个假设经常不成立。某些输入下的预测天然更不确定比如传感器在低温段噪声明显大于高温段。这时候如果你还用固定 σ² 的 MSE模型会把所有样本的误差一视同仁结果是在噪声大的区域过度拟合、在噪声小的区域欠拟合。解决办法是把 σ² 也作为模型输出的一部分这正好是极大似然可以自然扩展的地方。4.3 用 torch 实现带不确定性输出的回归头下面的模型输出两个值预测均值 μ 和 log σ²然后用极大似然形式构造损失。这样模型不仅能给出预测还能给出每个样本的不确定性估计在异常检测、主动学习和融合任务里非常实用。import torch import torch.nn as nn class UncertaintyHead(nn.Module): 同时输出 mu 和 log_var 的回归头 def __init__(self, feat_dim): super().__init__() self.fc_mu nn.Linear(feat_dim, 1) self.fc_log_var nn.Linear(feat_dim, 1) def forward(self, x): mu self.fc_mu(x) log_var self.fc_log_var(x) # 输出 log 方差保证 var 0 return mu, log_var def gaussian_nll_loss(mu, log_var, y): 高斯负对数似然损失。 相比 MSE额外考虑了每个样本自身的方差方差大的样本贡献自动降低。 var torch.exp(log_var) # 0.5 * (log(2pi) log(var)) (y - mu)^2 / (2 * var) return (0.5 * (log_var torch.log(torch.tensor(2 * torch.pi)) (y - mu)**2 / var)).mean() # 使用示例 model UncertaintyHead(feat_dim8) x torch.randn(16, 8) y torch.randn(16, 1) mu, log_var model(x) loss gaussian_nll_loss(mu, log_var, y)这里最关键的参数是log_var不用加约束因为它是方差的对数天然可以取到整个实数轴。训练时如果发现log_var一路猛降到负值说明某个样本的残差特别大模型正在试图通过降低方差来惩罚它——这是正常行为但如果你不想要这种自适应权重可以给log_var加一个最小值下限比如log_var log_var.clamp(min-5)。实际使用中还有一个常见坑初始化fc_log_var的偏置为 0等价于初始方差为 1损失起步会比较正常。如果不做任何初始化方差从随机值出发前期 loss 可能出现剧烈抖动看起来像发散了其实只是方差在重新校准。5. 极大似然实战避坑五个反复踩的现象与根治办法5.1 log(0) 导致的 NaN加 epsilon 与 log-sum-exp现象计算对数似然时出现 NaN尤其在预测概率接近 0 或 1 的时候梯度直接变成 NaN训练中断。原因似然函数是大量概率的连乘取对数后每一项是单个概率的对数。只要有一个样本的概率被算成 0整个对数似然就是负无穷。在神经网络里softmax 输出经过指数运算后可能下溢为 0log(0) 就是 NaN。解决优先用 log-sum-exp 技巧框架里对应log_softmax自己写代码时给对数输入加一个微小下界比如torch.log(prob 1e-8)。但注意 epsilon 不能太大否则会引入偏差通常 1e-7 到 1e-6 够用。如果是手写高斯 NLL注意让 σ 的预测不要走到 0log σ 表示法天然规避了这个问题。5.2 参数不可辩识似然面存在平谷现象多次拟合得到完全不同的参数但 NLL 几乎一样或者优化器报告成功收敛参数的置信区间却宽到离谱。原因某些参数组合对似然函数的影响是等价的。最典型的是混合高斯模型里两个成分交换标签似然不变还有线性模型里两个特征高度共线时参数的组合无法唯一确定。解决先观察似然面是否是平的。做法是固定其他参数扫描目标参数的范围把 NLL 画出来。如果是一条平线说明数据里根本没有足够信息识别这个参数。此时要么增加数据要么固定其中一个参数要么改用带先验的 MAP 估计给参数一个合理的取值范围。5.3 局部最优与初始值敏感多起点策略现象同一份数据换一组初始值拟合结果差距很大而且其中一组明显不合理。原因NLL 函数非凸存在多个局部极小值。优化器从某个初始点出发很可能掉进最近的局部坑里出不来。解决不要只跑一次。常见做法是随机生成 20~50 组初始值分别优化取 NLL 最小的结果。每组的最终 NLL 也可以打印出来看看分布如果几十组结果离散程度很大说明似然面的地形很复杂单次优化结果不可信。这个策略在混合模型上几乎是必须的算力成本可控但收益非常直接。5.4 数值下溢与精度丢失scaling 与 dtype现象数据数值很大时比如温度几万、时间戳十几位NLL 动辄变成 1e10 的量级优化器收敛很慢甚至直接不收敛。原因目标函数值太大梯度也相应放大有限差分近似时小扰动会被整数精度吃掉梯度失效。解决先把数据标准化让均值为 0、方差为 1拟合完之后把参数再反变换回原始尺度。另外检查 dtypenumpy 默认 float64但某些库会默认 float32精度不够时对数值梯度的影响非常大。用data.astype(np.float64)避免隐式降精度。5.5 把 MLE 当万能过拟合与退化问题现象模型参数很多、数据量不大时MLE 拟合出的参数极度过拟合甚至出现「方差趋近 0、似然无穷大」的退化情况。原因MLE 没有任何正则化机制它只负责最大化训练数据上的似然。参数越多模型越容易钻到训练样本的细节里。高斯混合模型里某个成分收缩到单个样本点就是这种退化。解决给目标函数加正则项或者直接改用 MAP——把参数的先验分布加到对数似然上。从数学形式上加一个高斯先验等价于 L2 正则这个视角能帮你理解为什么正则项系数不该拍脑袋定而是和先验的方差直接挂钩。遇到数据量小的拟合任务我一般会先把模型复杂度降下来实在降不下来就用 MAP。6. 进阶验证技巧用似然比检验和 BIC 判断模型值不值得留拟合完参数只是第一步更难的问题是你同时拟合了两个模型一个是简单模型一个是加了参数的复杂模型该留哪个只看 NLL 的话复杂模型永远不会输因为参数多了、灵活度高了。这时候用似然比检验和 BIC 做筛选比拍脑袋可靠。似然比检验的思路是复杂模型相比简单模型多出的参数是否显著提升了似然。统计量 D 2*(ℓ_full - ℓ_reduced)在大样本下近似服从自由度为参数差数的卡方分布。import numpy as np from scipy.stats import chi2 # 假设两个模型已经拟合完成得到对数似然值 ll_full -1234.5 # 复杂模型比如二次项回归的对数似然 ll_reduced -1240.1 # 简单模型比如线性回归的对数似然 df 1 # 参数个数之差 D 2.0 * (ll_full - ll_reduced) p_value 1.0 - chi2.cdf(D, dfdf) print(f似然比统计量 {D:.3f}, p 值 {p_value:.4f}) # BIC越小越好 n_samples 200 k_full 4 # 复杂模型参数个数 bic_full -2 * ll_full k_full * np.log(n_samples) print(fBIC_full {bic_full:.2f})p 值小于 0.05 时可以认为复杂模型带来了显著提升BIC 则更严格它对参数个数的惩罚随样本量增大而增大。我个人的习惯是样本量大时看 p 值样本量小时看 BIC因为小样本下似然比检验容易过度拒绝原假设。在应用似然比检验前先确认两个模型是嵌套关系——也就是复杂模型在特定参数取值下能退化成简单模型。如果模型不嵌套比如一个是高斯、一个是威布尔比较 NLL 没有统计意义这时候只能靠 BIC 或交叉验证。这套验证流程给我最大的教训是很多拟合项目做完之后大家只关心「参数估计出来是多少」很少有人追问「这个模型到底值不值得留」。前几年我处理一批时间序列数据花了很久拟合一个带季节项的复杂模型NLL 确实好看但一算 BIC比简单模型还差等于白忙一场。后来所有项目我都把似然比检验和 BIC 写进默认流程参数估计完顺手就验证一轮选型决策快了很多。希望帮到你。本文还有配套的精品资源点击获取