知行合一

极大似然估计(MLE)和最大后验概率估计(MAP)

MLE 与 MAP 分别对应频率学派与贝叶斯学派的参数估计框架:MLE 通过最大化似然获得估计,MAP 在此基础上引入先验并最大化后验。

频率学派与贝叶斯学派


关于频率学派与贝叶斯学派的标准解释网上有很多,下面是个人的一些理解。

频率学派

频率学派认为世界是确定的,也即某个事件发生的概率是确定的。在一次试验中无论该事件发生还是不发生都只是其发生概率$p$的一次体现而已。所以当进行多次重复试验时,事件发生的次数占总试验次数的比例就会逐渐接近这个事件的发生概率$p$。

在上述的前提下,如果想要知道这个$p$,那么只需要进行多次试验,然后解出$p$的最可能的取值即可。这就是从最初开始学习统计学时就接触到的极大似然估计(MLE)。很显然,在频率学派的方法下,当进行足够多的试验后即积累足够的数据量时我们可以得到一个很接近$p$的估计。

贝叶斯学派

与频率学派不同,贝叶斯学派认为世界是不确定的,它是随着我们对世界的认知改变的(就感觉很主观,但是又主观得很有道理…)。因此,在贝叶斯学派下,事件发生的概率$p$并不是一个确定的数字,根据不同的认知(即先验假设)对$p$应该作出不一样的估计,即最大后验概率估计(MAP)

随着试验次数的增加,我们的先验假设的作用被逐渐淡化,数据中体现的信息将会在估计中占据主导作用,因此当数据量足够大时,会发现两种估计方法会得到相同的结论。

极大似然估计(MLE)


以丢硬币为例,丢硬币10次,其中正面朝上7次,反面朝上3次,使用极大似然估计求解正面朝上的概率$p$。

有似然函数:

$$ f(X, p) = P(X|p) = p^7(1-p)^3 $$

取对数:

$$ \ln f(X, p) = 7\ln p + 3\ln(1-p) $$

求导得到:

$$ (\ln f(X, p))' = {7p + 3p - 7 \over p(p-1)} = {10p - 7 \over p^2-p} $$

令上式等于0可以得到$p$等于0.7。即当丢一个硬币10次,其中正面朝上7次,反面朝上3次时,极大似然估计认为正面朝上的概率最可能是0.7。

我们知道一个均匀的硬币正面朝上的概率应该是0.5,所以当数据量较小时,极大似然估计的结果与我们的先验知识就产生了偏差(在数据量大的情况下并不会有这样的问题,因为如果丢一个硬币10000次,正面朝上7000次,反面朝上3000次,那么应该做的是好好查查这个硬币是不是真的是均匀的了)。

最大后验概率估计(MAP)


在极大似然估计中我们寻找使得$P(X|p)$取到最大值的$p$,而在最大后验概率估计中,我们寻找的是使得$P(X|p)P(p)$取到最大值的$p$。

根据贝叶斯公式有:

$$ P(X|p)P(p) = P(p|X)P(X) $$

这里的$P(X)$是一个已知的值(7次正面,3次反面),所以我们要求解的实际是使得$P(p|X)$(即后验概率,故称为最大后验概率估计)取到最大值的$p$。

根据过去多年的丢硬币的体验,我们的先验知识是硬币正面朝上的概率是0.5,因此我们可以将$p$服从一个均值为0.5,标准差为$\sigma$的分布作为先验假设。

$$ f(X, p) = P(X|p)P(p) = p^7(1-p)^3{1 \over \sqrt{2\pi}\sigma} \exp \left( -{(p-0.5)^2 \over 2\sigma^2} \right) $$

取对数并求导:

$$ (\ln f(X, p))' = {7 \over p} + {3 \over p-1} - {p-0.5 \over \sigma^2} $$

求令上式等于0的$p$,到这里我已经不会了,所以我们用一个比较笨的办法吧…

先定义上方需要求解的方程

ln_map <- function(x, sigma) {
  item1 <- 7 / x
  item2 <- 3 / (x -1)
  item3 <- (x - 0.5) / sigma^2
  return(item1 + item2 - item3)
}

然后定义一个二分法查找方程解的函数

binary_search <- function(func, ...) {
  tolerantion <- 1e-5
  min <- 0 + tolerantion
  max <- 1 - tolerantion
  while (TRUE) {
    x <- (max + min) / 2
    LN_MAP <- do.call(func, args = c(list(x = x), ...))
    if (LN_MAP > 0) {
      min <- x
    } else if (LN_MAP < 0) {
      max <- x
    } else {
      break
    }
    if (abs(LN_MAP) < tolerantion) {
      break
    }
  }
  return(x)
}

对于不同的标准差$\sigma$使用二分法查找方程的解

for (sigma in seq(0.01, 0.1, 0.01)) {
  cat("标准差为", sigma, "时,")
  solution <- binary_search(ln_map, sigma = sigma)
  cat("p的最大后验概率估计为", solution, "\n")
}

可以看到随着$p$所服从的正态分布的标准差增大,$p$的最大后验概率逐渐偏离0.5。可以理解为,我们对硬币正面朝上的概率是0.5这个先验假设越没有信心,$p$的最大后验概率估计越远离0.5。当我们对硬币正面朝上的概率是0.5这个概率完全失去信心时,$p$的分布退化为均匀分布,$(\ln f(X, p))'$回到了与最大似然估计同样的形式,于是对$p$的最大后验概率估计也为0.7。

现在设想另一种情况,我们对硬币正面朝上的概率是0.5这个先验假设非常有信心(即标准差非常小),但是在试验中使用的这个硬币真的是个不均匀的硬币,那我们要如何通过最大后验概率估计去得出一个正确的$p$的估计呢。这时候需要做的就是增加试验次数,比如,将试验次数增加到10000次,如果是一个不均匀的硬币,那么正面朝上和反面朝上的次数的比例应该约为7:3,于是可以修改我们的ln_map()函数:

ln_map <- function(x, sigma) {
  item1 <- 7000 / x
  item2 <- 3000 / (x -1)
  item3 <- (x - 0.5) / sigma^2
  return(item1 + item2 - item3)
}

此时即使给一个非常小的$\sigma$,我们也能得到一个接近0.7的最大后验概率估计。

binary_search(ln_map, sigma = 0.01)