relliptical (version 1.3.0)

rtelliptical: Sampling Random Numbers from Truncated Multivariate Elliptical Distributions

Description

This function generates observations from a truncated multivariate elliptical distribution with location parameter mu, scale matrix Sigma, lower and upper truncation points lower and upper via Slice Sampling algorithm with Gibbs sampler steps.

Usage

rtelliptical(n = 10000, mu = rep(0, length(lower)),
  Sigma = diag(length(lower)), lower, upper = rep(Inf, length(lower)),
  dist = "Normal", nu = NULL, expr = NULL, gFun = NULL,
  ginvFun = NULL, burn.in = 0, thinning = 1)

Value

It returns a matrix of dimensions \(n\)x\(p\) with the random points sampled.

Arguments

n

number of observations to generate. Must be an integer >= 1.

mu

numeric vector of length \(p\) representing the location parameter.

Sigma

numeric positive definite matrix with dimension \(p\)x\(p\) representing the scale parameter.

lower

vector of lower truncation points of length \(p\).

upper

vector of upper truncation points of length \(p\).

dist

represents the truncated distribution to be used. The values are 'Normal', 't', 'Laplace', 'PE', 'PVII', 'Slash', and 'CN' for the truncated Normal, Student-t, Laplace, Power Exponential, Pearson VII, Slash, and Contaminated Normal distribution, respectively.

nu

additional parameter or vector of parameters depending on the density generating function. See Details.

expr

a character with the density generating function. See Details.

gFun

an R function with the density generating function. See Details.

ginvFun

an R function with the inverse of the density generating function defined in gFun. See Details.

burn.in

number of samples to be discarded as a burn-in phase.

thinning

factor for reducing the autocorrelation of random points.

Author

Katherine L. Valeriano, Christian E. Galarza and Larissa A. Matos

Details

The dist argument represents the truncated distribution to be used. The values are Normal, t, 't', PE, PVII, Slash, and CN, for the truncated Normal, Student-t, Laplace, Power Exponential, Pearson VII, Slash, and Contaminated Normal distribution, respectively.

The argument nu is a parameter or vector of parameters depending on the density generating function (DGF). For the truncated Student-t, Power Exponential, and Slash distribution, nu is a positive number. For the truncated Pearson VII, nu is a vector with the first element greater than \(p/2\) and the second element a positive number. For the truncated Contaminated Normal distribution, nu is a vector of length 2 assuming values between 0 and 1.

This function also allows generating random numbers from other truncated elliptical distributions not specified in the dist argument, by supplying the density generating function (DGF) through arguments either expr or gFun. The DGF must be a non-negative and strictly decreasing function on (0, Inf). The easiest way is to provide the DGF expression to argument expr as a character. The notation used in expr needs to be understood by package Ryacas0, and the environment of R. For instance, for the DGF \(g(t)=e^{-t}\), the user must provide expr = "exp(1)^(-t)". See that the function must depend only on variable \(t\), and any additional parameter must be passed as a fixed value. For this case, when a character expression is provided to expr, the algorithm tries to compute a closed-form expression for the inverse function of \(g(t)\), however, this is not always possible (a warning message is returned). See example 2.

If it was no possible to generate random samples by passing a character expression to expr, the user may provide a custom R function to the gFun argument. By default, its inverse function is approximated numerically, however, the user may also provide its inverse to the ginvFun argument to gain some computational time. When gFun is provided, arguments dist and expr are ignored.

References

fang2018symmetricrelliptical ho2012somerelliptical neal2003slicerelliptical robert2010introducingrelliptical valeriano2023momentsrelliptical

See Also

mvtelliptical

Examples

library(ggplot2)
library(ggExtra)
library(gridExtra)

# Example 1: Sampling from the Truncated Normal distribution
set.seed(1234)
mu  = c(0, 1)
Sigma = matrix(c(1,0.70,0.70,3), 2, 2)
lower = c(-2, -3)
upper = c(3, 3)
sample1 = rtelliptical(5e4, mu, Sigma, lower, upper, dist="Normal")

# Histogram and density for variable 1
ggplot(data.frame(sample1), aes(x=X1)) +
   geom_histogram(aes(y=after_stat(density)), colour="black", fill="grey", bins=15) +
   geom_density(color="red") + labs(x=bquote(X[1]), y="Density") + theme_bw()

# Histogram and density for variable 2
ggplot(data.frame(sample1), aes(x=X2)) +
   geom_histogram(aes(y=after_stat(density)), colour="black", fill="grey", bins=15) +
   geom_density(color="red") + labs(x=bquote(X[2]), y="Density") + theme_bw()

# \donttest{
# Example 2: Sampling from the Truncated Logistic distribution

# Function for plotting the sample autocorrelation using ggplot2
acf.plot = function(samples){
 p = ncol(samples); n = nrow(samples); q1 = qnorm(0.975)/sqrt(n); acf1 = list(p)
 for (i in 1:p){
   bacfdf = with(acf(samples[,i], plot=FALSE), data.frame(lag, acf))
   acf1[[i]] = ggplot(data=bacfdf, aes(x=lag,y=acf)) + geom_hline(aes(yintercept=0)) +
     geom_segment(aes(xend=lag, yend=0)) + labs(x="Lag", y="ACF", subtitle=bquote(X[.(i)])) +
     geom_hline(yintercept=c(q1,-q1), color="red", linetype="twodash") + theme_bw()
 }
 return (acf1)
}

set.seed(5678)
mu  = c(0, 0)
Sigma = matrix(c(1,0.70,0.70,1), 2, 2)
lower = c(-2, -2)
upper = c(3, 2)
# Sample autocorrelation with no thinning
sample2 = rtelliptical(10000, mu, Sigma, lower, upper, expr="exp(1)^(-t)/(1+exp(1)^(-t))^2")
grid.arrange(grobs=acf.plot(sample2), top="Logistic distribution with no thinning", nrow=1)

# Sample autocorrelation with thinning = 3
sample3 = rtelliptical(10000, mu, Sigma, lower, upper, expr="exp(1)^(-t)/(1+exp(1)^(-t))^2",
                       thinning=3)
grid.arrange(grobs=acf.plot(sample3), top="Logistic distribution with thinning = 3", nrow=1)
# }

# Example 3: Sampling from the Truncated Kotz-type distribution
set.seed(5678)
mu  = c(0, 0)
Sigma = matrix(c(1,-0.5,-0.5,1), 2, 2)
lower = c(-2, -2)
upper = c(3, 2)
sample4 = rtelliptical(2000, mu, Sigma, lower, upper, gFun=function(t){t^(-1/2)*exp(-2*t^(1/4))})
f1 = ggplot(data.frame(sample4), aes(x=X1,y=X2)) + geom_point(size=0.50) +
     labs(x=expression(X[1]), y=expression(X[2]), subtitle="Kotz(2,1/4,1/2)") + theme_bw()
ggMarginal(f1, type="histogram", fill="grey")