Karl Ove Hufthammer
2009-Jun-22 14:09 UTC
[R] The gradient of a multivariate normal density with respect to its parameters
Does anybody know of a function that implements the derivative (gradient) of the multivariate normal density with respect to the *parameters*? It?s easy enough to implement myself, but I?d like to avoid reinventing the wheel (with some bugs) if possible. Here?s a simple example of the result I?d like, using numerical differentiation: library(mvtnorm) library(numDeriv) f=function(pars, xx, yy) { mu=pars[1:2] sig1=pars[3] sig2=pars[4] rho=pars[5] sig=matrix(c(sig1^2,rho*sig1*sig2,rho*sig1*sig2,sig2^2),2) dmvnorm(cbind(x,y),mu,sig) } mu1=1 mu2=2 sig1=3 sig2=4 rho=.5 x=2 # or a x vector y=3 # or a y vector jacobian(f,c(mu1,mu2,sig1,sig2,rho),xx=x,yy=y) # (Can replace ?jacobian? with ?grad? if x and y have length 1.) -- Karl Ove Hufthammer
Ravi Varadhan
2009-Jun-22 14:48 UTC
[R] The gradient of a multivariate normal density with respect to its parameters
Hi Karl, I am not aware of any. May I ask what your purpose is? You don't really need this if you are going to use it in optimization, since most optimizers use a simple finite-difference approximation if you don't provide the gradient. Using the numerical approximation from "numDeriv" will be quite time-consuming in an optimization routine, since numDeriv uses a high-order Richardosn extrapolation to compute an accurate approximation of the gradient. Regardless of your purpose, there is a small bug in your function. You should change `dmvnorm(cbind(x,y),mu,sig)' to `dmvnorm(cbind(xx,yy),mu,sig)'. Also, taking the derivative of log-density might be better. Ravi. ---------------------------------------------------------------------------- ------- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: rvaradhan at jhmi.edu Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h tml ---------------------------------------------------------------------------- -------- -----Original Message----- From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On Behalf Of Karl Ove Hufthammer Sent: Monday, June 22, 2009 10:10 AM To: r-help at stat.math.ethz.ch Subject: [R] The gradient of a multivariate normal density with respect to its parameters Does anybody know of a function that implements the derivative (gradient) of the multivariate normal density with respect to the *parameters*? It's easy enough to implement myself, but I'd like to avoid reinventing the wheel (with some bugs) if possible. Here's a simple example of the result I'd like, using numerical differentiation: library(mvtnorm) library(numDeriv) f=function(pars, xx, yy) { mu=pars[1:2] sig1=pars[3] sig2=pars[4] rho=pars[5] sig=matrix(c(sig1^2,rho*sig1*sig2,rho*sig1*sig2,sig2^2),2) dmvnorm(cbind(x,y),mu,sig) } mu1=1 mu2=2 sig1=3 sig2=4 rho=.5 x=2 # or a x vector y=3 # or a y vector jacobian(f,c(mu1,mu2,sig1,sig2,rho),xx=x,yy=y) # (Can replace 'jacobian' with 'grad' if x and y have length 1.) -- Karl Ove Hufthammer ______________________________________________ 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.