Home > chronux_2_00 > spikesort > utility > datatools > histxy.m

histxy

PURPOSE ^

HISTXY 2-Dimensional Density Histogram.

SYNOPSIS ^

function [counts,x_inds,y_inds] = histxy(x, y, varargin)

DESCRIPTION ^

HISTXY            2-Dimensional Density Histogram.
   COUNTS = HISTXY(X,Y,D), where X and Y are matrices of the same size,
   returns a D x D matrix COUNTS containing the numbers of points in a
   binned scatter plot of the entries of X vs the corresponding entries
   of Y.  The elements of X and Y are each independently discretized
   into D evenly-spaced bin centers, offset such that the extrema are
   are bin centers.  If D is not specified (or is empty) it defaults
   to 100.

   COUNTS = HISTXY(X,Y,[D1,D2]) produces a D2 x D1 matrix COUNTS, with
   the corresponding binning of X and Y.

   [COUNTS,X_INDS,Y_INDS] = HISTXY(X,Y,D) also returns the bin centers.
   The density is then visualized with IMAGESC(X_INDS,Y_INDS,COUNTS).

   [...] = HISTXT(X,Y,D,BANDWIDTH), for scalar BANDWIDTH, convolves
   COUNTS with an isotropic 2-D Gaussian kernel of standard deviation
   BANDWIDTH bins.  If BANDWIDTH is a two element vector, an anisotropic
   Gaussian is used with column standard deviation BANDWIDTH(1) bins and
   row standard deviation BANDWIDTH(2) bins.

   [...] = HISTXY(...,'log') uses the log of the counts (0's yield -Inf).

   HISTXY(...) without output arguments produces an image of the density.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [counts,x_inds,y_inds] = histxy(x, y, varargin)
0002 %HISTXY            2-Dimensional Density Histogram.
0003 %   COUNTS = HISTXY(X,Y,D), where X and Y are matrices of the same size,
0004 %   returns a D x D matrix COUNTS containing the numbers of points in a
0005 %   binned scatter plot of the entries of X vs the corresponding entries
0006 %   of Y.  The elements of X and Y are each independently discretized
0007 %   into D evenly-spaced bin centers, offset such that the extrema are
0008 %   are bin centers.  If D is not specified (or is empty) it defaults
0009 %   to 100.
0010 %
0011 %   COUNTS = HISTXY(X,Y,[D1,D2]) produces a D2 x D1 matrix COUNTS, with
0012 %   the corresponding binning of X and Y.
0013 %
0014 %   [COUNTS,X_INDS,Y_INDS] = HISTXY(X,Y,D) also returns the bin centers.
0015 %   The density is then visualized with IMAGESC(X_INDS,Y_INDS,COUNTS).
0016 %
0017 %   [...] = HISTXT(X,Y,D,BANDWIDTH), for scalar BANDWIDTH, convolves
0018 %   COUNTS with an isotropic 2-D Gaussian kernel of standard deviation
0019 %   BANDWIDTH bins.  If BANDWIDTH is a two element vector, an anisotropic
0020 %   Gaussian is used with column standard deviation BANDWIDTH(1) bins and
0021 %   row standard deviation BANDWIDTH(2) bins.
0022 %
0023 %   [...] = HISTXY(...,'log') uses the log of the counts (0's yield -Inf).
0024 %
0025 %   HISTXY(...) without output arguments produces an image of the density.
0026 
0027 %%%%%%%%%%%%%%%%%%%%%%%%%%%% Parse Inputs %%%%%%%%%%%%%%%%%%%%%%%%%%%%
0028 D = [100,100];  logflag = 0;   blur = 0;      % defaults
0029 
0030 x = x(:);  M = length(x);
0031 y = y(:);  N = length(y);
0032 if ((M~=N) || ~isnumeric(x) || ~isnumeric(y) || any(isinf(x(:))) || any(isinf(y(:))))
0033     error('First 2 input arguments must be numeric matrices of equal size with no +/- Inf elements.');
0034 end
0035 mask = (isnan(x) | isnan(y));   nanflag = any(mask(:));
0036 if (nanflag),  warning('HISTXY:NaNflag', 'NaN elements will be ignored.');  end
0037 
0038 if (length(varargin) > 0)
0039     tail = varargin{end};
0040     if (ischar(tail) && strcmpi(tail,'log'))      % If the last arg was 'log' ...
0041         varargin = varargin(1:(end-1));           % ... chomp it and set a flag.
0042         logflag = 1;
0043     end
0044     if (length(varargin) > 1)                     % If two args left, ...
0045         tail = varargin{end};                     % ... try to chomp bandwidth.
0046         if (isnumeric(tail) && length(tail)<=2)
0047             varargin = varargin(1:(end-1));
0048             if (length(tail) == 2),  sigma = tail;
0049             elseif (length(tail) == 1), sigma = [tail, tail];
0050             else   error('Bandwidth can not be empty.');
0051             end
0052             blur = 1;
0053             sigma(sigma == 0) = 0.01;  % equiv to 0 b/c 100 stds till next bin => below realmin
0054         end
0055     end
0056     if (length(varargin) > 0)
0057         tail = varargin{end};
0058         if (isnumeric(tail) && length(tail)<=2)   % If remaining arg was a len<=2 vector, ...
0059             varargin = varargin(1:end-1);         % ... chomp it and set the bin count.
0060             if (length(tail) == 2),  D = tail;
0061             elseif (length(tail) == 1),  D = [tail, tail];
0062             end  % if tail was empty, use default ...
0063         end
0064     end
0065     if (length(varargin) > 0),  error('Unknown syntax.');  end;
0066 end
0067 
0068 if (~isequal(D,round(D))), error('Number of bins D must be integer-valued.');  end;
0069 if (any(D==1)),  error('Number of bins D must be greater than 1.');  end;
0070 
0071 %%%%%%%%%%%%%%%%%%%%%%%%%%%% Rescale Data %%%%%%%%%%%%%%%%%%%%%%%%%%%%
0072 % Separately scale the x and y data
0073 [x,min1,max1] = rescale(x,1,D(1));   x = round(x);
0074 [y,min2,max2] = rescale(y,1,D(2));   y = round(y);
0075 
0076 % Bin centers are equally spaced over the range of the unscaled data
0077 x_inds = linspace(min1,max1,D(1));
0078 y_inds = linspace(min2,max2,D(2));
0079 
0080 % Mask NaNs
0081 if (nanflag),  D=D+1;  x(mask)=D(1);  y(mask)=D(2);  end;
0082 
0083 
0084 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Histogram %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0085 counts = CORE_histxy(x,y,D(1),D(2));
0086 
0087 
0088 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Clean Up %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0089 if (nanflag),  counts(end,:) = []; counts(:,end) = [];   D = D-1;  end;
0090 
0091 if (blur)   % use separable kernels for convolution with least odd length  >= num bins
0092     x_kernel = exp(-((-floor(D(1)/2):floor(D(1)/2)).^2)./(2*sigma(1).^2));
0093     x_kernel = x_kernel ./ sum(x_kernel);
0094     y_kernel = exp(-((-floor(D(2)/2):floor(D(2)/2)).^2)./(2*sigma(2).^2))';
0095     y_kernel = y_kernel ./ sum(y_kernel);
0096     
0097     counts = conv2(y_kernel, x_kernel, counts, 'same');
0098 end
0099 
0100 if (logflag), o=warning('MATLAB:log:logOfZero', 'off');  counts=log(counts);  warning(o);  end;
0101 
0102 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Graphics %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0103 if (nargout > 0),  return;  end;
0104 imagesc(x_inds, y_inds, counts);   axis xy;
0105 clear counts x_inds y_inds  % clear these so nothing is dumped to output

Generated on Fri 15-Aug-2008 11:35:42 by m2html © 2003