[R] [FORGED] Problem running mvrnorm

[peter dalgaard]

I have a sense of a tiny bell ringing faintly from the distant past... There may be an issue with some versions of BLAS/LAPACK on some systems showing up in eigendecompositions. 

Checking... hmm, there have been several issues over time but I don't see anything that would be directly related to the present issue. There was a really stupid bug in R's reference BLAS, but that was found and fixed two years ago, and besides, it involved complex eigenvalues which I don't see creeping into mvrnorm. Also PR#15211 located an issue that could cause a hang (not segfault) which ended in a WONTFIX situation because it really needed to be fixed in LAPACK.

Anyways, some details about versions, hardware (you may be using a library optimized for the wrong CPU), would be needed. Even better if you can set up to run under Valgrind or a debugger to pinpoint the cause of the trouble. Also, check that you aren't simply running out of memory.

Pragmatically speaking, couldn't you get away with using a Choleski decomposition of Sigma? 

-pd 


> On 18 Nov 2015, at 21:56 , Rolf Turner <r.turner at auckland.ac.nz> wrote:
> 
> 
> I cobbled together a 190 x 190 positive definite matrix Sigma and ran your example.  I got a result instantaneously, with no error message. (I'm running Linux; an ancient Fedora 17 system.)
> 
> So the problem is peculiar to your particular Sigma.
> 
> As the error message tells you, the problem comes from doing an eigendecomposition of Sigma.  So start your investigation by doing
> 
>    E <- eigen(Sigma,symmetric=TRUE)
> 
> Presumably that will lead to the same error.  How to get around this error is beyond the scope of my capabilities.
> 
> You *might* get somewhere by using the singular value decomposition
> (equivalent for a positive definite matrix) rather than the eigendecomposition.  I have the vague notion that the svd is more numerically robust than eigendecomposition.  However I might well have that wrong.
> 
> Doing anything in 190 dimensions is bound to be fraught with numeric
> peril.
> 
> cheers,
> 
> Rolf Turner
> 
> -- 
> Technical Editor ANZJS
> Department of Statistics
> University of Auckland
> Phone: +64-9-373-7599 ext. 88276
> 
> On 19/11/15 08:28, Gang Chen wrote:
>>  I’m running R 3.2.2 on a Linux server (Redhat 4.4.7-16), and having the
>> following problem.
>> 
>> It works fine with the following:
>> 
>> require('MASS’)
>> var(mvrnorm(n = 1000, rep(0, 2), Sigma=matrix(c(10,3,3,2),2,2)))
>> 
>> However, when running the following in a loop with simulated data (Sigma):
>> 
>> # Sigma defined somewhere else
>> mvrnorm(n=1000, rep(0, 190), Sigma)
>> 
>> I get this opaque message:
>> 
>>  *** caught illegal operation ***
>> address 0x7fe78f8693d2, cause 'illegal operand'
>> 
>> Traceback:
>>  1: eigen(Sigma, symmetric = TRUE)
>>  2: mvrnorm(n = nr, rep(0, NN), Sigma)
>> 
>> Possible actions:
>> 1: abort (with core dump, if enabled)
>> 2: normal R exit
>> 3: exit R without saving workspace
>> 4: exit R saving workspace
>> 
>> I tried to do core dump (option 1), but it didn’t go anywhere (hanging
>> there forever). I also ran the same code on a Mac, and there was no problem
>> at all. What is causing the problem on the Linux server? In case the
>> variance-covariance matrix ‘Sigma’ is needed, I can provide its definition
>> later.
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

-- 
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd.mes at cbs.dk  Priv: PDalgd at gmail.com