pca.m

来自「这是盲信号的代码 都已经通过编译了 做这方面的同仁可以参考一下 我觉得蛮惯用的」· M 代码 · 共 42 行

function y=pca(x)
%主分量分析（观测信号大于源信号）
%创建时间：2005.7.8
%作者： lucky zhang
%y=pac(x);
[Dim, NumOfSampl] = size(x);
fprintf('Number of signals: %d\n', Dim);
fprintf('Number of samples: %d\n', NumOfSampl);
%remean
newVectors = zeros (size (x));
meanValue = mean (x')';
newVectors = x - meanValue * ones (1,size (x, 2));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Calculate the covariance Matrix
covarianceMatrix = cov(x', 1);
% Calculate the eigenvalues and eigenvectors of covariance matrix.
[E, D] = eig (covarianceMatrix);
%reduce dim
rankTolerance = 1e-7;
maxLastEig = sum (diag (D) > rankTolerance);
eigenvalues = flipud(sort(diag(D)));
firstEig = 1;
lastEig = maxLastEig;
fprintf('Dimension reduced to %d due to the singularity of covariance matrix\n',lastEig-firstEig+1);
oldDimension = size (x, 1);
if lastEig < oldDimension
  lowerLimitValue = (eigenvalues(lastEig) + eigenvalues(lastEig + 1)) / 2;
else
  lowerLimitValue = eigenvalues(oldDimension) - 1;
end

lowerColumns = diag(D) > lowerLimitValue;
selectedColumns = lowerColumns;
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));
E = selcol(E, selectedColumns);
D = selcol(selcol(D, selectedColumns)', selectedColumns);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
whiteningMatrix = inv (sqrt (D)) * E';
dewhiteningMatrix = E * sqrt (D);
y =  whiteningMatrix * x;