Re: Principal Components Analysis

[Rafael Laboissiere]

>>>>> "JWE" == John W Eaton <address@hidden> writes:

    JWE> Perhaps the following code will be of some help.  It will take a
    JWE> 2x2 matrix and generate a set of points for plotting an ellipse,
    JWE> its major and minor axes, and the bounding box that just encloses
    JWE> the ellipse.

    JWE> [...]

Hi John,

When I try to use your function to generate an ellipse of dispersion
correponding to the covariance matrix of a randomly generated set of
points, it does not seem to work properly.  (I am not sure whether your
function was intended to do that.)  Here is my code:

##################################################################
t = [rand()+randn(100,1), rand()+randn(100,1)];
level = 3;
npts = 181;
[x, y] = ellipse (cov(t), level, npts);

gset size ratio 1;
gset nokey

plot (t(:,1), t(:,2), "*", x+mean(t(:,1)), y+mean(t(:,2)), "-b");
##################################################################

In order to have a correct plot, I had to modify your ellipse function as
bellow :

diff -u /home/rafael/octave/ellipse.m /home/rafael/octave/ellipse-orig.m
--- ellipse-orig.m      Thu Jan 28 14:35:47 1999
+++ ellipse.m   Thu Jan 28 14:28:45 1999
@@ -12,7 +12,7 @@
 
   if (nargin > 1)
 
-    [v, l] = eig (amat / level);
+    [v, l] = eig (amat);
 
     dl = diag(l);
     if (any (imag (dl)) || any (dl <= 0))
@@ -21,8 +21,8 @@
 
     ## Generate contour data.
 
-    a = l / sqrt (l(1,1));
-    b = l / sqrt (l(2,2));
+    a = level * sqrt (l(1,1));
+    b = level * sqrt (l(2,2));
 
     t = linspace (0, 2*pi, n)';
  

--
Rafael Laboissiere
Institut de la Communication Parlee | Email: address@hidden
UPRESS A CNRS 5009 / INPG           | Voice: +33 4.76.57.48.49
46, av. Felix Viallet               |   Fax: +33 4.76.57.47.10
F-38031 Grenoble CEDEX 1 France     |   URL: http://www.icp.inpg.fr/~rafael