On Oct 18, 2013, at 8:26 AM, Steven LeBlanc wrote:
> On Oct 17, 2013, at 11:37 PM, David Winsemius <dwinsem...@comcast.net> wrote:
>
>>
>> On Oct 17, 2013, at 9:11 PM, Steven LeBlanc wrote:
>>
>>> 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()
>>
>> There are at least 5 different version of dmvnorm. None of them are in the
>> default packages.
>>
>>> within the likelihood function. Particular details are provided below.
>>
>> Complete? Except for the name of the package that has `dmvnorm`.
>>
>
> Package: ‘mvtnorm’ version 0.9-9992
> complete was the name of the data set.
I was clear that "complete" was the name of the dataset:
library(mvtnorm)
# First five rows of your complete:
complete <-
structure(c(0.84761637, 0.91487059, 0.84527267, 2.53821358, 1.16646209,
3.994261, 4.952595, 4.521837, 8.37488, 6.255022), .Dim = c(5L,
2L))
dmvnorm( complete, mean=c(1.267198, 5.475045),
sigma= matrix( c( 0.6461647, 2.2289513, 2.228951,
5.697834), 2) )
Error in dmvnorm(complete, mean = c(1.267198, 5.475045), sigma =
matrix(c(0.6461647, :
sigma must be a symmetric matrix
So trimming the covariance elements to be exactly equal:
> matrix( c( 0.6461647, 2.2289513, 2.228951, 5.697834), 2)
[,1] [,2]
[1,] 0.6461647 2.228951
[2,] 2.2289513 5.697834
> dmvnorm( complete, mean=c(1.267198, 5.475045),
sigma= matrix( c( 0.6461647, 2.228951, 2.228951,
5.697834), 2) )
[1] NaN NaN NaN NaN NaN
Warning message:
In log(eigen(sigma, symmetric = TRUE, only.values = TRUE)$values) :
NaNs produced
> dmvnorm(x=c(0,0))
[1] 0.1591549
> dmvnorm( complete, mean=c(1.267198, 5.475045) )
[1] 0.048690952 0.130494869 0.092440480 0.001059309 0.116818598
> eigen(sigma, symmetric = TRUE, only.values = TRUE)$values
[1] 6.5406882 -0.1966895
So the specified variance covariance matrix is not invertible. This can happen
if you use sample statistics:
http://stats.stackexchange.com/questions/49826/what-to-do-when-sample-covariance-matrix-is-not-invertible
I'm not sure what mixtools::dvnorm is doing that avoids the problem that
mvtnorm::dvnorm is identifying. Perhaps a pseudo-inverse if being constructed
and use as a substitute for sigma.
>
> Best Regards,
> Steven
David Winsemius
Alameda, CA, USA
______________________________________________
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.
