From ben@eno.Princeton.EDU  Tue Jul  1 17:33:12 1997
From: ben@eno.Princeton.EDU (Ben Bolker)
Date: Tue, 1 Jul 1997 12:33:12 -0400 (EDT)
Subject: R-alpha: bug in eigen()
Message-ID: <Pine.SUN.3.95.970701122527.2221A-100000@paca>


  I don't have time to track this down right now -- if no-one does it in
the next few days before I finish the current piece of work I will try to
find it myself.  
  R gets the wrong answer for the eigensystem of a (particular) 4x4
matrix (SunOS 4.1.4).

Running these commands:

x _ matrix(
     c(-10,   0    ,0     ,0,
         0  ,-0.3  ,0     ,0,
         1  , 0.03 ,-0.594,0.0020,
         1.2, 0.036,0.0792,-0.0176),byrow=T,ncol=4)
eigen(x)


R gives:

$values
[1] -10.2425681  -0.5515844  -0.3023825   0.1849350
 
$vectors
              [,1]        [,2]         [,3]      [,4]
[1,] -0.9881475178  0.06506841 -0.008135724 0.1387962
[2,] -0.0007232329 -0.07331061  0.993467182 0.0874528
[3,]  0.1014706558  0.95415639  0.046029849 0.2777959
[4,]  0.1151853550 -0.28280247 -0.104105822 0.9465290
 
S-PLUS 3.1 gives:

$values:
[1] -10.00000000  -0.59427468  -0.30000000  -0.01732532
 
$vectors:
           [,1]        [,2]       [,3]         [,4] 
[1,]  1.0000000  0.00000000  0.0000000  0.000000000
[2,]  0.0000000  0.00000000  1.0000000  0.000000000
[3,] -0.1062897 -0.06249770  0.1009810 -0.003466636
[4,] -0.1193683  0.00858338 -0.1557992 -0.999560734


Note that the structure of the matrix actually gives us (-10.0) and (-0.3)
for free, so the S-PLUS answer looks right.

  Even worse:

eigen(x[3:4,3:4])

  gives (R)

$values
[1] -0.604684392 -0.006915608
 
$vectors
           [,1]      [,2]
[1,]  0.9910228 0.1336929
[2,] -0.1336929 0.9910228


 (S-PLUS)
$values:
[1] -0.59427468 -0.01732532
 
$vectors:
            [,1]         [,2] 
[1,] -0.06249770 -0.003466636
[2,]  0.00858338 -0.999560734

  so it's not just some funky pivoting problem.

  Or am I doing something *really* dumb?

Benjamin Bolker          Dep't of Ecology and Evolutionary Biology,
ben@eno.princeton.edu    Princeton University
tel: (609) 258-6886,     fax: (609) 258-1334

=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-devel-request@stat.math.ethz.ch
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-


R-alpha: bug in eigen()
Ben Bolker