| rtelliptical {relliptical} | R Documentation |
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)
Arguments
n |
number of observations to generate. Must be an integer >= 1. |
mu |
numeric vector of length |
Sigma |
numeric positive definite matrix with dimension |
lower |
vector of lower truncation points of length |
upper |
vector of upper truncation points of length |
dist |
represents the truncated distribution to be used. The values are |
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
|
burn.in |
number of samples to be discarded as a burn-in phase. |
thinning |
factor for reducing the autocorrelation of random points. |
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.
Value
It returns a matrix of dimensions nxp with the random points sampled.
Note
The Normal distribution is a particular case of the Power Exponential distribution when nu = 1.
The Student-t distribution with \nu degrees of freedom results from the Pearson VII
distribution when nu = ((\nu+p)/2, \nu).
Author(s)
Katherine L. Valeriano, Christian E. Galarza and Larissa A. Matos
References
Fang KT, Kotz S, Ng KW (2018). Symmetric multivariate and related distributions. Chapman and Hall/CRC.
Ho HJ, Lin TI, Chen HY, Wang WL (2012). “Some results on the truncated multivariate t distribution.” Journal of Statistical Planning and Inference, 142(1), 25–40. doi:10.1016/j.jspi.2011.06.006.
Neal RM (2003). “Slice sampling.” Annals of statistics, 705–741.
Robert CP, Casella G (2010). Introducing Monte Carlo Methods with R, volume 18. New York: Springer.
Valeriano KA, Galarza CE, Matos LA (2023). “Moments and random number generation for the truncated elliptical family of distributions.” Statistics and Computing, 33(1), 32.
See Also
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()
# 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")