Learn R Programming

tlrmvnmvt (version 1.1.2)

pmvn: Quasi-Monte Carlo method for multivariate normal probabilities

Description

Compute multivariate normal probabilities with the dense-matrix based Quasi-Monte Carlo method and the tile-low-rank-matrix based Quasi-Monte Carlo method.

Usage

pmvn(lower = -Inf, upper = Inf, mean = 0, sigma = NULL, 
uselog2 = FALSE, algorithm = GenzBretz(), ...)

Arguments

lower

lower integration limits, a numeric vector of length n

upper

upper integration limits, a numeric vector of length n

mean

the mean parameter, a numeric vector of length n

sigma

the covariance matrix of dimension n

uselog2

whether return the result as the logarithm to the base 2

algorithm

an object of class GenzBretz or TLRQMC defining the hyper parameters of this algorithm

...

additional parameters used to construct `sigma` when it is not given:

  • geom a matrix of dimension n-by-2, specifying n spatial locations in the 2D unit square

  • kernelType the name of the covariance kernel, a string. Currently, only the Matern covariance function, e.g., "matern", is supported. Not case-sensitive. It should be given when `sigma` is not given

  • para the parameter for the covariance kernel, a numeric vector. When `type` is "matern", the length of `para` should be 4, representing the scale, range, smoothness, and nugget parameters of the covariance function. It should be given when `sigma` is not given

Value

When `uselog2` is set FALSE, the function returns the estimated probability with one attribute of the estimation error. When `uselog2` is set TRUE, the function only returns the estimated log-probability to the base 2. This is useful when the estimated probability is smaller than the machine precision.

Details

When `algorithm` is of the class `GenzBretz`, the Quasi-Monte Carlo sampling described in Genz, A. (1992) is used. When `algorithm` is of the class `TLRQMC`, the Quasi-Monte Carlo sampling with the tile-low-rank representation of the covariance matrix, described in Cao et al. (2020), is used. When `sigma`, is given, `geom`, `kernelType`, and `para` are not used. Otherwise, a covariance matrix is created with the information from `geom`, `kernelType`, and `para`.

References

Genz, A. (1992), "Numerical computation of multivariate normal probabilities," Journal of Computational and Graphical Statistics, 1, 141-149. Cao, J., Genton, M. G., Keyes, D. E., & Turkiyyah, G. M. (2022), "tlrmvnmvt: Computing High-Dimensional Multivariate Normal and Student-t Probabilities with Low-Rank Methods in R," Journal of Statistical Software, 101.4, 1-25.

Examples

Run this code
# NOT RUN {
  n = 225
  set.seed(0)
  a = rep(-10, n)
  b = rnorm(n, 3, 2)
  m = 15
  epsl = 1e-4
  vec1 = 1 : m
  vec2 = rep(1, m)
  geom = cbind(kronecker(vec1, vec2), kronecker(vec2, vec1))
  geom = geom / m
  beta = 0.3
  idx = zorder(geom)
  geom = geom[idx, ]
  a = a[idx]
  b = b[idx]
  distM = as.matrix(dist(geom))
  covM = exp(-distM / beta)
  pmvn(lower = a, upper = b, mean = 2, sigma = covM, uselog2 = FALSE, 
       algorithm = GenzBretz(N = 521))
  pmvn(lower = a, upper = b, mean = 2, uselog2 = TRUE, geom = geom, 
       kernelType = "matern", para = c(1.0, 0.3, 0.5, 0.0))
  pmvn(lower = a, upper = b, mean = 2, sigma = covM, uselog2 = FALSE, 
       algorithm = TLRQMC(N = 521, m = m, epsl = epsl))
  pmvn(lower = a, upper = b, mean = 2, uselog2 = TRUE, geom = geom, 
       algorithm = TLRQMC(N = 521, m = m, epsl = epsl),
       kernelType = "matern", para = c(1.0, 0.3, 0.5, 0.0))
# }

Run the code above in your browser using DataLab