pca2.m

该源码可用于分析气候多因子的主分量分析

%%%---------------------------------
%% R-PCA主成分分析
%------------------------------------
%X=[243.87 16521 6.46 34.57 149.85;
 %   240.31 8210 9.89 16.92 55.89;
 %   211.15 15349 10.09 29.77 80.13;
 %   413.18 16760 7.67 24.14 105.36;
%    349.60 7721 6.47 16.27 99.41;
 %   205.47 8123 12.33 18.47 46.18;
  %  298.11 13308 5.05 27.35 138.76;
   % 414.94 13781 4.10 16.65 98.20;
  %  287.25 14043 4.29 17.67 58.35;
 %   303.93 11126 7.63 18.39 74.23;
 %   608.40 22392 2.94 24.56 233.37;
 %   433.92 12508 0.69 20.06 118.70;
 %   572.63 12102 2.76 12.08 110.43;
  %  533.78 11990 3.80 11.59 75.55;
  %  545.70 9678 3.55 9.46 61.19;
  %  284.61 6513 6.41 12.83 48.15;
%    572.07 18664 2.31 17.76 162.11;
  %  409.86 7329 5.89 12.23 76.68;
  %  564.02 14911 4.93 28.50 233.58;
  %  221.20 6443 14.08 30.25 80.48]
function [Z,R,LL]=PCA(X)
%% X 为样本观测值
% Step1 样本标准化 Y
[n,p]=size(X);
Y=zeros(n,p);
for j=1:p
    XX(j)=mean(X(:,j));
    S(j)=std(X(:,j));
    Y(:,j)=(X(:,j)-XX(j))./S(j);
end
%%%-----------Step2  计算相关系数据矩阵R
R=zeros(p,p);
save yy.txt Y -ASCII -TABS;
for k=1:p
    for j=1:k
      for i=1:n
        R(k,j)=R(k,j)+(Y(i,j)*Y(i,k))/(n-1);
     end
R(k,k)=1;
R(j,k)=R(k,j);
  end
end 
%%% Step3 求R的特征值与特征向量
[L,D]=eig(R);
%% V为特征向量,D为特征值
%% 特征值排序以及选择相应的正则化特征向量
%% 选择排序
for i=1:p
    temp=i;
    for j=i+1:p
        if D(temp,temp)<D(j,j)
            temp=j;
            x=D(temp,temp);
            D(temp,temp)=D(i,i);
            D(i,i)=x;
            
            xx=L(:,temp);
            L(:,temp)=L(:,i);
            L(:,i)=xx;
        end
    end
end
%D
%L
%% Step4 建立主成分
 S=0.0;
 kk=1;
while (S<=0.85)
S=S+D(kk,kk)/p;
kk=kk+1;
end
kk=kk-1;
%% Step 5 计算前kk个主成分的样本值
Z=zeros(n,kk);
for j=1:kk
   LL(:,j)=L(:,j) ;
   for i=1:p
   Z(:,j)=Z(:,j)+L(i,j).*Y(:,i);
   end
end
%Z
%R
%L 
%clear i,j;
%% Step 5 计算前kk个主成分的样本值
%for i=1:n
 %   for j=1:kk
 %       for t=1:p
%        Z(i,j)=Z(i,j)+Y(i,t)*L(t,j);
%        end
 %   end
%end

%R
%L