概率统计随机过程之最大似然估计拓展

Apr 16, 2021 · 概率统计随机过程  ·
分享到:

概率统计随机过程之最大似然估计拓展

最大似然估计是参数点估计中非常常用且有效的估计方法，本文对其相关性质进行拓展介绍，包括渐进正态性，一致最小方差无偏估计，C-R界以及存在隐变量时采用的EM算法。

本文阅读需建立在已经理解最大似然估计基础上，首先我们精简的概述最大似然估计。

最大似然估计概述

核心思想：概率大的事件比概率小的事件更容易发生。要估计的参数能够使产生这个样本的概率最大。

步骤：

  1. 总体的概率/密度函数。
  2. 写出似然估计函数\(L(\theta)\)，其中\(\theta\)是待估计参数。
  3. 两边取\(\ln\)，即为\(l(\theta)=\ln(L(\theta))\)。(注意，对数似然函数是小写\(l\))
  4. 对\(l(\theta)\)的\(\theta\)求导（多参数估计就是偏导），令导数为0。
  5. 求出使导数为0的\(\theta\)即为最大似然估计参数\(\hat\theta_{MLE}\)。

TIPS：对于有些分布的极大似然估计没法直接求，比如均匀分布。

似然函数取对数的原因:

  1. 减少计算量。乘法变成加法，从而减少了计算量；同时，如果概率中含有指数项，如高斯分布，能把指数项也化为求和形式，进一步减少计算量；另外，在对联合概率求导时，和的形式会比积的形式更方便。
  2. 计算时更准确。为概率值都在[0,1]之间，因此，概率的连乘将会变成一个很小的值，可能会引起浮点数下溢，尤其是当数据集很大的时候，联合概率会趋向于0，非常不利于之后的计算。
  3. 取对数后，可以是一个上凸函数，更有利于求取最大值。

需要指出的是：取对数不影响单调性。 \[ p(x|\theta_1)>p(x|\theta_2)\Leftrightarrow \ln(p(x|\theta_1))>\ln(p(x|\theta_2)) \] 因为相同的单调性，它确保了概率的最大对数值出现在与原始概率函数相同的点上。因此，可以用更简单的对数似然来代替原来的似然。

最大似然估计的不变原理

我们介绍一个致使最大似然估计得到广泛应用的定理：

不变定理：设\(X\sim p(x;\theta),\;\theta\in\mathcal{\Theta}\)，若\(\theta\)的最大似然估计为\(\hat\theta\)，则对任意函数\(\gamma=g(\theta)\)，\(\gamma\)的最大似然估计为\(\hat\gamma=g(\hat\theta)\)。

最大似然估计服从渐进正态的条件

最大似然估计除了不变原理，还有另一个非常好的性质就是在一定条件下满足渐进正态性，从而保证了样本增加时的参数估计的收敛速度（通常为\(\frac{1}{\sqrt{n}}\)）。

渐近正态性

渐进正态性（大样本性质）：估计量的渐进正态性来源于中心极限定理，若统计量在样本容量\(n\rightarrow \infty\)时，统计量的分布也渐近于正态分布，称为渐进正态性。具体可定义为：如果存在一序列\(\{\sigma_n^2\}\)，满足\((\hat\theta_n-\theta)/\sigma_n(\theta)\overset{L}{\rightarrow}N(0,1)\)，则称\(\hat\theta_n\)是\(\theta\)的渐进正态估计，\(\sigma_n^2\)称为\(\hat\theta_n\)的渐进方差。

严格定义如下：

渐近正态性：设\(\hat\theta_n=\hat\theta(x_1,x_2,\dotsb,x_n)\)是\(\theta\)的一个相合估计序列，若存在一个趋于零的正数列\(\sigma_n(\theta)\)，使得规范变量\(y_n=\frac{\hat\theta_n-\theta}{\sigma_n(\theta)}\)的概率分布函数\(F_n(y)\)收敛于标准正态分布函数\(\varPhi(y)\)，即： \[ F_n(y)=P(\frac{\hat\theta_n-\theta}{\sigma_n(\theta)}\leq y)\rightarrow \varPhi(y)(n\rightarrow \infty) \] 如果用依分布收敛符号\(L\)可表示为： \[ \frac{\hat\theta_n-\theta}{\sigma_n(\theta)}\overset{L}{\rightarrow} N(0,1)(n\rightarrow\infty) \] 则称\(\hat\theta_n\)是\(\theta\)的渐进正态估计，或称\(\hat\theta_n\)具有渐进正态性，即： \[ \hat\theta_n\sim AN(\theta,\sigma^2_n(\theta)) \] 其中，\(\sigma^2_n(\theta)\)称为\(\hat\theta_n\)的渐进方差。

从\((\hat\theta_n-\theta)/\sigma_n(\theta)\)来看，分子项依概率收敛于\(\theta\)的速度与分母项\(\sigma_n(\theta)\)趋近于0的速度相同时，其比值才会稳定与正态分布。因此，\(\hat\theta_n\)收敛速度与渐近方差直接相关，渐近方差越小，收敛越快。

还需要指出，渐进方差并不是唯一的，如果存在另一\(\tau_n(\theta)\)有\(\frac{\sigma_n(\theta)}{\tau_n(\theta)}\rightarrow 1(n\rightarrow\infty)\)，根据依概率收敛定义可知\(\tau_n(\theta)\)也是\(\theta\)的渐近方差。

渐进正态性和相合性的关系类似于中心极限定理和大数定律。相合性是对估计的一种较低要求，它只要求估计序列\(\{\hat\theta_n\}\)在样本数量\(n\)增加的时候也趋近于\(\theta\)，但是并没有指出趋近的速度（例如是\(1/n,1/\sqrt{n}\)或\(1/\ln n\)）。而渐进正态性补充了这一点，收敛速度与渐进方差相关。经验来看，大多数渐进正态估计都是以\(1/\sqrt{n}\)的速度收敛于被估参数的。

最大似然估计具有正态性定理

在一定条件下，最大似然估计具有渐进正态性。我们将通过如下定理阐释。需要指出的是，定理是以连续分布的形式给出，但是对于离散场景也是适用的。

设\(p(x;\theta)\)是某密度函数，其参数空间\(\varTheta=\{\theta\}\)是直线上的非退化区间（即不是一个点），假如：

  1. 对一切\(\theta\in\varTheta\)，\(p(x;\theta)\)对\(\theta\)如下偏导都存在：\(\frac{\partial\ln p}{\partial\theta},\frac{\partial^2\ln p}{\partial\theta^2},\frac{\partial^3\ln p}{\partial\theta^3}\)
  2. 对一切\(\theta\in\varTheta\)，有\(|\frac{\partial\ln p}{\partial\theta}|<F_1(x),|\frac{\partial^2\ln p}{\partial\theta^2}|<F_2(x),\frac{\partial^3\ln p}{\partial\theta^3}<H(x)\)成立，其中\(F_1(x)\)与\(F_2(x)\)在实数轴上可积，而\(H(x)\)满足：\(\int_{-\infty}^\infty H(x)p(x;\theta)<M\)，这里\(M\)与\(\theta\)无关。
  3. 对一切\(\theta\in\varTheta\)，有\(0<I(\theta)=E[(\frac{\partial\ln p}{\partial \theta})^2]<+\infty\)

则在参数真值\(\theta\)为参数空间\(\varTheta\)内点的情况下，其似然方程有一个解存在，且此解\(\hat\theta_n=\theta(x_1,x_2,\dotsb,x_n)\)依概率收敛于\(\theta\)，且： \[ \hat\theta_n\sim AN(\theta,[nI(\theta)]^{-1}) \] 其中，\(I(\theta)\)为分布\(p(x;\theta)\)中含有\(\theta\)的费雪信息量，简称信息量。

这个定理的意义在于给定了最大似然分布有渐进正态性的条件，其中渐进方差（体现大样本效率）完全由样本数量\(n\)和分布的费雪信息量\(I(\theta)\)决定，且费雪信息量越大（分布中含有\(\theta\)）的信息越多，渐进方差在同等样本数量下越小，从而最大似然估计效果越好。

我们在这个定理中引入了费雪信息量这一概念，因此还要注意费雪信息量是否存在。对此，我们给出一个简单的概要：Cramer-Rao正则分布族中的分布费雪信息都是存在的。Cramer-Rao正则分布族定义如下：

Cramer-Rao正则分布族：分布\(p(x;\theta)\)，\(\theta\in\varTheta\)属于Cramer-Rao正则分布族，则需要满足以下五个条件：

  1. 参数空间\(\varTheta\)是直线上的开区间；
  2. \(\frac{\partial\ln p}{\partial\theta}\)对所有\(\theta\in\varTheta\)都存在；
  3. 分布的支撑\(\{x:p(x;\theta)>0\}\)与\(\theta\)无关；
  4. \(p(x;\theta)\)的微分与积分运算可交换；
  5. 对所有\(\theta\in\varTheta\)，期望\(0<E[(\frac{\partial\ln p(x;\theta)}{\partial\theta})^2]<+\infty\)。

常用的分布大多数都属于Cramer-Rao正则分布族，但是均匀分布\(U(0,\theta)\)不是，因为其分布的支撑与\(\theta\)相关。

Cramer-Rao正则分布族下费雪信息量的另一种表示

若\(p(x;\theta)\)为Cramer-Rao正则分布族，若其二节偏导数\(\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}\)对一切\(\theta\)存在，其费雪信息量还可以写为： \[ I(\theta)=-E[\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}] \] 这种方式可以帮助我们简便计算费雪信息量。此外，这是理解费雪信息量的另一个角度。我们知道二阶导数衡量一个函数的“凸程度”，绝对值越大，越陡峭，说明真值集中性越强。意味着，费雪信息量越大，我们更可能在凸峰出取到真实值附近的值。下面我们证明这两种表示是相等的。

证明：由于\(p(x;\theta)\)为Cramer-Rao正则分布族，且其二阶偏导数总是存在，对于定积分 \[ 1=\int_{-\infty}^{+\infty}p(x;\theta)\mathrm{d}x \] 两侧对\(\theta\)求偏导，由于定积分是个常数，因此导数必为0： \[ 0=\frac{\partial}{\partial\theta}\int_{-\infty}^{+\infty}p(x;\theta)\mathrm{d}x \] 根据Cramer-Rao正则分布族的第4条定义，我们交换上式中积分与微分的顺序： \[ 0=\frac{\partial}{\partial\theta}\int_{-\infty}^{+\infty}p(x;\theta)\mathrm{d}x=\int_{-\infty}^{+\infty}\frac{\partial p(x;\theta)}{\partial\theta}\mathrm{d}x \] 下面用到了一个计算技巧：\(\frac{\partial p(x;\theta)}{\partial\theta}=\frac{\partial\ln p(x;\theta)}{\partial\theta}p(x;\theta)\)，代入其中有： \[ 0=\int_{-\infty}^{+\infty}\frac{\partial\ln p(x;\theta)}{\partial\theta}p(x;\theta)\mathrm{d}x=E[\frac{\partial\ln p(x;\theta)}{\partial\theta}] \] 我们得到一个副结论：Log似然函数的一阶导期望为0。但我们还没有得到的结论，结论中需要二阶导，所以我们再对\(\int_{-\infty}^{+\infty}\frac{\partial\ln p(x;\theta)}{\partial\theta}p(x;\theta)\mathrm{d}x\)求\(\theta\)的偏导有： \[ \begin{aligned} 0&=\frac{\partial}{\partial\theta}\int_{-\infty}^{+\infty}\frac{\partial\ln p(x;\theta)}{\partial\theta}p(x;\theta)\mathrm{d}x(交互积分微分顺序)\\ &=\int_{-\infty}^{+\infty}[\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}p(x;\theta)+\frac{\partial\ln p(x;\theta)}{\partial\theta}\frac{\partial p(x;\theta)}{\partial\theta}]\mathrm{d}x\\ &\because \frac{\partial p(x;\theta)}{\partial\theta}=\frac{\partial\ln p(x;\theta)}{\partial\theta}p(x;\theta)\\ 0&=\int_{-\infty}^{+\infty}[\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}p(x;\theta)+(\frac{\partial\ln p(x;\theta)}{\partial\theta})^2p(x;\theta)]\mathrm{d}x\\ 0&=\int_{-\infty}^{+\infty}[\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}p(x;\theta)\mathrm{d}x+\underbrace{\int_{-\infty}^{+\infty}(\frac{\partial\ln p(x;\theta)}{\partial\theta})^2p(x;\theta)]\mathrm{d}x}_{I(\theta)=E[(\frac{\partial\ln p(x;\theta)}{\partial\theta})^2]}\\ \end{aligned} \] 从而有： \[ \begin{aligned} I(\theta)&=0-\int_{-\infty}^{+\infty}\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}p(x;\theta)\mathrm{d}x\\ &=-E[\frac{\partial^2\ln p(x;\theta)}{\partial\theta^2}] \end{aligned} \] 得证。

关于副结论：Log似然函数的一阶导期望为0。我们有时候也把\(p\)的Log似然函数称为分数函数，即分数函数的期望为0，而费雪信息量正好是分数函数的二阶原点矩，同时由于其期望为0，二阶原点矩等于二阶中心矩，即费雪信息量也正好是分数函数的方差。

Cramer-Rao不等式与界

从最大似然估计的渐进正态性，我们知道估计量始终会是一个随机变量，有自己的概率分布，而不是一个准确的值。Cramer-Rao除了给出了Cramer-Rao正则分布族这种费雪信息的存在条件，还有另一个更重要的贡献：C-R不等式，可以说给了统计学理论上的绝望。

C-R不等式，其实就是在说：统计，对真实的概率分布参数估计能力是有限的。举个不太恰当的类比，有点像量子理论中的测不准原理 （二者证明有相似之处哦）。C-R不等式告诉我们，无论我们如何抽样充足，无论我们统计方法如何科学，我们对参数的估计值，永远不可能无限逼近是逻辑上的真实值！它有自己估计的上限！

最大似然估计与相对熵（KL散度）、交叉熵的等价性

在机器学习、广义线性模型的应用场景中，常常使用最小（大）化交叉熵或相对熵，又叫KL散度来替换一些最大似然估计的操作，从理论上来看它们是具有等价性的。关于相对熵的概念，可以先看笔记《无线通信之信息熵回忆总结》。最大似然估计的目标是最大化似然函数，我们将其写为： \[ \argmax_{\theta} \ln(L(\theta))=\argmax_{\theta} \ln(\prod_{i=1}^n p(x_i;\theta)) \] 显然，似然函数\(L(\theta)\)是由抽样出来的已知样本\(x_i\)的概率密度函数（pdf）连乘得来的。当参数\(\theta\)已知时，\(\prod_{i=1}^n p(x_i;\theta)\)就是\(x=(x_1,x_2,\dotsb,x_n)\)的联合概率密度值；如果\(\theta\)未知，那么就会得到一个关于\(\theta\)的函数\(L(\theta)\)，而找到让函数\(L(\theta)\)最大的\(\theta^\star\)就是最大似然估计的核心。

如果将\(\ln\)符号和\(\prod\)符号交换，上式可改写成： \[ \argmax_{\theta} \ln(\prod_{i=1}^n p(x_i;\theta))=\argmax_{\theta}\sum_{i=1}^n \ln(p(x_i;\theta)) \] 简单地，我 们可以对似然函数加个负号，然后取最大值就变成了取最小值： \[ \argmax_{\theta}\sum_{i=1}^n \ln(p(x_i;\theta))=\argmin_{\theta}\sum_{i=1}^n -\ln(p(x_i;\theta)) \] 由于\(n\)是一个固定的数，因此在上式前面乘以一个常数\(\frac{1}{n}\)不改变\(\argmin\limits_{\theta}\)的结果，即\(\argmin\limits_{\theta} \frac{1}{n}\sum\limits_{i=1}^n -\ln(p(x_i;\theta))\)。

接下来是体现抽样与概率关系的一步。我们知道\(x_i\)都是从总体抽出来的样本，概率密度大的地方抽样出来的样本多，概率密度小的地方抽出来的样本就少。那么，想想我们如何使用蒙特卡洛法求一个随机变量函数的期望\(E[f(Y)]\)的？其中\(Y\)是一个随机变量。理论上来说： \[ E[f(Y)]=\int_{-\infty}^\infty p(y)f(y)\mathrm{d}y \] 但是，我们常常不知道\(Y\)精确的概率密度函数\(p(y)\)，因此我们会从总体\(Y\)中抽取一定量的\(y_i\)，然后用 \[ \tilde{E}[f(Y)]=\frac{1}{n}\sum_{i=1}^n f(y_i) \] 来近似计算期望，我们抽样的数量越多，这个值越准确。当抽样的数量趋近于无穷时，\(\tilde{E}[f(y)]{\rightarrow}E[f(y)](P=1)\)。这个表述是不是看得眼熟，没错这就是强大数定律！由此我们可以得到一个结论：

结论1：当样本数足够大时，可以使用\(\frac{1}{n}\sum_{i=1}^n f(y_i)\)替代期望\(E[f(Y)]\)。

回到原来\(\argmin\limits_{\theta} \frac{1}{n}\sum\limits_{i=1}^n -\ln(p(x_i;\theta))\)的场景，我们令\(f(x)=-\ln(p(x_i;\theta))\)，并且\(x_i\)都是从总体的抽样出来的，那么当样本数量较大时，我们同样可以用求期望的方式去近似函数均值，即： \[ \argmin\limits_{\theta} \frac{1}{n}\sum\limits_{i=1}^n -\ln(p(x;\theta))\approx\begin{cases} \argmin\limits_{\theta} \int\limits_x p(x) [(-\ln(p(x;\theta))]\mathrm{d}x,\quad x\text{为连续随机变量}\\ \argmin\limits_{\theta} \sum\limits_x p(x) [(-\ln(p(x;\theta))],\quad x\text{为离散随机变量} \end{cases}\\ =\argmin_\theta E[-\ln(p(x;\theta))]=\argmin_\theta H(x;\theta) \] 需要注意的是，作为概率的\(p(x)\)是与\(\theta\)无关的事实概率密度函数。因为抽样是从真实分布中抽样的，抽样结果只和\(x\)的真实分布有关，或者说抽样结果已经是包含参数\(\theta\)信息之后，体现出来的结果。到这里，我们又发现了一个重要的结论：

结论2：最大似然估计的效果等同于调节参数让信息熵最小。

这个直观上很好理解，最大似然估计是找出让事件发生概率最大的参数\(\theta\)，而熵是衡量事件不确定性的度量，越可能发生的事情，其信息含量越小，熵也越小。因此，最大似然估计会使得\(\argmin\limits_\theta H(x;\theta)\)。

到这里，我们距离交叉熵与相对熵（KL散度）貌似还有一些距离。因为上式中，我们使用的是一般的熵，而不是交叉熵、相对熵，那么区别在哪里呢？

答案就是似然函数\(-\ln(p(x;\theta))\)。我们在上文中都默认使用了真实的概率分布\(p(x;\theta)\)，然而我们实际上，我们并不知道！并且真实分布是客观存在的，实际上也不存在参数\(\theta\)。我们只是选取了一种带有参数\(\theta\)概率分布族\(q(x;\theta)\)的去近似真实的概率\(p(x)\)；或者假定数据\(x_i\)服从某一带有参数\(\theta\)概率分布族\(q(x;\theta)\)。然而，实际情况是\(q(x;\theta)\)与真实分布\(p(x)\)是有差距的。因此，我们写出实际情况: \[ \argmin\limits_{\theta} \frac{1}{n}\sum\limits_{i=1}^n -\ln(p(x;\theta))\approx\begin{cases} \argmin\limits_{\theta} \int\limits_x p(x) [(-\ln(q(x;\theta))]\mathrm{d}x,\quad x\text{为连续随机变量}\\ \argmin\limits_{\theta} \sum\limits_x p(x) [(-\ln(q(x;\theta))],\quad x\text{为离散随机变量} \end{cases}\\ =\argmin_\theta E[-\ln(q(x;\theta))] \] 注意只是将\(p(x;\theta)\)换成了\(q(x;\theta)\)，抽样的\(x_i\)反应的是真实的概率分布，因此就是约等于\(p(x)\)，且根据大数定律，样本数越大，越接近真实分布。根据笔记《无线通信之信息熵回忆总结》的定义，我们知道 \[ E[-\ln(q(x;\theta))]=\sum_{x}p(x)(-\ln q(x;\theta))=H(p,q) \] 其中，\(H(p,q)\)成为交叉熵，只有我们设定的参数分布族\(q(x;\theta)\)才与\(\theta\)相关，真实分布作为客观存在与\(\theta\)无关。

根据交叉熵和相对熵的关系 \[ D_{KL}(p||q)=H(p,q)-H(p) \] 可知，相对熵和交叉熵相差一个真实熵\(H(p)\)。因为\(p(x)\)是客观存在的，因此\(H(p)\)也是一个客观存在的常数，虽然我们不知道具体是多少，但是它确实只是个常数，与\(\theta\)无关，因此 \[ \argmin_\theta E[-\ln(q(x;\theta))]=\argmin_\theta H(p,q)=\argmin_\theta D_{KL}(p||q) \] 其中，\(q\)是与\(\theta\)相关的分布族，\(p\)与\(\theta\)无关。

综上可知，最大似然估计等同于最小化交叉熵或相对熵（KL散度）。

最大似然估计不具备的性质

  1. 最大似然估计不一定是无偏的。
  2. 当待估计参数与自变量\(x\)相关时，无法用求导的方式求解。我们需要回到最大似然估计的基本定义，既让似然函数概率最大，来考虑问题。