PCA的MATLAB实现

% ========================================================================= % power_ratio is for power-saving % main_feature is principal value % main_vectors is responding principal vector function [main_feature,main_vector]=PCA_Nicolas(data,power_ratio) % ========================================================================= [M,N]=size(data); % gain size of data mean_value = mean(data); % gain mean value % Get max value for ploting graphic Width = max(max(data)); x=-Width:0.1:Width; yy=-Width:0.1:Width; % plot for 2d data % high dimension will dismiss % -------------- Preparing for PCA ------------------------- % -------------- Preparing for PCA ------------------------- data_sub = zeros(M,N); %-------------------------------------------------------------- %-------------------------------------------------------------- switch N case 2 subplot(2,2,1) plot(data(:,1),data(:,2),'r*') hold on plot(x,0,'b') plot(0,yy,'b') title('Original Dataset Distributions') hold off % Original Dataset Distributions Plot %--------------------------------------------------------------- for i=1:M for j=1:N data_sub(i,j)=data(i,j)-mean_value(1,j); end end clear data % Gain the dataset which subtract the mean value and % perform a distribution with mean value zero subplot(2,2,2) plot(data_sub(:,1),data_sub(:,2),'r*') hold on plot(x,0,'b'); plot(0,yy,'b') title('Subtracted Dataset Distributions') hold off save data_sub % Subtracted Dataset Distributions Plot %------------------------------------------------------------------ % calculate covariance matrix,then eigenvectors,eigenvalues % and sort the eigenvalue and its eigenvector %%% Key section begins Cov_matrix = data_sub*data_sub'; % covariance matrix clear data_sub % release the memory [vector,feature]=eig(Cov_matrix);% covariance matrix decomposition clear Cov_matrix % release memory feature=diag(feature); % get eigenvalues on the diagonal feature_sum=sum(feature); % sum the energy of the subtracted dataset [junk,rindices]=sort(-1*feature);%《A Tutorial on Principal Component % Analysis》2005Y second edition feature=feature(rindices); vector=vector(:,rindices); %select the principal components main_feature_sum=0; MM=0; for i=1:M main_feature_sum = main_feature_sum+feature(i,1); if main_feature_sum > feature_sum*power_ratio break end end % read in the main feature MM=i; % 记录主要特征值的个数 clear junk % 释放内存 clear rindices % 释放内存 main_feature=zeros(MM,1);% create a array to save main feature for i=1:MM main_feature(i,1)=feature(i,1); end main_feature=diag([main_feature]); % diagonalize the main feature % and gain the main feature matrix % read in the responding vectors main_vector=zeros(M,MM); for i=1:MM main_vector(:,i)=vector(:,i); end clear vector % Gain the responding main vectors %%% Key section over % ----------------------------------------------------------------- % derive the origibal dataset with mean value 0 derived_data=main_vector*sqrt(main_feature); subplot(2,2,3) plot(derived_data(:,1),derived_data(:,2),'r*') hold on plot(x,0,'b'); plot(0,yy,'b') title('Derived Dataset Distributions') hold off save derived_data %----------------------------------------------------------------- data_derived=zeros(M,N); for i=1:M for j=1:N data_derived(i,j)=derived_data(i,j)+mean_value(1,j); end end clear derived_data save data_derived % derive the original dataset %------------------------------------------------------------------ subplot(2,2,4) plot(data_derived(:,1),data_derived(:,2),'r*') hold on plot(x,0,'b'); plot(0,yy,'b') title('Derived Dataset Distributions') hold off main_feature main_vector %%% Ups is the 2-dimension case %%% ----------------------------------------------------------------------- %%% ----------------------------------------------------------------------- case 3 subplot(2,2,1) plot3(data(:,1),data(:,2),data(:,3),'r*') title('Original Dataset Distributions') grid on % Original Dataset Distributions Plot %------------------------------------------------------------------ for i=1:M for j=1:N data_sub(i,j)=data(i,j)-mean_value(1,j); end end clear data % Gain the dataset which subtract the mean value and % perform a distribution with mean value zero subplot(2,2,2) plot3(data_sub(:,1),data_sub(:,2),data_sub(:,3),'r*') title('Subtracted Dataset Distributions') grid on save data_sub % Subtracted Dataset Distributions Plot %------------------------------------------------------------------ % calculate covariance matrix,then eigenvectors,eigenvalues % and sort the eigenvalue and its eigenvector %%% Key section begins Cov_matrix = data_sub*data_sub'; % covariance matrix clear data_sub % release the memory [vector,feature]=eig(Cov_matrix);% covariance matrix decomposition clear Cov_matrix % release memory feature=diag(feature); % get eigenvalues on the diagonal feature_sum=sum(feature); % sum the energy of the subtracted dataset [junk,rindices]=sort(-1*feature);%《A Tutorial on Principal Component % Analysis》2005Y second edition feature=feature(rindices); vector=vector(:,rindices); %select the principal components main_feature_sum=0; MM=0; for i=1:M main_feature_sum = main_feature_sum+feature(i,1); if main_feature_sum > feature_sum*power_ratio break; end end % read in the main feature MM=i; % 记录主要特征值的个数 clear junk % 释放内存 clear rindices % 释放内存 main_feature=zeros(MM,1);% create a array to save main feature for i=1:MM main_feature(i,1)=feature(i,1); end main_feature=diag([main_feature]); % diagonalize the main feature %%%% Gain the main feature matrix % read in the responding vectors main_vector=zeros(M,MM); for i=1:MM main_vector(:,i)=vector(