经典机器学习系列(十)【变分推断】

  贝叶斯推断的优势在于可以结合一些我们已知的先验信息。对于复杂的一些问题建模也非常灵活，很适合用于应用统计方面的工作。但是贝叶斯模型有一些问题，如果用传统的MCMC来求的话速度会非常慢，而现今问题的数据量和维度都比较大，在这样的问题上做贝叶斯推理的话就会有很大的挑战。

贝叶斯与近似贝叶斯

  贝叶斯推断(Bayesian Inference)，在贝叶斯推断中我们有观测数据 x = { x 1 , x 2 , ⋯   , x n } x=\{x_{1},x_{2},\cdots, x_{n}\} x={x1​,x2​,⋯,xn​}，将模型的参数设置为$z$，它具有$m$个参数 z = { z 1 , z 2 , ⋯   , z m } z=\{z_{1},z_{2},\cdots,z_{m}\} z={z1​,z2​,⋯,zm​}。给定一个先验(prior)分布$p(z)$，在我们观测到数据$x$之后，我们依据这个数据对$z$的分布有一个更新。这种更新我们通常称之为计算后验分布(posterior distribution) $p(z|x)$。

$$\begin{gathered} p(z|x)=\frac{p(x,z)}{p(x)} \\ =\frac{p(x|z)p(z)}{p(x)} \end{gathered}$$

  其中$p(z|x)$是后验分布；$p(x|z)$是likelihood，相当于是一个方程，将数据$x$和我们所关心的参数$z$联系在了一起；$p(z)$是参数的先验prior信息；而$p(x)$是我们对于数据的信息，通常称之为evidence。

  一般来说，贝叶斯建模需要走以下几步：

1. 选一个prior，然后去选一个给定的likelihood。一般我们需要有一些假设，比如prior或者likelihood是服从某种分布的。
2. 基于观测数据，去计算后验分布。
3. 如果后验分布变量$z$的维度非常高的话，我们就没办法将其全部表示出来，我们一般用后验分布的一些统计量，来代表整个后验分布。比如常见的用posterior mean and variances。

  但是有一个问题，就是要去计算后验分布$p(z|x)$的话，我们需要去计算$p(x)={\int }_{z}p(x,z)dz$这个积分，如果是高维度上的积分，数值上计算是非常困难的。并且有些情况是可能不存在解析解的。常用的方法有Markov Chain Monte Carlo(MCMC)，但是计算比较慢。

  另外一种思路，把贝叶斯模型的计算问题转化为一个优化问题。之后就可以通过一些优化上的手段来进行加速。

KL散度

  而用一个分布去拟合另一个分布通常需要衡量这两个分布之间的相似性，通常采用KL散度，当然还有其他的一些方法，像JS散度这种。下面介绍KL散度：

  机器学习中比较重要的一个概念—相对熵(relative entropy)。相对熵又被称为KL散度(Kullback–Leibler divergence) 或信息散度 (information divergence)，是两个概率分布间差异的非对称性度量 。在信息论中，相对熵等价于两个概率分布的信息熵的差值，若其中一个概率分布为真实分布，另一个为理论（拟合）分布，则此时相对熵等于交叉熵与真实分布的信息熵之差，表示使用理论分布拟合真实分布时产生的信息损耗 。其公式如下：

$$D_{KL}(p\Vert q)=\sum _{i=1}^N\left[p\left(x_i\right)\log p\left(x_i\right)-p\left(x_i\right)\log q\left(x_i\right)\right]$$

  合并之后表示为：

$$D_{KL}(p\Vert q)=\sum _{i=1}^{N}p\left(x_i\right)\log \left(\frac{p\left(x_i\right)}{q\left(x_i\right)}\right)$$

  假设理论拟合出来的事件概率分布$q(x)$跟真实的分布$p(x)$一模一样，即$p(x)=q(x)$，那么$p\left(x_i\right)\log q\left(x_i\right)$就等于真实事件的信息熵，这一点显而易见。在理论拟合出来的事件概率分布跟真实的一模一样的时候，相对熵等于0。而拟合出来不太一样的时候，相对熵大于0。其证明如下：

$$\sum _{i=1}^{N}p\left(x_i\right)\log \frac{q\left(x_i\right)}{p\left(x_i\right)}\le \sum _{i=1}^{N}p\left(x_i\right)\left(\frac{q\left(x_i\right)}{p\left(x_i\right)}-1\right)=\sum _{i=1}^N\left[p\left(x_i\right)-q\left(x_i\right)\right]=0$$

  其中第一个不等式是由$ln(x)\le x-1$推导出来的，只在$p(x_i)=q(x_i)$时取到等号。

  这个性质很关键，因为它正是深度学习梯度下降法需要的特性。假设神经网络拟合完美了，那么它就不再梯度下降，而不完美则因为它大于0而继续下降。

  但它有不好的地方，就是它是不对称的。也就是用$P$ 来拟合$Q$ 和用$Q$ 来拟合 $P$ 的相对熵居然不一样，而他们的距离是一样的。这也就是说，相对熵的大小并不跟距离有一一对应的关系。

变分贝叶斯

  我们经常利用贝叶斯公式求posterior distribution $P(Z|X)$

$$P(Z|X)=\frac{p(X,Z)}{{\int }_{z}p(X,Z=z)dz}$$

  但posterior distribution $P(Z|X)$求解用贝叶斯的方法是比较困难的，因为我们需要去计算${\int }_{z}p(X=x,Z=z)dz$，而$Z$通常会是一个高维的随机变量，这个积分计算起来就非常困难。在贝叶斯统计中，所有的对于未知量的推断(inference)问题可以看做是对后验概率(posterior)的计算。因此提出了Variational Inference来计算posterior distribution。

  那Variational Inference怎么做的呢？其核心思想主要包括两步：

1. 假设一个分布$q(z;\lambda )$ (这个分布是我们搞得定的，搞不定的就没意义了)
2. 通过改变分布的参数$\lambda$，使$q(z;\lambda )$ 靠近 $p(z|x)$。

  总结称一句话就是，为真实的后验分布引入了一个参数话的模型。 即：用一个简单的分布$q(z;\lambda )$ 拟合复杂的分布$p(z|x)$。

  这种策略将计算 $p(z|x)$ 的问题转化成优化问题了

$${\lambda }^{\ast }=\arg \underset{\lambda }{\min }\mathrm{divergence}(p(z|x),q(z;\lambda ))$$

  收敛后，就可以用 $q(z;\lambda )$ 来代替 $p(z|x)$了。

  本文的目的还是来求这个变分推理，不要走偏了。下面涉及一些公式等价转换：

$$\begin{array}{rl}\log P(x)& =\log P(x,z)-\log P(z|x)\\ & =\log \frac{P(x,z)}{Q(z;\lambda )}-\log \frac{P(z|x)}{Q(z;\lambda )}\end{array}$$

  等式两边同时对$Q(z)$求期望，得：

$$\begin{array}{rl}{\mathbb{E}}_{q(z;\lambda )}\log P(x)& ={\mathbb{E}}_{q(z;\lambda )}\log P(x,z)-{\mathbb{E}}_{q(z;\lambda )}\log P(z|x)\\ \log P(x)& ={\mathbb{E}}_{q(z;\lambda )}\log \frac{p(x,z)}{q(z;\lambda )}-{\mathbb{E}}_{q(z;\lambda )}\log \frac{p(z|x)}{q(z;\lambda )}\\ & =KL(q(z;\lambda )\Vert p(z|x))+{\mathbb{E}}_{q(z;\lambda )}\log \frac{p(x,z)}{q(z;\lambda )}\\ \log P(x)& =KL(q(z;\lambda )\Vert p(z|x))+{\mathbb{E}}_{q(z;\lambda )}\log \frac{p(x,z)}{q(z;\lambda )}\end{array}$$

  到这里我们需要回顾一下我们的问题，从另一个角度再来思考一遍，寻找一个近似后验$q^{\ast }(z)\in Q$去近似$p(z|x)$。

$$q^{\ast }(z)=\underset{q(z)\in Q}{\mathrm{argmin}}L(q(z),p(z|x))$$

  这里有两个注意点。1. 如何来选取这个$Q$。2. $L$表示它们之间的距离度量，通常用KL-divergence来表示：$L(q(z),p(z|x))=KL(q(z)||p(z|x))$。用KL的原因在于它满足一些凸性，可导等性质。另一方面在于它具备一些性质：

$$\begin{gathered} q^{\ast }(z)=\underset{q(z)\in Q}{\mathrm{argmin}}KL(q(z),p(z|x)) \\ =\underset{q(z)\in Q}{\mathrm{argmin}}-{\int }_{z}q(z)log[\frac{p(z|x)}{q(z)}] \end{gathered}$$

  $p(z|x)$是我们的target，$q(z)$是我们提出的，需要去优化的分布。要去算这个优化的目标函数，我们需要已知$p(z|x)$这个后验分布，但是我们不知道这个后验分布是什么。因此，上述这个公式我们是无法直接计算的。这就变成了一个循环的问题，为了要去算$p(z|x)$，要去计算$q(z)$，而要计算$q(z)$又需要计算$p(z|x)$。在没有其它约束的情况下，当$q^{\ast }=p(z|x)$的时候能够取到理论的最优值，但是$p(z|x)$并不知道。

  我们对KL Divergence最小化的问题做一些形式上的转变，期望去找到一个等价的问题：

$$\begin{array}{rl}\mathrm{KL}(q(z)\Vert p(z\mid x))& =-{\int }_{z}q(z)\log \left[\frac{p(z\mid x)}{q(z)}\right]dz\\ & ={\int }_{z}q(z)\log q(z)dz-{\int }_{z}q(z)\log p(z\mid x)dz\\ & ={\mathbb{E}}_q[\log q(z)]-{\mathbb{E}}_q[\log p(z\mid x)]\\ & ={\mathbb{E}}_q[\log q(z)]-{\mathbb{E}}_q\left[\log \left[\frac{p(x,z)}{p(x)}\right]\right]\\ & ={\mathbb{E}}_q[\log q(z)]-{\mathbb{E}}_q[\log p(x,z)]+{\mathbb{E}}_q[\log p(x)]\\ & =\underset{-\text{ELBO }}{\underbrace{{\mathbb{E}}_q[\log q(z)]-{\mathbb{E}}_q[\log p(x,z)]}}+\log p(x)\end{array}$$

  上述公式经过一定的变型之后，ELBO就是log-likelihood的一个下界函数：

$$\begin{array}{rl}\log p(x)& =\mathrm{ELBO}(q)+\mathrm{KL}(q(z)\Vert p(z\mid x))\\ & \ge \mathrm{ELBO}(q)\end{array}$$

  因为$\log p(x)$是一个与$q$无关的常数，那么对$\mathrm{KL}(q(z)\Vert p(z\mid x))$的Minimizing问题就可以等价于Maximizing ELBO($q$)。

  我们为什么要去做这样的一个替换呢？之前我们并不知道$q^{\ast }(z)={\mathrm{argmin}}_{q(z)\in Q}KL(q(z),p(z|x))$中的后验分布$p(z|x)$。因为想要计算它的话，需要计算一个积分。而ELBO的表达式如下：

$$\begin{array}{rl}\mathrm{ELBO}(q)& ={\mathbb{E}}_q[\log p(x,z)]-{\mathbb{E}}_q[\log q(z)]\\ & ={\int }_{z}q(z)\log [\stackrel{\text{known }}{\overbrace{\frac{p(x,z)}{q(z)}}}]dz\end{array}$$

  其中$p(x,z)=p(x|z)p(z)$，是一个先验分布与likelihood的乘积。因此可以对ELBO做优化求解。现在，variational inference的目标变成：

$$\underset{\lambda }{\max }{\mathbb{E}}_{q(z;\lambda )}[\text{log}p(x,z)-\text{log}q(z;\lambda )]$$

  ${\mathbb{E}}_{q(z;\lambda )}[\text{log}p(x,z)-\text{log}q(z;\lambda )]$称为Evidence Lower Bound(ELBO)。$p(x)$一般被称之为evidence，又因为 $KL(q||p)>=0$, 所以 $p(x)>=E_{q(z;\lambda )}[\text{log}p(x,z)-\text{log}q(z;\lambda )]$, 这就是为什么被称为ELBO。

  实际上EM算法（Expectation-Maximization)就是利用了这一特征，它分为交替进行的两步：E step假设模型参数不变，$q(z)=p(z|x)$，计算对数似然率，在M step再做ELBO相对于模型参数的优化。与变分法比较，EM算法假设了当模型参数固定时， $p(z|x)$ 是易计算的形式，而变分法并无这一限制，对于条件概率难于计算的情况，变分法仍然有效。

  那如何来求解上述公式呢？下面介绍平均场（mean-field）、蒙特卡洛、和黑盒变分推断 (Black Box Variational Inference) 的方法。

平均场变分族(mean-field variational family)

  之前我们说我们选择一族合适的近似概率分布 $q(Z;\lambda )$，那么实际问题中，我们可以选择什么形式的 $q(Z;\lambda )$ 呢？

  一个简单而有效的变分族为平均场变分族(mean-field variational family)。它假设了隐藏变量间是相互独立的：

$$q(Z;\lambda )=\prod _{k=1}^K{q}_k(Z_k;{\lambda }_k)$$

  这个假设看起来似乎比较强，但实际应用范围还是比较广泛，我们可以将其延展为将有实际相互关联的隐藏变量分组，而化为各组联合分布的乘积形式即可。

  利用ELBO和平均场假设，我们就可以利用coordinate ascent variational inference（简称CAVI)方法来处理：

• 利用条件概率分布的链式法则有 ：

$$p\left(z_{1:m},x_{1:n}\right)=p\left(x_{1:n}\right)\prod _{j=1}^{m}p\left(z_j|z_{1:(j-1)},x_{1:n}\right)$$

• 变分分布的期望为：

$$E\left[\log q\left(z_{1:m}\right)\right]=\sum _{j=1}^m{E}_j\left[\log q\left(z_j\right)\right]$$

  将其代入ELBO的定义得到：

$$ELBO=\mathrm{logp}\left(x_{1:n}\right)+\sum _{j=1}^{m}E\left[\log p\left(z_j|z_{1:(j-1)},x_{1:n}\right)\right]-E_j\left[\log q\left(z_j\right)\right]$$

  将其对 $z_k$ 求导并令导数为零有：

$$\frac{dELBO}{dq\left(z_k\right)}=E_{-k}\left[\log p\left(z_k|z_{-k},x\right)\right]-\log q\left(z_k\right)-1=0$$

  由此得到coordinate ascent 的更新法则为：

$$q^{\ast }\left(z_k\right)\propto \exp E_{-k}\left[\log p\left(z_k,z_{-k},x\right)\right]$$

  我们可以利用这一法则不断的固定其他的$z$ 的坐标来更新当前的坐标对应的$z$值，这与Gibbs Sampling过程类似，不过Gibbs Sampling是不断的从条件概率中采样，而CAVI算法中是不断的用如下形式更新：

$$q^{\ast }\left(z_k\right)\propto \exp E[\log (\text{conditional})]$$

Proof：

$$\begin{array}{rl}\mathrm{ELBO}(q)=& {\mathbb{E}}_q[\log p(x,z)]-{\mathbb{E}}_q[\log q(z)]\\ =& {\int }_{z}q(z)\log p(z,x)dz-\left[{\mathbb{E}}_{q_j}\left[\log q_j\left(z_j\right)\right]+\underset{\text{const }}{\underbrace{\sum _{i\ne j}{\mathbb{E}}_{q_i}\left[\log q_i\left(z_i\right)\right]}}\right]\\ =& {\int }_{z_j}{q}_j\left(z_j\right)\left[{\int }_{z_{-j}}{q}_{-j}\left(z_{-j}\right)\log p(z,x)d{z}_{-j}\right]d{z}_j\\ & -{\mathbb{E}}_{q_j}\left[\log q_j\left(z_j\right)\right]-\text{ const }\\ =& {\int }_{z_j}{q}_j\left(z_j\right){\mathbb{E}}_{q-j}[\log p(z,x)]d{z}_j-{\int }_{z_j}{q}_j\left(z_j\right)\log q_j\left(z_j\right)d{z}_j-\text{ const }\end{array}$$

  定义一个新的分布：

$$\log {\tilde{p}}_j\left(z_j,x\right)={\mathbb{E}}_{q-j}[\log p(z,x)]+\text{ const }$$

  此时的ELBO就可以写成如下形式：

$$\begin{array}{rl}\mathrm{ELBO}(q)& ={\int }_{z_j}{q}_j\left(z_j\right)\log {\tilde{p}}_j\left(z_j,x\right)d{z}_j-{\int }_{z_j}{q}_j\left(z_j\right)\log q_j\left(z_j\right)d{z}_j+\text{ const }\\ & ={\int }_{z_j}{q}_j\left(z_j\right)\log \frac{{\tilde{p}}_j\left(z_j,x\right)}{q_j\left(z_j\right)}d{z}_j+\text{ const }\\ & =-\mathrm{KL}\left(q_j\left(z_j\right)\Vert {\tilde{p}}_j\left(z_j,x\right)\right)+\text{ const }\end{array}$$

  此时的KL divergence达到最小值时，当：

q j ∗ ( z j ) = p ~ j ( z j , x ) ∝ exp ⁡ { E q − j [ log ⁡ p ( z , x ) ] } (  normalize  ) = exp ⁡ { E q − j [ log ⁡ p ( z , x ) ] } ∫ z j exp ⁡ { E q − j [ log ⁡ p ( z , x ) ] } d z j \begin{aligned} q_{j}^{*}\left(z_{j}\right) &=\tilde{p}_{j}\left(z_{j}, x\right) \\ & \propto \exp \left\{\mathbb{E}_{q-j}[\log p(z, x)]\right\} \\ (\text { normalize }) &=\frac{\exp \left\{\mathbb{E}_{q-j}[\log p(z, x)]\right\}}{\int_{z_{j}} \exp \left\{\mathbb{E}_{q_{-j}}[\log p(z, x)]\right\} d z_{j}} \end{aligned} qj∗​(zj​)( normalize )​=p~​j​(zj​,x)∝exp{Eq−j​[logp(z,x)]}=∫zj​​exp{Eq−j​​[logp(z,x)]}dzj​exp{Eq−j​[logp(z,x)]}​​

  其完整算法如下所示：

MCMC

  MCMC方法是利用马尔科夫链取样来近似后验概率，变分法是利用优化结果来近似后验概率，那么我们什么时候用MCMC，什么时候用变分法呢？

  首先，MCMC相较于变分法计算上消耗更大，但是它可以保证取得与目标分布相同的样本，而变分法没有这个保证：它只能寻找到近似于目标分布一个密度分布，但同时变分法计算上更快，由于我们将其转化为了优化问题，所以可以利用诸如随机优化(stochastic optimization)或分布优化(distributed optimization)等方法快速的得到结果。所以当数据量较小时，我们可以用MCMC方法消耗更多的计算力但得到更精确的样本。当数据量较大时，我们用变分法处理比较合适。

  另一方面，后验概率的分布形式也影响着我们的选择。比如对于有多个峰值的混合模型，MCMC可能只注重其中的一个峰而不能很好的描述其他峰值，而变分法对于此类问题即使样本量较小也可能优于MCMC方法。

黑盒变分推断(BBVI)

  ELBO公式表达为：

$${\mathbb{E}}_{q(z;\lambda )}[\text{log}p(x,z)-\text{log}q(z;\lambda )]$$

  对用参数$\theta$替代$\lambda$，并对其求导：

$${\nabla }_{\theta }\mathrm{ELBO}(\theta )={\nabla }_{\theta }{\mathbb{E}}_q\left(\log p(x,z)-\log q_{\theta }(z)\right)$$

  直接展开计算如下：

$$\begin{array}{rl}& \frac{\partial }{\partial \theta }\int q_{\theta }(z)\left(\log p(x,z)-\log q_{\theta }(z)\right)dz\\ =& \int \frac{\partial }{\partial \theta }\left[q_{\theta }(z)\left(\log p(x,z)-\log q_{\theta }(z)\right)\right]dz\\ =& \int \frac{\partial }{\partial \theta }\left(q_{\theta }(z)\log p(x,z)\right)-\frac{\partial }{\partial \theta }\left(q_{\theta }(z)\log q_{\theta }(z)\right)dz\\ =& \int \frac{\partial q_{\theta }(z)}{\partial \theta }\log p(x,z)-\frac{\partial q_{\theta }(z)}{\partial \theta }\log q_{\theta }(z)-\frac{\partial q_{\theta }(z)}{\partial \theta }dz\end{array}$$

  由于：

$$\int \frac{\partial q_{\theta }(z)}{\partial \theta }dz=\frac{\partial }{\partial \theta }\int q_{\theta }(z)dz=\frac{\partial }{\partial \theta }1=0$$

  因此：

$$\begin{array}{rl}{\nabla }_{\theta }\mathrm{ELBO}(\theta )& =\int \frac{\partial q_{\theta }(z)}{\partial \theta }\left(\log p(x,z)-\log q_{\theta }(z)\right)dz\\ & =\int q_{\theta }(z)\frac{\partial \log q_{\theta }(z)}{\partial \theta }\left(\log p(x,z)-\log q_{\theta }(z)\right)dz\\ & =\int q_{\theta }(z){\nabla }_{\theta }\log q_{\theta }(z)\left(\log p(x,z)-\log q_{\theta }(z)\right)dz\\ & ={\mathbb{E}}_q\left[{\nabla }_{\theta }\log q_{\theta }(z)\left(\log p(x,z)-\log q_{\theta }(z)\right)\right]\end{array}$$

  然后写成 SGD，就是所谓 Black Box Variational Inference (BBVI)。

$$\begin{array}{rl}& {\mathbb{E}}_{z\sim q_{\theta }(z)}\left[{\nabla }_{\theta }\log q_{\theta }(z)\left(\log p(x,z)-\log q_{\theta }(z)\right)\right]\\ \approx & \frac{1}{N}\sum _{i=1}^N{\nabla }_{\theta }\log q_{\theta }\left(z_i\right)\left(\log p\left(x,z_i\right)-\log q_{\theta }\left(z_i\right)\right)\end{array}$$

  其中$z_i\sim q_{\theta }(z)$。

参考

• https://www.jianshu.com/p/b8ecabf9af07
• https://blog.csdn.net/weixinhum/article/details/85064685
• https://zhuanlan.zhihu.com/p/51567052
• https://zhuanlan.zhihu.com/p/49401976
• https://space.bilibili.com/491707363/channel/detail?cid=162191

我的微信公众号名称：小小何先生
公众号介绍：主要研究分享深度学习、机器博弈、强化学习等相关内容！期待您的关注，欢迎一起学习交流进步！