|
wolffd@0
|
1 function h=plotgauss2d(mu, Sigma)
|
|
wolffd@0
|
2 % PLOTGAUSS2D Plot a 2D Gaussian as an ellipse with optional cross hairs
|
|
wolffd@0
|
3 % h=plotgauss2(mu, Sigma)
|
|
wolffd@0
|
4 %
|
|
wolffd@0
|
5
|
|
wolffd@0
|
6 h = plotcov2(mu, Sigma);
|
|
wolffd@0
|
7 return;
|
|
wolffd@0
|
8
|
|
wolffd@0
|
9 %%%%%%%%%%%%%%%%%%%%%%%%
|
|
wolffd@0
|
10
|
|
wolffd@0
|
11 % PLOTCOV2 - Plots a covariance ellipse with major and minor axes
|
|
wolffd@0
|
12 % for a bivariate Gaussian distribution.
|
|
wolffd@0
|
13 %
|
|
wolffd@0
|
14 % Usage:
|
|
wolffd@0
|
15 % h = plotcov2(mu, Sigma[, OPTIONS]);
|
|
wolffd@0
|
16 %
|
|
wolffd@0
|
17 % Inputs:
|
|
wolffd@0
|
18 % mu - a 2 x 1 vector giving the mean of the distribution.
|
|
wolffd@0
|
19 % Sigma - a 2 x 2 symmetric positive semi-definite matrix giving
|
|
wolffd@0
|
20 % the covariance of the distribution (or the zero matrix).
|
|
wolffd@0
|
21 %
|
|
wolffd@0
|
22 % Options:
|
|
wolffd@0
|
23 % 'conf' - a scalar between 0 and 1 giving the confidence
|
|
wolffd@0
|
24 % interval (i.e., the fraction of probability mass to
|
|
wolffd@0
|
25 % be enclosed by the ellipse); default is 0.9.
|
|
wolffd@0
|
26 % 'num-pts' - the number of points to be used to plot the
|
|
wolffd@0
|
27 % ellipse; default is 100.
|
|
wolffd@0
|
28 %
|
|
wolffd@0
|
29 % This function also accepts options for PLOT.
|
|
wolffd@0
|
30 %
|
|
wolffd@0
|
31 % Outputs:
|
|
wolffd@0
|
32 % h - a vector of figure handles to the ellipse boundary and
|
|
wolffd@0
|
33 % its major and minor axes
|
|
wolffd@0
|
34 %
|
|
wolffd@0
|
35 % See also: PLOTCOV3
|
|
wolffd@0
|
36
|
|
wolffd@0
|
37 % Copyright (C) 2002 Mark A. Paskin
|
|
wolffd@0
|
38
|
|
wolffd@0
|
39 function h = plotcov2(mu, Sigma, varargin)
|
|
wolffd@0
|
40
|
|
wolffd@0
|
41 if size(Sigma) ~= [2 2], error('Sigma must be a 2 by 2 matrix'); end
|
|
wolffd@0
|
42 if length(mu) ~= 2, error('mu must be a 2 by 1 vector'); end
|
|
wolffd@0
|
43
|
|
wolffd@0
|
44 [p, ...
|
|
wolffd@0
|
45 n, ...
|
|
wolffd@0
|
46 plot_opts] = process_options(varargin, 'conf', 0.9, ...
|
|
wolffd@0
|
47 'num-pts', 100);
|
|
wolffd@0
|
48 h = [];
|
|
wolffd@0
|
49 holding = ishold;
|
|
wolffd@0
|
50 if (Sigma == zeros(2, 2))
|
|
wolffd@0
|
51 z = mu;
|
|
wolffd@0
|
52 else
|
|
wolffd@0
|
53 % Compute the Mahalanobis radius of the ellipsoid that encloses
|
|
wolffd@0
|
54 % the desired probability mass.
|
|
wolffd@0
|
55 k = conf2mahal(p, 2);
|
|
wolffd@0
|
56 % The major and minor axes of the covariance ellipse are given by
|
|
wolffd@0
|
57 % the eigenvectors of the covariance matrix. Their lengths (for
|
|
wolffd@0
|
58 % the ellipse with unit Mahalanobis radius) are given by the
|
|
wolffd@0
|
59 % square roots of the corresponding eigenvalues.
|
|
wolffd@0
|
60 if (issparse(Sigma))
|
|
wolffd@0
|
61 [V, D] = eigs(Sigma);
|
|
wolffd@0
|
62 else
|
|
wolffd@0
|
63 [V, D] = eig(Sigma);
|
|
wolffd@0
|
64 end
|
|
wolffd@0
|
65 % Compute the points on the surface of the ellipse.
|
|
wolffd@0
|
66 t = linspace(0, 2*pi, n);
|
|
wolffd@0
|
67 u = [cos(t); sin(t)];
|
|
wolffd@0
|
68 w = (k * V * sqrt(D)) * u;
|
|
wolffd@0
|
69 z = repmat(mu, [1 n]) + w;
|
|
wolffd@0
|
70 % Plot the major and minor axes.
|
|
wolffd@0
|
71 L = k * sqrt(diag(D));
|
|
wolffd@0
|
72 h = plot([mu(1); mu(1) + L(1) * V(1, 1)], ...
|
|
wolffd@0
|
73 [mu(2); mu(2) + L(1) * V(2, 1)], plot_opts{:});
|
|
wolffd@0
|
74 hold on;
|
|
wolffd@0
|
75 h = [h; plot([mu(1); mu(1) + L(2) * V(1, 2)], ...
|
|
wolffd@0
|
76 [mu(2); mu(2) + L(2) * V(2, 2)], plot_opts{:})];
|
|
wolffd@0
|
77 end
|
|
wolffd@0
|
78
|
|
wolffd@0
|
79 h = [h; plot(z(1, :), z(2, :), plot_opts{:})];
|
|
wolffd@0
|
80 if (~holding) hold off; end
|
|
wolffd@0
|
81
|
|
wolffd@0
|
82 %%%%%%%%%%%%
|
|
wolffd@0
|
83
|
|
wolffd@0
|
84 % CONF2MAHAL - Translates a confidence interval to a Mahalanobis
|
|
wolffd@0
|
85 % distance. Consider a multivariate Gaussian
|
|
wolffd@0
|
86 % distribution of the form
|
|
wolffd@0
|
87 %
|
|
wolffd@0
|
88 % p(x) = 1/sqrt((2 * pi)^d * det(C)) * exp((-1/2) * MD(x, m, inv(C)))
|
|
wolffd@0
|
89 %
|
|
wolffd@0
|
90 % where MD(x, m, P) is the Mahalanobis distance from x
|
|
wolffd@0
|
91 % to m under P:
|
|
wolffd@0
|
92 %
|
|
wolffd@0
|
93 % MD(x, m, P) = (x - m) * P * (x - m)'
|
|
wolffd@0
|
94 %
|
|
wolffd@0
|
95 % A particular Mahalanobis distance k identifies an
|
|
wolffd@0
|
96 % ellipsoid centered at the mean of the distribution.
|
|
wolffd@0
|
97 % The confidence interval associated with this ellipsoid
|
|
wolffd@0
|
98 % is the probability mass enclosed by it. Similarly,
|
|
wolffd@0
|
99 % a particular confidence interval uniquely determines
|
|
wolffd@0
|
100 % an ellipsoid with a fixed Mahalanobis distance.
|
|
wolffd@0
|
101 %
|
|
wolffd@0
|
102 % If X is an d dimensional Gaussian-distributed vector,
|
|
wolffd@0
|
103 % then the Mahalanobis distance of X is distributed
|
|
wolffd@0
|
104 % according to the Chi-squared distribution with d
|
|
wolffd@0
|
105 % degrees of freedom. Thus, the Mahalanobis distance is
|
|
wolffd@0
|
106 % determined by evaluating the inverse cumulative
|
|
wolffd@0
|
107 % distribution function of the chi squared distribution
|
|
wolffd@0
|
108 % up to the confidence value.
|
|
wolffd@0
|
109 %
|
|
wolffd@0
|
110 % Usage:
|
|
wolffd@0
|
111 %
|
|
wolffd@0
|
112 % m = conf2mahal(c, d);
|
|
wolffd@0
|
113 %
|
|
wolffd@0
|
114 % Inputs:
|
|
wolffd@0
|
115 %
|
|
wolffd@0
|
116 % c - the confidence interval
|
|
wolffd@0
|
117 % d - the number of dimensions of the Gaussian distribution
|
|
wolffd@0
|
118 %
|
|
wolffd@0
|
119 % Outputs:
|
|
wolffd@0
|
120 %
|
|
wolffd@0
|
121 % m - the Mahalanobis radius of the ellipsoid enclosing the
|
|
wolffd@0
|
122 % fraction c of the distribution's probability mass
|
|
wolffd@0
|
123 %
|
|
wolffd@0
|
124 % See also: MAHAL2CONF
|
|
wolffd@0
|
125
|
|
wolffd@0
|
126 % Copyright (C) 2002 Mark A. Paskin
|
|
wolffd@0
|
127
|
|
wolffd@0
|
128 function m = conf2mahal(c, d)
|
|
wolffd@0
|
129
|
|
wolffd@0
|
130 m = chi2inv(c, d); % matlab stats toolbox
|
