主成分分析的原理及Matlab源程序

1、主成分分析的原理及 Matlab源程序2010-07-03 20:59:02 作者：huright 来源： 浏览次数：82网友评论0 条主成分分析也称主分量分析，旨在利用降维的思想，把多指标转化为少数几个综合指标。在实证问题研究中，为了全面、系统地分析问题，我们必须考虑众多影响因素。这些涉及 的因素一般称为指标，在多元统计分析中也称为变量。因为每个变量都在不同程度上反映 了所研究问题的某些信息，并且指标之间彼此有一定的相关性，因而所得的统计数据反映 的信息在一定程度上有重叠。在用统计方法研究多变量问题时，变量太多会增加计算量和增加分析问题的复杂性，人们希望在进行定量分析的过程中，涉及的变量较少

2、，得到的信 息量较多。主成分分析正是适应这一要求产生的。主成分分析法是一种数学变换的方法，它把给定的一组相关变量通过线性变换转成另一组不相关的变量，这些新的变量按照方差依次递减的顺序排列。在数学变换中保持变量的 总方差不变，使第一变量具有最大的方差，称为第一主成分，第二变量的方差次大，并且 和第一变量不相关，称为第二主成分。依次类推，I个变量就有I个主成分。其中Li为p维正交化向量（Li衣Li = 1） , Zi之间互不相关且按照方差由大到小排列， 则称Zi为X的第I个主成分。设 X的协方差矩阵为 2,则工必为半正定对称矩阵，求特征 值入i（按从大到小排序）及其特征向量，可以证明，入所对应的正

3、交化特征向量，即为第I个主成分Zi所数据标准化；分析步骤求相关系数矩阵；一系列正交变换，使非对角线上的数置0,加到主对角上；得特征根xi （即相应那个主成分引起变异的方差）,并按照从大到小的顺序把特征根排列求各个特征根对应的特征向量；用下式计算每个特征根的贡献率Vi;Vi=xi/（x1+x2+）贡献率定义为入入通常要求提取的主成分的数量k满足工入k工入j>0.85% function y = pca(mixedsig)%程序说明： y = pca(mixedsig) ，程序中 mixedsig 为 n*T 阶混合数据矩阵， n 为信号个数，T为采样点数% y 为 m*T 阶主分量矩阵。%

4、 n 是维数， T 是样本数。if nargin = 0error('You must supply the mixed data as input argument.');endif length(size(mixedsig)>2error('Input data can not have more than two dimensions. ');endif any(any(isnan(mixedsig)error('Input data contains NaN''s.');end%去均值meanValue = mean

5、(mixedsig')'当数据本身维数很大时m,n = size(mixedsig);%mixedsig = mixedsig - meanValue*ones(1,size(meanValue); %容易出现 Out of memoryfor s = 1:mfor t = 1:nmixedsig(s,t) = mixedsig(s,t) - meanValue(s);endendDim,NumofSampl = size(mixedsig);oldDimension = Dim;fprintf('Number of signals: %d ',Dim);fpr

6、intf('Number of samples: %d ',NumofSampl);fprintf('Calculate PCA.');firstEig = 1;lastEig = Dim;计算协方差矩阵计算协方差矩阵的特征值和特征向量covarianceMatrix = corrcoef(mixedsig'); %E,D = eig(covarianceMatrix); %计算协方差矩阵的特征值大于阈值的个数lastEig%rankTolerance = 1;%maxLastEig = sum(diag(D) >= rankTolerance);

7、%lastEig = maxLastEig;lastEig = 10;%降序排列特征值 eigenvalues = flipud(sort(diag(D);%选择相应的特征值和特征向量if lastEig < oldDimensionlowerLimitValue = (eigenvalues(lastEig) + eigenvalues(lastEig + 1)/2;elselowerLimitValue = eigenvalues(oldDimension) - 1;endlowerColumns = diag(D) > lowerLimitValue;%去掉较大的特征值 (一

8、般没有这一步 ) if firstEig > 1higherLimitValue = (eigenvalues(firstEig - 1) + eigenvalues(firstEig)/2;elsehigherLimitValue = eigenvalues(1) + 1;endhigherColumns = diag(D) < higherLimitValue;%合并选择的特征值 selectedColumns =lowerColumns & higherColumns;%输出处理的结果信息 fprintf('Selected %d dimensions. &#

9、39;,sum(selectedColumns);fprintf('Smallest remaining (non-zero) eigenvalue %g ',eigenvalues(lastEig); fprintf('Largest remaining (non-zero) eigenvalue %g ',eigenvalues(firstEig); fprintf('Sum of removed eigenvalue %g ',sum(diag(D) .* (selectedColumns);E = selcol(E,selectedCol

10、umns);D = selcol(selcol(D,selectedColumns)',selectedColumns);%计算白化矩阵 whiteningMatrix = inv(sqrt(D) * E' dewhiteningMatrix = E * sqrt(D);%提取主分量 y = whiteningMatrix * mixedsig;%行选择子程序function newMatrix = selcol(oldMatrix,maskVector)if size(maskVector,1)= size(oldMatrix,2)error('The mask vector and matrix are of uncompatible size.'); end