mnf

 function [phi,psi] = mnf(X)
 Extracts the maximum noise fraction components of the data in X.
 and projects points onto these components

 In:
    X - an NxP data matrix  (P is number of samples)

 Out:
    phi - the basis vectors (NxP)
    psi - the mixing coefficients (NxN)


 Error: This code will return null values for the two output parameters
        if there is an error during execution.

 This follows Hundley, Kirby and Anderle (not anymore)

0001 function [phi,psi] = mnf(X)
0002 % function [phi,psi] = mnf(X)
0003 % Extracts the maximum noise fraction components of the data in X.
0004 % and projects points onto these components
0005 %
0006 % In:
0007 %    X - an NxP data matrix  (P is number of samples)
0008 %
0009 % Out:
0010 %    phi - the basis vectors (NxP)
0011 %    psi - the mixing coefficients (NxN)
0012 %
0013 %
0014 % Error: This code will return null values for the two output parameters
0015 %        if there is an error during execution.
0016 %
0017 % This follows Hundley, Kirby and Anderle (not anymore)
0018 %
0019 
0020 X = X';
0021 
0022 [numSamples dim] = size(X);
0023 
0024 %mean center the columns
0025 %X = X - repmat(mean(X),numSamples,1);
0026 
0027 try
0028     [Vx,Lx] = sorteig(X'*X);
0029 catch
0030     phi = [];
0031     psi = [];
0032     fprintf('An error occured while computing the eigen decomposition of the covariance matrix.\n')
0033     lasterr
0034     return;
0035 end
0036 
0037 DSx = sqrt(diag(Lx));
0038 INVDSx = diag(1./(DSx));
0039 
0040 %estimate covariance matrix of Z = Q'*Q
0041 SX(1,:) = X(2,:);
0042 SX(2:size(X,1),:) = X(1:size(X,1)-1,:);
0043 Z1 = SX'*X;
0044 s1 = norm(X-SX,'fro')/norm(X,'fro');
0045 Z2 = X'*SX;
0046 s2=norm(Z2);
0047 Z3 = X'*X;
0048 Z4 = SX'*SX;
0049 Z = 0.5*(Z3+Z4-Z1-Z2);
0050 ZHAT = INVDSx*Vx'*Z*Vx*INVDSx;
0051 
0052 try
0053     [PSIHAT,DPSIHAT] = sorteig(ZHAT);
0054 catch
0055     phi = [];
0056     psi = [];
0057     fprintf('An error occured while computing the eigen decomposition of the ZHAT.\n')
0058     fprintf('%s\n\n',lasterr)
0059     return;
0060 end
0061 
0062 psi = Vx*INVDSx*PSIHAT;
0063 phi = X*psi;
0064 
0065 phi = phi';
0066 psi = psi';
0067 
0068 % this will make all signals have positive activation
0069 %phi(find(sum(psi)<0),:) = -phi(find(sum(psi)<0),:);
0070