nt_find_outlier_trials

[idx,d,mn,idx_unsorted]=nt_find_outlier_trials(x,criterion,plot,regress_flag) - find outlier trials

  idx: indices of trials to keep
  d: relative deviations from mean
  
  x: data (time * channels * trials)
  criterion: keep trials less than criterion from mean
  disp: if true plot trial deviations before and after 
  regress_flag: if true regress out mean, rather than subtract

  For example criterion=2 rejects trials that deviate from the mean by
  more than twice the average deviation from the mean.

0001 function [idx,d]=nt_find_outlier_trials(x,criterion,disp_flag,regress_flag)
0002 %[idx,d,mn,idx_unsorted]=nt_find_outlier_trials(x,criterion,plot,regress_flag) - find outlier trials
0003 %
0004 %  idx: indices of trials to keep
0005 %  d: relative deviations from mean
0006 %
0007 %  x: data (time * channels * trials)
0008 %  criterion: keep trials less than criterion from mean
0009 %  disp: if true plot trial deviations before and after
0010 %  regress_flag: if true regress out mean, rather than subtract
0011 %
0012 %  For example criterion=2 rejects trials that deviate from the mean by
0013 %  more than twice the average deviation from the mean.
0014 %
0015 
0016 if nargin<2; criterion=inf; end
0017 if nargin<3; disp_flag=1; end
0018 if nargin<4; regress_flag=0; end
0019 if ndims(x)>3; error('x should be 2D or 3D'); end
0020 
0021 if ndims(x)==3;
0022     [m,n,o]=size(x);
0023     x=reshape(x,m*n,o);
0024 else
0025     [~,o]=size(x);
0026 end
0027 
0028 if regress_flag; error('need to take nt_tsregress.m out of DISUSE'); end
0029 
0030 mn=mean(x,2);
0031 if regress_flag
0032     mn=nt_tsregress(x,mean(x,2));  % regression
0033 else
0034     mn=repmat(mn(:),1,o);       % mean
0035 end
0036 d=x-mn; % difference from mean
0037 dd=zeros(1,size(d,2));
0038 for k=1:size(d,2); dd(k)=sum(d(:,k).^2); end
0039 d=dd; clear dd;
0040 d=d/(sum(x(:).^2)/o);
0041 
0042 idx=find(d<criterion(1));
0043 
0044 if nargout==0;
0045     % just plot deviations
0046     plot(d,'.-');
0047     xlabel('trial'); ylabel('normalized deviation from mean'); 
0048     clear idx d mn idx_unsorted
0049 else
0050     if disp_flag
0051         % plot deviations before & after outlier removal
0052         figure(100); clf
0053         nt_banner('outlier trials');
0054         
0055         subplot 121; 
0056         plot(d,'.-'); hold on; 
0057         plot(setdiff(1:o,idx), d(setdiff(1:o,idx)), '.r');
0058         xlabel('trial'); ylabel('normalized deviation from mean'); title(['before, ',num2str(numel(d))]);
0059         drawnow
0060         
0061         subplot 122; 
0062         [~,dd]=nt_find_outlier_trials(x(:,idx),0,[]);
0063         plot(dd,'.-');
0064         xlabel('trial'); ylabel('normalized deviation from mean'); title(['after, ',num2str(numel(idx))]);
0065         drawnow
0066     end
0067     
0068 end
0069      
0070 criterion=criterion(2:end);
0071 if ~isempty(criterion)
0072     idx=nt_find_outlier_trials(x(:,idx),criterion,disp_flag,regress_flag);
0073     idx = idx(idx2); % otherwise, idx doesn?t correspond to original matrix anymore
0074 end