ojafirst.m

matlab实现主成分分析PAC神经网络的程序

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Copy Right:  the Computational Intelligence Laboratory
%% System name: PCA
%% File name:   OjaFirstEx.m
%% Description: This is a algorithm implement about Oja's PCA.
%%
%% Author: 
%% Date:   4/14/2004
%% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear all

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  The original image is digitized into the 64x4096 matrix.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

originalImage = imread('lina.bmp');
preserveOrigianl = originalImage;

[rr,cc] = size(originalImage);

Img = originalImage;
mc=256;

r = 8; c = 8 ; p = r*c ;
digitalImg = blkM2vc(Img, [r c]) ;
[p N] = size(digitalImg) ; %digitalImage is a 64x4096 matrix.

Xm = mean(digitalImg')' ;
Xmax = max(max(abs(digitalImg)));
digitalImg = digitalImg - Xm(:, ones(1, N)) ;
digitalImg = digitalImg/max(max(abs(digitalImg))) ;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  The original image is showed.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

tempXX = digitalImg*Xmax;
tempXX = tempXX+Xm(:, ones(1, N)) ;
Imr = vc2blkM(tempXX, r, rr) ;
Imr = round(mc*Imr/max(max(Imr))) ;
figure(1), set(1,'Name','the original picture'),image(Imr), colormap(gray(mc))



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% the single neuron network is trained until it converge. 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

m = 1 ; % number of neurons
%for m=1:16
time(m)=rem(Now,1);
W = 0.6*(rand(m, p)-0.5); % weight initialisation

er = .001; % the length convergence error
e0=digitalImg;
eta = 0.01
for j=1:m
    

    for s=1:16
    [rs rn] = sort(rand(1, N));
    R=0;
    for n=1:N

        x = e0(:, rn(n)) ; % randomised selection of patterns
        y = W*x ;
        dW = (eta*y)*(x' - W*y);
        W = W + dW;
        
       if sqrt(sum(dW.^2))<er
           break;           
       end
       
       sqrt(sum(dW.^2))
    end
    end
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Extracting the first principal component.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
yy = W*digitalImg ;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% The first principal component is showed.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[rTem,cTem] = size(yy);
tV = rr/8;
for j = 1:rTem
  rVc = 256*yy(j,:);
  Imr = eigvcToMask(rVc,tV);
  Imr = round(mc*Imr/max(max(Imr))) ;
  Imr = round(Imr) ;
  figure(2),set(2,'Name','the Principal Component'), image(Imr), colormap(gray(mc))

end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% The weight vector mask is showed.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[rTem,cTem] = size(W);
tV = 64/8;
for j = 1:rTem
  rVc = 256*W(j,:);
  Imr = eigvcToMask(rVc,tV);
  Imr = round(mc*Imr/max(max(Imr))) ;
  Imr = round(Imr) ;
  figure(3),set(3,'Name','the weight vector mask'), image(Imr), colormap(gray(mc))

end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% The image is reconstructed by the first principal component.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

recoverImg = W'*yy;

recoverImg=recoverImg*Xmax;
recoverImg = recoverImg + Xm(:, ones(1, N)) ;
Imr = vc2blkM(recoverImg, r, rr) ;
Imr = round(mc*Imr/max(max(Imr))) ;
Imr = round(Imr) ;
figure(4),set(4,'Name','the reconstructed image'), image(Imr), colormap(gray(mc))

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% this is end.