Normal Distribution

正态分布介绍

神说,要有正态分布,就有了正态分布。

神看正态分布是好的,就让随机误差就服从了正态分布. 然而,对正太分布推动最大的是中心极限定理.设随机变量 \(X_n(n=1,2,\cdots)\) 服从参数为 \(p\) 的二项分布,则对任意的 \(x\), 恒有

$$\lim_{n\rightarrow\infty}P\{ \frac{X_n – np}{\sqrt{np(1-p)}} \le x \} \int_{-\infty}^x \frac{1}{\sqrt{2\pi}} e^{\frac{-t^2}{2}}dt$$

中心极限定理后续又被其它数学家们推广到了其它任意分布的情形,而不限于二项分布。后续的统计学家发现,一系列的重要统计量,在样本量 N 趋于无穷的时候, 其极限分布都有正态的形式, 这构成了数理统计学中大样本理论的基础 .除了中心界限定理的推动之外,高斯分布自己的一些特殊性质也使得该分布非常方便使用。

正态分布与熵

熵是混乱程度或者不确定度的度量。对于离散型分布来说,我们用\(H_n(p_1,p_2,\cdots,p_n)\)来表示该离散分布的熵的数值度量。对于这个熵函数,我们要求该函数对于任意的\(p_i\)是连续的。否则\(p_i\)的微小变化会引起不确定度的剧烈改变,这与直觉是相违背的。同时在直觉的限定下,我们还有下式

\begin{equation}H_3(p_1,p_2,p_3)=H_2(p_1,p_2+p_3)+(p_2+p_3)H_2(\frac{p_2}{p_2+p_3},\frac{p_3}{p_2+p_3})\end{equation}

这个等式相当于决策树的判定过程,此时我们把\(p_2,p_3\)分在一组作为顶层判断,然后再去细分。由这个公式,我们可以推演到下面这个公式

\begin{equation}H(p_1,\cdots,p_n)=H(w_1,\cdots,w_r)+w_1H(\frac{p_1}{w_1},\cdots,\frac{p_k}{w_1})+w_2H(\frac{p_{k+1}}{w_2},\cdots,\frac{p_{k+m}}{w_2})+\cdots\end{equation}

其中,每一个\(w_i\)可以看作一组\(p_j\)的聚合。当\(n\)固定时,等概率分布不确定度最大,我们定义

\begin{equation}h(n)=H_n(\frac{1}{n},\frac{1}{n},\cdots,\frac{1}{n})\end{equation}

又由于在离散分布的情况下\(p_i=\frac{n_i}{\sum{n_i}}\),此时我们可以根据上式来得到

\begin{equation}h(\sum{n_i})=H(p_1,\cdots,p_n)+\sum{p_ih(n_i)}\end{equation}

根据之前的熵计算公式,我们考虑这种情况\(n_i=m\),在这种情况下有下面的结果

\begin{equation} h(m*n)=h(m)+h(n)\end{equation}

很显然,我们可以得出下面这个结论

\begin{equation}h(n)=K ln(n)\end{equation}

其中\(K\)为任意的常数。由于我们只需要知道相对度量,而K的取值并不影响相对度量。此时熵函数也可以解析的表示出来

$$H(p_1,\cdots,p_n)=-K\sum_{i=1}^{n}{p_iln(p_i)}$$

由于我们只需要知道相对度量,而\(K\)的取值并不影响相对度量。我们可以把\(K\)指定为1,此时得到

\begin{equation} H(p_1,\cdots,p_n)=-\sum_{i=1}^{n}{p_iln(p_i)}\end{equation}

这就是脍炙人口的熵计算公式。

虽然我们当前把熵定义在离散型分布上,但是其实熵公式对连续型分布也同样适用,此时熵公式变形为

\begin{equation}H[x]=-\int{p(x)ln(p(x))} dx\end{equation}

这里顺带提一下相对熵。相对熵(relative entropy),也叫做Kullback-Leibler divergence。这个量表征的是我们用概率分布函数\(p(x)\)去近似概率分布函数\(q(x)\)时所引起的偏差,其计算方法如下

\begin{equation}-\int{p(x)ln(\frac{q(x)}{p(x)}})dx\end{equation}

这个值在\(p(x)=q(x)\)时得到最小值0.在PRML中,对于这个最小值的证明,使用了\(−ln(x)\)是凸函数这一性质。对于凸函数\(f\),我们有

\begin{equation} f(E[x])\leq E[f(x)]\end{equation}

对于连续型变量,我们有

\begin{equation}f(\int{p(x)xdx})\leq\int{f(p(x))dx}\end{equation}

所以

\begin{equation}-\int{p(x)ln(\frac{q(x)}{p(x)}}dx)\geq -ln \int{q(x)dx}=0\end{equation}

其实这里还有一个更简单的证明,我们可以注意到有\(ln(x)\leq x-1\),等于号只有在\(x=1\)的时候成立。所以,我们有

\begin{equation}-\int{p(x)ln(\frac{q(x)}{p(x)}}dx)\geq -\int{p(x)(\frac{q(x)}{p(x)}-1)dx}=\int{p(x)dx}-\int{q(x)dx}=0\end{equation}

等式成立条件跟之前一样,即完全拟合。

由熵到正态分布

对于某一\((-\infty ,+\infty)\),我们知道其期望\(\mu\)和方差\(\delta ^2\),则其概率密度函数\(p(x)\)满足下面的式子

\begin{equation}\begin{split}\int_{-\infty}^{+\infty}{p(x)dx}&=1\\ \int_{-\infty}^{+\infty}{p(x)xdx}&=\mu\\ \int_{-\infty}^{+\infty}{p(x){(x-\mu)}^2dx}&={\sigma}^2\end{split}\end{equation}

此时我们要求符合上述条件且使熵最大的\(p(x)\),根据拉格朗日乘子,原问题等价于使下式最大化

\begin{equation}-\int_{-\infty}^{+\infty}{p(x)ln(p(x))dx} +{\Lambda}_1(\int_{-\infty}^{+\infty}{p(x)dx}-1)+{\Lambda}_2(\int_{-\infty}^{+\infty}{p(x)xdx}-u)+{\Lambda}_3(\int_{-\infty}^{+\infty}{p(x){(x-\mu)}^2dx}-{\sigma}^2)\end{equation}

使导数等于0,我们可以得到

\begin{equation} p(x)=e^{-1+{\Lambda}_1+{\Lambda}_2x+{\Lambda}_3(x-u)^2}\end{equation}

将这个结果带入原来对\(p(x)\)的三个积分限制,最终得到了高斯分布方程

\begin{equation}p(x)=\frac{1}{\sqrt{2\pi}\sigma}e^{-\frac{(x-\mu)^2}{2{\sigma}^2}}\end{equation}

就这样,我们从最大熵原理变出了高斯分布公式。 推导出高斯分布有好几种方法,其中最直观的推导是假设期望的最大似然估计为样本均值推导出高斯分布。关于高斯分布的前世今生,想了解更多的话,可以参考此链接。事实上,这个链接所说的内容绝大部分都是在《Probility Theory The Logic of Science》中的直接翻译,这本书是一本非常好的书,对于概率论有非常深刻的了解,非常值得一看.

如何生成高斯分布

计算机中我们只能生成伪均匀分布的伪随机数,因为我们无法生成真随机数。Any one who considers arithmetical methods of producing random digits is ,of course,in a state of sin.(参见The Art of Computer Programming第二卷第一页)但是我们经常会有生成随机数的需求,生成高斯分布是其中比较常见的需求。最作死的方法莫过于\(k\)次重复模拟\((0,1)\)均匀分布,然后求均值。此均值在极限情况下服从\(N(0,1/{k^2})\)的高斯分布。另外一种方法就是Box和muller在1958年给出的由均匀分布的随机变量生成高斯分布的随机变量的算法。假设\(u,v\)是两个独立的随机变量,且都符合\((0,1)\)均匀分布。则下面两个变量都是服从\(N(0,1)\)分布的随机变量,且相互独立。

\begin{equation}\begin{split} x_1&=\sqrt{(-2ln(u))}*cos(2\pi v) \\ x_2&=\sqrt{(-2ln(v))}*cos(2\pi u)\end{split}\end{equation}

多元高斯分布

前面的高斯分布公式只是一元高斯分布,现在我们来看多元高斯分布。多元高斯分布的方程如下

\begin{equation}\mathcal{N}(x|\mu,{\tiny\sum})=\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)}\end{equation}

其中\(D\)为维数。现在我们来验证这个方程的正确性,即验证下面的等式

\begin{equation}\begin{split}\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}dx}=1\\\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}xdx}=\mu\\\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}(x-\mu)^T(x-\mu)dx}=|{\small\sum}|\end{split}\end{equation}

对于这三个等式,PRML中已经证明过了,不过我还是想再去证明一下,以加强记忆。首先对\(\sum\)进行谱分解。由于\(\sum\)为协方差矩阵,所以是对称正定的,有\(D\)个线性无关的特征向量\(u_i\),相应的特征值为\(\Lambda _i\)。所以我们可以完成该协方差矩阵的谱分解

\begin{equation}{\tiny\sum}=\sum_{i=1}^{D}{{\Lambda}_i{u_i}{u_i}^T}=diag(\Lambda)[u_1,\cdots,u_D][u_1,\cdots,u_D]^T\end{equation}

由于最后的三个矩阵都是方阵且可逆,所以我们可以非常简单的得到下式

\begin{equation}\begin{aligned}{{\tiny \sum}}^{-1}=[u_1,\cdots,u_D][u_1,\cdots,u_D]^Tdiag(\frac{1}{\Lambda})&=diag(\frac{1}{\Lambda})[u_1,\cdots,u_D][u_1,\cdots,u_D]^T\\&=[u_1,\cdots,u_D]diag(\frac{1}{\Lambda})[u_1,\cdots,u_D]^T\end{aligned}\end{equation}

我们令\(U=[u_1,\cdots,u_D]\),所以

\begin{equation}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)=(x-\mu)^TUdiag(\frac{1}{\Lambda})U^T(x-\mu)=(U^T(x-\mu))^Tdiag(\frac{1}{\Lambda})(U^T(x-\mu))\end{equation}

此时我们令\(y=U^T(x-\mu)\)\(dx=Udy\),此时\(y\)\(x\)的Jacobian矩阵为\(Jacobian(x)=U^T\)。根据多元积分的替代规则

\begin{equation}\int{f(x)dx}=\int{f(y)|Jacobian(x)|dy}\end{equation}

又因为

\begin{equation}det(|{\small \sum}|)=\prod_{i=1}^{D}{{\Lambda}_i}\end{equation}

所以

\begin{equation}p(y)=p(x)|U|=\prod_{i=1}^{D}{\frac{1}{{{({2\pi}{\Lambda}_i})}^{1/2}}exp\{-\frac{{y_i}^2}{2{\Lambda}_i}\}}\end{equation}

所以我们最终得到了

\begin{equation}\begin{aligned}&\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}dx}\\=&\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}y^Tdiag(\frac{1}{\Lambda})y\}dy}=1\end{aligned}\end{equation}

同理

\begin{equation}\begin{aligned}&\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}xdx}\\=&\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}y^Tdiag(\frac{1}{\Lambda})y\}(Uy+\mu)dy}\end{aligned}\end{equation}

\(y\)的对称性可得

\begin{equation}\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}y^Tdiag(\frac{1}{\Lambda})y\}Uydy}=0\end{equation}

同样,我们可以证明

\begin{equation}\begin{split}&\quad\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}(x-\mu)^T(x-\mu)dx}\\&=\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}y^Tdiag(\Lambda)y\}y^2dy}\\&=\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}\int{exp\{-\frac{1}{2}y^Tdiag(\Lambda)y\}y^2dy}\\&=|{\small\sum}|\end{split}\end{equation}

综上,我们证明了\(\mathcal{N}(x|\mu,{\tiny\sum})\)的确是一个概率分布,且期望为\(\mu\),协方差为\(|\sum|\)

高斯分布的边际分布

高斯分布的边际分布形式为\(p(x_1)\),其中\(x_1\)为一个或多个变量,且是总体变量的一部分。根据边际分布的计算公式,我们有

\begin{equation}p(x_1)=\int{p(x_1,x_2)dx_2}\end{equation}

其中\(x_2\)为一个或多个变量,\(x_1+x_2=x\)。由于

\begin{equation}p(x_1,x_2)=\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}\end{equation}

,所以这个积分直接计算比较困难,我们考虑分解\((x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\).令 \(x=[x_1,x_2]\),此时可以把期望分解为\(\mu=[{\mu}_1,{\mu}_2]\),同时可以把协方差分解为

\begin{equation}{\small \sum}=\left[\begin{matrix} {\tiny \sum}_{11}&{\tiny \sum}_{12}\\{\tiny \sum}_{21}&{\tiny \sum}_{22}\end{matrix}\right]\end{equation}

同时我们令

\begin{equation}\Lambda={\small \sum}^{-1}\end{equation}

此时

\begin{equation}\begin{split}(x-\mu)^T{\tiny \sum}^{-1}(x-\mu)&=(x-\mu)^T\Lambda(x-\mu)\\&=(x_1-{\mu}_1)^T{\Lambda}_{11}(x_1-{\mu}_1)+(x_1-{\mu}_1)^T{\Lambda}_{12}(x_2-{\mu}_2)\\&\quad +(x_2-{\mu}_2)^T{\Lambda}_{21}(x_1-{\mu}_1)+(x_2-{\mu}_2)^T{\Lambda}_{22}(x_2-{\mu}_2)\end{split}\end{equation}

由于

\begin{equation}-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)=-\frac{1}{2}x^T{{\tiny\sum}}^{-1}x+x^T{{\tiny\sum}}^{-1}\mu+-\frac{1}{2}{\mu}^T{{\tiny\sum}}^{-1}{\mu}\end{equation}

所以我们考虑把前一条式子转换为这种形式,把与\(x_2\)相关的都提取出来,类似于配方法,这样利用高斯分布的一些性质来简化积分。不过,这里我们先换元\(v_1=x_1-{\mu}_1,v_2=x_2-{\mu}_2\),所以

\begin{equation}\begin{split}&\quad(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\\&={v_1}^T{\Lambda}_{11}{v_1}+{v_1}^T{\Lambda}_{12}{v_2}+{v_2}^T{\Lambda}_{21}{v_1}+{v_2}^T{\Lambda}_{22}{v_2}\\&={v_2}^T{\Lambda}_{22}{v_2}+2{v_2}^T{\Lambda}_{21}{v_1}+{v_1}^T{\Lambda}_{11}{v_1}\\&={v_2}^T{\Lambda}_{22}{v_2}+2{v_2}^T{\Lambda}_{22}{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1+{v_1}^T{\Lambda}_{11}{v_1}\\&=(v_2+{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)^T{\Lambda}_{22}(v_2+{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)+const\end{split}\end{equation}

其中

\begin{equation}const={v_1}^T{\Lambda}_{11}{v_1}-({{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)^T{\Lambda}_{22}({{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)\end{equation}

此时假设\(v_1,v_2\)分别为\(a,b\)维的向量,所以

\begin{equation}\begin{split}p(v_1)&=\int{\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{-\frac{1}{2}(x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)\}dx_2}\\&=\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}exp\{const\}\int{exp\{-\frac{1}{2}(v_2+{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)^T{\Lambda}_{22}(v_2+{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)\}dv_2}\\&=\frac{1}{{(2\pi)}^{D/2}}\frac{1}{{|{\small\sum}|}^{1/2}}(2\pi)^{b/2}{|{\tiny \sum}_{22}|}^{1/2}exp\{-\frac{1}{2}const\}\end{split}\end{equation}

现在我们终于把与\(v_2\)相关的项弄消失了,下一步就是重新表示\(const\)

\begin{equation}\begin{split}const&={v_1}^T{\Lambda}_{11}{v_1}-({{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)^T{\Lambda}_{22}({{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1)\\&={v_1}^T{\Lambda}_{11}{v_1}-{v_1}^T{{\Lambda}_{21}}^T{{{\Lambda}_{22}}^{-1}}{\Lambda}_{22}{{\Lambda}_{22}}^{-1}{\Lambda}_{21}v_1\\&={v_1}^T{\Lambda}_{11}{v_1}-{v_1}^T{\Lambda}_{12}{{\Lambda}_{22}}^{-1}{\Lambda}_{21}{v_1}={v_1}^T({\Lambda}_{11}-{\Lambda}_{12}{{\Lambda}_{22}}^{-1}{\Lambda}_{21}){v_1}\end{split}\end{equation}

现在我们再用分块矩阵求逆来表示各个\(Λ\),由于

\begin{equation}\left[ \begin{matrix} A&B\\C&D\end{matrix}\right]=\left[ \begin{matrix}M&-MBD^{-1}\\-D^{-1}CM&D^{-1}+D^{-1}CMBD^{-1}\end{matrix}\right]\end{equation}

其中\(M\)是传说中的舒尔补

\begin{equation}M={(A-BD^{-1}C)}^{-1}\end{equation}

所以我们可以得到

\begin{equation}{\Lambda}_{11}={({\tiny \sum}_{11}-{\tiny \sum}_{12}{{\tiny \sum}_{22}}^{-1}{\tiny \sum}_{21})}^{-1}\end{equation}

将这个式子中的\(\Lambda,\sum\)相替换,我们可以得到

\begin{equation}const={v_1}^T{{\tiny \sum}_{11}}^{-1}{v_1}\end{equation}

。利用谱分解和Cholesky,我们可以得到

\begin{equation}|{\small \sum}|=|{\small \sum}_{11}|*|{\small \sum}_{22}|\end{equation}

我操,不想再写了,反正最后可以得到\(\mathcal{N}(0,{\small \sum}_{11})\),不要问我为什么。。。。所以最终\(x_1\)服从\(\mathcal{N}({\mu}_1,{\small \sum}_{11})\)

正态分布的条件分布

条件分布的形式为\(p(x_1|x_2)\),根据贝叶斯定理

\begin{equation}p(x_1|x_2)=\frac{p(x_1,x_2)}{\int{p(x_2)dx_2}}\end{equation}

这里我们只讨论\(x_1,x_2\)组成全集的情况,其他情况可以通过缩小全集的形式推广。根据我们之前求的边际分布公式

\begin{equation}\int{p(x_2)dx_2}=\frac{1}{{(2\pi)}^{b/2}}\frac{1}{{|{\small\sum}_{22}|}^{1/2}}exp\{-\frac{1}{2}(x_2-{\mu}_2)^T{{\tiny\sum}_{22}}^{-1}(x_2-{\mu}_2)\}\end{equation}

所以我们有

\begin{equation}p(x_1|x_2)=\frac{1}{{(2\pi)}^{a/2}}\frac{1}{{|{\small\sum}_{11}|}^{1/2}}exp\{-\frac{1}{2}((x-\mu)^T{{\tiny\sum}}^{-1}(x-\mu)-(x_2-{\mu}_2)^T{{\tiny\sum}_{22}}^{-1}(x_2-{\mu}_2))\}\end{equation}

再次利用之前的分块矩阵公式,并使用配方方法。这里我就不在去详细的算了,直接给出结果。此时的条件分布仍然是一个高斯分布,期望和条件方差为

\begin{equation}{\mu}_{1|2}={\mu}_1+{\tiny \sum}_{12}{{\tiny \sum}_{22}}^{-1}(x_2-{\mu}_2)\end{equation}
\begin{equation}{\tiny \sum}_{1|2}={\tiny \sum}_{11}-{\tiny \sum}_{12}{{\tiny \sum}_{22}}^{-1}{{\tiny \sum}_{21}}\end{equation}
Published:
2014-07-17 20:31
Category:
Tag: