参数估计之极大似然估计、极大后验估计和贝叶斯参数估计

本文来源1：Parameter estimation for text analysis， Gregor Heinrich技术报告的翻译

在统计推断中，我们一般面临两项推断任务：1）估计一组分布参数集合$$\theta$$使其能最好地解释观测量$$\mathcal{X}$$; 2）在给定之前观测值的情况下，计算新观测量的概率，即计算$$p(\bar{x}|\mathcal{X})$$。我们把前者称作估计问题，后者称为预测或者回归问题。

假设，我们有数据集$$\mathcal{X} \cong \{ x_{i} \} _{i=1}^{|\mathcal{X}|} $$，是随机变量（$$\bold{r.v.}$$）$$X$$的一组独立同分布（$$\bold{i.i.d}$$）的样本。参数$$\vartheta$$是依赖分布的，比如对于高斯分布，$$\vartheta = \{\mu,\sigma^2\}$$。

本文接下来将分别依次说明极大似然估计、极大后验估计和贝叶斯估计的参数方法。我们将以N重伯努利试验（Bernoulli）实验为案例说明这些参数估计的实现。

1、极大似然估计（Maximum likelihood estimation）

极大似然估计是寻找参数以最大化如下似然的估计方法：

L(\vartheta|\mathcal{X}) \cong p(\mathcal{X}|\vartheta) = \bigcap_{x\in \mathcal{X}} \{X=x|\vartheta \} = \prod _{x \in X} p(x|\vartheta) \ \ \ \ \ \ \ (1)

上式也就是X产生数据$$\mathcal{X}$$的联合事件概率。这里的出发点是我们寻找一个参数，使其最有可能导致当前这些事件的发生，那么这个参数就是分布的参数了。比如，扔一个硬币100次，有90次正面以上，那么我们就寻找一个正面朝上概率，使得100次扔硬币90次朝上的概率最大就好了（显然这个硬币不是均匀的）。上式是连乘形式，计算不够方便，因此我们通常求解log。这个在之前的博客贝叶斯分析推断的一些基础知识中已经说明了：

\hat{\vartheta}_{ \mathrm{ML} } = \mathrm{argmax}_{\vartheta} \mathcal{L}(\vartheta|\mathcal{X}) = \mathrm{argmax}_{\vartheta} = \sum_{x\in \mathcal{X}}\log p(x|\vartheta)  \ \ \ \ \ \ \ (2)

公式（4）求解很简单。只需要计算其在$$\theta_k$$处的一阶偏导等于0的点即可：

\frac { \partial \mathcal{L} } {\partial \vartheta{k}} = 0 \ \ \ \ \forall \vartheta_k \in \vartheta \ \ \ \ \ \ \ (3)

在给定数据$$\mathcal{X}$$的情况下，对于新的观测量$$\tilde{x}$$的预测概率可以近似等于（公式（4）的推导过程参见先前的博客贝叶斯分析推断的一些基础知识）：

\begin{aligned}

p(\tilde{x}|\mathcal{X}) &= \int _{\vartheta \in \Theta} p(\tilde{x}|\vartheta)p(\vartheta|\mathcal{X}) \mathrm{d} \vartheta \ \ \ \ \ \ \ (4)\\

& \approx \int _{\vartheta \in \Theta} p(\tilde{x}|\hat{\vartheta}_{\mathrm{ML}})p(\vartheta|\mathcal{X}) \mathrm{d} \vartheta = p(\tilde{x}|\hat{\vartheta}_{\mathrm{ML}}) \ \ \ \ \ \ \ (5)

\end{aligned}

也就是说，下一个样本被认为是服从原始分布，其参数为$$\hat{\theta}_{\mathrm{ML}}$$。

2、极大后验估计（Maximum a posteriori estimation）

极大后验估计与极大似然估计非常相似，但是引入了一些参数的先验知识$$p(\vartheta)$$。名字来源于它要求给定数据情况下参数的极大后验：

\hat{\theta}_\mathrm{MAP} = \mathrm{argmax}_{\vartheta} p(\vartheta|\mathcal{X}) \ \ \ \ \ \ \ (6)

通过使用贝叶斯公式，上式可以写成：

\begin{aligned}

\hat{\vartheta}_{\mathrm{MAP}} &= \ \mathrm{argmax}_{\vartheta} \frac{ p(\mathcal{X})p(\theta) }{p(\mathcal{X})} \\
&= \ \mathrm{argmax}_{\vartheta} p(\mathcal{X}|\vartheta)p(\vartheta) \\
&= \ \mathrm{argmax}_{\vartheta} \{\mathcal{L}(\vartheta | \mathcal{X}) + \log p(\vartheta)\} \\
&= \ \mathrm{argmax}_{\vartheta} \{ \sum_{x\in \mathcal{X}} \log p(x|\vartheta) + \log p(\vartheta) \}
 \ \ \ \ \ \ \ (7)
\end{aligned}

和公式（2）比较我们可以看出极大后验估计是先验加上似然的结果。在实际中，先验$$p(\vartheta)$$是用来组织额外的知识来避免过拟合的结果。这也称作是奥卡姆剃刀(Occam's razor)

在引入先验$$p(\vartheta)$$之后，极大后验估计是服从贝叶斯的数据建模方法的，其中$$\vartheta$$可以理解为随机变量的统计量，其先验就是以它们自己作为参数的，即$$p(\vartheta) := p(\vartheta|\alpha)$$，$$\alpha$$是超参数。

极大后验估计可以用类似公式（3）来求解，我们只需要保持$$\mathcal{L}(\vartheta|\mathcal{X}+\log p(\vartheta))$$最大即可。与上面的内容相似，对于新的观测值，其近似等于：

p(\tilde{x} | \mathcal{X}) \approx \int _{\theta\in \vartheta} p(\tilde{x}|\hat{\vartheta}_{\mathrm{MAP}})p(\vartheta|\mathcal{X}) \mathrm{d} \vartheta = p(\tilde{x}|\hat{\vartheta}_{\mathrm{MAP}})
 \ \ \ \ \ \ \ (8)

3、贝叶斯推断（Bayesian inference）

贝叶斯推断扩展了极大后验方法，它把参数集合$$\vartheta$$当作是一个分布，而不是直接估计$$\vartheta$$的值。它不仅考虑了产生这个极大值的数据参数，也将期望、方差作为另外的参数引入用来衡量估计的质量和可信度。这里的意思就是它把参数当作分布而不是某个具体值来对待。比如假设投硬币正面向上的次数服从概率为p的二项式分布。那么贝叶斯推断认为这个p应该服从某个一分布，而不是具体值，我们可以算出p的分布，然后用p的均值或者其他统计量来做预测或者运算。贝叶斯推断的主要步骤是，首先要写出后验;

p(\vartheta|\mathcal{X}) = \frac { P(\mathcal{X}|\vartheta) \cdot p(\vartheta)}{p(\mathcal{X})}
 \ \ \ \ \ \ \ (9)

由于我们并没有限制一定要找最大值，所以我们需要计算公式（8）的分母：

p(\mathcal{X}) = \int_{\vartheta \in \Theta} p(\mathcal{X}|\vartheta)p(\vartheta) \mathrm{d} \vartheta
 \ \ \ \ \ \ \ (10)

对于预测问题，贝叶斯推断也扩展了极大后验推断&你好：

\begin{aligned}
p(\tilde{x}|\mathcal{X}) &= \int_{\vartheta\in \Theta}p(\vartheta|\mathcal{X}) \mathrm{d} \vartheta
 \ \ \ \ \ \ \ (11) \\
&= \int _{\vartheta \in \Theta} p(\tilde{x}|\vartheta)\frac{P(\mathcal{X}|\vartheta)p(\vartheta)}{p(\mathcal{X})}  \ \ \ \ \ \ \ (12)

\end{aligned}

这里的后验$$p(\vartheta|\mathcal{X})$$替换了之前具体的$$\vartheta$$积分掉之后，先验就自动被引入到预测中了，它就是$$\tilde{x}$$的分布，可以通过信度被再一次分析利用。值。

4、极大似然估计、极大后验估计和贝叶斯估计的实例

假设我们扔一个不均匀的硬币，其正面向上是事件c=1，概率是p，未知。那么正面向下的事件为c=0，其概率为1-p。假设我们实验了20次，12次正面向上，问这个硬币正面向上事件发生的概率是多少？

首先，一次试验，随机变量C（正面向上或者正面向下）的密度函数为：

p(C=c|p) = p^c (1-p)^{1-c} \cong \mathrm{Bern}(c|p)  \ \ \ \ \ \ \ \ (13)

假设我们有N次试验结果。使用极大似然估计，其似然函数的log为：

\begin{aligned}
\mathcal{L} &=  \log \prod_{i=1}^{N}p(C=c_i | p)=\sum_{i=1}^N \log p(C=c_i |p)  \ \ \ \ \ \ \ \ (14)\\
&= n^{(1)}\log p(C=1|p) + n^{(0)}\log p(C=0|p)\\
&= n^{(1)} \log p + n^{(0)}\log (1-p)  \ \ \ \ \ \ \ \ (15)
\end{aligned}

其中$$n^{(c)}$$表示一个伯努利实验中时间c产生的次数。求参数p的微分，我们有：

\frac{ \partial \mathcal{L}} { \partial p} = \frac{n^{(1)}}{p} - \frac{n^{(0)}}{1-p}=0 \Leftrightarrow \hat{p}_{ML} = \frac{n^{(1)}}{n^{(1)}+n^{(0)}} = \frac{n^{(1)}}{N}   \ \ \ \ \ \ \ \ (16)

也就是说，极大似然估计的结果就是某个事件发生的次数除以事件总数。

让我们再使用极大后验估计来思考这个问题。假设我们认为这个硬币是均匀分布的。那么也就是说，我们认为p的先验大约是0.5。假设我们取beta作为p的先验：

p(p|\alpha,\beta) = \frac{1}{\mathrm{\bold{B}}(\alpha,\beta)} p^{\alpha-1} (1-p)^{(\beta-1)} \cong \mathrm{\bold{Beta}}(p|\alpha,\beta) \ \ \ \ \ \ \ \ (17)

其中，$$\mathrm{\bold{Beta}}(\alpha,\beta)=\frac{\Gamma{\alpha}\Gamma{\beta}}{\Gamma(\alpha+\beta)}$$，$$\Gamma(x)$$是Gamma函数，这个可以理解为阶乘在实数上的推广。Beta分布的区间在[0,1]之间，可以用来产生正规化概率值。

在我们的例子中，我们认为硬币是均匀的，因此取$$\alpha=\beta=5$$，使得该分布的众数为0.5，进而p可以取0.5。根据前面的极大后验估计公式，我们得到估计公式为：

\begin{aligned}
\frac{\partial}{\partial p}{\mathcal{L}}+\log p(p) &= \frac{n^{(1)}}{p} - \frac{n^{(0)}}{1-p}+\frac{\alpha-1}{p}-\frac{\beta-1}{1-p}=0  \ \ \ \ \ \ \ \ (18) \\
&\Leftrightarrow \hat{p}_{\mathrm{MAP}} = \frac{n^{(1)}+\alpha-1}{n^{(1)}+n^{(0)}+\alpha+\beta-2} = \frac{n^{(1)}+4}{(n^{1})+n^{(0)}+8} \ \ \ \ \ \ \ \ (19)
\end{aligned}

这里的结果有两个很有意思的地方。第一个是关于$$n^{(c)}$$的估计结果发生了改变：由于分子和分母都加了一个值（4和8），使得其结果趋向于0.5。超参数$$\alpha$$和$$\beta$$的值越高，要修正这个先验知识就需要更多的实际观测值。第二点是，计数是来自实际观测值的推导还是超变量先验的推导是不相关的，也就是说，超参数的$$\alpha$$和$$\beta$$也可以理解成物理计数结果。这也就是为什么$$\alpha$$和$$\beta$$被称为伪计数的原因。根据上面试验结果，我们得到了正面向上的概率是0.571小于极大似然估计（因为结果得到了先验的修正）。这就是先验调整的结果。