[R] eigen() may fail for some symmetric matrices, affects mvrnorm()

From: Jouni Kerman <jpk2005_at_columbia.edu>
Date: Sun 01 May 2005 - 12:19:30 EST

Hi all,

Recently our statistics students noticed that their Gibbs samplers were crashing due to some NaNs in some parameters. The NaNs came from mvrnorm (Ripley & Venables' MASS package multivariate normal sampling function) and with some more investigation it turned out that they were generated by function eigen, the eigenvalue computing function. The problem did not seem to happen when the option "symmetric" was set to FALSE; at first I thought that only matrices that were not completely symmetric (such as ones with small errors of the magnitude 1e-17 or so) were affected but then we encountered matrices that were perfectly symmetric but eigen still generated NaNs when symmetric=TRUE.

I wrote a trivial patch to mvrnorm to detect this problem and to recompute the eigenvalues with symmetric=TRUE and EISPACK=TRUE. This seems to work. symmetric=FALSE also seems to work, but is somewhat slower.

If you use mvrnorm then you may want to try this. Details and the patch are on my web page, along with one example of a problematic symmetric matrix (to be loaded with the 'load' function). The NaNs should appear in the eigenvectors. This particular matrix may not produce an error on all platforms.

      http://www.stat.columbia.edu/~kerman/Software/rprog.html

Thanks to Radford Neal who suggested on our research web log (link below) that setting EISPACK=TRUE may help.

     http://www.stat.columbia.edu/~cook/movabletype/mlm/

Jouni Kerman
Department of Statistics
Columbia University
New York



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sun May 01 12:24:58 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:31:31 EST