princomp.m

%princomp源程序
function P = princomp(X, q) 
%PRINCOMP Obtain principal-component vectors and related quantities. 
%   P = PRINCOMP(X, Q) Computes the principal-component vectors of 
%   the vector population contained in the rows of X, a matrix of 
%   size K-by-n where K is the number of vectors and n is their 
%   dimensionality. Q, with values in the range [0, n], is the number 
%   of eigenvectors used in constructing the principal-components 
%   transformation matrix. P is a structure with the following 
%   fields: 
% 
%     P.Y      K-by-Q matrix whose columns are the principal- 
%              component vectors. 
%     P.A      Q-by-n principal components transformation matrix 
%              whose rows are the Q eigenvectors of Cx corresponding 
%              to the Q largest eigenvalues.  
%     P.X      K-by-n matrix whose rows are the vectors reconstructed 
%              from the principal-component vectors. P.X and P.Y are 
%              identical if Q = n.  
%     P.ems    The mean square error incurred in using only the Q 
%              eigenvectors corresponding to the largest 
%              eigenvalues. P.ems is 0 if Q = n.   
%     P.Cx     The n-by-n covariance matrix of the population in X. 
%     P.mx     The n-by-1 mean vector of the population in X. 
%     P.Cy     The Q-by-Q covariance matrix of the population in 
%              Y. The main diagonal contains the eigenvalues (in 
%              descending order) corresponding to the Q eigenvectors. 
 
%   Copyright 2002-2004 R. C. Gonzalez, R. E. Woods, & S. L. Eddins 
%   Digital Image Processing Using MATLAB, Prentice-Hall, 2004 
%   $Revision: 1.4 $  $Date: 2003/10/26 23:16:16 $ 
 
[K, n] = size(X); 
X = double(X); 
 
% Obtain the mean vector and covariance matrix of the vectors in X. 
[P.Cx, P.mx] = covmatrix(X); 
P.mx = P.mx'; % Convert mean vector to a row vector. 
 
% Obtain the eigenvectors and corresponding eigenvalues of Cx.  The 
% eigenvectors are the columns of n-by-n matrix V.  D is an n-by-n 
% diagonal matrix whose elements along the main diagonal are the 
% eigenvalues corresponding to the eigenvectors in V, so that X*V = 
% D*V. 
[V, D] = eig(P.Cx); 
 
% Sort the eigenvalues in decreasing order.  Rearrange the 
% eigenvectors to match.  
d = diag(D); 
[d, idx] = sort(d); 
d = flipud(d); 
idx = flipud(idx); 
D = diag(d); 
V = V(:, idx); 
 
% Now form the q rows of A from first q columns of V. 
P.A = V(:, 1:q)'; 
 
% Compute the principal component vectors. 
Mx = repmat(P.mx, K, 1); % M-by-n matrix.  Each row = P.mx. 
P.Y = P.A*(X - Mx)'; % q-by-K matrix. 
 
% Obtain the reconstructed vectors. 
P.X = (P.A'*P.Y)' + Mx; 
 
% Convert P.Y to K-by-q array and P.mx to n-by-1 vector. 
P.Y = P.Y'; 
P.mx = P.mx'; 
 
% The mean square error is given by the sum of all the  
% eigenvalues minus the sum of the q largest eigenvalues. 
d = diag(D); 
P.ems = sum(d(q + 1:end)); 
 
% Covariance matrix of the Y's: 
P.Cy = P.A*P.Cx*P.A';