本篇简要介绍深度信念网络DBN，是一个不太常见的神经网络。

首先对深度信念网络（DBN）进行简要介绍，然后对其组成原件RBM的结构、原理和训练过程进行介绍。接着对DBN的训练过程进行介绍。

目录

一、DBN概述

二、首先玻尔兹曼机（RBM）

1. RBM基本结构

2. 神经元之间的条件独立性

3. 使用RBM的过程

三、深度信念网络DBN

1. DBN训练过程

 2. 调优过程

四、将DBN用于分类任务

1. 用于BP网络分类

2. DBN的拓展

总结：

---

一、DBN概述

深度信念网络（Deep Belief Nets，DBN），既可用于非监督学习，类似于一个自编码机；也可以用于监督学习，作为分类器使用。

从非监督学习来讲，其目的是尽可能地保留原始特征的特点，同时降低特征的维度。从监督学习来讲，提高分类的准确率。二者都是对特征进行学习的过程，即如何得到更好的特征表达。

作为神经网络，DBN由若干神经元组成，组成原件是受限玻尔兹曼机（RBM）。训练DBN的过程是一层一层进行的，每一层中，用数据向量推断隐层，再把这一隐层当作下一层的数据向量。

二、首先玻尔兹曼机（RBM）

1. RBM基本结构

RBM有两层神经元，一层叫做显层（visible layer），用于输入训练数据；一层叫做隐层（hidden layer），用于提取特征，也称特征检测器。

上图为RBM的结构。每一层都可以用一个向量表示，每一维表示一个神经元。层间是双向（兑成）全连接，层内不存在连接。隐层单元被训练去捕捉在可视层表现出来的高阶数据的相关性。

2. 神经元之间的条件独立性

由于显层和隐层内部神经元没有互联，只有层间是对称全连接，这样的好处是：

（1）在给定所有显元的值的情况下，每一个隐元取值是互不相关的。同理，在给定所有隐元下，显元取值不相关；

（2）这样，计算每个神经元取值就可以并行计算。

3. 使用RBM的过程

假设我们已经得到一个训练好的RBM，隐元和显元之间的权重用矩阵W表示：

 M代表显元个数，N代表隐元个数。

当我们将新来的数据赋给显层后，RBM依照权值W决定开启或关闭隐元，具体操作如下：

（1）首先，（并行）计算出每个隐元的激励值；’

（2）将每个隐元的激励值用S形函数进行标准化，变成它们处于开启状（用1表示）的概率值；

（3）由此我们得到每个隐元的开启和关闭概率，此时我们将开启的概率与从一个[0，1]均匀分布中的随机抽取值进行比较，决定每个神经元开启或者关闭的概率；

（4）给定隐层，显层的计算概率一样的。

这时我们就要进行训练RBM。它实际上是求一个最能产生训练样本的概率分布，也就是说，要求一个分布，在这个分布里，训练样本的概率最大。

 G.Hitton提出了名为对比散度的学习算法，简述如下：

首先，对于训练集的每一条记录x，将其赋给显层，计算隐层神经元被开启的概率，然后从（隐层）计算出的概率中抽取一个样本，用来重构显层，再次使用重构过的显层计算隐层神经元概率，按公式更新权重矩阵W。

上述反复向上和向下的操作称为吉布斯采样（Gibbs）。隐层激活单元和可视层输入之间的相关性差别就作为权值更新的主要依据。

如此训练之后的RBM就能较为准确的提取显层的特征了，或者根据隐层所代表的特征还原显层。

三、深度信念网络DBN

1. DBN训练过程

DBN是由多层RBM组成的一个神经网络，它既可以被看作一个生成模型，也可以看作一个判别模型（最后一层加label进行BP算法输出训练）。

其训练过程是，使用非监督贪婪逐层方法预训练获得权值。训练过程如下：

（1）首先充分训练第一个RBM；

（2）固定第一个RBM的权重和偏移量，然后使用其隐性神经元的状态，作为第二个RBM的输入向量；

（3）充分训练第二个RBM后，将第二个RBM堆叠在第一个RBM的上方；

（4）重复以上三个步骤任意多次；

（5）如果训练集中的数据有标签，那么在顶层的RBM训练时，这个RBM的显层中除了显性神经元，还需要有代表分类标签的神经元，一起进行训练：

          ——假设顶层 RBM 的显层有 500 个显性神经元，训练数据的分类一共分成了 10 类；

          ——那么顶层 RBM 的显层有 510 个显性神经元，对每一训练训练数据，相应的标签神经元被打开设为1，而其他的则被关闭设为0。

训练结果如下图：

 2. 调优过程

生成模型使用Contrastive Wake-Sleep算法进行调优，其算法过程是：

（1）除了顶层RBM，其他层RBM的权重被分成向上的认知权重和向下的生成权重；

（2）Wake阶段：认知过程，通过外界的特征和向上的权重 (认知权重) 产生每一层的抽象表示 (结点状态) ，并且使用梯度下降修改层间的下行权重 (生成权重) ；

（3）Sleep阶段：生成过程，通过顶层表示 (醒时学得的概念) 和向下权重，生成底层的状态，同时修改层间向上的权重。

使用过程 ：

（1）使用随机隐性神经元状态值，在顶层 RBM 中进行足够多次的吉布斯抽样；

（2）向下传播，得到每层的状态。

四、将DBN用于分类任务

1. 用于BP网络分类

RBM的训练特点导致训练时间显著减少，因为只需要单个步骤就可以接近最大似然学习。网络的每一层都可会改进训练数据的对数概率。

在预训练后，DBN可以通过利用带标签数据用BP算法去对判别性能做调整。在这里，一个标签集将被附加到顶层（推广联想记忆），通过一个自下向上的，学习到的识别权值获得一个网络的分类面。

这个性能会比单纯的BP算法训练的网络好。这可以很直观的解释，DBNs的BP算法只需要对权值参数空间进行一个局部的搜索，这相比前向神经网络来说，训练是要快的，而且收敛的时间也少。

2. DBN的拓展

（1）用于二维数据的卷积DBNs

DBNs并没有考虑到图像的2维结构信息，因为输入是简单的从一个图像矩阵一维向量化的。而CDBNs就是考虑到了这个问题，它利用邻域像素的空域关系，通过一个称为卷积RBMs的模型区达到生成模型的变换不变性，而且可以容易得变换到高维图像。

（2）用于时间序列的堆叠时间RBMs

（3）DBN中RBM的替换

目前，和DBNs有关的研究包括堆叠自动编码器，它是通过用堆叠自动编码器来替换传统DBNs里面的RBMs。

这就使得可以通过同样的规则来训练产生深度多层神经网络架构，但它缺少层的参数化的严格要求。与DBNs不同，自动编码器使用判别模型，这样这个结构就很难采样输入采样空间，这就使得网络更难捕捉它的内部表达。

但是，降噪自动编码器却能很好的避免这个问题，并且比传统的DBNs更优。它通过在训练过程添加随机的污染并堆叠产生场泛化性能。训练单一的降噪自动编码器的过程和RBMs训练生成模型的过程一样。

总结：

（1）注意RBM层间的双向全连接和层内无连接，由此有概率独立性和并行计算的特点；

（2）熟悉单个RBM的训练算法流程，即对比散度算法；简单来说就是根据显层的初始值计算隐层，用隐层的样本重构显层，再用显层计算一遍隐层，将两遍计算隐层的结果结合起来，得到权值矩阵W；

（3）同时，重构后的显层也可以向下传播，被随机选择后去重构原始的输入信号（因为它既是本层的显层，也是下层的隐层）。上面这种反复过程称为吉布斯采样（Gibbs）；

（3）了解DBN的训练过程，即逐层训练；简单来说就是训练好一个RBM，将这个隐层作为下一个RBM的显层，训练好第二个，如此循环，一直训练好整个网络；

（4）如果要做分类器使用的话，就在顶层的RBM显层加上分类概率神经元即可；

（5）了解一些DBN的拓展。

三、DBN

3.1 生成模型

深度信念网络是一个生成模型，用来生成符合特定分布的样本。隐变量用来描述在可观测变量之间的高阶相关性。假如加入服从分布 𝑝(𝑣)的训练数据，通过训练得到一个深度信念网络。

生成样本时，先在最顶两层进行足够多的吉布斯采样，在达到热平衡时生成样本ℎ^((𝐿−1))，然后依次计算下一层隐变量的分布。因为在给定上一层变量取值时，下一层的变量是条件独立的，故可独立采样。这样，从第𝐿−1层开始，自顶向下进行逐层采样，最终得到可观测层的样本。

3.2 参数学习

深度信念网络最直接的训练方式是最大化可观测变量的边际分布𝑝(𝑣)在训练集上的似然 。但是在深度信念网络中，隐变量ℎ之间的关系十分复杂，由于“贡献度分配问题”，很难直接学习。即使对于简单的单层Sigmoid信念网络：

$$p(v=1\mid h)=\sigma \left(b+{\omega }^{T}h\right)$$

在已知可观测变量时，其隐变量的联合后验概率𝑝(ℎ|𝑣)不再相互独立，因此很难估计所有隐变量的后验概率，早期深度信念网络的后验概率一般通过蒙特卡洛方法或变分方法来近似估计，但效率低，从而导致其参数学习比较困难。

为了有效训练深度信念网络，我们将每一层的Sigmoid信念网络转换为受限玻尔兹曼机，这样做的好处是隐变量的后验概率事相互独立的，从而可容易进行采样。这样，深度信念网络可看做多个受限玻尔兹曼机从下到上进行堆叠，第𝑙层受限玻尔兹曼机的隐层作为第𝑙+1受限玻尔兹曼机的可观测层。进一步，深度信念网络可采用逐层训练的方式来快速训练，即从最底层开始，每次只训练一层，直到最后一层。

深度信念网络的训练过程可分为逐层预训练和精调两个阶段，先通过逐层预训练将模型的参数初始化为较优的值，然后通过传统机器学习方法对参数进行精调。

3.2 逐层预训练

采用逐层训练的方式，将深度信念网络的训练简化为对多个受限玻尔兹曼机的训练。具体的逐层训练过程为自下而上依次训练每一层的首先玻尔兹曼机。假设已训练好前𝑙−1层的受限玻尔兹曼机，可计算隐变量自下而上的条件概率：

$$p\left(h^{(i)}\mid h^{(i-1)}\right)=\sigma \left(b^i+W^{(i)}{h}^{(i-1)}\right),1\le i\le (l-1)$$

这样可按照$v=h^{(0)}\to \cdots \to h^{(l-1)}$的顺序生成一组$h^{l-1}$的样本，记为 $H^{(l-1)}=h^{(l,1)},...,h^{(l,M)}$。然后将$h^{(l-1)}$ 和$h^{(l)}$组成一个受限玻尔兹曼机，用$H^{(l-1)}$作为训练集充分训练第𝑙层的受限玻尔兹曼机.

大量实践表明，逐层预训练可以产生非常好的参数初始值，从而极大地降低了模型的学习难度。

3.3 精调

经过预训练，再结合具体的任务（监督或无监督学习），通过传统的全局学习算法对网络进行精调（fine-tuning），使模型收敛到更好的局部最优点。

作为生成模型的精调：除了顶层的受限玻尔兹曼机，其它层之间的权重可以被分为向上的认知权重（Recognition Weight）$W^{\ast }$和向下的生成权重（Generative Weight）𝑊。认知权重用来计算后验概率，生成权重用来定义模型，认知权重初始值 𝑊^{(‘𝑙)=𝑊}(𝑙)𝑇

深度信念网络一般采用 Contrastive Wake-Sleep算法进行精调：

1.Wake阶段：认知过程，通过外界输入（可观测变量）和向上的认知权重，计算每一层隐变量的后验概率并采样。修改下行的生成权重使得下一层的变量的后验概率最大。

2.Sleep阶段：生成过程，通过顶层的采样和向下的生成权重，逐层计算每一层的后验概率并采样。然后，修改向上的认知权重使得上一层变量的后验概率最大。

3.交替进行Wake和Sleep过程，直到收敛

作为判别模型的精调 ：深度信念网络的一个应用是作为深度神经网络的预训练模型，提供神经网络的初始权重，这时只需要向上的认知权重，作为判别模型使用：

四 代码演示

import sys
import numpy


numpy.seterr(all='ignore')
 
def sigmoid(x):
    return 1. / (1 + numpy.exp(-x))
 
def softmax(x):
    e = numpy.exp(x - numpy.max(x))  # prevent overflow
    if e.ndim == 1:
        return e / numpy.sum(e, axis=0)
    else:  
        return e / numpy.array([numpy.sum(e, axis=1)]).T  # ndim = 2


class DBN(object):
    def __init__(self, input=None, label=None,\
                 n_ins=2, hidden_layer_sizes=[3, 3], n_outs=2,\
                 numpy_rng=None):
        
        self.x = input
        self.y = label

        self.sigmoid_layers = []
        self.rbm_layers = []
        self.n_layers = len(hidden_layer_sizes)  # = len(self.rbm_layers)

        if numpy_rng is None:
            numpy_rng = numpy.random.RandomState(1234)

        
        assert self.n_layers > 0


        # construct multi-layer
        for i in range(self.n_layers):
            # layer_size
            if i == 0:
                input_size = n_ins
            else:
                input_size = hidden_layer_sizes[i - 1]

            # layer_input
            if i == 0:
                layer_input = self.x
            else:
                layer_input = self.sigmoid_layers[-1].sample_h_given_v()
                # print('=============')
                
            # construct sigmoid_layer
            sigmoid_layer = HiddenLayer(input=layer_input,
                                        n_in=input_size,
                                        n_out=hidden_layer_sizes[i],
                                        numpy_rng=numpy_rng,
                                        activation=sigmoid)
            self.sigmoid_layers.append(sigmoid_layer)


            # construct rbm_layer
            rbm_layer = RBM(input=layer_input,
                            n_visible=input_size,
                            n_hidden=hidden_layer_sizes[i],
                            W=sigmoid_layer.W,     # W, b are shared
                            hbias=sigmoid_layer.b)
            self.rbm_layers.append(rbm_layer)


        # layer for output using Logistic Regression
        self.log_layer = LogisticRegression(input=self.sigmoid_layers[-1].sample_h_given_v(),
                                            label=self.y,
                                            n_in=hidden_layer_sizes[-1],
                                            n_out=n_outs)

        # finetune cost: the negative log likelihood of the logistic regression layer
        self.finetune_cost = self.log_layer.negative_log_likelihood()



    def pretrain(self, lr=0.1, k=1, epochs=100):
        # pre-train layer-wise
        # print('pre-training')
        for i in range(self.n_layers):
            if i == 0:
                layer_input = self.x
            else:
                layer_input = self.sigmoid_layers[i-1].sample_h_given_v(layer_input)
                # print(layer_input)
            rbm = self.rbm_layers[i]
            
            for epoch in range(epochs):
                rbm.contrastive_divergence(lr=lr, k=k, input=layer_input)
                # cost = rbm.get_reconstruction_cross_entropy()
                # print >> sys.stderr, \
                #        'Pre-training layer %d, epoch %d, cost ' %(i, epoch), cost

    # def pretrain(self, lr=0.1, k=1, epochs=100):
    #     # pre-train layer-wise
    #     for i in range(self.n_layers):
    #         rbm = self.rbm_layers[i]
            
    #         for epoch in range(epochs):
    #             layer_input = self.x
    #             for j in range(i):
    #                 layer_input = self.sigmoid_layers[j].sample_h_given_v(layer_input)
            
    #             rbm.contrastive_divergence(lr=lr, k=k, input=layer_input)
    #             # cost = rbm.get_reconstruction_cross_entropy()
    #             # print >> sys.stderr, \
    #             #        'Pre-training layer %d, epoch %d, cost ' %(i, epoch), cost


    def finetune(self, lr=0.1, epochs=100):
        # print('finetune')
        layer_input = self.sigmoid_layers[-1].sample_h_given_v()

        # train log_layer
        epoch = 0
        done_looping = False
        while (epoch < epochs) and (not done_looping):
            self.log_layer.train(lr=lr, input=layer_input)
            # self.finetune_cost = self.log_layer.negative_log_likelihood()
            # print >> sys.stderr, 'Training epoch %d, cost is ' % epoch, self.finetune_cost
            
            lr *= 0.95
            epoch += 1


    def predict(self, x):
        layer_input = x
        
        for i in range(self.n_layers):
            sigmoid_layer = self.sigmoid_layers[i]
            # rbm_layer = self.rbm_layers[i]
            layer_input = sigmoid_layer.output(input=layer_input)

        out = self.log_layer.predict(layer_input)
        return out



class HiddenLayer(object):
    def __init__(self, input, n_in, n_out,\
                 W=None, b=None, numpy_rng=None, activation=numpy.tanh):
        
        if numpy_rng is None:
            numpy_rng = numpy.random.RandomState(1234)

        if W is None:
            a = 1. / n_in
            initial_W = numpy.array(numpy_rng.uniform(  # initialize W uniformly
                low=-a,
                high=a,
                size=(n_in, n_out)))

            W = initial_W

        if b is None:
            b = numpy.zeros(n_out)  # initialize bias 0


        self.numpy_rng = numpy_rng
        self.input = input
        self.W = W
        self.b = b

        self.activation = activation

    def output(self, input=None):
        if input is not None:
            self.input = input
        
        linear_output = numpy.dot(self.input, self.W) + self.b

        return (linear_output if self.activation is None
                else self.activation(linear_output))

    def sample_h_given_v(self, input=None):
        if input is not None:
            self.input = input

        v_mean = self.output()
        # print('v_mean:\n',v_mean)
        h_sample = self.numpy_rng.binomial(size=v_mean.shape,
                                           n=1,
                                           p=v_mean)
        # print('h_sample:\n',h_sample)
        return h_sample



class RBM(object):
    def __init__(self, input=None, n_visible=2, n_hidden=3, \
        W=None, hbias=None, vbias=None, numpy_rng=None):
        
        self.n_visible = n_visible  # num of units in visible (input) layer
        self.n_hidden = n_hidden    # num of units in hidden layer

        if numpy_rng is None:
            numpy_rng = numpy.random.RandomState(1234) # 使用RandomState获得随机数生成器。

        if W is None:
            a = 1. / n_visible
            initial_W = numpy.array(numpy_rng.uniform(  # initialize W uniformly
                low=-a,
                high=a,
                size=(n_visible, n_hidden)))  # 生成权重矩阵

            W = initial_W

        if hbias is None:
            hbias = numpy.zeros(n_hidden)  # initialize h bias 0

        if vbias is None:
            vbias = numpy.zeros(n_visible)  # initialize v bias 0


        self.numpy_rng = numpy_rng
        self.input = input
        self.W = W
        self.hbias = hbias
        self.vbias = vbias
        # self.params = [self.W, self.hbias, self.vbias]


    def contrastive_divergence(self, lr=0.1, k=1, input=None):
        if input is not None:
            self.input = input
        
        ''' CD-k '''
        ph_mean, ph_sample = self.sample_h_given_v(self.input)

        chain_start = ph_sample

        for step in range(k):
            if step == 0:
                nv_means, nv_samples,\
                nh_means, nh_samples = self.gibbs_hvh(chain_start)
            else:
                nv_means, nv_samples,\
                nh_means, nh_samples = self.gibbs_hvh(nh_samples)
        # nv_means, nv_samples, nh_means, nh_samples = self.gibbs_hvh(chain_start)
        # for _ in range(1, k):
        #         nv_means, nv_samples, nh_means, nh_samples = self.gibbs_hvh(nh_samples)

        # chain_end = nv_samples

        self.W += lr * (numpy.dot(self.input.T, ph_sample)
                        - numpy.dot(nv_samples.T, nh_means))
        self.vbias += lr * numpy.mean(self.input - nv_samples, axis=0)
        self.hbias += lr * numpy.mean(ph_sample - nh_means, axis=0)

        # cost = self.get_reconstruction_cross_entropy()
        # return cost


    def sample_h_given_v(self, v0_sample):
        h1_mean = self.propup(v0_sample)  # 
        h1_sample = self.numpy_rng.binomial(size=h1_mean.shape,   # discrete: binomial 二项分布
                                       n=1,
                                       p=h1_mean)

        return [h1_mean, h1_sample]


    def sample_v_given_h(self, h0_sample):
        v1_mean = self.propdown(h0_sample)
        v1_sample = self.numpy_rng.binomial(size=v1_mean.shape,   # discrete: binomial
                                            n=1,
                                            p=v1_mean)
        return [v1_mean, v1_sample]

    def propup(self, v):
        # 返回隐藏层被激活的概率，应该是一个矩阵，对应的是每一个神经元被激活的概率
        pre_sigmoid_activation = numpy.dot(v, self.W) + self.hbias
        return sigmoid(pre_sigmoid_activation)

    def propdown(self, h):
        pre_sigmoid_activation = numpy.dot(h, self.W.T) + self.vbias
        return sigmoid(pre_sigmoid_activation)

    def gibbs_hvh(self, h0_sample):
        v1_mean, v1_sample = self.sample_v_given_h(h0_sample)
        h1_mean, h1_sample = self.sample_h_given_v(v1_sample)
        return [v1_mean, v1_sample, h1_mean, h1_sample]
    

    def get_reconstruction_cross_entropy(self):
        pre_sigmoid_activation_h = numpy.dot(self.input, self.W) + self.hbias
        sigmoid_activation_h = sigmoid(pre_sigmoid_activation_h)
        
        pre_sigmoid_activation_v = numpy.dot(sigmoid_activation_h, self.W.T) + self.vbias
        sigmoid_activation_v = sigmoid(pre_sigmoid_activation_v)

        cross_entropy =  - numpy.mean(
            numpy.sum(self.input * numpy.log(sigmoid_activation_v) +
            (1 - self.input) * numpy.log(1 - sigmoid_activation_v),
                      axis=1))
        
        return cross_entropy

    def reconstruct(self, v):
        h = sigmoid(numpy.dot(v, self.W) + self.hbias)
        reconstructed_v = sigmoid(numpy.dot(h, self.W.T) + self.vbias)
        return reconstructed_v


class LogisticRegression(object):
    def __init__(self, input, label, n_in, n_out):
        self.x = input
        self.y = label
        self.W = numpy.zeros((n_in, n_out))  # initialize W 0
        self.b = numpy.zeros(n_out)          # initialize bias 0

    def train(self, lr=0.1, input=None, L2_reg=0.00):
        if input is not None:
            self.x = input

        p_y_given_x = softmax(numpy.dot(self.x, self.W) + self.b)
        d_y = self.y - p_y_given_x
        
        self.W += lr * numpy.dot(self.x.T, d_y) - lr * L2_reg * self.W
        self.b += lr * numpy.mean(d_y, axis=0)

    def negative_log_likelihood(self):
        sigmoid_activation = softmax(numpy.dot(self.x, self.W) + self.b)

        cross_entropy = - numpy.mean(
            numpy.sum(self.y * numpy.log(sigmoid_activation) +
            (1 - self.y) * numpy.log(1 - sigmoid_activation),
                      axis=1))

        return cross_entropy

    def predict(self, x):
        return softmax(numpy.dot(x, self.W) + self.b)




def test_dbn(pretrain_lr=0.1, pretraining_epochs=1000, k=1, \
             finetune_lr=0.1, finetune_epochs=200):

    x = numpy.array([[1,1,1,0,0,0],
                     [1,0,1,0,0,0],
                     [1,1,1,0,0,0],
                     [0,0,1,1,1,0],
                     [0,0,1,1,0,0],
                     [0,0,1,1,1,0]])
    y = numpy.array([[1, 0],
                     [1, 0],
                     [1, 0],
                     [0, 1],
                     [0, 1],
                     [0, 1]])

    
    rng = numpy.random.RandomState(123)

    # construct DBN
    # print('construct DBN')
    dbn = DBN(input=x, label=y, n_ins=6, hidden_layer_sizes=[3, 3], n_outs=2, numpy_rng=rng)

    # pre-training (TrainUnsupervisedDBN)
    # k是gibbs的次数
    
    dbn.pretrain(lr=pretrain_lr, k=1, epochs=pretraining_epochs)
    
    # fine-tuning (DBNSupervisedFineTuning)
    dbn.finetune(lr=finetune_lr, epochs=finetune_epochs)


    # test
    x = numpy.array([1, 1, 0, 0, 0, 0])
    print(dbn.predict(x))


if __name__ == "__main__":
    test_dbn()

[0.72344411 0.27655589]