bivgui.m

常用ROBUST STATISTICAL

function bivgui(arg)
% BIVGUI  Distribution Shapes - Bivariate GUI
%
% This GUI function allows one to explore the bivariate distributions in
% the data set. These are the distributions of pairs of columns of the data
% matrix X.
%
% This GUI provides the means to create bivariate polar smooths in
% scatterplots, hexagonal binning in scatterplots and bivariate histograms.
%
% This can be called from the edagui main GUI or as a stand-alone function
% from the command line. To call from the command line, use
%
%       bivgui
%
% NOTE: Coloring by groups does not affect the plots in this GUI.
%
%   Exploratory Data Analysis Toolbox, December 2006
%   Martinez and Martinez, Exploratory Data Analysis with MATLAB
%   CRC Press

% First set up the layout if it does not exist.
flg = findobj('tag','bivgui');
if isempty(flg)
    % then create the gui
    dsblayout
elseif nargin == 0
    % bring it forward
    figure(flg)
end

if nargin == 0
    arg = ' ';
end
if strcmp(arg,'polarplot')
    % Display all possible pairwise scatterplots with polar smooths.
    disppolar
elseif strcmp(arg,'disphist')
    % Display all possible pairwise bivariate histograms.
    disphist
elseif strcmp(arg,'hexplot')
    % Display scatterplots with hexagonal smoothing.
    disphex
elseif strcmp(arg,'close')
    % in other gui's we will do some housekeeping. With this gui, we do not
    % have to do so. Obviously, the user will want to keep the data from
    % the loadgui for other applications. That is the purpose.
    tg = findobj('tag','bivgui');
    H = get(tg,'userdata');
    if ~isempty(H.plots)
        button = questdlg('Closing this GUI will close the associated plot windows.',...
            'Closing GUI Warning','OK','Cancel','Cancel');
        if strcmp(button,'Cancel')
            return
        else
            close(H.plots)
        end
    end
    delete(tg)
end

%%%%%%%%%%%%%%%%%%%%%%%%%  SUB FUNCTIONs %%%%%%%%%%%%%%%%%%%%%%%%%%%%
function disppolar    
% This function displays all possible pairwise scatterplots.
% Adds polar smoothing.
% Get the data matrix. Get the GUI info.
ud = get(0,'userdata');
if isempty(ud.X)
    errordlg('You must load up some data first.')
    return
end
[n,p] = size(ud.X);
tg = findobj('tag','bivgui');
H = get(tg,'userdata');
% Get the value of the radio button indicating robust loess or not.
rflag = get(H.robflag,'value');
% Get the value of the smoothing parameter
lam = str2double(get(H.smooth,'string'));
if ~(lam > 0 & lam < 1)
    errordlg('Smoothing parameter must be between 0 and 1.','Data Entry Error')
    return
end
% Get the degree of the polynomial
deg = str2double(get(H.degree,'string'));
if deg ~= 1 & deg ~= 2
    errordlg('Degree of polynomial must be 1 (linear) or 2 (quadratic).')
    return
end
% Get the dimensions to plot.
dims = get(H.dimpol,'string');
try
    if strcmp('all',dims)
        % Wants to plot all dimensions
        dim = 1:p;
    else
        % Wants just a subset - convert to numbers.
        eval(['dim = [' dims '];'])
    end
catch
    errordlg('Data entry error in edit box.','Data Enty Error')
    return
end
if length(dims) == 1
    errordlg('Must have more than 1 dimension specified in edit box.','Data Entry Error')
    return
end
hf = figure;
set(hf,'numbertitle','off','name','EDA: Pairwise Scatterplots with Polar Smooths')
% Upon figure close, this should delete from the array.
set(hf,'CloseRequestFcn',...
    'tg = findobj(''tag'',''bivgui''); H = get(tg,''userdata''); H.plots(find(H.plots == gcf)) = []; set(tg,''userdata'',H); delete(gcf)')
H.plots = [H.plots, hf];
set(tg,'userdata',H)
plotpolar(ud.X(:,dim),lam,deg,rflag,dim)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
function disphist
% This function displays the histograms.
  
% Get the data matrix. Get the GUI info.
ud = get(0,'userdata');
if isempty(ud.X)
    errordlg('You must load up some data first.')
    return
end
[n,p] = size(ud.X);
tg = findobj('tag','bivgui');
H = get(tg,'userdata');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
% This whole section was changed on August 19. The bandwidths do not work
% correctly and it is difficult for data enty. So, change to allow the user
% to enter the number of bins instead.
% % Get the bandwidths
% bs = get(H.band,'string');
% if strcmp(bs,'Default')
%     % Use the default bandwidths
%     bw = 3.5*std(ud.X)*n^(-1/(p+1));
%     bws = num2str(bw,2);
%     set(H.band,'string',bws)
% else
%     % Get what the user enetered.
%     eval(['bw = [' bs '];'])
% end
% % Get the dimensions to plot.
% dims = get(H.dimhist,'string');
% try
%     if strcmp('all',lower(dims))
%         % Wants to plot all dimensions
%         dim = 1:p;
%     else
%         % Wants just a subset - convert to numbers.
%         eval(['dim = [' dims '];'])
%     end
% catch
%     errordlg('Data entry error in edit box.','Data Enty Error')
%     return
% end
% if length(dim) == 1
%     errordlg('Must have more than 1 dimension specified in edit box.','Data Entry Error')
%     return
% end
% if length(bw) ~= length(dim)
%     errordlg('You must enter one bandwidth per dimension.','Date Entry Error')
%     return
% end
% % check the number of bins. If too small, then the user should enter
% % the h values. Provide an error.
%  % This is just a check.
% rng = max(ud.X(dims)) - min(ud.X(dims));
% nb = rng./bw;
% if any(nb < 4)
%     errordlg('Some of the bandwidths are too large, yielding too few bins. Make them smaller.')
%     return
% end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Get the number of bins
% Default number of bins is 10.
nb = round(str2double(get(H.band,'string')));
% Get the dimensions to plot.
dims = get(H.dimhist,'string');
try
    if strcmp('all',lower(dims))
        % Wants to plot all dimensions
        dim = 1:p;
    else
        % Wants just a subset - convert to numbers.
        eval(['dim = [' dims '];'])
    end
catch
    errordlg('Data entry error in edit box.','Data Enty Error')
    return
end
if length(dim) == 1
    errordlg('Must have more than 1 dimension specified in edit box.','Data Entry Error')
    return
end

% check the number of bins. 
if nb < 5
    errordlg('Should have 5 or more bins.')
    return
end

rng = max(ud.X(:,dim)) - min(ud.X(:,dim));
bw = rng/nb;



hf = figure;
set(hf,'numbertitle','off','name','EDA: Pairwise Bivariate Histograms')
% Upon figure close, this should delete from the array.
set(hf,'CloseRequestFcn',...
    'tg = findobj(''tag'',''bivgui''); H = get(tg,''userdata''); H.plots(find(H.plots == gcf)) = []; set(tg,''userdata'',H); delete(gcf)')
H.plots = [H.plots, hf];
set(tg,'userdata',H)
plothist(ud.X(:,dim),bw,dim)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
function disphex
% This function displays the histograms.
  
% Get the data matrix. Get the GUI info.
ud = get(0,'userdata');
if isempty(ud.X)
    errordlg('You must load up some data first.')
    return
end
[n,p] = size(ud.X);
tg = findobj('tag','bivgui');
H = get(tg,'userdata');
% Get the bandwidths
Nb = round(str2double(get(H.bins,'string')));
% Get the dimensions to plot.
dims = get(H.dimhex,'string');
try
    if strcmp('all',lower(dims))
        % Wants to plot all dimensions
        dim = 1:p;
    else
        % Wants just a subset - convert to numbers.
        eval(['dim = [' dims '];'])
    end
catch
    errordlg('Data entry error in edit box.','Data Enty Error')
    return
end
if length(dim) == 1
    errordlg('Must have more than 1 dimension specified in edit box.','Data Entry Error')
    return
end


hf = figure;
set(hf,'numbertitle','off','name','EDA: Pairwise Scatterplots with Hexagonal Binning')
% Upon figure close, this should delete from the array.
set(hf,'CloseRequestFcn',...
    'tg = findobj(''tag'',''bivgui''); H = get(tg,''userdata''); H.plots(find(H.plots == gcf)) = []; set(tg,''userdata'',H); delete(gcf)')
H.plots = [H.plots, hf];
set(tg,'userdata',H)
plothex(ud.X(:,dim),Nb,dim)


%%%%%%%%%%%%%%%%%%%%   HELPER FUNCTIONS %%%%%%%%%%%%%%%%%%%%%%%

function plothex(X,Nb,dim)
% does the scatterplot matrix, but adds polar smooths.
% lam is the smoothing parameter, deg is the degree of the polynomial and
% rflag is 1 means to use robust loess, and a 0 means to use regular.
ud = get(0,'userdata'); % for variable names

Hfig = findobj('name','EDA: Pairwise Scatterplots with Hexagonal Binning');
set(Hfig, 'units','normalized',...
    'position',  [0 0.0365 0.9678 0.8750],...
    'toolbar','none',...
    'menubar','none');
[n,p] = size(X);
% Order of axes is left-right, top-bottom
pp = 0;
I = 0; J = 0;
for j = (p-1):-1:0
    I = I + 1;
    for i = 0:(p-1)
        J = J + 1;
        box on
        pp = pp + 1;
        subplot(p,p,pp)
        if I~=J
            % The following is the column index (to data) for the X and Y
            % variables.
            set(gca,'yticklabel','','xticklabel','','ticklength',[0 0])
            % Do the bivariate histogram.
            Xt = [X(:,J),X(:,I)];
            hexplot(Xt,Nb)
            drawnow
        else
            axis off
            % This is a center axes - plot the variable name.
            if isnumeric(ud.varlab)
                text(0.35,0.45, num2str(ud.varlab(dim(I))))
            elseif iscell(ud.varlab)
                text(0.35,0.45,ud.varlab{dim(I)})
            end
        end  % if stmt
    end   % for j loop
    J = 0;
end   % for i loop


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5

function plothist(X,bw,dim)
% does the scatterplot matrix, but adds polar smooths.
% lam is the smoothing parameter, deg is the degree of the polynomial and
% rflag is 1 means to use robust loess, and a 0 means to use regular.
ud = get(0,'userdata'); % for variable names

Hfig = findobj('name','EDA: Pairwise Bivariate Histograms');
set(Hfig, 'units','normalized',...
    'position',  [0 0.0365 0.9678 0.8750],...
    'toolbar','none',...
    'menubar','none');
[n,p] = size(X);
% Order of axes is left-right, top-bottom
pp = 0;
I = 0; J = 0;
for j = (p-1):-1:0
    I = I + 1;
    for i = 0:(p-1)
        J = J + 1;
        box on
        pp = pp + 1;
        subplot(p,p,pp)
        if I~=J
            % The following is the column index (to data) for the X and Y
            % variables.
            set(gca,'yticklabel','','xticklabel','','ticklength',[0 0])
            % Do the bivariate histogram.
            Xt = [X(:,J),X(:,I)];
            Z = cshist2d(Xt,1,[bw(J),bw(I)]);
            box off
          
            drawnow
        else
            axis off
            % This is a center axes - plot the variable name.
            if isnumeric(ud.varlab)
                text(0.35,0.45, num2str(ud.varlab(dim(I))))
            elseif iscell(ud.varlab)
                text(0.35,0.45,ud.varlab{dim(I)})
            end
        end  % if stmt
    end   % for j loop
    J = 0;
end   % for i loop

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function  plotpolar(X,lam,deg,rflag,dim)
% does the scatterplot matrix, but adds polar smooths.
% lam is the smoothing parameter, deg is the degree of the polynomial and
% rflag is 1 means to use robust loess, and a 0 means to use regular.
ud = get(0,'userdata'); % for variable names

Hfig = findobj('name','EDA: Pairwise Scatterplots with Polar Smooths');
set(Hfig, 'units','normalized',...
    'position',  [0 0.0365 0.9678 0.8750],...
    'toolbar','none',...
    'menubar','none');
[n,p] = size(X);
minx = min(X);
maxx = max(X);
rngx = range(X);
% set up the axes
H.IndX = zeros(p,p);    % X dim for data
H.IndY = zeros(p,p);    % Y dim for data
H.AxesLims = cell(p,p); % Axes limits.
H.Haxes = zeros(p,p);   % Axes handles.
H.HlineReg = zeros(p,p);    % Line handles to non-highlighted data.

% Order of axes is left-right, top-bottom
% Take up the entire figure area with axes.
I = 0; J = 0;
for j = (p-1):-1:0
    I = I + 1;
    for i = 0:(p-1)
        J = J + 1;
        pos = [i/p j/p 1/p 1/p];
        pos = floor(pos*100)/100;
        H.Haxes(I,J) = axes('pos',[i/p j/p 1/p 1/p]);
        box on
        if I~=J
            % The following is the column index (to data) for the X and Y
            % variables.
            H.IndX(I,J) = J;
            H.IndY(I,J) = I;
            set(gca,'yticklabel','','xticklabel','','ticklength',[0 0])
            % Do the scatterplot.
            H.HlineReg(I,J) = line('xdata',X(:,J),'ydata',X(:,I),...
                'markersize',3,'marker','o','linestyle','none');
            % Do the polar smooth and add to plot.
            hold on
            [xhat2,yhat2] =  polarloess(X(:,J),X(:,I),lam,deg,rflag);
            plot(xhat2,yhat2)