% OjaPCA.m
%          5 May 2005

clear
xS = [4  3  5  2  1  6  3  3  2
      1  7  5  4  6  3  2  4  7] ;
nn = max(size(xS)) ;
mnX = mean(xS,2) ;
xx = xS - mnX(:, ones(1,nn)) ;

% Oja's rule -- multiply epochs
kk = 5 ; % number of epochs
alf = 3e-2 ;
w = zeros(kk*nn+1,2) ;
w(1,:) = 0.4*rand(1,2) ;
k = 1 ;
for ep = 1:kk
  for n = 1:nn
    x = xx(:,n) ;
    y = w(k,:)*x ;
    w(k+1,:) = w(k,:) + alf*y*(x' - y*w(k,:)) ;
    k = k+1 ;
  end
end
ws = w(end,:)

load p3bW
 myfigA(1,17,8)
subplot(1,2,1)
plot(w(:,1), w(:,2), '-', w(1,1), w(1,2), '>', ws(1), ws(2), '*'), grid on
xlabel('w_1'), ylabel('w_2'),  axis([-0.6 0.6 0 1])
subplot(1,2,2)
plot(0:nn*kk, sqrt(sum(w'.^2))), grid on, axis([0 kk*nn 0 1.2])
xlabel('learning steps, n'), ylabel('||w||')
text(-30,1.25,'Evolution of the weight vector and its magnitude')
 print -f1 -depsc2 p3LOjaD

  myfigA(2, 8,7)
 plot(xx(1,:), xx(2,:), 'd', ...
    4*ws(1)*[-1 1], 4*ws(2)*[-1 1],'-', ...
    ws(1), ws(2),'*', 0, 0, '+k'), 
 grid on, axis([-4 4 -4 4])
 title('The principal direction of a 2-D shape')
 text(ws(1)+0.2*abs(ws(1)), ws(2), 'w_s')
 xlabel('x_1'), ylabel('x_2'), 
 print -f2 -depsc2 p3ShapeD
