Greets,
I'm using nlminb() to estimate the parameters of a multivariate normal random
sample with missing values and ran into an unexpected result from my call to
dmvnorm() within the likelihood function. Particular details are provided
below. It appears that dmvnorm() makes a call to log(eigen(sigma)). Whereas
eigen(sigma) is returning a negative number, I understand log()'s complaint.
However, it is a mystery to me why this data set should produce such a result.
Any suggestions?
Best Regards,
Steven
> complete
[,1] [,2]
[1,] 0.84761637 3.994261
[2,] 0.91487059 4.952595
[3,] 0.84527267 4.521837
[4,] 2.53821358 8.374880
[5,] 1.16646209 6.255022
[6,] 0.94706527 4.169510
[7,] 0.48813564 3.349230
[8,] 3.71828469 9.441518
[9,] 0.08953357 1.651497
[10,] 0.68530515 5.498403
[11,] 1.52771645 8.484671
[12,] 1.55710697 5.231272
[13,] 1.89091603 4.152658
[14,] 1.08483541 5.401544
[15,] 0.58125385 5.340141
[16,] 0.24473250 2.965046
[17,] 1.59954401 8.095561
[18,] 1.57656436 5.335744
[19,] 2.73976992 8.572871
[20,] 0.87720252 6.067468
[21,] 1.18403087 3.526790
[22,] -1.03145244 1.776478
[23,] 2.88197343 7.720838
[24,] 0.60705218 4.406073
[25,] 0.58083464 3.374075
[26,] 0.87913427 5.247637
[27,] 1.10832692 3.534508
[28,] 2.92698371 8.682130
[29,] 4.04115277 11.827360
[30,] -0.57913297 1.476586
[31,] 0.84804365 7.009075
[32,] 0.79497940 3.671164
[33,] 1.58837762 5.535409
[34,] 0.63412821 3.932767
[35,] 3.14032433 9.271014
[36,] -0.18183869 1.666647
[37,] 0.57535770 6.881830
[38,] 3.21417723 10.901636
[39,] 0.29207932 4.120408
[40,] 0.65938218 5.209301
> u
[1] 1.267198 5.475045
> sigma
[,1] [,2]
[1,] 0.6461647 2.228951
[2,] 2.2289513 5.697834
> dmvnorm(x=complete,mean=u,sigma=sigma)
[1] NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
[30] NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
Warning message:
In log(eigen(sigma, symmetric = TRUE, only.values = TRUE)$values) :
NaNs produced
______________________________________________
R-help@r-project.org mailing list
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.
