wolffd@0: function [P,V,me,l] = pcaproj(D,arg1,arg2)
wolffd@0: 
wolffd@0: %PCAPROJ Projects data vectors using Principal Component Analysis.
wolffd@0: %
wolffd@0: % [P,V,me,l] = pcaproj(D, odim)
wolffd@0: % P =          pcaproj(D, V, me)
wolffd@0: %
wolffd@0: %  Input and output arguments ([]'s are optional)
wolffd@0: %   D      (matrix) size dlen x dim, the data matrix
wolffd@0: %          (struct) data or map struct            
wolffd@0: %   odim   (scalar) how many principal vectors are used
wolffd@0: %  
wolffd@0: %   P      (matrix) size dlen x odim, the projections
wolffd@0: %   V      (matrix) size dim x odim, principal eigenvectors (unit length)
wolffd@0: %   me     (vector) size 1 x dim, center point of D
wolffd@0: %   l      (vector) size 1 x odim, the corresponding eigenvalues, 
wolffd@0: %                   relative to total sum of eigenvalues
wolffd@0: %                   
wolffd@0: % See also SAMMON, CCA.
wolffd@0: 
wolffd@0: % Contributed to SOM Toolbox 2.0, February 2nd, 2000 by Juha Vesanto
wolffd@0: % Copyright (c) by Juha Vesanto
wolffd@0: % http://www.cis.hut.fi/projects/somtoolbox/
wolffd@0: 
wolffd@0: % juuso 191297 070200
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0: 
wolffd@0: error(nargchk(2, 3, nargin)); % check the number of input arguments
wolffd@0: 
wolffd@0: % the data
wolffd@0: if isstruct(D), 
wolffd@0:   if strcmp(D.type,'som_map'), D=D.codebook; else D=D.data; end
wolffd@0: end
wolffd@0: [dlen dim] = size(D);
wolffd@0: 
wolffd@0: if nargin==2, 
wolffd@0: 
wolffd@0:   odim = arg1;
wolffd@0:     
wolffd@0:   % autocorrelation matrix
wolffd@0:   A = zeros(dim);
wolffd@0:   me = zeros(1,dim);
wolffd@0:   for i=1:dim, 
wolffd@0:     me(i) = mean(D(isfinite(D(:,i)),i)); 
wolffd@0:     D(:,i) = D(:,i) - me(i); 
wolffd@0:   end  
wolffd@0:   for i=1:dim, 
wolffd@0:     for j=i:dim, 
wolffd@0:       c = D(:,i).*D(:,j); c = c(isfinite(c));
wolffd@0:       A(i,j) = sum(c)/length(c); A(j,i) = A(i,j); 
wolffd@0:     end
wolffd@0:   end
wolffd@0:   
wolffd@0:   % eigenvectors, sort them according to eigenvalues, and normalize
wolffd@0:   [V,S]   = eig(A);
wolffd@0:   eigval  = diag(S);
wolffd@0:   [y,ind] = sort(abs(eigval)); 
wolffd@0:   eigval  = eigval(flipud(ind));
wolffd@0:   V       = V(:,flipud(ind)); 
wolffd@0:   for i=1:odim, V(:,i) = (V(:,i) / norm(V(:,i))); end
wolffd@0:   
wolffd@0:   % take only odim first eigenvectors
wolffd@0:   V = V(:,1:odim);
wolffd@0:   l = abs(eigval)/sum(abs(eigval));
wolffd@0:   l = l(1:odim); 
wolffd@0: 
wolffd@0: else % nargin==3, 
wolffd@0: 
wolffd@0:   V = arg1;
wolffd@0:   me = arg2;
wolffd@0:   odim = size(V,2);    
wolffd@0:   D = D-me(ones(dlen,1),:);
wolffd@0:   
wolffd@0: end
wolffd@0:   
wolffd@0: % project the data using odim first eigenvectors
wolffd@0: P = D*V;
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%