最大似然估计：从样本到参数的统计推断核心方法

做数据分析这几年有一类对话反复出现业务方跑过来扔给我几百条用户留存数据问我到底该把留存率估成多少训练完一个回归模型问我这些系数背后的含义再或者做 A/B 测试时两组样本差了几个百分点问我能不能下结论。表面上看是三件事内核却共用同一把钥匙我们手里有一批观测样本想知道背后生成这批样本的分布参数是什么。最大似然估计MLE就是回答这类问题的第一把钥匙。很多人第一次接触最大似然估计是在概率论课上背下一个“求导、令其为零、解方程”的三段式流程转头就忘。我真正理解它也是工作以后的事——难的不是求导而是搞清楚为什么要这样求以及求出来的数值在什么条件下才可信。这篇文章按我自己的理解路径来写先讲清楚估计问题到底长什么样再把似然函数的构造、解析解推导、数值优化、统计性质和工程落地时常踩的坑逐个过一遍。适合正在学统计推断的学生也适合一直把 MLE 当黑盒在调用的建模工程师。1. 最大似然估计到底在最大化什么先修正那个“反直觉”的设定1.1 我们面对的是“有样本、没参数”的问题统计学里有一类问题叫参数估计。所谓参数就是决定分布形状的那些未知常量。比如一个二项分布参数 p 决定了成功的概率一个正态分布参数 μ 和 σ² 决定了中心位置和离散程度。现实世界中参数永远不知道我们手里只有有限个样本。我举个自己常遇到的例子。产品方说新用户第二天留存率可能做到 40%给我一个月的数据一共看了 320 个新用户其中 121 个第二天又回来了。那么问题来了留存率到底是多少在没有其他信息的情况下最朴素的办法就是拿 121/320 ≈ 37.8% 去交差。这个 37.8% 其实就是一个最大似然估计的结果只不过它简单到我们往往没意识到。模型化以后就是假设每个用户的次日留存服从一个伯努利分布参数是 p。每一位用户是否留存在给定 p 的前提下都是独立事件。我们要做的就是从这一串 0 和 1 的观测里反推 p。所有参数估计方法本质上都是“给定样本给参数打个分选一个最合理的”。1.2 “似然”这个词的方向和概率是相反的初学最让人困惑的一点就是“似然”和“概率”看起来长得一模一样方向却截然相反。概率是已知参数去预测数据。假设我知道 p0.37我可以算出 320 个人里正好出现 121 个留存者的概率有多大。似然是反过来数据已经固定下来了我去评价不同的 p 值在多大程度上“配得上”这批数据。写出来就是概率P(数据 | 参数) ——参数在前数据是随机的。似然L(参数 | 数据) ——数据在前参数是可变的。这也是为什么英文里管它叫 likelihood而不是 probability就是刻意想绕开“概率”这个已经被占用理解方向的概念。你甚至可以把它理解成对一批已经发生的事实我们追问“在某个参数假设下这个事实显得多自然”。越自然的参数得分越高。1.3 那个 arg max 究竟在做什么最大似然估计的定义其实一句话就能说完在所有可能的参数取值里找出一个让当前观测数据出现概率最大的参数。写成数学符号是θ̂ argmax L(θ | data) argmax P(data | θ)关键是“argmax”它不是问你最大值是多少而是问“最大值是在哪个参数那取到的”。很多人把这个弄混写程序求出一个很大的分数却忘了把对应的参数取出来。沿用留存的例子。p 取 0.1 的时候出现“320 人里 121 人留存”的概率非常小p 取 0.378 的时候这个概率在所有 p 里最大p 取 0.9 的时候概率又变得极小。MLE 做的就是把整条 p 从 0 扫到 1找到最高点所在的 p。它可能不是真实 p但在“当前这批数据”这个证据之下它是最合理的解释。2. 似然函数的两种构造离散分布和连续分布各有写法2.1 离散情形伯努利和二项分布的标准构造先从最简单的离散情况做起。假设观测样本是 x₁, x₂, …, xₙ每个 xᵢ 只有 0 和 1 两种取值概率分别是 1-p 和 p。由于独立同分布的假设多个事件同时发生的联合概率就是每个事件概率的乘积L(p) p^(x₁) (1-p)^(1-x₁) × p^(x₂) (1-p)^(1-x₂) × … × p^(xₙ) (1-p)^(1-xₙ)如果记 K Σxᵢ 是“成功次数”乘积就能压缩成L(p) p^K (1-p)^(n-K)你看它本质上是一个二项分布的概率质量函数只不过把 x 固定住剩下的自变量是 p。这里很容易产生一个误解既然 L(p) 和某个概率分布的表达式相同那 L 是不是也代表“概率”呢从数值上看它确实是一个概率值但它回答的问题已经变了——我们不再问“抽到 K 次成功的概率”而是问“给定 K哪个 p 最可能”。2.2 连续情形为什么是概率密度的连乘积连续分布的情况稍微抽象一点。以正态分布为例单次观测 xᵢ 的“概率”在严格意义上是 0因为连续变量取到某个精确值的概率测度是 0。那为什么 MLE 还要用概率密度函数 f(xᵢ|θ) 去相乘这里的关键是“局部概率的近似”。密度函数 f(x) 的意义是x 落在非常小的区间 [xᵢ, xᵢdx] 里的概率约等于 f(xᵢ)dx。如果 dx 足够小每个观测点附近小邻域的概率近似就是密度乘以宽度。把 n 个观测的邻域事件求交联合概率近似为∏ f(xᵢ|θ) · (dx)ⁿ其中那个 (dx)ⁿ 和参数无关在比较不同 θ 的时候被约掉了。所以我们干脆直接最大化 ∏ f(xᵢ|θ)也就是概率密度的乘积。这就是为什么 MLE 在连续分布下照样能用的底层原因。不过要敲个警钟密度乘积不再是一个真正的概率它可以大于 1。尤其在样本量很小的分布中密度函数的数值可能比 1 大不少。所以解读它的时候关注相对大小就够了别把它当成“概率”去对外汇报。2.3 为什么总要对似然取对数对数似然的三重价值标准教科书几乎都会给你下一步操作取对数。初学者往往觉得这只是为了求导方便其实它有三层意义。第一层是数值稳定性。样本量一大几百个小于 1 的小数相乘很快就下溢成 0。比如 0.1 的 500 次方在计算机里直接就变成 0连比较都无从谈起。取对数之后乘法变加法数值范围从指数级落到一个可以安全处理的量纲。第二层是数学便利。指数族分布的密度函数里几乎都有 exp(...)取对数把指数摘掉求导时不需要反复用乘法法则。以正态分布为例不取对数时对 μ 求导还得被 exp 里的二次项牵连取对数之后目标函数变成二次函数求导轻快很多。第三层是单调性保证。对数函数在正数区间是严格递增的所以 ln L(θ) 和 L(θ) 的最大值点完全一致。也就是说你求对数后的峰值对应的 θ和原始似然的峰值对应的 θ 是同一个。这就给了我们一个理论许可随便切换。3. 从手推到代码解析解和数值优化两条路都得懂3.1 正态分布的均值与方差最经典的完整推导我们把最大似然估计手推一遍正态分布这对后面理解线性回归太重要了。假设 x₁, …, xₙ 独立同分布且都来自 N(μ, σ²)单点的密度是f(xᵢ) 1 / (σ√(2π)) · exp(-(xᵢ-μ)² / (2σ²))似然函数是 n 个这样的东西相乘取对数之后变成ℓ(μ, σ²) -n/2 · ln(2π) - n · ln σ - 1/(2σ²) · Σ(xᵢ-μ)²先固定 σ²对 μ 求偏导并令其为零。Σ(xᵢ-μ)² 对 μ 求导等于 -2Σ(xᵢ-μ)于是得到 Σ(xᵢ-μ) 0解出来就是 μ̂ (1/n)Σxᵢ也就是样本均值。再看 σ²。对 σ² 求导之前先令 u σ²把 -ln(σ) 写成 -1/2 ln u。对 u 求导得到∂ℓ/∂u -n/(2u) 1/(2u²)Σ(xᵢ-μ)²令其为零解得 u (1/n)Σ(xᵢ-μ̂)²。这里有个很多初学者忽略的细节分母是 n 而不是 n-1。也就是说方差的最大似然估计是有偏的它比常用的样本方差略小一点。后面我会专门说这个“偏差换有效性”的取舍。你在实操里用 Python 的统计库或者手写 numpy 时经常看到默认除以 n-1那就是为了无偏性而做的一次纠偏不是 MLE 本来的输出。3.2 线性回归里的最小二乘原来就是 MLE线性回归的最小二乘法学的时候有点“横空出世”的感觉为什么非要最小化平方误差用 MLE 一推就明白了。假设 y wᵀx ε误差 ε 是独立同分布的正态噪声 N(0, σ²)。于是给定 w 和 σ² 时y 的条件分布是 N(wᵀx, σ²)。对全部 n 个样本构造对数似然ℓ(w, σ²) -n/2 · ln(2π) - n · ln σ - 1/(2σ²) · Σ(yᵢ - wᵀxᵢ)²注意到最后一项前面的系数与 w 无关最大化 ℓ 关于 w等价于最小化 Σ(yᵢ - wᵀxᵢ)²。这就是最小二乘。换句话说最小二乘不是一个独立发明出来的“经验损失”它是在“高斯噪声 独立同分布”假设下自然冒出来的等价写法。这个等价关系也解释了为什么线性回归对极端离群点那么敏感——因为正态分布的尾巴很薄偶尔一个离谱的观测会极大拉高平方误差项整个 w 的解都会被拽走。怀里没有这个假设就谈不上最小二乘的合理性。这也是为什么后来出现 Huber 损失、分位数回归这些稳健方法本质就是在放松“正态噪声”这个前提。3.3 当解析解不存在时梯度法、牛顿法与拟牛顿法怎么挑并不是每个模型都能像正态分布那样求出闭式解。逻辑回归的似然函数里带着 sigmoid对参数求导出来的方程没有解析解一些定制分布的 MLE 更是没法手推闭式解。这时候就要动用数值优化。最稳妥的起点是梯度上升或者对负对数似然做梯度下降。每次迭代沿着对数似然的梯度方向走一小步θ_new θ_old η · ∇ℓ(θ_old)这里的 η 是学习率。操作上我一般先看一轮迭代后对数似然的增量如果增量变成负的说明步子迈大了把 η 调小一半重来。梯度法实现简单对问题结构依赖少缺点也很明显收敛慢尤其在参数维度较高或者目标函数“峡谷”形状明显的时候。如果需要更快收敛就用牛顿法。牛顿法把目标函数在当前点做二阶泰勒展开移动步长需要用到 Hessian 矩阵 Hθ_new θ_old - H⁻¹∇ℓ(θ_old)这里 H 是 ℓ 对 θ 的二阶偏导矩阵。牛顿法在目标函数接近二次函数时收敛极快但现实中 Hessian 矩阵继续解和存储成本都很高。对 MLE 场景一个常见让步是用拟牛顿法比如 BFGS 或 L-BFGS用梯度信息近似 Hessian 的逆矩阵既保留快速收敛又大幅省内存。重复一句不管用什么优化器最终目标都是让负对数似然降下来。所以实践上一定在训练循环里打印逐轮的负对数似然值用它判断有没有走过头、有没有卡住而不是只看参数本身的变化。4. 同一条线索上的三个分支交叉熵、矩估计和贝叶斯估计4.1 逻辑回归的交叉熵损失就是从负对数似然来的很多做机器学习的人熟悉“交叉熵损失”这个名字却不清楚它怎么来的。逻辑回归里对二分类问题模型输出的是正类概率p σ(wᵀx) 1/(1e^(-wᵀx))这个 p 被看作伯努利分布的参数。对单个样本 (xᵢ, yᵢ)它的“贡献”是 p^(yᵢ) · (1-p)^(1-yᵢ)。将所有样本的似然乘起来再取负对数得到NLL(w) -Σ [ yᵢ ln pᵢ (1-yᵢ) ln(1-pᵢ) ]这个式子就是交叉熵损失。换句话说训练逻辑回归时你真正在做的事情就是最大似然估计在给定数据下找一组 w 使伯努利模型的似然最高。你看到的一切“损失函数设计”说到底是在选一个概率模型分布然后用负对数似然当优化目标。4.2 矩估计和 MLE 的对比谁更简单谁更可靠矩估计是另一条老牌推断路线。它的思路很朴素让样本矩等于理论矩。一阶矩让样本均值等于总体期望二阶矩让样本二阶原点矩等于总体二阶原点矩。这么做最大的好处是“无脑”——不需要搞懂似然函数不需要优化直接联立方程解参数。但可靠性常常不如 MLE。举个经典例子均匀分布 U(0, θ)真实 θ 未知。矩估计的方法是令 E[X] θ/2 等于样本均值 x̄解出 θ̂_moment 2x̄。MLE 的方法是最大化 L(θ) 1/θⁿ要求所有观测 ≤ θ解出 θ̂_mle max(x₁,…,xₙ)。你看当样本里最大值明明已经达到了 9.8你硬把它估成 7.5显然不合理。MLE 给出的是观测到的最大值它把所有观测的“证据”用到了极限。在正确设定模型的情况下MLE 的渐近方差往往比矩估计更小。所以说矩估计适合作为快速初筛MLE 才是拿来出正式结论的东西。4.3 贝叶斯估计只是比 MLE 多了一个先验项贝叶斯估计经常被拿来和 MLE 对比。它的目标是最大化后验概率p(θ|data) ∝ P(data|θ) · p(θ)其中 p(θ) 是参数的先验分布。取对数之后后验对数等于似然对数加先验对数ln p(θ|data) ℓ(θ) ln p(θ) 常数所以你看到的“最大后验估计” MAP本质上就是“最大似然估计加了一个先验的正则项”。当先验是均匀分布时ln p(θ) 是个常数MAP 和 MLE 完全等价。当先验是正态分布时它对应 L2 正则化当先验是拉普拉斯分布时它对应 L1 正则化。这个视角特别有用——每当你给 MLE 加正则项你其实是在默认参数更靠近某个位置。原本“只信数据”的估计变成了“带着经验看数据”。5. MLE 的理论保证到底在保证什么一致性、无偏性与边界风险5.1 一致性是所有保证里最值得信赖的一条MLE 最经典的理论性质叫一致性当样本量趋于无穷时估计值会依概率收敛到真实参数。直觉很好理解样本越多信息越多我们推断的“最优解释”就应该越来越接近真相。这正是我们日常建模最依赖的性质。手头有 30 个样本时MLE 可能偏得离谱收集到 3000 个样本后你会看到估计值慢慢稳定下来不再因为增加几条数据就跳出一个大弯。实际操作中如果你想快速估一个参数先跑一次“样本量逐渐增大时参数估计轨迹”的曲线假如曲线在某一处开始水平说明信息已经足够假如还在陡峭爬升那别急着把当前估计值写进报告。5.2 无偏性与 Cramér-Rao 下界效率的天花板MLE 并不总是无偏的。前面提到的方差估计分母是 n 而非 n-1就是有偏的典型。但是随着样本量增加偏差会逐步缩小渐近无偏这一条腿可以被满足。更有意思的是方差下界。在一定的正则性条件下任何无偏估计量的方差都下不了某个值Var(θ̂) ≥ 1 / I(θ)这里 I(θ) 是 Fisher 信息量等于 -E[∂²ℓ/∂θ²]。Fisher 信息量刻画的是“这个分布对参数有多敏感”。直观理解就是似然函数在真实参数附近的弯曲程度越大我们越有把握锁定真实参数对方差下界的要求就越低。MLE 的优秀之处在于它通常能渐近达到这个下界因此被称作“渐近有效”。这给实际操作提供了一个检查办法如果你用别的方法做推断可以计算一下估计值的经验方差和 CRB 的比值如果比值远远大于 1说明你浪费了数据里的信息这时换个估计方法往往就有收益。5.3 MLE 什么时候会翻车理论性质再好MLE 也有无能为力的时刻。下面这张表是我自己整理的常见风险情形每一次都是切切实实踩过之后才记住的。问题类型现场表现应对思路模型不可识别多个参数组合给出一模一样的似然值约束参数空间或改参数化方式多峰似然数值优化结果强烈依赖初始值多起点、网格扫描、EM 算法非正则似然真实参数在参数空间边界上导数不存在换约束优化或用变换后参数分布误设模型并非真实分布参数出现系统性偏差做残差诊断换稳健估计或更灵活分布以非正则似然为例如果真实参数落在边界上那在你求导的那个内点附近梯度根本不会为零。这时候“求导令其为零”的套路直接失效必须在边界约束下另想办法。多峰问题在混合模型里尤为常见天天有人跑出不同的结果又说不清谁对——答案通常是谁都不对你要么用多随机初值挑似然最大的那个要么用 EM 等专门算法去逼近。6. 工程落地时绕不开的三个细节数值、初值和正则化6.1 概率连乘的下溢问题log 域是底线不是技巧我在前面已经提过下溢问题但实际操作中的危险程度远超很多人的想象。计算机里 double 能表示的最小正数大概是 10^(-324)几个很小的密度乘上几百次结果直接落 0。一旦变成 0后面取对数只能得到负无穷梯度也无法计算。所以只要模型样本量超过几百一律先在日志里算。一般我会直接把损失函数改成负对数似然把模型内部所有涉及连乘的地方都换成求和。如果你需要计算一组归一化权重比如 softmax 类型公式最稳的办法是先把每个分量减去最大值再做 exp这就是 log-sum-exp 技巧。别嫌它啰嗦数值稳定性是所有做概率模型的人都必须过的第一道坎。6.2 初值敏感性和多峰问题永远不要只跑一次当似然函数有多个局部极大值时MLE 的数值解会跟着初始值走。类似爬一座有多座山峰的山从不同山脚出发终点可能落在不同山头且中间山峰之间还被谷地隔开简单梯度法翻不过去。我的习惯性做法是随机生成 5 到 20 组初始值分别跑完整轮优化记录每组参数对应的对数似然。最终选对数似然最大的那一组作为结果。如果发现好几组初值落在完全相同或非常接近的峰上那基本可以认为找到了主峰如果结果散得很开那就得回头审视模型本身——很可能是参数化导致的不可识别。6.3 样本量小与过拟合惩罚似然是自然的兜底MLE 在样本量不足时有一个明显毛病会把模型拟合得过于极端。比如高斯混合模型如果某个子类只有两三个样本MLE 得出的方差可能趋近于 0似然函数的值则冲向无穷大。这种情况用“正规化”的说法来兜底就是给对数似然加上惩罚项。在优化目标里直接加 λ‖θ‖² 或 λ‖θ‖₁逻辑上很清楚等价于 MAP 中选择了不同的先验。工程上有一些专门针对小样本的调整比如在估计协方差矩阵时加一个小对角阵或者用缩水估计。别看这些做法在数学上不起眼实际项目里它们对小样本的稳定性影响远比想象的大。我自己的体会是碰到这种场景先别急着调模型结构看一眼似然面和 Hessian 的特征值。如果 Hessian 有接近零的负特征值说明目标函数在某个方向上几乎是平面的这时参数本身就没有足够证据支撑加正则项往往是更合理的自我约束。做概率推断这么久最大似然估计给我的感觉不是一条公式而是一种思维方式。你先给数据定义一个生成模型然后相信“能解释当前证据的才是值得选择的”。它看起来朴素却把参数估计、回归、分类、正则化这些内容都串到了同一条线上。真到了项目里请记得MLE 的结果只在你对分布假设有信心的前提下才是可信的输出模型之前留时间看看残差、跑几次多初始值、检查一遍对数似然曲线这些细活比任何炫技都有用。