Essayez d'implémenter RBM avec chainer.

0. Introduction

RBM est une méthode utilisée comme pré-apprentissage lorsque plusieurs couches sont utilisées dans le Deep Learning. Le code source est publié sur theano. http://deeplearning.net/tutorial/rbm.html Cette fois, j'ai porté ce code sur le chainer.

Dans le cas de RBM, la fonction d'erreur est une fonction de l'énergie libre et est compliquée (?). En raison de la propagation d'erreur inverse, si vous différenciez avec la matrice de poids et les variables de biais, vous obtiendrez une belle formule de calcul. Il peut être calculé à partir de la valeur du traitement de la divergence contrastive (CD). Lors de la mise en œuvre avec chainer, plutôt que de créer une fonction d'erreur (perte) et d'effectuer une propagation d'erreur inverse (loss.backward) Il serait plus simple de simplement l'implémenter avec cupy (numpy + GPU), mais j'ai osé l'implémenter en utilisant la fonction d'erreur. Le code source de cupy (numpy + GPU) ne sera publié que lorsque l'opportunité se présentera. Cependant, le code source suivant utilise également cupy pour accélérer le traitement des CD. .. ..

# -*- coding: utf-8 -*-

import os,sys
import argparse
import time

import numpy as np

import chainer
from chainer import computational_graph
from chainer import cuda, Variable

import chainer.functions as F
import chainer.links as L
import chainer.optimizers as O


import pickle
import gzip

parser = argparse.ArgumentParser()
parser.add_argument('--gpu', '-g', default=-1, type=int,
                    help='GPU ID (negative value indicates CPU)')
parser.add_argument('--pcd', '-p', default=1, type=int,
                    help='pcd_flag')
parser.add_argument('--kcd', '-k', default=1, type=int,
                    help='cd-k')


'''Lors du calcul avec GPU, cupy= numpy +Processus avec GPU.'''
args = parser.parse_args()
if args.gpu >= 0:
    cuda.check_cuda_available()
xp = cuda.cupy if args.gpu >= 0 else np

def sigmoid(x):
    return 1. / (1 + xp.exp(-x))

class RBM(chainer.Chain):
    '''
    n_visbile:Dimension de couche visible
    n_hidden:Dimension du calque masqué
    l.W:Matrice de poids
    l.b:hbias
    l.a:vbaias
    real:La couche visible est 0,Les observations sont divisées selon qu'il s'agit de 1 ou d'un nombre réel. La couche initiale suppose un nombre réel, mais la couche intermédiaire est 0,Je ne pense qu'à un seul.
h en chaîne= l(v)La matrice de poids qui apparaît en général RBM est transposée de manière à se multiplier.
    '''
    def __init__(self,n_visible,n_hidden,real=0,k=1,pcd_flag=0):
        super(RBM,self).__init__(
            l=L.Linear(n_visible,n_hidden),
            # chainer 1.Cette utilisation était obsolète dans 5 manuels, alors ajoutez_Utilisez param.
#            a=L.Parameter(xp.zeros(n_visible,dtype=xp.float32)),
        )
        self.l.add_param("a",(n_visible),dtype=xp.float32)

        ''' l.N'apprend-il pas à moins que a ne soit initialisé et initialisé à l'aide du chainer?'''
        self.l.a.data.fill(0)

        self.n_visible = n_visible
        self.n_hidden = n_hidden
        self.real = real
        self.k = k
        self.pcd_flag = pcd_flag

    def __call__(self,v_data,v_prev_data=None):

        batch_size = v_data.shape[0]

        if self.pcd_flag == 0:
            v_prev_data = v_data
        elif self.pcd_flag ==1 and v_prev_data is None:
            v_prev_data = v_data
        vh_data = self.constrastive_divergence(v_prev_data,self.k)
        v = Variable(v_data)
        vh = Variable(vh_data.astype(np.float32))

        '''
        http://deeplearning.net/tutorial/rbm.expression html(5)De
Fonction d'erreur(loss)Est(Énergie libre des données de couche visible)Quand(CD de données de couche visible-k énergie libre)の差Quandなる。
        loss = (self.free_energy(v) - self.free_energy(vh)) / batch_size
        v->Le mappage vers vh n'est pas inclus dans l'apprentissage, donc le CD-Variable après que k traitement est terminé
        '''

        loss = (self.free_energy(v) - self.free_energy(vh)) / batch_size
        return loss

    def free_energy(self,v):
        ''' Function to compute the free energy '''
        '''
La valeur d'entrée v doit être variable
À l'origine, après avoir pris SUM en unités de ligne, unités de colonne(Nombre de lots)C'est un processus qui devrait prendre SUM
Après tout, je vais prendre SUM, donc SUM tout à la fois
        '''
        batch_size = v.data.shape[0]
        n_visible = self.n_visible
        real = self.real
        if real == 0:
            '''
La couche visible[0,1]Quand
            vbias_term = -1 * SUM((a(i) * v(i))
            '''
            vbias_term = F.sum(F.matmul(v,self.l.a))
        else:
            '''
Lorsque la couche visible est un nombre réel
            vbias_term = -0.5 * SUM((v(i)-a(i)) * (v(i)-a(i)))
Nombre de lots dans le chainer*Pour gérer dans la couche d'entrée, pour soustraire le biais de chaque ligne
Utiliser de force m * n
            '''
            m = Variable(xp.ones((batch_size,1),dtype=xp.float32))
            n = F.reshape(self.l.a,(1,n_visible))
            v_ = v - F.matmul(m,n)
            vbias_term = -F.sum(0.5 * v_ * v_)

        wx_b = self.l(v)
        hidden_term = F.sum(F.log(1+F.exp(wx_b)))
        return -vbias_term-hidden_term

    def propup(self,vis):
        '''
Calculer la distribution de probabilité des couches cachées à partir des données de couches visibles
La couche cachée[0,1]Donc, il renvoie la probabilité d'être 1.
La valeur d'entrée vis est toujours xp(np)
        '''
        pre_sigmoid_activation = xp.dot(vis,self.l.W.data.T) + self.l.b.data
        return sigmoid(pre_sigmoid_activation)

    def propdown(self,hid):
        '''
Calculer la distribution de probabilité de la couche visible à partir des données de la couche masquée
La couche visible[0,1]Quand est, la probabilité de devenir 1 est renvoyée.
Lorsque la couche visible est un nombre réel, elle renvoie la moyenne de la distribution normale. La distribution est fixée à 1.
La valeur d'entrée masquée est toujours xp(np)
        '''
        real = self.real
        if real == 0:
            pre_sigmoid_activation = xp.dot(hid,self.l.W.data) + self.l.a.data
            v_mean = sigmoid(pre_sigmoid_activation)
        else:
            v_mean = xp.dot(hid,self.l.W.data) + self.l.a.data
        return v_mean

    def sample_h_given_v(self,v0_sample):
        '''
        v0_sample → h1_Échantillon d'échantillonnage Gibbs
        [0,1]Créez une valeur aléatoire de et comparez-la à la probabilité de 1 calculée par pop up
        h1_mean[i] >  p[i]Puis h1_sample[i] = 1
        h1_mean[i] <= p[i]Puis h1_sample[i] = 0
        '''
        h1_mean = self.propup(v0_sample)
        h1_sample = xp.random.binomial(size=h1_mean.shape,n=1,p=h1_mean)
        return h1_mean,h1_sample

    def sample_v_given_h(self,h0_sample):
        '''
        h0_sample → v1_Échantillon d'échantillonnage Gibbs
La couche visible[0,1]dans le cas de
        [0,1]La valeur de est une valeur aléatoire p[i]Et comparez-le à la probabilité de 1 calculée par popdown
        v1_mean[i] >  p[i]Puis v1_sample[i] = 1
        v1_mean[i] <= p[i]Puis v1_sample[i] = 0
Lorsque la couche visible est un nombre réel
        v1_sample[i]Est moyenne v1_mean[i], Parce que c'est une distribution normale avec variance 1
Valeur aléatoire u de la distribution normale avec moyenne 0 et variance 1[i]Et ajoutez la valeur calculée par popdown
        v1_sample[i] = v1_mean[i] + u[i]
        '''
        v1_mean = self.propdown(h0_sample)
        if self.real == 0:
            v1_sample = xp.random.binomial(size=v1_mean.shape,n=1,p=v1_mean)
        else:
            batch_number = h0_sample.shape[0]
            v1_sample = v1_mean + xp.random.randn(batch_number,self.n_visible)
        return v1_mean,v1_sample

    def gibbs_hvh(self,h0_sample):
        ''' h->v->échantillonnage gibbs h'''
        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 gibbs_vhv(self,v0_sample):
        ''' v->h->échantillonnage v gibbs'''
        h1_mean,h1_sample = self.sample_h_given_v(v0_sample)
        v1_mean,v1_sample = self.sample_v_given_h(h1_sample)
        return h1_mean,h1_sample,v1_mean,v1_sample

    def constrastive_divergence(self,v0_sample,k=1):
        ''' CD-le traitement de k, cupy est nécessaire pour en faire un GPU'''
        vh_sample = v0_sample
        for step in range(k):
            ph_mean,ph_sample,vh_mean,vh_sample = self.gibbs_vhv(vh_sample)
        return vh_sample

    def reconstruct(self, v):
        h = sigmoid(xp.dot(v,self.l.W.data.T) + self.l.b.data)
        reconstructed_v = sigmoid(xp.dot(h,self.l.W.data) + self.l.a.data)
        return reconstructed_v

Dans le cas de chainer, vous pouvez utiliser chainer.links.Linear pour définir la matrice de poids (W) et le biais (b) de la cartographie de chaque couche. Dans le cas de la GAR, un biais supplémentaire est nécessaire. l.add_param("a",(n_visible),dtype=xp.float32) J'ai ajouté un paramètre dans. l.W: matrice de poids l.b:hbias l.a:vbaias Ce sera. Pour gérer le mappage de la couche visible à la couche cachée avec chainer.links.Linear Il s'agit d'une transposition de la matrice de poids décrite dans un livre général RBM. (Wij → Wji)

1. Expérimentez avec MNIST

J'ai expérimenté le code source suivant avec des données MNIST. Puisqu'il n'y avait pas de serveur GPU, il a été formé par le CPU, donc le nombre de formations est de 50 fois. Cette fois, la couche masquée est définie sur 500 dimensions et les résultats avec un coefficient d'apprentissage de 0,01, une impulsion de 0,5 et un PCD-10 sont utilisés.

Voici les données originales. En fait, il y a 60000 cas, mais seulement 150 cas sont affichés

Ce qui suit est une image de l'image ci-dessus reconstruite. Il est presque reproduit. Seuls les éléments sont affichés.

Vient ensuite une image qui visualise la matrice de poids. Le calque caché est de 500 dimensions par rapport à la base de chaque vecteur Il est converti en 28 * 28 images.