function y = pca(mixedsig,q)
%程序说明：y = pca(mixedsig)，程序中mixedsig为 n*T 阶混合数据矩阵，n为信号个数，T为采样点数
% y为 m*T 阶主分量矩阵。
% n是维数，T是样本数。
if nargin == 0
error('You must supply the mixed data as input argument.');
end
if length(size(mixedsig))>2
error('Input data can not have more than two dimensions. ');
end
if any(any(isnan(mixedsig)))
error('Input data contains NaN''s.');
end
%——————————————去均值————————————
meanValue = mean(mixedsig')';
%%%%128*k 转置后为k*128 变成1*128 再转置 128*1 每行是k个点的均值
[m,n] = size(mixedsig);
%mixedsig = mixedsig - meanValue*ones(1,size(meanValue)); %当数据本身维数很大时容易出现Out of memory
for s = 1:m
for t = 1:n
mixedsig(s,t) = mixedsig(s,t) - meanValue(s);
end
end
%[m,n]=size(descr1)
%for j=1:n
% for i=1:m
% if descr1(i,j)==NaN
% descr1(:,j)=0;
% end
% end
%end
%J=isnan(descr1)
% num=0;
%for h=1:m
%for k=1:n
% num=num+J(h,k);
%end
%end
[Dim,NumofSampl] = size(mixedsig);
oldDimension = Dim;
fprintf('Number of signals: %d\n',Dim);
fprintf('Number of samples: %d\n',NumofSampl);
fprintf('----------------------------------------------------------------------------\n')
fprintf('Calculate PCA...\n');
firstEig = 1;
lastEig = Dim;
covarianceMatrix = corrcoef(mixedsig'); %计算协方差矩阵
[E,D] = eig(covarianceMatrix); %计算协方差矩阵的特征值和特征向量
%———计算协方差矩阵的特征值大于阈值的个数lastEig———
%rankTolerance = 1;
%maxLastEig = sum(diag(D) >= rankTolerance);
%lastEig = maxLastEig;
lastEig = q;
%——————————降序排列特征值——————————
eigenvalues = flipud(sort(diag(D)));
%—————————去掉较小的特征值——————————
if lastEig < oldDimension
lowerLimitValue = (eigenvalues(lastEig) + eigenvalues(lastEig + 1))/2;
else
lowerLimitValue = eigenvalues(oldDimension) - 1;
end
lowerColumns = diag(D) > lowerLimitValue;
%—————去掉较大的特征值(一般没有这一步)——————
%if firstEig > 1
% higherLimitValue = (eigenvalues(firstEig - 1) + eigenvalues(firstEig))/2;
%else
% higherLimitValue = eigenvalues(1) + 1;
%end
%higherColumns = diag(D) < higherLimitValue;
%—————————合并选择的特征值——————————
selectedColumns =lowerColumns ;%& higherColumns;
%—————————输出处理的结果信息—————————
fprintf('Selected [%d] dimensions.\n',sum(selectedColumns));
%fprintf('Smallest remaining (non-zero) eigenvalue[ %g ]\n',eigenvalues(lastEig));
%fprintf('Largest remaining (non-zero) eigenvalue[ %g ]\n',eigenvalues(firstEig));
%fprintf('Sum of removed eigenvalue[ %g ]\n',sum(diag(D) .* (~selectedColumns))/sum(diag(D)));
fprintf('Percents of selected eigenvalue[ %g ]\n',sum(diag(D) .* (selectedColumns))/sum(diag(D)));
%fprintf('Sum of eigenvalue[ %g ]\n',sum(diag(D) ));
%———————选择相应的特征值和特征向量———————
E = selcol(E,selectedColumns);%特征向量
D = selcol(selcol(D,selectedColumns)',selectedColumns);
%fprintf('D=%d rows and %d lines\n',b,v);
%——————————计算白化矩阵———————————
%whiteningMatrix = inv(sqrt(D)) * E';
%dewhiteningMatrix = E * sqrt(D);
%为什么要计算白化矩阵呢？
%——————————提取主分量————————————
y = E' * mixedsig;
[b,v]=size(y);
%y = whiteningMatrix * mixedsig;
%——————————行选择子程序———————————
%fprintf('%d rows and %d lines\n',b,v);
评论10