Home > NoiseTools > nt_statmatrix.m

nt_statmatrix

PURPOSE ^

stats=nt_statMatrix(x,plot_params) - calculate statistics arrays for each dim of matrix

SYNOPSIS ^

function stats=nt_statMatrix(x,plot_params)

DESCRIPTION ^

stats=nt_statMatrix(x,plot_params) - calculate statistics arrays for each dim of matrix

  stats: array of statistics arrays

  x: data to stat
  plot_params: parameters

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function stats=nt_statMatrix(x,plot_params)
0002 %stats=nt_statMatrix(x,plot_params) - calculate statistics arrays for each dim of matrix
0003 %
0004 %  stats: array of statistics arrays
0005 %
0006 %  x: data to stat
0007 %  plot_params: parameters
0008 
0009 if nargin<1; error('!'); end
0010 if nargin<2;
0011     plot_params.bottom=0.1;
0012     plot_params.height=0.8;
0013 end
0014 
0015 % count true number of dimensions
0016 nDims=ndims(x);
0017 sizeX=size(x);
0018 if nDims==2 && (sizeX(1)==1 || sizeX(2)==1)
0019     nDims=1;
0020 end
0021 
0022 % if we're plotting, concatenate dimensions beyond 4th
0023 if nargout==0 && ndims(x)>4; x=x(:,:,:,:); nDims=4; end
0024 
0025 stats={}; 
0026 if 0
0027     % computationally expensive
0028     for iDim=1:nDims
0029         sz=size(x);
0030         x=reshape(x,sz(1),[]);
0031         stats{iDim}.iDim=iDim;
0032         stats{iDim}.n=size(x,2);
0033         stats{iDim}.min=min(x')';
0034         stats{iDim}.max=max(x')';
0035         stats{iDim}.mean=nanmean(x')';
0036         stats{iDim}.var=nanvar(x')';
0037         x=reshape(x,sz);
0038         x=shiftdim(x,1);
0039     end
0040 else
0041     % cheaper
0042     sz=size(x);
0043     for iDim=1:nDims
0044         stats{iDim}.iDim=iDim;
0045         stats{iDim}.n=size(x,iDim);
0046         
0047         % min
0048         y=x;
0049         for k=1:iDim-1;
0050              y=squeeze(min(y,[],1));
0051         end
0052         if iDim<nDims; y=min(y(:,:),[],2); else y=min(y(:,:),[],1); end
0053         stats{iDim}.min=y(:); 
0054 
0055         % max
0056         y=x;
0057         for k=1:iDim-1;
0058             y=squeeze(max(y,[],1));
0059         end
0060         if iDim<nDims; y=max(y(:,:),[],2); else y=max(y(:,:),[],1); end
0061         stats{iDim}.max=y(:); 
0062 
0063         % mean
0064         y=x;
0065         for k=1:iDim-1;
0066             y=squeeze(nanmean(y,1));
0067         end
0068         if iDim<nDims; y=nanmean(y(:,:),2); else y=nanmean(y(:,:),1); end
0069         stats{iDim}.mean=y(:); 
0070       
0071         % var
0072         y=x.^2;
0073         for k=1:iDim-1;
0074             y=squeeze(nanmean(y,1));
0075         end
0076         if iDim<nDims; y=nanmean(y(:,:),2); else y=nanmean(y(:,:),1); end
0077         stats{iDim}.var=y(:)-stats{iDim}.mean.^2; 
0078         %nt_whoss
0079         
0080     end
0081 end
0082 
0083 
0084 % plot the statistics
0085 if nargout==0;
0086     if isempty(get(0,'currentfigure'));
0087         figH=figure;
0088         set(gcf,'color',[1 1 1]);
0089     end
0090     figH=figure(gcf); %clf
0091     %set(gcf,'name',inputname(1));
0092 
0093     fontsize=12;
0094     
0095     % plot in 1,2,3 or 4 panels depending on number of dimensions
0096     switch nDims
0097         case 1
0098             axes('position',[0.05, plot_params.bottom, 0.93, plot_params.height]);  set(gca,'box','on','fontsize',fontsize);
0099             plot(x, 'k'); title(['n = ', num2str(sizeX(1))]);
0100             xlabel('samples');
0101         case 2
0102             axes('position',[0.05, plot_params.bottom, 0.45, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0103             nt_plotstats(stats{1});
0104             xlabel('samples'); title(['n = ', num2str(sizeX(1))]);
0105             axes('position',[0.53, plot_params.bottom, 0.45, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0106             nt_plotstats(stats{2});
0107             xlabel('samples'); title(['n = ', num2str(sizeX(2))]); set(gca,'ytick',[]);
0108         case 3
0109             axes('position',[0.05, plot_params.bottom, 0.3, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0110             nt_plotstats(stats{1});
0111             xlabel('samples'); title(['n = ', num2str(sizeX(1))]);
0112             axes('position',[0.37, plot_params.bottom, 0.3, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0113             nt_plotstats(stats{2});
0114             xlabel('samples'); title(['n = ', num2str(sizeX(2))]); set(gca,'ytick',[]);
0115             axes('position',[0.69, plot_params.bottom, 0.3, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0116             nt_plotstats(stats{3});
0117             xlabel('samples'); title(['n = ', num2str(sizeX(3))]); set(gca,'ytick',[]);
0118         otherwise % limit to 4 panels (last dims are concatenated)
0119             axes('position',[0.05, plot_params.bottom, 0.2, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0120             nt_plotstats(stats{1});
0121             xlabel('samples'); title(['n = ', num2str(sizeX(1))]);
0122             axes('position',[0.27, plot_params.bottom, 0.2, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0123             nt_plotstats(stats{2});
0124             xlabel('samples'); title(['n = ', num2str(sizeX(2))]); set(gca,'ytick',[]);
0125             axes('position',[0.49, plot_params.bottom, 0.2, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0126             nt_plotstats(stats{3});
0127             xlabel('samples'); title(['n = ', num2str(sizeX(3))]); set(gca,'ytick',[]);
0128             axes('position',[0.71, plot_params.bottom, 0.2, plot_params.height]); set(gca,'box','on','fontsize',fontsize);
0129             nt_plotstats(stats{4});
0130             xlabel('samples'); 
0131             title_string=['n = ', num2str(sizeX(4))];
0132             for k=5:numel(sizeX)
0133                 title_string=[title_string,'*',num2str(sizeX(k))];
0134             end
0135             title(title_string); 
0136             set(gca,'ytick',[]);
0137             
0138     end
0139 end
0140             
0141 
0142 function nt_plotstats(stats)
0143 %h=nt_plotstats(stats)
0144 %
0145 %  stats: stats stucture
0146 
0147 holdStatus=ishold;
0148 hold on
0149 if isfield(stats,'min')
0150     nt_plot2(1:size(stats.min,1), [stats.min,stats.max], [1 1 1]*0.9);
0151 end
0152 if isfield(stats,'var')
0153     sd=sqrt(stats.var);
0154     nt_plot2(1:size(stats.mean,1), [stats.mean+sd,stats.mean-sd], [1 1 1]* 0.5);
0155 end
0156 nt_plot2(1:size(stats.mean,1),[stats.mean,stats.mean], [0 0 0]);
0157 if holdStatus;
0158     hold on;
0159 else
0160     hold off
0161 end
0162 
0163 function h=nt_plot2(x,y,c)
0164 %h=nt_plot2(x,y,c) - color region between two plots
0165 %
0166 %
0167 %  x: abscissa
0168 %  y: ordinate (1 or 2 columns)
0169 %  c: color (see 'fill')
0170 %
0171 %  h: vector of handles to patches
0172 %
0173 %  nt_plot2(x
0174 
0175 if nargin<1; error('!'); end
0176 
0177 % process parameters
0178 if nargin==1;
0179     y=x;
0180     x=(1:size(y,1))';
0181     c=[0 0 0];
0182 elseif nargin==2;
0183     c=[0,0,0];
0184 elseif nargin==3;
0185     ;
0186 else
0187     error('!');
0188 end
0189 
0190 % format data
0191 if size(y,2)==1
0192     y=[y,zeros(size(y))];
0193 elseif size(y,2)>2
0194     error('!');
0195 else
0196     ;
0197 end
0198 x=x(:);
0199 
0200 if 1
0201 % make sure that y(:,1)<y(:,2);
0202     yy=y;
0203     yy(:,1)=min(y,[],2);
0204     yy(:,2)=max(y,[],2);
0205     y=yy;
0206 
0207     % downsample if data array is too large
0208     TARGET_N=2000;
0209     if size(x,1)>TARGET_N
0210         DSR=ceil(size(x,1)/TARGET_N);
0211         n=floor(size(x,1)/DSR);
0212         x=x(1:DSR:end);
0213         x=x(1:n);
0214         y_extra=y(DSR*n:end,:);
0215         y=y(1:DSR*n,:);
0216         a=min(reshape(y(:,1),DSR,size(y,1)/DSR));
0217         b=max(reshape(y(:,2),DSR,size(y,1)/DSR));
0218         a(end)=min(a(end),min(y_extra(:,1)));
0219         b(end)=max(b(end),max(y_extra(:,2)));
0220         y=[a',b'];
0221     end
0222 end
0223 
0224 % draw plot
0225 h=fill([x;flipud(x)],[y(:,1);flipud(y(:,2))],c, 'edgecolor', c);
0226 if x(end)>x(1); 
0227     xlim([x(1)-1,x(end)+1]);
0228 end
0229 set(gca,'box','on');
0230 
0231 % from stats toolbox
0232 function m = nanmean(x,dim)
0233 %%
0234 
0235 % Find NaNs and set them to zero
0236 nans = isnan(x);
0237 x(nans) = 0;
0238 
0239 if nargin == 1 % let sum deal with figuring out which dimension to use
0240     % Count up non-NaNs.
0241     n = sum(~nans);
0242     n(n==0) = NaN; % prevent divideByZero warnings
0243     % Sum up non-NaNs, and divide by the number of non-NaNs.
0244     m = sum(x) ./ n;
0245 else
0246     % Count up non-NaNs.
0247     n = sum(~nans,dim);
0248     n(n==0) = NaN; % prevent divideByZero warnings
0249     % Sum up non-NaNs, and divide by the number of non-NaNs.
0250     m = sum(x,dim) ./ n;
0251 end
0252

Generated on Mon 10-Nov-2014 14:40:42 by m2html © 2005