认识Beta/Dirichlet分布

本文主要是对rickjin的《LDA数学八卦》以及PRML一书中关于Beta分布和Dirichlet分布的知识整理。

我们知道Gamma函数

$\Gamma(x) = \int_0^{\infty}t^{x-1}e^{-t}dt$

是阶乘运算在实数集上延伸。它具有如下性质

$\Gamma(x+1) = x\Gamma(x)$

所以，我们其实有 $\Gamma(n)=(n-1)!$

是不是感觉有些奇怪？为什么不是 $\Gamma(n)=n!$ 而是 $\Gamma(n)=(n-1)!$ ？
今天读了rickjin的《LDA数学八卦》才知道，原来欧拉是研究了Beta函数

$B(m,n) = \int_0^1x^{m-1}(1-x)^{n-1}dx$

之后，他发现，如果Gamma函数的定义选取满足 $\Gamma(n) = (n-1)!$ ,那么Beta函数会有一个很漂亮的对称形式

$B(m,n) = \frac{\Gamma(m)\Gamma(n)}{\Gamma(m+n)}$

而如果选取 $\Gamma(n) = n!$ 的定义，则有

$B'(m,n)= \frac{\Gamma(m)\Gamma(n)}{\Gamma(m+n+1)}$

这个形式显然不如B(m,n)优美，而数学家总是很在乎美感的。

讲完了Gamma函数，我们再来扒一扒Beta函数，他又有什么物理意义呢？

Beta分布

我之前找工作面试的时候，曾经被面过这么一个题目

$X \sim Uniform(0,1)$ ;
随机生成10个数，把这10个数排序后得到的顺序统计量是 $X_1,X_2,...,X_n$ ;
问第7大的数的概率分布?

那时候我是不知道Beta分布，否则肯定不会被虐的这么惨::>_<::

我们先将之一般化，对于一般的情况 $X_k$ 的概率密度是什么呢？下面，我们尝试计算一下 $X_k$ 落在一个区间 $[x,x+\Delta x]$ 的概率值

$P(x \le X_k \le x + \Delta x) = ?$

beta_distribution

如上图所示，我们把[0,1]区间分成三段 $[0,x),[x,x+\Delta x], (x+\Delta x,1]$ 三段。我们假定， $\Delta x$ 足够小，只能够容纳一个点,则由排列组合理论可得

$P(x \le X_k \le x + \Delta x) = {n\choose 1}\Delta x{n-1 \choose k-1}x^{k-1}(1-x-\Delta x)^{n-k}$

所以我们可以得到 $X_k$ 的概率密度函数为

$\begin{align}f(x) &= \lim_{x \to 0} \frac{P(x\le X_k \le x+\Delta x)}{\Delta x}\\ &= {n \choose 1} {n-1 \choose k-1}x^{k-1}(1-x)^{n-k}\\ &= \frac{n!}{(k-1)!(n-k)!}x^{k-1}(1-x)^{n-k}\\ &= \frac{\Gamma(n+1)}{\Gamma(k)\Gamma(n-k+1)}x^{k-1}(1-x)^{n-k}\end{align}$

我们取 $\alpha = k, \beta = n-k+1$ ,于是

$f(x) = \frac{\Gamma(\alpha + \beta)}{\Gamma(\alpha)\Gamma(\beta)}x^{\alpha -1}(1-x)^{\beta -1}$

这就是Beta分布！
回到上面那个面试题，把 $n=10,k=7$ 带入其中，得到密度函数

$f(x) = \frac{10!}{6! \times 3!}x^6(1-x)^3,x \in [0,1]$

Beta-Binomial共轭

上边的面试题还有几个衍化版本，我们先看第一个衍化版本：

$X \sim Uniform(0,1)$ ;
随机生成n个数,由小到大排序后为 $X_1,X_2,...,X_n$ ,我们要猜测第k大的数 $p=X_k$ ;
我们再生成m个数， $Y_1,Y_2,...,Y_m \sim Uniform(0,1)$ , 其中有 $m_i$ 个数比p小， $m_2$ 比p大；
求 $P(p\|Y_1,Y_2,...,Y_m)$ 的分布是什么。

容易看出，我们一共生成了 $m+n$ 个数，而 $p=X_k$ 在最终生成的m+n个数中，是第 $k+m_1$ 大的。按照我们之前讲过的Beta分布的逻辑，p其实应该服从 $\alpha = k+m_1,\beta = n-k+1+m_2$ 的Beta分布。我们知道贝叶斯学派进行参数估计的基本过程是

先验分布 + 后验数据 = 后验分布

对应到Beta分布，后验数据其实相当于是做了m次Bernoulli实验，其中 $m_1$ 次比p小， $m_2$ 次比p大，相当于

$Beta(p\|\alpha,\beta) + BinomCount(m_1,m_2)=Beta(p\|\alpha + m_1,\beta + m_2)$

上面这个式子描述的就是 $Beta-Bonomial共轭$ .

共轭的意思是，参数的先验分布和后验分布都能保持Beta分布的形式，这样的好处是，我们能够在先验分布中赋予参数明确的物理意义，并且这个物理意义可以通过后验数据，延续到后验分布中进行解释。

由上边的解释可知，Beta分布重的参数 $\alpha,\beta$ 其实都可以理解为物理技术，这两个参数也经常被称为伪计数（pseudo-count)。所以， $Beta(\alpha,\beta)$ 也可以理解为

$Beta(\alpha,\beta) = Beta(1,1) + BinomCount(\alpha -1,\beta-1)$

其中

$Beta(1,1)=\frac{\Gamma(1+1)}{\Gamma(1)\Gamma(1)}x^{1-1}(1-x)^{1-1}=1$

这恰好就是均匀分布 $Uniform(0,1)$ 。

贝叶斯学派和频率学派的不同
>假设有一个不均匀的硬币，抛出正面的概率为p，抛掷 $m$ 次后，出现正面和翻面的次数分别为 $m_1$ 和 $m_2$ ，那么按照传统频率学派的观点，p的估计值应该为 $\hat p=\frac{m_1}{m}$ ,而从贝叶斯学派的观点来看，开始对硬币不均匀性一无所知，所以应该假设p服从均匀分布 $Uniform(0,1)$ ,也就是 $Beta(1,1)$ ,于是在有了后验数据之后，我们得出p其实应该服从 $Beta(p|m_1 +1,m_2 +1)$ .

百变星君Beta分布

beta

Beta分布的概率密度如上图， $\alpha,\beta$ 的不同，他可以是凹的、凸的、单调上升的、单调下降的，可以是曲线也可以是直线；而且，如前所述，均匀分布也是Beta分布的一种特殊形式。正是由于Beta分布能够你和如何之多的形状，因此他经常被贝叶斯学派用作先验分布。

Dirichlet-Multinomial共轭：Beta分布的高维推广

更一步的问题

$X \sim Uniform(0,1)$ ;
随机生成n个数，排序后为 $X_1,X_2,...,X_n$ ;
求 $(X_{k_1},X_{k_2})$ 的联合分布。

同推导Beta分布类似,我们取 $\Delta x$ 足够小，只能容纳一个点.

Dirichlet

由于 $\Delta x$ 足够小，我们有 $x_1+x_2+x_3=1$ .

$\begin{align}&P(X_{k_1} \in (x_1,x_1+\Delta x),X_{k_1+k_2}\in(x_2,x_2+\Delta x))\\&={n\choose 1}{n-1 \choose 1}{n-2 \choose {k_1 -1,k_2 -1}}x_1^{k_1 -1}x_2^{k_2 -1}x_3^{n-k_1 - k_2}(\Delta x)^2\end{align}$

于是我们得到 $(X_{k_1},X_{k_2})$ 的联合分布为

$\begin{align}f(x_1,x_2,x_3) &= \frac{n!}{(k_1 -1)!(k_2 -1)!(n-k_1 -k_2)!}x_1^{k_1 -1}x_2{k_2 -1}x_3^{n-k_1 -k_2}\\&=\frac{Gamma(n+1)}{\Gamma(k_1)\Gamma(k_2)\Gamma(n-k_1 -k_2 +1)}x_1^{k_1 -1}x_2^{k_2 -1}x_3^{n-k_1 -k_2}\end{align}$

令 $\alpha_1=k_1,\alpha_2=k_2,\alpha_3=n-k_1-k_2+1$ ,我们得到

$f(x_1,x_2,x_3)=\frac{\Gamma(\alpha_1 +\alpha_2 +\alpha _3)}{\Gamma(\alpha_1)\Gamma(\alpha_2)\Gamma(\alpha_3)}x_1^{\alpha_1 -1}x_2^{\alpha_2 -1}x_3^{\alpha_3 -1}$

上边这个分布其实就是一个三维形式的Dirichlet分布 $Dir(\alpha_1,\alpha_2,\alpha_3)$ .同Beta分布类似，Dirichlet分布也是一个百变星君，下图为不同 $\alpha$ 值时Dirichlet分布的图像。

Dirichlet_distribution_pic

一般形式的Dirichlet分布定义如下

$Dir(\vec p\|\vec \alpha)=\frac{\Gamma(\sum_{k=1}^K\alpha_k)}{\prod_{k=1}K\Gamma(\alpha_k)}\prod_{k=1}^Kp_k^{\alpha_k -1}$

Dirichlet分布也是Binomial共轭的

$Dir( \vec p \| \vec \alpha) +MultCount(\vec m)=Dir(\vec p\| \vec \alpha + \vec m)$

我们同样也有

$Dir(\vec p\|\vec \alpha)=Dir(\vec p\|\vec 1)+MultCount(\vec m - \vec 1)$

Beta分布和Dirichlet分布的性质

如果 $p \sim Beta(t\|\alpha,\beta)$ ,则

$\begin{align}E(p)&=\int_0^1 t*Beta(t\|\alpha,\beta)dt\\&=\int_0^1 t*\frac{\Gamma(\alpha + \beta)}{\Gamma(\alpha)\Gamma(\beta)} t^{\alpha-1}(1-t)^{\beta -1}dt\\&=\frac{\Gamma(\alpha + \beta)}{\Gamma(\alpha)\Gamma(\beta)}\int_0^1 t^\alpha(1-t)^{\beta -1}dt\end{align}$

上式右边的积分对应到概率分布 $Beta(t\|\alpha +1,\beta)$

$Beta(t\|\alpha +1,\beta)=\int_0^1 t*\frac{\Gamma(\alpha + \beta +1 )}{\Gamma(\alpha +1)\Gamma(\beta)} t^\alpha(1-t)^{\beta -1}dt$

而且我们有

$\int_0^1Beta(t\|\alpha +1,\beta)dt=1$

所以我们有

$\int_0^1 t^\alpha(1-t)^{\beta -1}dt=\frac{\Gamma(\alpha+1)\Gamma(\beta)}{\Gamma(\alpha+\beta+1)}$

把上式带入E(p)中得

$E(p)=\frac{\Gamma(\alpha + \beta)}{\Gamma(\alpha)\Gamma(\beta)}\cdot\frac{\Gamma(\alpha+1)\Gamma(\beta)}{\Gamma(\alpha+\beta+1)}=\frac{\alpha}{\alpha + \beta}$

同样的，对于Dirichlet分布我们可以得到

$E(\vec p)=(\frac{\alpha_1}{\sum_{i=1}{K}\alpha_i},\frac{\alpha_2}{\sum_{i=1}{K}\alpha_i},...,\frac{\alpha_K}{\sum_{i=1}{K}\alpha_i})$