skewstudent {nvmix}

Functionalities for the skew-t distribution and copula

Description

Sampling and density evaluation for the multivariate skew-t distribution and copula.

Usage

rskewt(n, loc = rep(0, d), scale = diag(2), factor = NULL, gamma = rep(0, d),
       df = Inf, method = c("PRNG", "sobol", "ghalton"), skip = 0)
dskewt(x, loc = rep(0, d), scale = diag(2), gamma = rep(0, d), df,
       log = FALSE, scale.inv, ldet)

rskewtcopula(n, scale = diag(2), factor = NULL, gamma = rep(0, d), df = Inf,
             pseudo = TRUE, method = c("PRNG", "sobol", "ghalton"), skip = 0)
dskewtcopula(u, scale = diag(2), gamma = rep(0, d), df, log = FALSE,
             scale.inv, ldet)

Arguments

u

(n, d)-matrix of evaluation points or data; Have to be in (0,1).

x

(n, d)-matrix of evaluation points or data

n

sample size n (positive integer).

df

positive degress of freedom; can also be Inf in which case the copula is interpreted as the Gaussian copula.

loc

location of length d.

gamma

Skewness-vector of dimension d; if all(gamma == 0), the classical t distribution or copula results.

scale

scale matrix (a covariance matrix entering the distribution as a parameter) of dimension (d, d) (defaults to d = 2). Note that scale must be positive definite,

factor

(d, d)-matrix such that factor %*% t(factor) equals scale. If not provided, factor is internally determined via chol()

scale.inv

inverse of scale; if not provided, computed via pd.solve(scale).

ldet

log(det(scale)); if not provided, computed via pd.solve(scale).

log

logical indicating whether the logarithmic density is to be computed.

pseudo

logical; if TRUE, copula samples are computed via pobs() from a multivariate skew-t sample. If FALSE, the univariate skew t distribution functions are internally approximated via integrate(); see details below.

method

see rnvmix().

skip

see rnvmix().

Details

Functionalities for sampling from the multivariate skew-t distribution and copula; the former has stochastic representation \mu + W\gamma + \sqrt{W}AZ where AA^T=scale, W follows an inverse-gamma distrubution with parameters df/2 and is independent of the d-dimensional vector Z following a standard multivariate normal distribution. When gamma is the null-vector, the distribution becomes the multivariate t distribution.

A major computational challenge when working with the skew t copula is the lack of an available distribution and quantile function of the univariate skew t distribution. These are required in rskewtcopula(, pobs = FALSE) and in dskewtcopula(). The unviarate skew t distribution and quantile functions are currently implemented as described Yoshiba, T. (2018).

The functions described here are currently being further developed to improve stability, accuracy and speed, so that arguments may change in subsequent versions of nvmix.

Value

n-vector of (log-)density values and (n, d)-matrix of samples, respectively.

Author(s)

Erik Hintz, Marius Hofert and Christiane Lemieux

References

Hintz, E., Hofert, M. and Lemieux, C. (2020), Grouped Normal Variance Mixtures. Risks 8(4), 103.

Hintz, E., Hofert, M. and Lemieux, C. (2021), Normal variance mixtures: Distribution, density and parameter estimation. Computational Statistics and Data Analysis 157C, 107175.

Hintz, E., Hofert, M. and Lemieux, C. (2022), Multivariate Normal Variance Mixtures in R: The R Package nvmix. Journal of Statistical Software, doi:10.18637/jss.v102.i02.

McNeil, A. J., Frey, R. and Embrechts, P. (2015). Quantitative Risk Management: Concepts, Techniques, Tools. Princeton University Press.

Yoshiba, T. (2018). Maximum Likelihood Estimation of Skew-t Copulas with Its Applications to Stock Returns. Journal of Statistical Computation and Simulation 88 (13): 2489–2506.

See Also

rStudent(), dStudent(), rStudentcopula(), dStudentcopula()

Examples

## Sampling from the skew-t copula

n <- 100 # sample size
d <- 10 # dimension
rho <- 0.5
scale <- matrix(rho, ncol = d, nrow = d)
diag(scale) <- 1 # scale
gamma <- rep(1, d) # skewness
df <- 7 # degrees-of-freedom parameter
set.seed(1) # same random numbers for both runs
system.time(samplecop_pobs <- rskewtcopula(n, scale = scale, gamma = gamma,
                                           df = df, pseudo = TRUE))
set.seed(1)
system.time(samplecop_pskewt <- rskewtcopula(n, scale = scale, gamma = gamma,
                                             df = df, pseudo = FALSE))
## Plot first two coordinates
layout(rbind(1:2))
plot(samplecop_pobs, xlab = expression(U[1]), ylab = expression(U[2]))
mtext("pobs = TRUE")
plot(samplecop_pskewt, xlab = expression(U[1]), ylab = expression(U[2]))
mtext("pobs = FALSE")
layout(1)