[受限玻尔兹曼机] 原理、求解过程推导、深度信念网络

写在前面

本篇文章主要写受限玻尔兹曼机、Gibbs求解方法、CD对比散度求解方法和深度信念网络。并非全部都是原创，有部分来自于书籍和网络。

一、受限玻尔兹曼机

1、简介玻尔兹曼机结构

如下图所示（图片来源：知乎），玻尔兹曼机的网络是相互连接型网络，其有两种形式：
①网络节点全部由可见单元构成。
②网络节点由可见单元和隐藏单元共同组成。

这种网络结构的训练非常复杂，因此基于此网络，辛顿等人提出了受限玻尔兹曼机。

2、受限玻尔兹曼机

受限玻尔兹曼机（Restricted Boltzmann Machines, RBM）主要有伯努利-伯努利RBM、伯努利-高斯RBM、高斯-伯努利RBM，本文主要基于伯努利-伯努利RBM，即可观测变量服务伯努利分布，隐变量也服从伯努利分布。

在玻尔兹曼机的基础上，首先玻尔兹曼机增加了限制：层内单元之间无连接。如下图所示：

受限玻尔兹曼机是由可见层和隐藏层构成的两层结构。可见层由可见变量组成，隐藏层由隐藏变量组成。可见层和隐藏层之间相互连接，层内之间无连接。

受限玻尔兹曼机的目的： 通过求取连接权重和偏置项，使得重构的$\tilde{v}$与$v$误差最小。利用这些连接权重和偏置项我们可以进行分类或者进行生成。

受限玻尔兹曼机的能量函数如下：
$$E(v,h,\theta )=-\sum _{i=1}^n{b}_i{v}_i-\sum _{j=1}^m{c}_j{h}_j-\sum _{i=1}^n\sum _{j=1}^m{w}_{ij}{v}_i{h}_j$$
其中，$b_i$是可见变量的偏置，$c_j$是隐藏变量的偏置，$w_{ij}$是连接权重，$\theta$是表示所有连接权重和偏置的参数集合。状态$(v,h)$的联合概率分布如下所示：
p ( v , h ∣ θ ) = 1 Z e x p { − E ( v , h , θ ) } Z = ∑ v , h e x p { − E ( v , h , θ ) } p(\red{v},h|\theta)=\frac{1}{Z}exp\{-E(\red{v},h,\theta)\} \\ Z=\sum_{v,h}exp\{-E(v,h,\theta)\} p(v,h∣θ)=Z1​exp{−E(v,h,θ)}Z=v,h∑​exp{−E(v,h,θ)}
注意：标红的$v$表示是某一个$v$，而$Z$中未标红的$v$表示随机变量，不是某个固定值。

其对数似然函数为（关于此似然函数的原因，参见深度学习读书笔记之RBM（限制波尔兹曼机）1.3.2 从概率到极大似然）：
l o g L ( θ ∣ v ) = l o g 1 Z ∑ h e x p { − E ( v , h , θ ) } = l o g ∑ h e x p { − E ( v , h , θ ) } − l o g ∑ v , h e x p { − E ( v , h , θ ) } \begin{array}{lcl} logL(\theta|\red{v}) &=& log\frac{1}{Z}\sum_{h}exp\{-E(\red{v},h,\theta)\} \\ &=& log\sum_hexp\{-E(\red{v},h,\theta)\}-log\sum_{v,h}exp\{-E(v,h,\theta)\} \end{array} logL(θ∣v)​==​logZ1​∑h​exp{−E(v,h,θ)}log∑h​exp{−E(v,h,θ)}−log∑v,h​exp{−E(v,h,θ)}​
上述似然函数对$\theta$求导计算梯度，过程如下：
∂ l o g L ( θ ∣ v ) ∂ θ = ∂ ( l o g ∑ h e x p { − E ( v , h , θ ) } ) ∂ θ − ∂ ( l o g ∑ v , h e x p { − E ( v , h , θ ) } ) ∂ θ = − 1 ∑ h e x p { − E ( v , h , θ ) } ∑ h ( e x p { − E ( v , h , θ ) } ∂ E ( v , h , θ ) ∂ θ ) + 1 ∑ v , h e x p { − E ( v , h , θ ) } ∑ v , h ( e x p { − E ( v , h , θ ) } ∂ E ( v , h , θ ) ∂ θ ) = − ∑ h ( p ( h ∣ v ) ∂ E ( v , h , θ ) ∂ θ ) + ∑ v , h ( p ( v , h ) ∂ E ( v , h , θ ) ∂ θ ) \begin{array}{lcl} \frac{\partial logL(\theta|v)}{\partial \theta} &=& \frac{\partial (log\sum_h exp\{-E(\red{v},h,\theta)\})}{\partial \theta}-\frac{\partial (log\sum_{v,h}exp\{-E(v,h,\theta)\})}{\partial \theta} \\ &=& -\frac{1}{\sum_h exp\{-E(\red{v},h,\theta)\}}\sum_h (exp\{-E(\red{v},h,\theta)\}\frac{ \partial E(\red{v},h,\theta)}{\partial \theta}) \\ &+& \frac{1}{\sum_{v,h} exp\{-E(v,h,\theta)\}}\sum_{v,h} (exp\{-E(v,h,\theta)\}\frac{ \partial E(v,h,\theta)}{\partial \theta}) \\ &=& -\sum_h (p(h|\red{v})\frac{ \partial E(\red{v},h,\theta)}{\partial \theta})+\sum_{v,h} (p(v,h)\frac{ \partial E(v,h,\theta)}{\partial \theta}) \end{array} ∂θ∂logL(θ∣v)​​==+=​∂θ∂(log∑h​exp{−E(v,h,θ)})​−∂θ∂(log∑v,h​exp{−E(v,h,θ)})​−∑h​exp{−E(v,h,θ)}1​∑h​(exp{−E(v,h,θ)}∂θ∂E(v,h,θ)​)∑v,h​exp{−E(v,h,θ)}1​∑v,h​(exp{−E(v,h,θ)}∂θ∂E(v,h,θ)​)−∑h​(p(h∣v)∂θ∂E(v,h,θ)​)+∑v,h​(p(v,h)∂θ∂E(v,h,θ)​)​
上式中，$\theta$可以取$w_{ij}$、$b_i$和$c_j$，分别计算梯度如下：

（1）求对$w_{ij}$的梯度

$$\begin{array}{lcl}\frac{\partial logL(\theta |v)}{\partial w_{ij}}& =& -{\sum }_h(p(h|v)\frac{\partial E(v,h,\theta )}{\partial w_{ij}})+{\sum }_{v,h}p(v,h)(\frac{\partial E(v,h,\theta )}{\partial w_{ij}})\\ & =& {\sum }_{h_j}(p(h_j|v)v_i{h}_j)-{\sum }_{v_i,h_j}p(v_i,h_j)v_i{h}_j\\ & =& {\sum }_{h_j}(p(h_j|v)v_i{h}_j)-{\sum }_{v_i}(p(v_i){\sum }_{h_j}p(h_j|v_i)v_i{h}_j)\\ & =& p(h_j=0|v)v_i\cdot 0+p(h_j=1|v)v_i\cdot 1\\ & -& ({\sum }_{v_i}p(v_i)p(h_j=0|v_i)v_i\cdot 0+{\sum }_{v_i}p(v_i)p(h_j=1|v_i)v_i\cdot 1)\\ & =& p(h_j=1|v)v_i-{\sum }_{v_i}(p(v_i)p(h_j=1|v_i)v_i)\end{array}$$
注意： 在这里，因为求的是$w_{ij}$的梯度，因此求导之后参考其中的只有$h_j$和$v_i$。

（2）求对$b_i$的梯度

$$\begin{array}{lcl}\frac{\partial logL(\theta |v)}{\partial b_i}& =& -{\sum }_h(p(h|v)\frac{\partial E(v,h,\theta )}{\partial b_i})+{\sum }_{v,h}p(v,h)(\frac{\partial E(v,h,\theta )}{\partial b_i})\\ & =& {\sum }_h(p(h|v)v_i)-{\sum }_{v_i,h}p(v_i,h)v_i\\ & =& {\sum }_h(p(h|v)v_i)-{\sum }_{v_i}(p(v_i){\sum }_{h}p(h|v_i)v_i)\\ & =& v_i-{\sum }_{v_i}(p(v_i)v_i)\end{array}$$
注意： ${\sum }_{h}p(h|v_i)=1$，注意区分${\sum }_{v}p(h|v)=p(h)$。

（3）求对$c_j$的梯度

$$\begin{array}{lcl}\frac{\partial logL(\theta |v)}{\partial c_j}& =& -{\sum }_h(p(h|v)\frac{\partial E(v,h,\theta )}{\partial c_j})+{\sum }_{v,h}p(v,h)(\frac{\partial E(v,h,\theta )}{\partial c_j})\\ & =& {\sum }_{h_j}(p(h_j|v)h_j)-{\sum }_{v,h_j}p(v,h_j)h_j\\ & =& {\sum }_{h_j}(p(h_j|v)h_j)-{\sum }_v(p(v){\sum }_{h_j}p(h_j|v)h_j)\\ & =& p(h_j=0|v)\cdot 0+p(h_j=1|v)\cdot 1\\ & -& ({\sum }_{v}p(v)p(h_j=0|v)\cdot 0+{\sum }_{v}p(v)p(h_j=1|v)\cdot 1)\\ & =& p(h_j=1|v)-{\sum }_v(p(v)p(h_j=1|v))\end{array}$$

（4）更新方程

对$w_{ij}$、$b_i$和$c_j$求完梯度之后，可以按照如下方程更新参数：
$$\begin{array}{lcl}{w}_{ij}& \leftarrow & w_{ij}-\frac{\partial (logL(\theta |v))}{\partial w_{ij}}\\ b_i& \leftarrow & b_i-\frac{\partial (logL(\theta |v))}{\partial b_i}\\ c_j& \leftarrow & c_j-\frac{\partial (logL(\theta |v))}{\partial c_j}\end{array}$$
迭代更新即可。

3、求解存在的问题

受限玻尔兹曼机也存在效率的问题，上述公式可以看到，${\sum }_{v}p(v)$是需要计算的，而这个式子的计算需要所有可见变量的模式之和，其复杂度为$O(2^{n_v+n_h})$，这会使得效率非常低。因此可以使用Gibbs采样和CD对比散度算法来求解。

二、受限玻尔兹曼机的求解

在给定一个训练样本后，训练一个RBM的意义在于调整模型的参数，以拟合给定的训练样本，使得在该参数下RBM表示的可见层节点概率分布尽可能的与训练数据相符合。我们可以采用Gibbs采样或者CD算法求解参数，CD算法的效率更高。
Gibbs采样求解和CD算法求解都是近似求解模型，我们希望$p(v)=p_{train}(v)$（数据的真实、底层分布）。所以我们通过近似算法， 迭若干次之后，使得马尔科夫链拟合最终的$v$的分布。

1、Gibbs采样

我们希望$p(v)=p_{train}(v)$（数据的真实、底层分布），所以我们可以使用Gibss采样去收敛这个最终的分布$p$。
Gibbs采样的思想是虽然我们不知道样本数据$v_1,v_2,...,v_{n_v}$的联合概率，但是我们知道样本的条件概率$p(v_i|v_{-i})$，则我们可以先求出每个数据的条件概率值，得到$v$的任一状态$[v_1(0),v_2(0),...,v_{n_v}(0)]$，然后，用条件概率公式迭代对每一个数据求条件概率。最终，迭代$k$次的时候，$v$的某一个状态$[v_1(k),v_2(k),...,v_{n_v}(k)]$将收敛于$v$的联合概率分布$p(v)$，其过程图如下（图片来自参考博客）：

求解过程是： 假设给我一个训练样本$v_0$，根据公式$P(h_j=1|v)$求$h_0$中每一个节点的条件概率，再根据公式$p(v_i=1|h)$求$v_1$中每个节点的条件概率，然后依次迭代，直到K步，此时$p(v|h)$的概率收敛于$p(v)$的概率。
$$\begin{array}{l}{h}_0\sim p(h|v_0)\\ v_1\sim p(v|h_0)\\ h_1\sim p(h|v_1)\\ v_2\sim p(v|h_1)\\ ...\\ v_{k+1}\sim p(v|h_k)\end{array}$$

2、CD对比散度算法

上述使用Gibbs采样虽然已经相比较BP算法来说，代价减少了很多，但是Gibbs采样需要不断迭代直到收敛到平稳分布，所以迭代过程需要多次。基于此，hinton提出了对比散度算法，对比散度算法和Gibbs采样的区别在于：Gibbs采样$p(h|v)$中$v$的初始化是随机初始化，但是CD算法中$v$初始化为训练样本，其余和Gibbs采样相同。实验证明，CD算法只要迭代1次就可以有很好的效果。

（1）CD算法的算法流程图

（图片来源:受限玻尔兹曼机(RBM)+对比散度算法(CD-k)）

这是CD-1的算法流程，CD-k算法便是k次CD-1算法的迭代。

（2）条件分布$P(h|v)$和$P(v|h)$

因为在CD算法中需要求$p(h_j=1|v)$的取值，所以我们先来看看它的公式是什么。
CD算法的能量函数为：
$$E(v,h,\theta )=-\sum _{i=1}^n{b}_i{v}_i-\sum _{j=1}^m{c}_j{h}_j-\sum _{i=1}^n\sum _{j=1}^m{w}_{ij}{v}_i{h}_j$$
因此，定义给定$v,h$情况下RBM的状态概率分布为：
$$\begin{array}{l}P(v,h)=\frac{1}{Z}{e}^{-E(v,h)}\\ Z={\sum }_{v,h}{e}^{-E(v,h)}\end{array}$$
有了概率分布，我们来看条件分布$P(h|v)$：
P ( h ∣ v ) = P ( h , v ) P ( v ) = 1 P ( v ) 1 Z e x p { ∑ i = 1 n b i v i + ∑ j = 1 m c j h j + ∑ i = 1 n ∑ j = 1 m w i j v i h j } = 1 Z ′ e x p { ∑ j = 1 m c j h j + ∑ i = 1 n ∑ j = 1 m w i j v i h j } = 1 Z ′ ∏ j = 1 n h e x p { c j h j + w i j v i h j } \begin{array}{lcl} P(h|v) &=& \frac{P(h,v)}{P(v)} \\ &=& \frac{1}{P(v)}\frac{1}{Z}exp\{\sum_{i=1}^nb_iv_i+ \sum_{j=1}^mc_jh_j+\sum_{i=1}^n\sum_{j=1}^mw_{ij}v_ih_j\} \\ &=& \frac{1}{{Z}'}exp\{\sum_{j=1}^mc_jh_j+\sum_{i=1}^n\sum_{j=1}^mw_{ij}v_ih_j\} \\ &=& \frac{1}{{Z}'}\prod_{j=1}^{n_h}exp\{c_jh_j +w_{ij}v_ih_j \} \end{array} P(h∣v)​====​P(v)P(h,v)​P(v)1​Z1​exp{∑i=1n​bi​vi​+∑j=1m​cj​hj​+∑i=1n​∑j=1m​wij​vi​hj​}Z′1​exp{∑j=1m​cj​hj​+∑i=1n​∑j=1m​wij​vi​hj​}Z′1​∏j=1nh​​exp{cj​hj​+wij​vi​hj​}​
其中，$Z^{\prime }$为新的归一化常数：
1 Z ′ = 1 P ( v ) 1 Z e x p { ∑ i = 1 n b i v i } \frac{1}{{Z}'} = \frac{1}{P(v)}\frac{1}{Z}exp\{\sum_{i=1}^nb_iv_i\} Z′1​=P(v)1​Z1​exp{i=1∑n​bi​vi​}
下面以同样的方式求条件分布$P(v|h)$。

能量模型基于条件分布的激活函数如下，我们以$P(h_j=1|v)$为例推导如下：
P ( h j = 1 ∣ v ) = P ( h j = 1 ∣ v ) P ( h j = 0 ∣ v ) + P ( h j = 1 ∣ v ) = e x p { c j + w i j v i } 1 + e x p { c j + w i j v i } = 1 1 + e x p { − ( c j + w i j v i ) } = σ ( c j + w i j v i ) \begin{array}{lcl} P(h_j=1|v) &=& \frac{P(h_j=1|v) }{P(h_j=0|v) +P(h_j=1|v) } \\ &=& \frac{exp\{c_j +w_{ij}v_i\}}{1+exp\{c_j +w_{ij}v_i\}} \\ &=& \frac{1}{1+exp\{-(c_j +w_{ij}v_i)\}} \\ &=& \sigma (c_j +w_{ij}v_i) \end{array} P(hj​=1∣v)​====​P(hj​=0∣v)+P(hj​=1∣v)P(hj​=1∣v)​1+exp{cj​+wij​vi​}exp{cj​+wij​vi​}​1+exp{−(cj​+wij​vi​)}1​σ(cj​+wij​vi​)​
因此，$P(h_j=1|v)$实际上就就是sigmoid激活函数，同理，$P(v_j=1|h)$也是sigmoid激活函数，其公式如下：
$$P(v_i=1|h)=\sigma (b_i+w_{ij}{h}_i)$$
有了激活函数，我们就可以从可见层和参数推导出隐藏层的神经元的取值概率了。

（3）CD算法求解举例

如下图所示（图片来自图解深度学习），我们以此为例讲解对比散度算法（Contrastive Divergence，CD）的过程。

如何判断值为0还是1：若根据上面的公式求出$P(h_j=1|v)=0.7$，现在我们取一个0~1之间的随机数，如果随机数小于0.7，则我们认为$h_j=1$（条件概率求出$h_j=1$的概率是0.7，随机数小于0.7的概率也是0.7），否则$h_j=0$。

(A)如上图(a)所示，$v_1^{(0)},v_2^{(0)},v_i^{(0)},v_n^{(0)}$的取值分别为$1,0,1,1$。

(B) 我们根据公式$P(h_j=1|v)=\sigma (c_j+w_{ij}{v}_i)$可以求出$h_1^{(0)},h_j^{(0)},h_m^{(0)}$的概率分别为$0.8,0.5,0.9$。

(C ) $h_1^{(0)},h_j^{(0)},h_m^{(0)}$的概率分别为$0.8,0.5,0.9$，根据上述判断结果是0还是1的方法，$h_1^{(0)},h_j^{(0)},h_m^{(0)}$的取值分别为$0,1,1$。

(D) $h_1^{(0)},h_j^{(0)},h_m^{(0)}$的取值分别为$0,1,1$，根据$P(v_i=1|h)=\sigma (b_i+w_{ij}{h}_i)$我们求得$v_1^{(1)},v_2^{(1)},v_i^{(1)},v_n^{(1)}$的概率分别为$0.7,0.2,0.8,0.6$。

(E)$v_1^{(1)},v_2^{(1)},v_i^{(1)},v_n^{(1)}$的概率分别为$0.7,0.2,0.8,0.6$，根据上述判断结果是0还是1的方法，$v_1^{(1)},v_2^{(1)},v_i^{(1)},v_n^{(1)}$的取值分别为$1,1,0,1$。

(F)梯度的更新公式为：
$$\begin{array}{l}{w}_{ij}\leftarrow w_{ij}+\eta [p(h_j^{(0)}=1|v^{(0)})v_i^{(0)}-p(h_j^{(1)}=1|v^{(1)})v_i^{(1)}]\\ b_i\leftarrow b_i+\eta (v_i^{(0)}-v_i^{(1)})\\ c_j\leftarrow c_j+\eta [p(h_j^{(0)}=1|v^{(0)})-p(h_j^{(1)}=1|v^{(1)})]\end{array}$$

三、深度信念网络

深度信念网络（Deep Belief Network, DBN）是由受限玻尔兹曼机堆叠而成，深度信念网络和其他神经网络的不同是：

1）神经网络的权重是先使用前向传播求得网络目标结果，再根据反向传播算法，把误差传播到下一层，调整所有的连接权重和偏置。
2）深度信念网络是逐层来调整连接权重和偏置项。具体做法是受限调整输入层和隐藏层之间的参数，把训练后得到的参数作为下一层的输入，再调整该层与下一个隐藏层之间的参数。不断迭代。

下图是DBN的示意图（图片来源：Deep Belief Networks (DBNs)）：

四、应用

1、预训练
在预训练任务中，可以使用DBM逐层训练参数，因此当做其他任务的参数初始化结果。
2、生成模型和分类模型
DBN可以当做生成模型来使用，除此之外，也可以当做分类模型，在DBN的最顶层增加一层softmax层，即可变成分类模型。
3、降维
隐藏层可以看做是降维的结果
4、和自编码器关系
和自编码器关联极大，这点会在自编码器中介绍二者的异同。

参考文章：
[1] 深度学习读书笔记之RBM（限制波尔兹曼机）
[2] 受限玻尔兹曼机和深度置信网络
[3] 受限玻尔兹曼机（RBM）原理总结
[4] 受限玻尔兹曼机基础教程
[5] 深度信念网络
[6] Hinton G E, Osindero S. A Fast Learning Algorithm for Deep Belief Nets[J].
[7] 深度学习（书）
[8] 图解深度学习（书）
[9] 受限玻尔兹曼机（RBM）学习笔记（六）对比散度算法