dmnorm {mnormt} | R Documentation |
The probability density function, the distribution function and random number generation for the multivariate normal (Gaussian) probability distribution
dmnorm(x, mean = rep(0, d), varcov, log = FALSE) pmnorm(x, mean = rep(0, length(x)), varcov, ...) rmnorm(n = 1, mean = rep(0, d), varcov) sadmvn(lower, upper, mean, varcov, maxpts = 2000 * d, abseps = 1e-06, releps = 0)
x |
for dmnorm , this is either a vector of length d or
a matrix with d columns, where d=ncol(varcov) ,
giving the coordinates of the point(s) where the density must
be evaluated;
for pmnorm , only a vector of length d is allowed,
and d cannot exceed 20 |
mean |
a numeric vector representing the expected value of the
distribution; it must be of length d , as defined above |
varcov |
a positive definite matrix representing the
variance-covariance matrix of the distribution;
a vector of length 1 is also allowed
(in this case, d=1 is set) |
log |
a logical value; if TRUE ,
the logarithm of the density is computed
|
... |
parameters passed to sadmvn ,
among maxpts , abseps , releps |
n |
the number of random numbers to be generated |
lower |
a numeric vector of lower integration limits of
the density function; must be of maximal length 20;
+Inf and -Inf entries are allowed |
upper |
a numeric vector of upper integration limits
of the density function; must be of maximal length 20;
+Inf and -Inf entries are allowed |
maxpts |
the maximum number of function evaluations
(default value: 2000*d ) |
abseps |
absolute error tolerance (default value: 1e-6 ) |
releps |
relative error tolerance (default value: 0 ) |
The function pmnorm
works by making a suitable call to
sadmvn
if d>2
, or to biv.nt.prob
if d=2
.
Function sadmvn
is an interface to a Fortran-77 routine with
the same name written by Alan Genz, and available from his web page;
this makes uses of some auxiliary functions whose authors are
documented in the Fortran code. The routine uses an adaptive
integration method.
dmnorm
returns a vector of density values (possibly log-transformed);
pmnorm
and sadmvn
return a single probability with
attributes giving details on the achieved accuracy;
rmnorm
returns a matrix of n
rows of random vectors
The attributes error
and status
of the probability
returned by pmnorm
and sadmvn
indicate whether the function
had a normal termination, achieving the required accuracy. If
this is not the case, re-run the function with an higher value of
maxpts
Fortran code of SADMVN
and most auxiliary functions by Alan Genz,
some additional auxiliary functions by people referred to within his
program. Porting to R and additional R code by Adelchi Azzalini
Genz, A. (1992). Numerical Computation of Multivariate Normal Probabilities. J. Computational and Graphical Statist., 1, 141-149.
Genz, A. (1993). Comparison of methods for the computation of multivariate normal probabilities. Computing Science and Statistics, 25, 400-405.
Genz, A.: Fortran code available at http://www.math.wsu.edu/math/faculty/genz/software/mvn.f
x <- seq(-2,4,length=21) y <- 2*x+10 z <- x+cos(y) mu <- c(1,12,2) Sigma <- matrix(c(1,2,0,2,5,0.5,0,0.5,3), 3, 3) f <- dmnorm(cbind(x,y,z), mu, Sigma) p1 <- pmnorm(c(2,11,3), mu, Sigma) p2 <- pmnorm(c(2,11,3), mu, Sigma, maxpts=10000, abseps=1e-10) x <- rmnorm(10, mu, Sigma) p <- sadmvn(lower=c(2,11,3), upper=rep(Inf,3), mu, Sigma) # upper tail # p0 <- pmnorm(c(2,11), mu[1:2], Sigma[1:2,1:2]) p1 <- biv.nt.prob(0, lower=rep(-Inf,2), upper=c(2, 11), mu[1:2], Sigma[1:2,1:2]) p2 <- sadmvn(lower=rep(-Inf,2), upper=c(2, 11), mu[1:2], Sigma[1:2,1:2]) c(p0, p1, p2, p0-p1, p0-p2) # p1 <- pnorm(0, 1, 3) p2 <- pmnorm(0, 1, 3^2)