Home > NoiseTools > nt_mcca.m

nt_mcca

PURPOSE ^

[A,score,AA]=nt_mcca(C,N) - multi-set cca

SYNOPSIS ^

function [A,score,AA]=nt_mcca(C,N)

DESCRIPTION ^

[A,score,AA]=nt_mcca(C,N) - multi-set cca

  A: transform matrix
  score: commonality score (ranges from 1 to N)
  AA: array of subject-specific MCCA transform matrices
 
  C: covariance matrix of aggregated data sets
  N: number of channels of each data set

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [A,score,AA]=nt_mcca(C,N)
0002 %[A,score,AA]=nt_mcca(C,N) - multi-set cca
0003 %
0004 %  A: transform matrix
0005 %  score: commonality score (ranges from 1 to N)
0006 %  AA: array of subject-specific MCCA transform matrices
0007 %
0008 %  C: covariance matrix of aggregated data sets
0009 %  N: number of channels of each data set
0010 
0011 if nargin<2; error('!'); end
0012 if size(C,1) ~= size(C,2); error('!'); end
0013 if size(C,1) ~= round(size(C,1)/N)*N; error('!'); end
0014 
0015 % sphere by blocks
0016 nblocks=size(C,1)/N;
0017 for iBlock=1:nblocks
0018     idx=(iBlock-1)*N + (1:N);
0019     CC=C(idx,idx);
0020     [V, S] = eig(CC) ;  
0021     V=real(V); S=real(S);
0022     [E,idx2] = sort(diag(S)', 'descend');
0023     topcs=V(:,idx2);
0024     EXP=1-10^-12; 
0025     E=E.^EXP; % break symmetry when x and y perfectly correlated (otherwise cols of x*A and y*B are not orthogonal)
0026     EE=(1./E); EE(find(E<=0))=0;
0027     A(idx,idx)=topcs*diag(sqrt(EE));
0028 end
0029 C=A'*C*A;
0030 
0031 
0032 % final PCA
0033 [V, S] = eig(C) ;
0034 V=real(V); S=real(S);
0035 [E, idx] = sort(diag(S)', 'descend') ;
0036 topcs = V(:,idx);
0037 A=A*topcs;
0038 %A=A(:,1:N);
0039 
0040 C=topcs'*C*topcs;
0041 score=diag(C);
0042 
0043 
0044 if nargout>2;
0045     AA=[];
0046     for iBlock=1:nblocks
0047         AA{iBlock}=A(N*(iBlock-1)+(1:N),:);
0048     end
0049 end
0050 
0051 % if nargout>3;
0052 %     if isempty(Nkeep); error('must specify Nkeep'); end
0053 %     AAA=[];
0054 %     for iBlock=1:nblocks
0055 %         % covariance of subject's data
0056 %         idx=(iBlock-1)*N + (1:N);
0057 %         C11=C(idx,idx);
0058 %         % covariance of selected MCCA components
0059 %         tmp=A(:,1:Nkeep);
0060 %         C22=tmp'*C*tmp;
0061 %         % cross covariance between subject's data and transformed data
0062 %         C12=C(idx,:)*tmp; clear tmp
0063 %         C21=C12';
0064 %         % CCA:
0065 %         [tmp]=nt_cca([],[],[],[C11,C12;C21,C22],N);
0066 %         AAA{iBlock}=tmp;
0067 %     end
0068 % end
0069 
0070 
0071 return
0072 
0073 % test code
0074 
0075 % 3 uncorrelated data sets
0076 figure(1); clf;
0077 x1=randn(10000,10); x2=randn(10000,10); x3=randn(10000,10); 
0078 x=[x1,x2,x3];
0079 C=x'*x;
0080 [A,score,AA]=nt_mcca(C,10);
0081 subplot 131; nt_imagescc(A); title('mcca transform');
0082 subplot 132; nt_imagescc(A'*C*A); title('covariance of transformed data');
0083 subplot 133; nt_imagescc(x'*(x*A)); title('crosscorr between raw & transf'); xlabel('transformed'); ylabel('raw');
0084 z=x*A;
0085 figure(11); clf;
0086 plot(mean(z.^2));
0087 
0088 % 3 identical data sets
0089 figure(2); clf
0090 x1=randn(10000,10); x=[x1,x1,x1]; 
0091 C=x'*x; 
0092 [A,score,AA]=nt_mcca(C,10);
0093 subplot 131; nt_imagescc(A); title('mcca transform');
0094 subplot 132; nt_imagescc(A'*C*A); title('covariance of transformed data');
0095 subplot 133; nt_imagescc(x'*(x*A)); title('cross correlation'); xlabel('transformed'); ylabel('raw');
0096 
0097 % 3 data sets with shared parts
0098 figure(3); clf
0099 x1=randn(10000,5); x2=randn(10000,5); x3=randn(10000,5); x4=randn(10000,5); 
0100 x=[x2,x1,x3,x1,x4,x1];
0101 C=x'*x; 
0102 [A,score,AA]=nt_mcca(C,10);
0103 subplot 131; nt_imagescc(A); title('mcca transform');
0104 subplot 132; nt_imagescc(A'*C*A); title('covariance of transformed data');
0105 subplot 133; nt_imagescc(x'*(x*A)); title('cross correlation'); xlabel('transformed'); ylabel('raw');
0106 
0107 % 3 data sets with parts shared 2 by 2
0108 figure(4); clf
0109 x1=randn(10000,5); x2=randn(10000,5); x3=randn(10000,5); x4=randn(10000,5); 
0110 x=[x2,x1,x3,x1,x3,x4];
0111 C=x'*x; 
0112 [A,score,AA]=nt_mcca(C,10);
0113 subplot 131; nt_imagescc(A); title('mcca transform');
0114 subplot 132; nt_imagescc(A'*C*A); title('covariance of transformed data');
0115 subplot 133; nt_imagescc(x'*(x*A)); title('cross correlation'); xlabel('transformed'); ylabel('raw');
0116

Generated on Sat 29-Apr-2023 17:15:46 by m2html © 2005