[R] Fwd: dmnorm not meant for 1024-dimensional data?
Adelchi Azzalini
aa at tango.stat.unipd.it
Wed Apr 25 21:49:59 CEST 2007
On Wed, Apr 25, 2007 at 01:25:09PM -0500, Daniel Elliott wrote:
> Hello. Thank you for your mnormt package. Below is an email I sent to
> R-help regarding the function dmnorm.
>
> Any help you can provide would be greatly appreciated.
>
Hello.
You have not provided information on what specifically you
have tried as for parameters and evaluation points, so I have
chosen a particularly simple case, namely the 1024-dimensional
density with independent standard marginals, using this code
d<- 1024
x <- 1
const<- 0.5*log(2*pi)
log.pdf <-0
for (i in 1:d) log.pdf <- log.pdf - x^2/2 -const
#
S <- diag(d)
mu<- rep(0,d)
X <- rep(x,d)
log.pdf2 <- dmnorm(X,mu,S, log=TRUE)
cat(log.pdf, log.pdf2, abs(log.pdf-log.pdf2),"\n")
and the outcome was this one:
-1453 -1453 2.137e-11
which seems quite decent to me. Obviously, you must not take exp()
of this log-density, as otherwise a 0 is indeed produced, but
the reason is not in the package, rather in the possibility of
representing that number using standard floating-point
numerical devices.
best regards,
Adelchi Azzalini
> Thank you.
>
> - dan elliott
>
> ---------- Forwarded message ----------
> From: Daniel Elliott <danelliottster at gmail.com>
> Date: Apr 25, 2007 1:21 PM
> Subject: dmnorm not meant for 1024-dimensional data?
> To: r-help at stat.math.ethz.ch
>
> Hello,
>
> I have some data generated by a simple mixture of Gaussians (more like
> K-means) and (as a test) am using dmnorm to calculate the probability
> of each data point coming from each Gaussian. However, I get only
> zero probabilities.
>
> This code works in low dimensions (tried 2 and 3 already). I have run
> into many implementations that do not work in high dimension, but I
> thought that I was safe with dmnorm because it has an option to
> compute the log of the probability.
>
> So, is dmnorm not intended to be used with data of such high dimensionality?
>
> Thank you,
>
> dan elliott
--
Adelchi Azzalini <azzalini at stat.unipd.it>
Dipart.Scienze Statistiche, Università di Padova, Italia
tel. +39 049 8274147, http://azzalini.stat.unipd.it/
More information about the R-help
mailing list