R help - Jun 2011 - Multivariate normal density in C for R

Does anyone know of a package that uses C code to calculate a multivariate
normal density?

My goal is to find a faster way to calculate MVN densities and avoid R loops
or apply functions, such as when X and mu are N x K matrices, as opposed to
vectors, and in this particular case, speed really matters. I would like to
be able to use .C or .Call to pass X, mu, Sigma, and N to a C program and
have it return a vector of log densities to R.

I'm new to putting C in R, but am sure I'll figure it out. Thanks for
any
suggestions.

--
View this message in context:
http://r.789695.n4.nabble.com/Multivariate-normal-density-in-C-for-R-tp3624602p3624602.html
Sent from the R help mailing list archive at Nabble.com.

I use the following function which does not uses loops and seems to be 
pretty fast:

dmvnorm <- function (x, mu, Sigma, df, log = FALSE) {
     if (!is.matrix(x))
         x <- rbind(x)
     p <- nrow(Sigma)
     ed <- eigen(Sigma, symmetric = TRUE)
     ev <- ed$values
     if (!all(ev >= -1e-06 * abs(ev[1])))
         stop("'Sigma' is not positive definite")
     ss <- if (!is.matrix(mu)) {
         x - rep(mu, each = nrow(x))
     } else {
         x - mu
     }
     inv.Sigma <- ed$vectors %*% (t(ed$vectors)/ev)
     quad <- 0.5 * rowSums((ss %*% inv.Sigma) * ss)
     fact <- -0.5 * (p * log(2 * pi) + sum(log(ev)))
     if (log)
         as.vector(fact - quad)
     else
         as.vector(exp(fact - quad))
}


I hope it helps.

Best,
Dimitris


On 6/25/2011 3:58 PM, zerfetzen wrote:
> Does anyone know of a package that uses C code to calculate a multivariate
> normal density?
>
> My goal is to find a faster way to calculate MVN densities and avoid R
loops
> or apply functions, such as when X and mu are N x K matrices, as opposed to
> vectors, and in this particular case, speed really matters. I would like to
> be able to use .C or .Call to pass X, mu, Sigma, and N to a C program and
> have it return a vector of log densities to R.
>
> I'm new to putting C in R, but am sure I'll figure it out. Thanks
for any
> suggestions.
>
> --
> View this message in context:
http://r.789695.n4.nabble.com/Multivariate-normal-density-in-C-for-R-tp3624602p3624602.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at 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.
>
-- 
Dimitris Rizopoulos
Assistant Professor
Department of Biostatistics
Erasmus University Medical Center

Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands
Tel: +31/(0)10/7043478
Fax: +31/(0)10/7043014
Web: http://www.erasmusmc.nl/biostatistiek/

zerfetzen wrote:
> 
> Does anyone know of a package that uses C code to calculate a multivariate
> normal density?
> 
> My goal is to find a faster way to calculate MVN densities and avoid R
> loops or apply functions, such as when X and mu are N x K matrices, as
> opposed to vectors, and in this particular case, speed really matters. I
> would like to be able to use .C or .Call to pass X, mu, Sigma, and N to a
> C program and have it return a vector of log densities to R.
> 
> I'm new to putting C in R, but am sure I'll figure it out. Thanks
for any
> suggestions.
> 
Just a suggestion: package mvtnorm

Berend


--
View this message in context:
http://r.789695.n4.nabble.com/Multivariate-normal-density-in-C-for-R-tp3624602p3625954.html
Sent from the R help mailing list archive at Nabble.com.