[Deep Learning ]RBM trained with Contrastive Divergence

Boltzmann Machine은 [0,1]의 값을 갖는 binary unit들로 이루어진 network를 의미한다.

 위의 Figure1은 Boltzmann Machine을 안다면 누구나 한번쯤 봤을 그림이다. 먼저 왼쪽의 모형이 general BM이다. 이 BM의 특징은 full-connectivity에 있다. 그리고 오른쪽에 모형이 restricted BM이다. 이 모형은 visible node와 hidden node를 분리시켰다. 이것이 BM과 RBM의 차이이다. 이 간단한 차이로 RBM은 실제 구현이 가능하고, BM은 구현이 매우 어렵다. 

 앞서 설명하였듯이 BM에서 node는 0또는 1의 binary한 값을 갖는다. 그리고 각 node사이에는 symmetric하게 link가 있는데, 이 link에는 weight가 존재한다. 이 weight의 값은 굳이 양수일 필요 없이 모든 값을 가질 수 있다. RBM의 경우 각 node를 visible과 hidden으로 나눠 놓았고, 여기서 visible node는 우리의 data가 들어가는 곳을 의미하고, hidden의 경우 우리는 각 node가 1이 될 확률만을 알게된다. 

BM의 상태는 에너지를 통해서 설명될 수 있는데 엔트로피와 마찬가지로 에너지가 높을 수록 그 존재 확률이 낮아지게 된다. 먼저 특정 상태의 에너지는 다음과 같이 표시된다. 

 그리고 이때 해당 상태의 확률은 다음과 같다. 

 여기서 Z는 단순히 normalizing constant이지만 문제는 이 Z의 계산이 매우 어렵다. 그 이유는 뒤에서 나오는데 쉽게 풀어 말하자면 가능한 모든 경우에 대해서 계산을 해봐야하기 때문이다. 이렇게 구해진 확률의 log-likelihood를 계산해서 그 값이 최대한 크게 하는 gradient ascent방식을 사용하는 것이 일반적인 RBM을 학습시키는 방법이다. 이때 log-likelihood는 다음과 같다. 

 이 식에서 알 수 있듯이 우리가 알고자하는 link의 weight는 주어진 data로 구해진 expectation value에서 model 전체에 대한 expectation value를 빤 값이다. 보통 전자를 positive phase라고 부르고, 후자를 negative phase라고 부른다. 

 이러한 RBM을 train 시키는 다른 방법으로는 Contrastive Divergence(CD)라는 방법이 있다. 이 CD는 간략화 된 학습이라고 보면 되는데, positive phase를 주어진 data에 대한 expectation으로 보고, negative phase를 n번 gibbs sampling을 통해서 얻어진 distribution의 expectation으로 보는 것이다. 말은 어렵지만 gibbs sampling이 무엇인지 안다면 대충 이해가 갈 것이다. 

 

 *gibbs sampling에 대한 포스팅: http://enginius.tistory.com/296

 여튼 내가 구한 소스에는 1-CD를 사용했다. 무슨 말인고 하니 gibbs sampling을 한번만 한 것을 distribution으로 생각하겠다는 것이다.

 만약 n번 gibbs sampling을 할 경우 n-CD가 된다. 간단하다. 

 실험은 매트랩으로 하였고, dataset은 MNSIT digit을 사용하였다. MNIST는 10000개의 train image와 6000개의 test image를 갖는 dataset이다. dataset은 다음 사이트에서 구할 수 있다. 

 http://yann.lecun.com/exdb/mnist/

 여튼 전체적인 실험 과정은 이렇다. 이 dataset을 이용해서 RBM을 train시키고, 이렇게 train된 RBM에 노이즈가 섞인 data를 넣어주었을 때 어떤 결과를 보이는지를 확인하는 것이다. RBM의 경우 unsupervised learning이라 이 자체로 classification 문제를 풀긴 어렵다. 

<매트랩 코드>

 %% 초기화 하자

clc;

rand('state',100);          % make random distribution

%% 데이터를 가져오고, 한 epoch에 처리할 단위인 batch를 만들자//

if exist('batchdata', 'var') == 0

    clear all;

    

    fprintf(1,'1. Converting Raw files into Matlab format \n');

    converter;

    % make batchs from image data

    fprintf(1,'2. Make Batch Data \n');

    makebatches;

    [numcases numdims numbatches] = size(batchdata);

end

%% 실제 Restricted Boltzmann Machine을 training하자

numhid      = 100;

epoch       = 1;

maxepoch    = 10000;

%% 먼저 RBM에 필요한 변수들을 초기화 하자

epsilonW    = 0.05;     % Learning rate for weights

epsilonVb   = 0.05;     % Learning rate for biases of visible units

epsilonHb   = 0.05;     % Learning rate for biases of hidden units

CD          = 1;

weightCost  = 0.001;

initialMomentum  = 0.5;

finalMomentum    = 0.9;

[numCases numDims numBatches] = size(batchdata);

% Initializing symmetric weights and biases.

visHid      = 0.001*randn(numDims, numhid);

hidBiases   = zeros(1,numhid);

visBiases   = zeros(1,numDims);

posHidProbs = zeros(numCases,numhid);

negHidProbs = zeros(numCases,numhid);

posProds    = zeros(numDims,numhid);

negProds    = zeros(numDims,numhid);

visHidInc   = zeros(numDims,numhid);

hidBiasInc  = zeros(1,numhid);

visBiasInc  = zeros(1,numDims);

batchPosHidProbs = zeros(numCases, numhid, numBatches);

%% 실제로 RBM을 train하자

batchIndex = floor((numBatches-1)*rand())+1;

ii = 1;

figure(1);

dispims((batchdata(:, :, batch))',28,28); title(sprintf('[%dth epoch] %dth batch raw data ', epoch, batch)); 

            

for epoch = epoch:maxepoch

    errsum = 0;

    

    % 1. 전체 batch에 대해서

    %for batch = 1 : numBatches,            

    

    % 2. 토이 테스트를 위해서 한 batch에 대해서

    for batch = batchIndex : batchIndex     

        

        

        % 모든 batch에 대해서, (한 batch는 100개의 image로 이뤄져있다.)

        %fprintf(1,'epoch %  d batch %d \r\n', epoch, batch);

        

        % 하나의 batch를 한번에 처리하기 위해서 bias를 repmat

        visBias = repmat(visBiases, numCases, 1);

        hidBias = repmat(hidBiases, numCases, 1);

        %%%%%%%%% START POSITIVE PHASE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        % 1. batch 중에서 우리가 사용할 data를 읽어온다. 100개의 digit으로 이뤄져있다. 

        data = batchdata(:, :, batch);              

        % get this index's batch data(100개의 digit을 포함)

        % 이 data는 numCases, 즉 100개의 digit 정보를 가지고 있다. 이 값들은 0~1사이의 값이고,

        % RBM은 binary만 처리할 수 있어서 이러한 루틴을 넣는 것 같다.

        data = data > rand(numCases, numDims);      

        % 이 값을 0또는 1로 바꾼다.

        

        % 2. positive state에서 각 batch의 특정 hidden node가 1이 될 확률을 의미한다. 

        posHidProbs = 1./ (1 + exp(-data*(2*visHid) - hidBias));

        % posHidProbs는 R[numCases*numhid]의 공간을 갖고 있다.

        % i번째 행은 현재 batch의 i번째 item을 의미하고,

        % j번째 열을 해당 item의 j번째 hidden node가 1이 될 확률을 의미한다.

        % 3. 모든 batch에 대해서 positive state에서 hidden node가 1이 될 확률을 저장한다. 

        batchPosHidProbs(:, :, batch) = posHidProbs;

        % batchPosHidProbs는 R[nunCases*numhid*numBatches]의 크기를 갖는다.

        % 현재 이 loop는 각 batch마다 한번씩 들어오는데 이를 다 저장을 하고, batchPosHidProbs는 여기선 안쓰인다.

        % 4. data의 평균 

        posProds    = data' * posHidProbs;

        % posProds는 R[numDims, numhid]의 크기를 갖는다.

        % 현재 batch내의 모든 data의 각 dimension에 대한 expectation?

        % 뒤에서 numCases로 나눠져서 data의 평균을 구한다. 

        % 5. positve phase에서 hidden node의 기대값을 나타낸다. 

        poshidact   = sum(posHidProbs);

        % poshidact는 R[numhid]의 크기를 갖는다.

        % 각 item은 뒤에서 numCases로 나눠져서 현재 batch의 각 hidden node의 기대값을 나타낸다.

        % 6. positve phase에서 visible node의 기대값을 나타낸다. 

        posvisact   = sum(data);

        % poshidact는 R[numDims]의 크기를 갖는다.

        % 각 item은 뒤에서 numCases로 나눠져서 현재 batch의 각 visible node의 기대값을 나타낸다.

        %%%%%%%%% END OF POSITIVE PHASE  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        %%%%% START NEGATIVE PHASE  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        poshidstates = posHidProbs > rand(numCases,numhid);

        % posHidProbs 앞에서 구한 현재 batch의 각 case들의 각 hidden node가 1이 될 확률이다.

        % poshidstates는 RBM의 binary한 연산을 위해서 이를 1또는 0으로 바꾼다.

        negdata = 1./(1 + exp(-poshidstates*visHid' - visBias));

        % negdata는 앞서 positive phase에서 계산된 hid node들 (1 if 확률>rand)의 값에 지금 갖고 있는 

        % visHid, 즉 vis와 hid를 연결하는 weight를 곱해서 결국 neg phase에서 vis가 1이될 확률을 의미한다. 

        

        if rem(epoch, 100) == 1 && batch == batchIndex

            ii  = ii + 1;

            % raw data (batch)

            %figure(1);

            %dispims((batchdata(:, :, batch))',28,28); title(sprintf('[%dth epoch] %dth batch raw data ', epoch, batch)); drawnow; 

            

            % negative phase에서 계산된 image -

            figure(ii); 

            dispims(negdata',28,28); title(sprintf('[%dth epoch] %dth batch negative phase data', epoch, batch)); drawnow; 

            %pause();

        end

        

        negdata = negdata > rand(numCases,numDims);

        % RBM은 binary data밖에 쓸 수 없다. 즉 neg phase에서 vis data를 의미한다. 

        

        negHidProbs = 1./(1 + exp(-negdata*(2*visHid) - hidBias));

        % negHidProbs 는 현재 batch내에 각 data에 대해서 각 hidden node들의 확률을 나타낸다.

        negProds  = negdata'*negHidProbs;

        % negProds는 R[784*50]의 크기를 같는다.

        % 즉 현재 batch에 대해서 vis과 hid를 연결하는 weight를 의미? 

        

        neghidact = sum(negHidProbs);

        % neghidact R[numhid]의 크기를 갖는다.

        % 각 item은 뒤에서 numCases로 나눠져서 현재 batch의 각 hidden node의 기대값을 나타낸다.

        

        negvisact = sum(negdata);

        % negvisact R[numDims]의 크기를 갖는다.

        % 각 item은 뒤에서 numCases로 나눠져서 현재 batch의 각 visible node의 기대값을 나타낸다.

        %%%%%%%%% END OF NEGATIVE PHASE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        

        err= sum(sum( (data - negdata).^2 ));

        errsum = err + errsum; % 이 errsum은 batch마다 초기화 된다. 즉 현재 batch의 총 에러를 의미 

        if epoch > 5  

            momentum = finalMomentum;

        else

            momentum = initialMomentum;

        end;

        %%%%%%%%% UPDATE WEIGHTS AND BIASES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

        % inc를 계산한다. E[pos]-E[neg] 

        visHidInc = momentum*visHidInc + ...

            epsilonW*( (posProds-negProds)/numCases - weightCost*visHid ) ;

        visBiasInc = momentum*visBiasInc + (epsilonVb/numCases)*(posvisact - negvisact);

        hidBiasInc = momentum*hidBiasInc + (epsilonHb/numCases)*(poshidact - neghidact);

        

        % 실제 update 부분 

        visHid = visHid + visHidInc;

        visBiases = visBiases + visBiasInc;

        hidBiases = hidBiases + hidBiasInc;

        %%%%%%%%%%%%%%%% END OF UPDATES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    end

    fprintf(1, 'epoch %4i error %6.1f  \r\n', epoch, errsum);

end;

%% 구해진 RBM을 test해보자 

% 1. 먼저 테스트를 위해 random하게 data를 가져온다. 

%batchIndex = floor((numBatches-1)*rand())+1;

fprintf(' %d th batch\r\n', batchIndex);

data = batchdata(:, :, batchIndex); 

% 2. 이 data에 noise를 섞는다.

data_with_noise = data + 0.1*randn(numCases, numDims);

%. 3. RBM에서 hidden node로 갔다가 다시 visible로 오자. 

data2 = data_with_noise > rand(numCases, numDims);      % 이 값을 0또는 1로 바꾼다.

posHidProbs = 1./ (1 + exp(-data2*(2*visHid) - hidBias)); 

poshidstates = posHidProbs > rand(numCases,numhid);

negdata = 1./(1 + exp(-poshidstates*visHid' - visBias));

negdata = negdata > rand(numCases, numDims);

% 4. 결과를 plot한다. 

figure(1); dispims( (batchdata(:, :, batchIndex))', 28, 28); title(sprintf('1. raw data')); 

figure(2); dispims( data2', 28, 28); title(sprintf('2. data with noise')); 

figure(4); dispims( negdata', 28, 28); title('3. RBM result'); 

실험 결과

 - 이런식으로 epoch가 커질수록 원래 이미지에 더 가까이 가는 것을 알 수 있다. 이제 테스트를 위해 이 image에 gaussian noise를 입혔을 때 어떤 결과가 나오는지 확인해보자. 

1. 원래 이미지

2. 노이즈를 추가한 이미지 

3. RBM을 거친 후의 이미지 

 RBM을 거친 후의 이미지에서 노이즈가 줄어든 것을 알 수 있다. 

 * 참고로 DBN은 이 RBM을 여러 층 cascade한 것을 의미한다. 

 매트랩 프로젝트

8. Restricted Boltzmann Machine on MNIST.z01

8. Restricted Boltzmann Machine on MNIST.zip