受限玻尔兹曼机RBM

相关算法

python代码参考http://blog.csdn.net/zc02051126/article/details/9668439#（作少量修改与注释）

 #coding:utf8

 import matplotlib.pylab as plt

 import numpy as np

 import cPickle

 class RBM:

     def __init__(self,n_visul, n_hidden, max_epoch = 50, batch_size = 110, penalty = 2e-4):

         self.n_visible = n_visul

         self.n_hidden = n_hidden

         self.max_epoch = max_epoch

         self.batch_size = batch_size

         self.penalty = penalty

         self.w = np.random.random((self.n_visible, self.n_hidden)) * 0.1

         self.v_bias = np.zeros((1, self.n_visible))

         self.h_bias = np.zeros((1, self.n_hidden))

     def sigmoid(self, z):

         return 1.0 / (1.0 + np.exp( -z ))

     def forward(self, vis):

         return self.sigmoid(np.dot(vis.T, self.w) + self.h_bias)

     def backward(self, vis):

         return self.sigmoid(np.dot(vis, self.w.T) + self.v_bias)

     def batch(self):

         d, N = self.x.shape

         num_batchs = int(round(N / self.batch_size)) + 1

         groups = np.ravel(np.repeat([range(0, num_batchs)], self.batch_size, axis = 0))

         groups=groups[:N]

         np.random.shuffle(groups)

         batch_data = []

         for i in range(0, num_batchs):

             index = groups == i

             batch_data.append(self.x[:, index])

         return batch_data

     def rbmBB(self, x):

         self.x = x

         eta = 0.1

         momentum = 0.5  #动量项

         W = self.w

         b = self.h_bias

         c = self.v_bias

         Winc  = np.zeros((self.n_visible, self.n_hidden))

         binc = np.zeros(self.n_hidden)

         cinc = np.zeros(self.n_visible)

         batch_data = self.batch()

         num_batch = len(batch_data)

         errors = []

         for epoch in range(0, self.max_epoch):

             err_sum = 0.0

             for batch in range(0, num_batch):

                 num_dims, num_cases = batch_data[batch].shape

                 data = batch_data[batch]

                 # 已知可见层，采样出隐藏层

                 ph = self.forward(data)

                 ph_states = np.zeros((num_cases, self.n_hidden))

                 ph_states[ph > np.random.random((num_cases, self.n_hidden))] = 1

                 # 已知隐藏层，采样出可见层

                 neg_data = self.backward(ph_states)

                 neg_data_states = np.zeros((num_cases, num_dims))

                 neg_data_states[neg_data > np.random.random((num_cases, num_dims))] = 1

                 neg_data_states = neg_data_states.transpose()

                 nh = self.forward(neg_data_states)

                 # CD算法

                 dW = np.dot(data, ph) - np.dot(neg_data_states, nh)

                 dc = np.sum(data, axis = 1) - np.sum(neg_data_states, axis = 1)

                 db = np.sum(ph, axis = 0) - np.sum(nh, axis = 0)

                 # 刷新参数

                 Winc = momentum * Winc + eta * (dW / num_cases - self.penalty * W)

                 binc = momentum * binc + eta * (db / num_cases);

                 cinc = momentum * cinc + eta * (dc / num_cases);

                 W = W + Winc

                 b = b + binc

                 c = c + cinc

                 self.w = W

                 self.h_bais = b

                 self.v_bias = c

                 err = np.linalg.norm(data - neg_data.transpose())

                 err_sum += err

             print epoch, err_sum

             errors.append(err_sum)

         self.errors = errors

         self.hiden_value = self.forward(self.x)

         h_row, h_col = self.hiden_value.shape

         hiden_states = np.zeros((h_row, h_col))

         hiden_states[self.hiden_value > np.random.random((h_row, h_col))] = 1

         self.rebuild_value = self.backward(hiden_states)

     def visualize(self, X):  #可视化

         D, N = X.shape

         s = int(np.sqrt(D))

         num = int(np.ceil(np.sqrt(N)))

         a = np.zeros((num*s + num + 1, num * s + num + 1)) - 1.0

         x = 0

         y = 0

         for i in range(0, N):

             z = X[:,i]

             z = z.reshape(s,s,order='F')

             z = z.transpose()

             a[x*s+x:x*s+s+x , y*s+y:y*s+s+y] = z

             x = x + 1

             if(x >= num):

                 x = 0

                 y = y + 1

         return a

 def readData(path):

     data = []

     for line in open(path, 'r'):

         ele = line.split(' ')

         tmp = []

         for e in ele:

             if e != '':

                 tmp.append(float(e.strip(' ')))

         data.append(tmp)

     return data

 if __name__ == '__main__':

     f = open('mnist.pkl', 'rb')

     training_data, validation_data, test_data = cPickle.load(f)

     training_inputs = [np.reshape(x, 784) for x in training_data[0]]

     data =training_inputs[:5000]

     data = np.array(data)

     data = data.transpose()

     rbm = Rbm(784, 100,max_epoch = 50)

     rbm.rbmBB(data)

     a = rbm.visualize(data)  #(2060L, 2060L)

     fig = plt.figure(1)

     ax = fig.add_subplot(111)

     ax.imshow(a)

     plt.title('original data')

     rebuild_value = rbm.rebuild_value.transpose()

     b = rbm.visualize(rebuild_value)  #(2060L, 2060L)

     fig = plt.figure(2)

     ax = fig.add_subplot(111)

     ax.imshow(b)

     plt.title('rebuild data')

     hidden_value = rbm.hiden_value.transpose()

     c = rbm.visualize(hidden_value)  #(782L, 782L)

     fig = plt.figure(3)

     ax = fig.add_subplot(111)

     ax.imshow(c)

     plt.title('hidden data')

     w_value = rbm.w

     d = rbm.visualize(w_value)  #(291L, 291L)

     fig = plt.figure(4)

     ax = fig.add_subplot(111)

     ax.imshow(d)

     plt.title('weight value(w)')

     plt.show()