annotate _FullBNT/KPMtools/plotcov2New.m @ 9:4ea6619cb3f5 tip

removed log files
author matthiasm
date Fri, 11 Apr 2014 15:55:11 +0100
parents b5b38998ef3b
children
rev   line source
matthiasm@8 1 % PLOTCOV2 - Plots a covariance ellipsoid with axes for a bivariate
matthiasm@8 2 % Gaussian distribution.
matthiasm@8 3 %
matthiasm@8 4 % Usage:
matthiasm@8 5 % [h, s] = plotcov2(mu, Sigma[, OPTIONS]);
matthiasm@8 6 %
matthiasm@8 7 % Inputs:
matthiasm@8 8 % mu - a 2 x 1 vector giving the mean of the distribution.
matthiasm@8 9 % Sigma - a 2 x 2 symmetric positive semi-definite matrix giving
matthiasm@8 10 % the covariance of the distribution (or the zero matrix).
matthiasm@8 11 %
matthiasm@8 12 % Options:
matthiasm@8 13 % 'conf' - a scalar between 0 and 1 giving the confidence
matthiasm@8 14 % interval (i.e., the fraction of probability mass to
matthiasm@8 15 % be enclosed by the ellipse); default is 0.9.
matthiasm@8 16 % 'num-pts' - if the value supplied is n, then (n + 1)^2 points
matthiasm@8 17 % to be used to plot the ellipse; default is 20.
matthiasm@8 18 % 'label' - if non-empty, a string that will label the
matthiasm@8 19 % ellipsoid (default: [])
matthiasm@8 20 % 'plot-axes' - a 0/1 flag indicating if the ellipsoid's axes
matthiasm@8 21 % should be plotted (default: 1)
matthiasm@8 22 % 'plot-opts' - a cell vector of arguments to be handed to PLOT3
matthiasm@8 23 % to contol the appearance of the axes, e.g.,
matthiasm@8 24 % {'Color', 'g', 'LineWidth', 1}; the default is {}
matthiasm@8 25 % 'fill-color' - a color specifier; is this is not [], the
matthiasm@8 26 % covariance ellipse is filled with this color
matthiasm@8 27 % (default: [])
matthiasm@8 28 %
matthiasm@8 29 % Outputs:
matthiasm@8 30 % h - a vector of handles on the axis lines
matthiasm@8 31 %
matthiasm@8 32 % See also: PLOTCOV3
matthiasm@8 33
matthiasm@8 34 % Copyright (C) 2002 Mark A. Paskin
matthiasm@8 35 %
matthiasm@8 36 % This program is free software; you can redistribute it and/or modify
matthiasm@8 37 % it under the terms of the GNU General Public License as published by
matthiasm@8 38 % the Free Software Foundation; either version 2 of the License, or
matthiasm@8 39 % (at your option) any later version.
matthiasm@8 40 %
matthiasm@8 41 % This program is distributed in the hope that it will be useful, but
matthiasm@8 42 % WITHOUT ANY WARRANTY; without even the implied warranty of
matthiasm@8 43 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
matthiasm@8 44 % General Public License for more details.
matthiasm@8 45 %
matthiasm@8 46 % You should have received a copy of the GNU General Public License
matthiasm@8 47 % along with this program; if not, write to the Free Software
matthiasm@8 48 % Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
matthiasm@8 49 % USA.
matthiasm@8 50 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
matthiasm@8 51
matthiasm@8 52 function [h, s] = plotcov2New(mu, Sigma, varargin)
matthiasm@8 53
matthiasm@8 54 h = [];
matthiasm@8 55 s = [];
matthiasm@8 56
matthiasm@8 57 if size(Sigma) ~= [2 2], error('Sigma must be a 2 by 2 matrix'); end
matthiasm@8 58 if length(mu) ~= 2, error('mu must be a 2 by 1 vector'); end
matthiasm@8 59
matthiasm@8 60 Sigma = checkpsd(Sigma);
matthiasm@8 61
matthiasm@8 62 [p, ...
matthiasm@8 63 n, ...
matthiasm@8 64 label, ...
matthiasm@8 65 plot_axes, ...
matthiasm@8 66 plot_opts, ...
matthiasm@8 67 fill_color] = process_options(varargin, 'conf', 0.9, ...
matthiasm@8 68 'num-pts', 20, ...
matthiasm@8 69 'label', [], ...
matthiasm@8 70 'plot-axes', 1, ...
matthiasm@8 71 'plot-opts', {}, ...
matthiasm@8 72 'fill-color', []);
matthiasm@8 73 holding = ishold;
matthiasm@8 74 % Compute the Mahalanobis radius of the ellipsoid that encloses
matthiasm@8 75 % the desired probability mass.
matthiasm@8 76 k = conf2mahal(p, 2);
matthiasm@8 77 % Scale the covariance matrix so the confidence region has unit
matthiasm@8 78 % Mahalanobis distance.
matthiasm@8 79 Sigma = Sigma * k;
matthiasm@8 80 % The axes of the covariance ellipse are given by the eigenvectors of
matthiasm@8 81 % the covariance matrix. Their lengths (for the ellipse with unit
matthiasm@8 82 % Mahalanobis radius) are given by the square roots of the
matthiasm@8 83 % corresponding eigenvalues.
matthiasm@8 84 [V, D] = eig(full(Sigma));
matthiasm@8 85 V = real(V);
matthiasm@8 86 D = real(D);
matthiasm@8 87 D = abs(D);
matthiasm@8 88
matthiasm@8 89 % Compute the points on the boundary of the ellipsoid.
matthiasm@8 90 t = linspace(0, 2*pi, n);
matthiasm@8 91 u = [cos(t(:))'; sin(t(:))'];
matthiasm@8 92 w = (V * sqrt(D)) * u;
matthiasm@8 93 z = repmat(mu(:), [1 n]) + w;
matthiasm@8 94 h = [h; plot(z(1, :), z(2, :), plot_opts{:})];
matthiasm@8 95 if (~isempty(fill_color))
matthiasm@8 96 s = patch(z(1, :), z(2, :), fill_color);
matthiasm@8 97 end
matthiasm@8 98
matthiasm@8 99 % Plot the axes.
matthiasm@8 100 if (plot_axes)
matthiasm@8 101 hold on;
matthiasm@8 102 L = sqrt(diag(D));
matthiasm@8 103 h = plot([mu(1); mu(1) + L(1) * V(1, 1)], ...
matthiasm@8 104 [mu(2); mu(2) + L(1) * V(2, 1)], plot_opts{:});
matthiasm@8 105 h = [h; plot([mu(1); mu(1) + L(2) * V(1, 2)], ...
matthiasm@8 106 [mu(2); mu(2) + L(2) * V(2, 2)], plot_opts{:})];
matthiasm@8 107 end
matthiasm@8 108
matthiasm@8 109
matthiasm@8 110 if (~isempty(label))
matthiasm@8 111 th = text(mu(1), mu(2), label);
matthiasm@8 112 set(th, 'FontSize', 18);
matthiasm@8 113 set(th, 'FontName', 'Times');
matthiasm@8 114 set(th, 'FontWeight', 'bold');
matthiasm@8 115 set(th, 'FontAngle', 'italic');
matthiasm@8 116 set(th, 'HorizontalAlignment', 'center');
matthiasm@8 117 end
matthiasm@8 118
matthiasm@8 119 if (~holding & plot_axes) hold off; end