# [R] mnormt and mvtnorm

Shh shh102004 at yahoo.com
Sun Feb 15 19:56:29 CET 2009

```Hi R user,

I came across a question when I played around with mnormt package (mvtnorm as well), could any of you let me know where I made the mistake?  Thanks!

For instance, with a multivariate normal distribution with known mean and variance (to be simple, I tried the i.i.d case first)

u1 <- u2 <- u3 <- 1
test <- function(x) {2*x*dnorm(x,u1)* pmnorm(rep(x,1),c(u2),diag(1)) }
value1 <- integrate(test , lower=-Inf, upper=Inf)\$value
testA <- function(x) {2*x*dnorm(x,u1,1)*pnorm(x,u2)}
value1A <- integrate(testA , lower=-Inf, upper=Inf)\$value
value1
# 1.564190
value1A
# 1.564190

# when pmnorm handle the univariate case, it is the same as pnorm. However, when I extended to multivariate, I got a confusing result:

test2 <- function(x) {x*dnorm(x,u1)*pmnorm(rep(x,2),c(u2,u3),diag(2)) }
value2 <- integrate(test2 , lower=-Inf, upper=Inf)\$value
test2A <- function(x) {x*dnorm(x,u1)*pnorm(x,u2)*pnorm(x,u3)}
value2A <- integrate(test2A , lower=-Inf, upper=Inf)\$value
value2
#  0.6997612
value2A
#  0.6154281

# I assume value2=value2a since they are independent, am I right? Also when I used some acutual values, I got:

pmnorm(rep(0,2),c(u2,u3),diag(2))
#  0.02517149
pnorm(0,u2)*pnorm(0,u3)
#  0.02517149

test3  <- function(x) {x*dnorm(x,u1)*dmnorm(rep(x,2),c(u2,u3),diag(2)) }
value3 <- integrate(test3 , lower=-Inf, upper=Inf)\$value
test3A <- function(x) {x*dnorm(x,u1)*dnorm(x,u2)*dnorm(x,u3)}
value3A <- integrate(test3A , lower=-Inf, upper=Inf)\$value
value3
# 0
value3A
#   0.09188815

# I got they same with an actual value:

dmnorm(rep(0,2),c(u2,u3),diag(2))
# 0.05854983
dnorm(0,u2)*dnorm(0,u3)
#  0.05854983

Similar thing occurred in mvtnorm I got three different results with mvtnorm, mnormt and pnorm.

```