pmvnorm
Multivariate Normal Distribution
Computes the distribution function of the multivariate normal distribution for arbitrary limits and correlation matrices.
- Keywords
- distribution
Usage
pmvnorm(lower=-Inf, upper=Inf, mean=rep(0, length(lower)),
corr=NULL, sigma=NULL, algorithm = GenzBretz(), ...)
Arguments
- lower
the vector of lower limits of length n.
- upper
the vector of upper limits of length n.
- mean
the mean vector of length n.
- corr
the correlation matrix of dimension n.
- sigma
the covariance matrix of dimension n less than 1000. Either
corr
orsigma
can be specified. Ifsigma
is given, the problem is standardized. If neithercorr
norsigma
is given, the identity matrix is used forsigma
.- algorithm
an object of class
GenzBretz
,Miwa
orTVPACK
specifying both the algorithm to be used as well as the associated hyper parameters.- ...
additional parameters (currently given to
GenzBretz
for backward compatibility issues).
Details
This program involves the computation of multivariate normal probabilities with arbitrary correlation matrices. It involves both the computation of singular and nonsingular probabilities. The implemented methodology is described in Genz (1992, 1993) (for algorithm GenzBretz), in Miwa et al. (2003) for algorithm Miwa (useful up to dimension 20) and Genz (2004) for the TVPACK algorithm (which covers 2- and 3-dimensional problems for semi-infinite integration regions).
Note the default algorithm GenzBretz is randomized and hence slightly depends on
.Random.seed
and that both -Inf
and +Inf
may
be specified in lower
and upper
. For more details see
pmvt
.
The multivariate normal
case is treated as a special case of pmvt
with df=0
and
univariate problems are passed to pnorm
.
The multivariate normal density and random deviates are available using
dmvnorm
and rmvnorm
.
Value
The evaluated distribution function is returned with attributes
estimated absolute error and
status messages.
References
Genz, A. (1992). Numerical computation of multivariate normal probabilities. Journal of Computational and Graphical Statistics, 1, 141--150.
Genz, A. (1993). Comparison of methods for the computation of multivariate normal probabilities. Computing Science and Statistics, 25, 400--405.
Genz, A. (2004), Numerical computation of rectangular bivariate and trivariate normal and t-probabilities, Statistics and Computing, 14, 251--260.
Genz, A. and Bretz, F. (2009), Computation of Multivariate Normal and t Probabilities. Lecture Notes in Statistics, Vol. 195. Springer-Verlag, Heidelberg.
Miwa, A., Hayter J. and Kuriki, S. (2003). The evaluation of general non-centred orthant probabilities. Journal of the Royal Statistical Society, Ser. B, 65, 223--234.
See Also
Examples
# NOT RUN {
n <- 5
mean <- rep(0, 5)
lower <- rep(-1, 5)
upper <- rep(3, 5)
corr <- diag(5)
corr[lower.tri(corr)] <- 0.5
corr[upper.tri(corr)] <- 0.5
prob <- pmvnorm(lower, upper, mean, corr)
print(prob)
stopifnot(pmvnorm(lower=-Inf, upper=3, mean=0, sigma=1) == pnorm(3))
a <- pmvnorm(lower=-Inf,upper=c(.3,.5),mean=c(2,4),diag(2))
stopifnot(round(a,16) == round(prod(pnorm(c(.3,.5),c(2,4))),16))
a <- pmvnorm(lower=-Inf,upper=c(.3,.5,1),mean=c(2,4,1),diag(3))
stopifnot(round(a,16) == round(prod(pnorm(c(.3,.5,1),c(2,4,1))),16))
# Example from R News paper (original by Genz, 1992):
m <- 3
sigma <- diag(3)
sigma[2,1] <- 3/5
sigma[3,1] <- 1/3
sigma[3,2] <- 11/15
pmvnorm(lower=rep(-Inf, m), upper=c(1,4,2), mean=rep(0, m), corr=sigma)
# Correlation and Covariance
a <- pmvnorm(lower=-Inf, upper=c(2,2), sigma = diag(2)*2)
b <- pmvnorm(lower=-Inf, upper=c(2,2)/sqrt(2), corr=diag(2))
stopifnot(all.equal(round(a,5) , round(b, 5)))
# }