nt_subspace_prune5

[Y]=nt_subspace_prune(X,npass,thresh) - local cleaning matrices

  Y: denoised data
  
  X: data to denoise (nsamples X nchans X ntrials)
  npass: number of passes [default: 1]
  thresh: threshold power ratio between segment & all [default: 10]

0001 function [Y,scores]=nt_subspace_prune5(X,npass,thresh)
0002 %[Y]=nt_subspace_prune(X,npass,thresh) - local cleaning matrices
0003 %
0004 %  Y: denoised data
0005 %
0006 %  X: data to denoise (nsamples X nchans X ntrials)
0007 %  npass: number of passes [default: 1]
0008 %  thresh: threshold power ratio between segment & all [default: 10]
0009 
0010 if nargin<2||isempty(npass); npass=1; end
0011 if nargin<3||isempty(thresh); thresh=1; end
0012 
0013 if isnumeric(X)
0014     if ndims(X)<3; error('!'); end
0015     tmp={};
0016     for iTrial=1:size(X,3); 
0017         tmp{iTrial}=X(:,:,iTrial); 
0018     end
0019     X=tmp;
0020     [Y,scores]=nt_subspace_prune5(X,npass,thresh);
0021     tmp=zeros(size(Y{1},1),size(Y{2},2),numel(Y));
0022     for iTrial=1:numel(X) 
0023         tmp(:,:,iTrial)=Y{iTrial}; 
0024     end
0025     Y=tmp;
0026     return
0027 end
0028 
0029 ntrials=numel(X);
0030 nchans=size(X{1},2);
0031 
0032 original_power=nt_wpwr(X);
0033 scores=[]; D=[];
0034 for iPass=1:npass
0035     
0036     X=nt_demean2(X);
0037     
0038     % covariance of full data
0039     [C0,tw]=nt_cov(X); 
0040     C0=C0/tw;
0041     
0042     iBest=0; best_score=0; 
0043     CC1=zeros(nchans,nchans,ntrials);
0044     for iTrial=1:numel(X)
0045         a=X{iTrial};
0046         
0047         % covariance of this trial
0048         C1=nt_cov(a)/size(a,1);
0049         
0050         % contrast this trial with rest
0051         [todss,pwr0,pwr1]=nt_dss0(C0,C1,[],0);
0052         fromdss=pinv(todss);
0053         
0054         % denoising matrix for this trial
0055         D{iTrial}=todss(:,2:end)*fromdss(2:end,:);
0056         
0057         % is this trial the most excentric?
0058         if pwr1(1)/pwr0(1)>best_score
0059             iBest=iTrial;
0060             best_score=pwr1(1)/pwr0(1);
0061         end
0062         scores(iPass,iTrial)=pwr1(1)/pwr0(1);
0063     end
0064     
0065     % remove most excentric component of most excentric trials
0066     if best_score>thresh
0067         X{iBest}=nt_mmat(X{iBest},D{iBest});
0068     else
0069         break;
0070     end
0071 
0072     figure(10); clf; nt_imagescc(scores); colorbar;
0073     %disp(nt_wpwr(X)/original_power);
0074 end
0075     
0076 Y=X;