贝塔-伯努利模型(Beta-Binomial Model)

作者: 引线小白-本文永久链接:http://www.limoncc.com/机器学习/2017-03-07-机器学习笔记03/
知识共享许可协议: 本博客采用署名-非商业-禁止演绎4.0国际许可证

一、贝塔-伯努利模型(beta-binomial model)

伯努利分布:

$$\begin{align}
x\sim \mathrm{Ber}(\mu), x\in\{0,1\},\mu\in[0,1]
\end{align}$$

1、0-1分布(bernoulli distribution)概率质量函数可以表示为:

$$\begin{align}
\mathrm{PMF}:\mathrm{Ber}(x\mid \mu)=\mu^x(1-\mu)^{1-x}
\end{align}$$
其他表示方法:
$\displaystyle \mathrm{Ber}(x\mid \mu)=\mu^{\mathbb{I}(x=1)}(1-\mu)^{\mathbb{I}(x=0)} $

$\displaystyle \mathrm{Ber}(x\mid \mu)=\begin{cases}\mu&\text{if }x=1\\\nu=1-\mu&\text{if }x=0\end{cases} $

2、均值与方差:

知道 $\displaystyle \varphi_x(t)=\mu\mathrm{e}^{\mathrm{i}t}+\nu$

$\displaystyle \mathrm{E}[x]=(-\mathrm{i})^1\left.\frac{\partial{\varphi(t)}}{\partial{t}}\right| _{t=0}=\mu$

$\displaystyle \mathrm{var}[x]=(-\mathrm{i})^2\left.\frac{\partial^2{\varphi(t)}}{\partial{t}^2}\right| _{t=0}-\mathrm{E}^2[x]=\mu(1-\mu)=\mu\nu$

3、似然函数(likelihood)

有数据集: $\displaystyle \mathcal{D}=\{x_i\}_{i=1}^{n}$于是有似然函数

$$\begin{align}
\mathrm{L}(\mu)=p(\mathcal{D}\mid\mu)=p(\boldsymbol{X}\mid\mu)=\prod_{i=1}^{n}\mu^{x_i}(1-\mu)^{1-x_i}=\mu^{N_1}(1-\mu)^{N_0}
\end{align}$$

其中 $\displaystyle N_1=\sum_{i=1}^{n}\mathbb{I}(x_i=1)$, $\displaystyle N_0=\sum_{i=1}^{n}\mathbb{I}(x_i=0)$

4、对数似然函数

$$\begin{align}
\mathcal{L}(\mu)
&=\ln p(\mathcal{D}\mid\mu)=\ln \prod_{i=1}^{n}\mu^{x_i}(1-\mu)^{1-x_i}=\sum_{i=1}^{n}\left(x_i\ln\mu+(1-x_i)\ln(1-\mu)\right)\\
&=\boldsymbol{I}^\text{T}\boldsymbol{x}\ln\mu+\boldsymbol{I}^\text{T}(\boldsymbol{I}-\boldsymbol{x})\ln(1-\mu)
\end{align}$$

5、求极大似然估计:

$$ \frac{\partial{\mathcal{L}}}{\partial{\mu}}=\frac{1}{\mu}\boldsymbol{I}^\text{T}\boldsymbol{x}-\frac{1}{1-\mu}\boldsymbol{I}^\text{T}(\boldsymbol{I}-\boldsymbol{x})=0\\
\mu\boldsymbol{I}^\text{T}\boldsymbol{I}=\boldsymbol{I}^\text{T}\boldsymbol{x}
$$$\begin{align}
\mu_{MLE}=\frac{\boldsymbol{I}^\text{T}\boldsymbol{x}}{\boldsymbol{I}^\text{T}\boldsymbol{I}}=\frac{1}{n}\sum_{i=1}^{n}x_i=\bar{x}
\end{align}$
极大似然估计分析:
我们知道 $\displaystyle x\in\{0,1\}$,在0-1分布的n次试验中, $\displaystyle x=1$的次数为 $\displaystyle k$。于是
$$
\begin{align}
\mu_{MLE}=\frac{1}{n}\sum_{i=1}^{n}=\frac{k}{n}
\end{align}
$$
现在我们假设抛硬币3次,3次都是正面朝上,那么 $\displaystyle n=k=3$,那么 $\displaystyle \mu_{MLE}=1$。这种情况下,极大似然的结果预测所有未来的观察值都是正面向上。常识告诉我们这是不合理的。事实上,这就是极大似然中过拟合现象的一个极端例子。下面我们引入 $\displaystyle \mu$的先验分布,我们会得到一个更合理的结论。

6、二项式分布(binomial distribution)

我们现在把 $\displaystyle k$作为随机变量 $\displaystyle k=x_1+x_2+…+x_n$,于是我们有 $\displaystyle k\in\{N_1\}$,易知二项式分布和0-1分布的似然函数有相同的形式,是正比关系。最大化它们的 $\displaystyle \mu$都是 $\displaystyle \frac{k}{n}$。

$$\begin{align}
\mathrm{L}(\mu)\propto\mathrm{Bin}(k\mid \mu,n)=\mathrm{C}_{n}^k\,\mu^k(1-\mu)^{n-k}
\end{align}$$

为了与后面符号衔接:我们令 $\displaystyle k=N_1,n=N,N_0=n-k$于是又有:
$$\begin{align}
p(\mathcal{D}\mid\mu)=\mathrm{L}(\mu)\propto\mathrm{Bin}(N_1\mid \mu,N_1+N_0)=\mathrm{C}_{N}^{N_1}\,\mu^{N_1}(1-\mu)^{N_0}
\end{align}$$

我们知道 $\displaystyle \mathrm{C}_{n}^{k}=\frac{n!}{(n-k)!k!}$。
同时我们有离散随机变量特征函数$\displaystyle \varphi(t)=\mathrm{E}[\mathrm{e}^{\mathrm{i}t x_i}]=\sum_{i=1}^{\infty}p_i\mathrm{e}^{\mathrm{i}t x_i}$,注意虚数 $\displaystyle \mathrm{i}$和变量$\displaystyle i$的区别。我们令 $\displaystyle \nu=1-\mu$于是有:
$$\begin{align}
\varphi_k(t)=(\mu\mathrm{e}^{\mathrm{i}t}+\nu)^n
\end{align} $$
$\displaystyle \mathrm{E}[k]=\sum_{k=0}^{n}k\mathrm{Bin}(k\mid \mu,n)=(-\mathrm{i})^1\left.\frac{\partial{\varphi(t)}}{\partial{t}}\right| _{t=0}=n\mu$
$\displaystyle \mathrm{var}[k]=\sum_{k=0}^{n}(k-\mathrm{E}[k])^2\mathrm{Bin}(k\mid \mu,n)=(-\mathrm{i})^2\left.\frac{\partial^2{\varphi(t)}}{\partial{t}^2}\right| _{t=0}-\mathrm{E}^2[k]=n\mu(1-\mu)=n\mu\nu$

7、共轭先验分布

如果先验和似然函数有相同形式,那就非常方便:这样后验也有相同形式。这个时候我们称之为共轭先验(conjugate prior)。现在,我们把 $\displaystyle \mu$做为一个变量观察似然函数得到
$$\begin{align}
p(\mu)\propto\mu^{\rho_1}(1-\mu)^{\rho_2}
\end{align}$$
这样计算后验就很容易了:

$\displaystyle p(\mu\mid\mathcal{D})\propto p(\mathcal{D}\mid \mu)p(\mu)=\mu^{N_1}(1-\mu)^{N_0}\mu^{\rho_1}(1-\mu)^{\rho_2}=\mu^{N_1+\rho_1}(1-\mu)^{N_0+\rho_2}$

在给出先验分布之前,我们看看这个两个函数。之后我们会补充 $\displaystyle \Gamma(x),\mathrm{B}(a,b)$函数的性质的证明。
$\displaystyle \Gamma(x)=\int_0^\infty u^{x-1}\mathrm{e}^{-u}\mathrm{d}u, x>0$
$\displaystyle \mathrm{B}(a,b)=\int_{0}^{1}x^{a-1}(1-x)^{b-1}\mathrm{d}x,a>0,b>0$
$\displaystyle \Gamma(x+1)=x\Gamma(x)$
$\displaystyle \Gamma(n+1)=n!$
$\displaystyle \mathrm{B}(a,b)=\mathrm{B}(b,a)$
$\displaystyle \mathrm{B}(a,b)=\frac{b-1}{a+b-1}\mathrm{B}(a,b-1),a>0,b>1$

$\displaystyle \mathrm{B}(a,b)=\frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)}$

然后我们知道贝塔分布: $\displaystyle \mathrm{Beta}(\mu\mid a,b)=\frac{\Gamma(a+b)}{\Gamma{a}\Gamma(b)}\mu^{a-1}(1-\mu)^{b-1} $。这正是我们要的共轭先验。

$$\begin{align}
\mathrm{Beta}(\mu\mid a,b)=\frac{1}{\mathrm{B}(a,b)}\mu^{a-1}(1-\mu)^{b-1}
\end{align}$$利用 $\displaystyle \Gamma$函数性质易知:
$\displaystyle \mathrm{E}[\mu]=\frac{a}{a+b} $
$\displaystyle \mathrm{var}[\mu]=\frac{ab}{(a+b)^2(a+b+1)}$
$\displaystyle \mathrm{mode}[\mu]=\mathop{\mathrm{argmax}}_{\mu}\,\mathrm{Beta}(\mu\mid a,b)=\frac{a-1}{a+b-2},a>0,b>0$

8、后验分布

我们用似然函数乘以贝塔先验得到后验分布:
$$\begin{align}
p(\mu\mid\mathcal{D})
\propto\mathrm{Bin}(N_1\mid\mu, N_1+N_0)\mathrm{Beta}(\mu \mid a,b)\end{align}$$ 归一化得:
$$\begin{align}
p(\mu\mid\mathcal{D})=\mathrm{Beta}(\mu\mid N_1+a,N_0+b)
\end{align}$$

1、在线学习

容易证明 $\displaystyle \mathrm{Beta}$分布具有再生性质。于是我们可以假设我们有两个数据集 $\displaystyle \mathcal{D}_1,\mathcal{D}_2$。于是有:
1、当我们观察到 $\displaystyle \mathrm{D}_1$时,有
$$p(\mu\mid\mathcal{D}_1)=\mathrm{Beta}(\mu\mid N_1^1+a,N_0^1+b)$$
2、当我们观察到 $\displaystyle \mathcal{D}_2$时,有
$$p(\mu\mid\mathcal{D}_1,\mathcal{D}_2)\propto p(\mathcal{D}_2\mid\mu) p(\mu\mid\mathcal{D}_1)$$
易得:
$$\begin{align}
p(\mu\mid\mathcal{D}_1,\mathcal{D}_2)=\mathrm{Beta}(\mu\mid N_1^1+N_1^2+a,N_0^1+N_0^2+b)
\end{align}$$
这表明 该贝叶斯推断具有在线学习的良好性质。

2、后验均值与众数

最大后验估计:
$$\begin{align}
\mu_{MAP}=\mathrm{mode}[\mu]=\mathop{\mathrm{argmax}}_{\mu}\,p(\mu\mid\mathcal{D})=\frac{N_1+a-1}{N+a+b-2}=\frac{k+a-1}{n+a+b-2}
\end{align}$$
极大似然估计:

$$\begin{align}
\mu_{MLE}=\frac{N_1}{N}=\frac{k}{n}
\end{align}$$

后验均值:

$$\begin{align}
\mathrm{E}[\mu\mid\mathcal{D}]=\frac{N_1+a}{N+a+b}=\frac{k+a}{n+a+b}\end{align}$$
我们发现众数和均值不同。
我们还发现 后验均值是先验均值和最大似然估计的凸组合。下面我们下证明这一点。 令先验均值$\displaystyle \frac{a}{a+b}=\mu_0,a+b=c$。

$$\begin{align}
\mathrm{E}[\mu\mid\mathcal{D}]=\frac{k+a}{n+a+b}=\frac{c}{n+c}\mu_0+\frac{n}{n+c}\mu_{MLE}=\lambda\mu_0+(1-\lambda)\mu_{MLE}
\end{align}$$

我们发现 $\displaystyle c$可以理解为先验对于后验的等价样本大小。 $\displaystyle \lambda=\frac{c}{n+c}$是先验对于后验的等价样本大小比率。

同理可知最大后验估计是先验众数和最大似然估计的凸组合,而且更加倾向于最大似然估计:
$$\begin{align}
\mu_{MAP}=\frac{c-2}{n+c-2}\mathrm{mode}[\mu_0]+\frac{n}{n+c-2}\mu_{MLE}=\eta\mathrm{mode}[\mu_0]+(1-\eta)\mu_{MLE}
\end{align}$$

3、后验方差

$$\begin{align}
\mathrm{var}(\mu\mid\mathcal{D}]=\frac{(N_1+a)(N_0+b)}{(N_1+a+N_0+b)^2(N_1+a+N_0+b+1)}
\end{align}$$
当 $\displaystyle N$足够大时:
$$\begin{align}
\mathrm{var}(\mu\mid\mathcal{D}]\approx\frac{N_1N_0}{NNN}=\frac{\mu_{MLE}(1-\mu_{MLE})}{n}
\end{align}$$
后验标准差

$$\begin{align}
\sigma=\sqrt{\frac{\mu_{MLE}(1-\mu_{MLE})}{n}}
\end{align}$$
1、它随着我们数据的增加以 $\displaystyle \sqrt{\frac{1}{n}}$的速度下降。
2、 $\displaystyle \mu=\mathop{\mathrm{argmax}}_{\mu}\sigma(\mu)=0.5$
3、 $\displaystyle \mu=\mathop{\mathrm{argmin}}_{\mu}\sigma(\mu)=0\,or\,1$

4、后验预测分布(posterior predictive distribution)

上述,我们讨论了未知参数的推断,现在我们来讨论预测问题。在这之前:我们令后验分布 $\displaystyle p(\mu\mid\mathcal{D})=\mathrm{Beta}(a,b)$以简化符号 。考虑 $\displaystyle x_{n+1}$发生了,那么我们想预测 $\displaystyle x_{n+1}=0\,or\,1$的概率。为了简记,令 $\displaystyle \tilde{x}$为一随机变量,代表 $\displaystyle n$次观测后的值。同时考虑到 $\displaystyle \tilde{x}\in\{0,1\}$于是有:
$$\begin{align}
p(\tilde{x}\mid \mathcal{D})
&=\int_0^1p(x\mid \mu)p(\mu\mid\mathcal{D})\mathrm{d}\mu\\
&=\int_0^1\mu^\tilde{x}(1-\mu)^{1-\tilde{x}}\frac{\Gamma(a+b)}{\Gamma{a}\Gamma(b)}\mu^{a-1}(1-\mu)^{b-1}\mathrm{d}\mu\\
&=\frac{\Gamma(a+b)}{\Gamma{a}\Gamma(b)}\int_0^1\mu^{a+\tilde{x}-1}(1-\mu)^{b+1-\tilde{x}-1}\mathrm{d}\mu\\
&=\frac{\Gamma(a+b)}{\Gamma{a}\Gamma(b)}\frac{\Gamma(a+\tilde{x})\Gamma(b+1-\tilde{x})}{\Gamma(a+b+1)}\\
&=\frac{a^{\tilde{x}}b^{1-\tilde{x}}}{(a+b)^{\tilde{x}}(a+b)^{1-\tilde{x}}}\\
&=\left(\frac{a}{a+b}\right)^{\tilde{x}}\left(\frac{b}{a+b}\right)^{1-\tilde{x}}\\
&=\mathrm{Ber}\left(\tilde{x}\mid\mathrm{E}[\mu\mid\mathcal{D}]\right)
\end{align}$$

所以可以看到,后验预测分布等价于plug-in$\displaystyle \mathrm{E}[\mu\mid\mathcal{D}]$。这里的plug-in x是插值近似 plug-in approximation的意思。

5、过拟合与黑天鹅悖论

如果我们使用 $\displaystyle plug-in\,\mu_{MLE}$:
$$\begin{align}
p(\tilde{x}\mid \mathcal{D})\approx\mathrm{Ber}\left(\tilde{x}\mid\mu_{MLE}\right)
\end{align}$$
1、当数据集较小,例如 $\displaystyle n=3,k=0,then,\mu_{MLE}=\frac{0}{3}=0$。这样我们预测 $\displaystyle \tilde{x}=0$的概率为 $\displaystyle 0$。这叫零数问题(zero count problem),或者叫数据匮乏问题(sparse data problem)。这种问题在小样本情况,经常发生。当然有人就要说了,现在是大数据时代,干嘛关注这个问题?如果我们基于特定的目的,划分数据:例如特定的人做特别的事——个人购物推荐。这种情况下即使是在大数据时代,数据也不会很大。所以 $$即使在大数据当道的时代,贝叶斯方法依然是有用的![Jordan2011]$$
2、零数问题类似于黑天鹅事情。
3、我们来用贝叶斯方法来解决这个问题: 使用均匀分布先验 $\displaystyle \mathrm{Beat}(1,1)=\mathrm{U}[0,1]$,plug-$\displaystyle in\,\mathrm{E}[\mu\mid\mathcal{D}]$给出了拉普拉斯继承规则(laplace’s rule of succession):
$$\begin{align}
p(\tilde{x}=1\mid\mathcal{D})=\mathrm{E}[\mu\mid\mathcal{D}]=\frac{k+1}{n+2}
\end{align}$$
可以理解为:在实践中,通常在数据集中加一应该是正确的做法。称之为加一平滑(add-one smoothing)
4、注意到:plug-in$\displaystyle \mu_{MAP}$没有这个性质。

6、多试验预测

考虑这个问题,我们又观察到了 $\displaystyle m$个数据,那么发生 $\displaystyle \tilde{x}=1$的次数是 $\displaystyle s$次的概率。
$$\begin{align}
p(s\mid\mathcal{D},m)
&=\int_0^1\mathrm{Bin}(s\mid\mu,m)\mathrm{Beta}(\mu\mid a.b)\mathrm{d}\mu\\
&=\frac{\mathrm{C}_m^s}{\mathrm{B}(a.b)}\int_0^1\mu^{s+a-1}(1-\mu)^{m-s+b-1}\mathrm{d}\mu\\
&=\mathrm{C}_m^s\frac{\mathrm{B}(s+a,m-s+b)}{\mathrm{B}(a.b)}
\end{align}$$
我们称 $\displaystyle \mathrm{Bb}(s\mid a,b,m)=\mathrm{C}_m^s\frac{\mathrm{B}(s+a,m-s+b)}{\mathrm{B}(a.b)}$为贝塔-二项式分布(beta-binomial distribution)。易知(使用积分表达式,同时利用 $\displaystyle \mathrm{Bin}(s\mid \mu,m)$容易求得):
$\displaystyle \mathrm{E}[s\mid\mathcal{D},m]=m\frac{a}{a+b}$

$\displaystyle \mathrm{var}[s\mid\mathcal{D},m]=\frac{mab}{(a+b)^2}\frac{a+b+m}{a+b+1}$

$\displaystyle \mathrm{var}[s\mid\mu_{MAP}]=\frac{m(a-1)(b-1)}{(a+b-2)^2}$(插值求得的)
于是有:
$$\begin{align}
\mathrm{var}[s\mid\mathcal{D},m]\geqslant\mathrm{var}[s\mid\mu_{MAP}]
\end{align}$$
至于证明,可以分析这个式子: $\displaystyle f(a,b)=ab \left( a+b-2 \right) ^{2} \left( a+b+m \right) - \left( a-1
\right) \left( b-1 \right) \left( a+b \right) ^{2} \left( a+b+1
\right) \geqslant 0,(a>0,b>0)
$可以使用微积分一阶导数,二阶导数证明。所以 贝叶斯预测的概率密度更宽、拥有长尾,从而避免了过拟合和黑天鹅悖论$$

9、评述

通过贝塔-伯努利模型,我们熟悉了贝叶斯方法的基本概念、流程、特点。一旦我们熟悉了伯努利、二项、贝塔、均匀、贝塔-伯努利分布后,很多关键的东西就可以用文字表述清楚了。

1、上帝给 $\displaystyle x$选择了一个分布:$\displaystyle x\sim \mathrm{Ber}(\mu), x\in\{0,1\},\mu\in[0,1] =\mu^x(1-\mu)^{1-x} $ 这个分布由 $\displaystyle \mu$确定。

2、犹豫种种原因,人类也知晓了 $\displaystyle x$的分布类型,但是不知道 $\displaystyle \mu$。于是人类搞了假设空间 $\displaystyle \mathcal{H}=\{\mu_i\}$,为了找到上帝的那个 $\displaystyle \hat{\mu}$,于是人类的征程开始了。

3、人类的智者,选择了用概率来解决这个难题。因为我们可以根据我们在某个时刻,不断知晓的信息流 $\displaystyle \mathcal{D}_t,\mathcal{D}_{t+1}$,来给假设空间的每个 $\displaystyle \mu$更新概率。于是这个表达式横空出世 $\displaystyle p(\mu\mid \mathcal{D}_{t+1},\mathcal{D}_t)\propto p(\mathcal{D}_{t+1}\mid\mu)p(\mu\mid \mathcal{D}_{t})$。于是我们得到了:
$$p(\mu\mid\mathcal{D}_t)=\mathrm{Beta}(\mu\mid N_1^t+a,N_0^t+b) $$
这样我们就把假设空间的 $\displaystyle \mu$的都给了个概率。至于人类相信那一个。贝叶斯方法说,基于不同的目的,我们要重新开发一门学问,叫决策论。

4、上帝微微一笑, 人类继续开始征程:令 $\tilde{x}$为一随机变量,代表 n次观测后的值,于是得到:
$$ p(\tilde{x}\mid \mathcal{D}_t)=\int_0^1p(x\mid \mu)p(\mu\mid\mathcal{D}_t)\mathrm{d}\mu=\mathrm{Ber}\left(\tilde{x}\mid\mathrm{E}[\mu\mid\mathcal{D}_t]\right)$$
于是基于对假设空间赋概,我们对未来的可能值也赋概。至于人类相信那一个。贝叶斯方法说,基于不同的目的,我们要重新开发一门学问,叫决策论。

5、上帝开始感觉人类是个聪明的物种,甚为满意! 我们又观察到了 m个数据,那么发生 $\tilde{x}=1$的次数是$s$次的概率
$$\begin{align}
p(s\mid\mathcal{D}_t,m)
&=\int_0^1\mathrm{Bin}(s\mid\mu,m)p(\mu\mid\mathcal{D}_t)\mathrm{d}\mu\\
&=\int_0^1\mathrm{Bin}(s\mid\mu,m)\mathrm{Beta}(\mu\mid N_1^t+a,N_0^t+b)\mathrm{d}\mu\\
&=\mathrm{Bb}(s\mid N_1^t+a,N_0^t+b,m)
\end{align}$$
这样我们对 $\displaystyle s$的可能值也赋概了。至于人类相信那一个。贝叶斯方法说,基于不同的目的,我们要重新开发一门学问,叫决策论。

6、至此,当 $\displaystyle \lim_{t\to\infty}\mathcal{D}_t$人类基本可以看到上帝裸体了。

7、靠那个贝叶斯是谁,怎么总要出来嘟囔一下。决策论是啥?被冷落的贝叶斯微微一笑。


版权声明
引线小白创作并维护的柠檬CC博客采用署名-非商业-禁止演绎4.0国际许可证。
本文首发于柠檬CC [ http://www.limoncc.com ] , 版权所有、侵权必究。
本文永久链接http://www.limoncc.com/机器学习/2017-03-07-机器学习笔记03/

予汝玫瑰,渡人沃土。