Home > matutils > outlierDims.m

outlierDims

PURPOSE ^

function flag=outlier(pt, dim, nsig,upperCut, lowerCut, epsil)

SYNOPSIS ^

function flag=outlierDims(pt,dim, nsig, upperCut, lowerCut, epsil)

DESCRIPTION ^

 function flag=outlier(pt, dim, nsig,upperCut, lowerCut, epsil)

 flag points that is 'nsig' sigmas away from the mean.  the returned
 flag array is NOT bit oriented.

 if the std is smaller than epsil, nothing will get flagged.  epsil
 default value is 0.01

 dim is the dimension over which we want to flag.
 only works on 3D data..if you want more, contact stephen.
 Michael Loh

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function flag=outlierDims(pt,dim, nsig, upperCut, lowerCut, epsil)
0002 
0003 % function flag=outlier(pt, dim, nsig,upperCut, lowerCut, epsil)
0004 %
0005 % flag points that is 'nsig' sigmas away from the mean.  the returned
0006 % flag array is NOT bit oriented.
0007 %
0008 % if the std is smaller than epsil, nothing will get flagged.  epsil
0009 % default value is 0.01
0010 %
0011 % dim is the dimension over which we want to flag.
0012 % only works on 3D data..if you want more, contact stephen.
0013 % Michael Loh
0014 
0015 if (~exist('epsil'))
0016   epsil=0.01;
0017 end
0018 
0019 if (~exist('upperCut'))
0020   upperCut = inf;
0021 end
0022 
0023 if (~exist('lowerCut'))
0024   lowerCut = -inf;
0025 end
0026 
0027 
0028 if(ndims(pt)>3)
0029     error('only supports up to 3D data');
0030 end
0031 
0032 % let's permute everything so that we always operate on teh first
0033 % dimension.
0034 switch dim
0035   case 1
0036     data = pt;
0037     
0038   case 2
0039     data = permute(pt, [2 1 3]);
0040 
0041   case 3
0042     data = permute(pt, [3 1 2]);
0043 end
0044 
0045 % initialize flag array
0046 flag=zeros(size(data));
0047 
0048 
0049 % first we do the upperCut and lower cuts
0050 flag = ((data>upperCut) | (data < lowerCut))>0;
0051 
0052 
0053 if(size(flag,1)<500)
0054     % look for single points driving the rms WAY up
0055     % unnecessary if we have many many points, as the prob is unlikely.
0056     tic
0057     for i=1:size(pt,1)
0058         x=data;x(i,:,:)=NaN;
0059         r(i,:,:)=nanstd(x);
0060     end
0061    toc
0062    meanR = nanmean(r);
0063    meanR = repmat(meanR,  [size(r,1) 1 1]);
0064    f = find(r./meanR < (1 - 3/sqrt(size(r,1))));
0065    data(f)=nan; flag(f) = 1;
0066 end
0067 
0068 data(flag) = nan;
0069 % exit if the variance is tiny
0070 sig=nanstd(data);
0071 epsil = ones(size(sig))*epsil;
0072 
0073 noFlag = (sig<epsil);
0074 if(all(noFlag(:)))
0075     % no more need to flag
0076     return;
0077 end
0078 
0079 
0080 % flag
0081 meanVal = nanmean(data);
0082 meanVal = repmat(meanVal, [size(data,1) 1 1]);
0083 sigVal  = repmat(sig, [size(data,1) 1 1]);
0084 
0085 devs    = abs((pt-meanVal)./sigVal);
0086 devFlag = devs>nsig;
0087 
0088 if(any(noFlag(:)))
0089     [x y] = find(squeeze(noFlag));
0090     if(~isempty(x))
0091         for m=1:length(x)
0092             devFlag(:,x(m),y(m)) = 0;
0093         end
0094     end
0095 end
0096 
0097 flag = (flag + devFlag)>0;
0098 
0099 % now we permute back to original dimension
0100 switch dim
0101   case 1
0102     flag = flag;
0103 
0104   case 2
0105     flag = permute(flag, [2 1 3]);
0106 
0107   case 3
0108     flag = permute(flag, [3 1 2]);
0109 end
0110 
0111 return;

Generated on Sun 14-Jun-2015 17:12:45 by m2html © 2005