mvTProbQuasiMonteCarlo: Multivariate t distribution function

View source: R/mvTQMC.R

mvTProbQuasiMonteCarloR Documentation

Multivariate t distribution function

Description

Estimate the multivariate t distribution function with quasi-Monte Carlo method.

Usage

mvTProbQuasiMonteCarlo(p, upperBound, cov, nu, genVec, ...)

Arguments

p

Number of samples used for quasi-Monte Carlo estimation. Must be a prime number.

upperBound

Vector of probabilities, i.e., the upper bound of the integral.

cov

Covariance matrix of the multivariate normal distribution. Must be positive semi-definite. WARNING: for performance in high-dimensions, no check is done to ensure positive-definiteness of the covariance matrix. It is the user responsibility to ensure that this property is verified.

nu

Degrees of freedom of the t distribution.

genVec

Generating vector for the quasi-Monte Carlo procedure. Can be computed using genVecQMC.

...

Additional arguments passed to Cpp routine.

Details

The function uses a quasi-Monte Carlo procedure based on randomly shifted lattice rules to estimate the distribution function a multivariate normal distribution as described in Genz and Bretz (2009) on page 50.

For compatibility reasons, the function handles the univariate case, which is passed on to pt.

Value

A named vector with components estimate estimate of the distribution function along error, 3 times the empirical Monte Carlo standard error over the nrep replications.

Author(s)

Raphael de Fondeville

Source

Adapted from the QSILATMVTV Matlab routine by A. Genz (2013)

References

Genz, A. and Bretz, F. (2009). Computations of Multivariate Normal and t Probabilities, volume 105. Springer: Dordrecht.

Examples


#Define locations
loc <- expand.grid(1:4, 1:4)
ref <- sample.int(16, 1)

#Define degrees of freedom
nu <- 3

#Compute variogram matrix
variogramMatrix <- ((sqrt((outer(loc[,1],loc[,1],"-"))^2 +
(outer(loc[,2],loc[,2],"-"))^2)) / 2)^(1.5)

#Define an upper bound
upperBound <- variogramMatrix[-ref,ref]

#Compute covariance matrix
cov <-  (variogramMatrix[-ref,ref]%*%t(matrix(1, (nrow(loc) - 1), 1)) +
t(variogramMatrix[ref,-ref]%*%t(matrix(1, (nrow(loc) - 1), 1))) -
variogramMatrix[-ref,-ref])

#Compute generating vector
p <- 499
latticeRule <- genVecQMC(p, (nrow(loc) - 1))

#Estimate the multivariate distribution function
mvTProbQuasiMonteCarlo(latticeRule$primeP, upperBound, cov, nu, latticeRule$genVec)

mvPot documentation built on April 4, 2025, 5:28 a.m.