# rtmvt: Sampling Random Numbers From The Truncated Multivariate... In tmvtnorm: Truncated Multivariate Normal and Student t Distribution

## Description

This function generates random numbers from the truncated multivariate Student-t distribution with mean equal to `mean` and covariance matrix `sigma`, lower and upper truncation points `lower` and `upper` with either rejection sampling or Gibbs sampling.

## Usage

 ```1 2 3 4``` ```rtmvt(n, mean = rep(0, nrow(sigma)), sigma = diag(length(mean)), df = 1, lower = rep(-Inf, length = length(mean)), upper = rep(Inf, length = length(mean)), algorithm=c("rejection", "gibbs"), ...) ```

## Arguments

 `n` Number of random points to be sampled. Must be an integer >= 1. `mean` Mean vector, default is `rep(0, length = ncol(x))`. `sigma` Covariance matrix, default is `diag(ncol(x))`. `df` Degrees of freedom parameter (positive, may be non-integer) `lower` Vector of lower truncation points,\ default is `rep(-Inf, length = length(mean))`. `upper` Vector of upper truncation points,\ default is `rep( Inf, length = length(mean))`. `algorithm` Method used, possible methods are rejection sampling ("rejection", default) and the R Gibbs sampler ("gibbs"). `...` additional parameters for Gibbs sampling, given to the internal method `rtmvt.gibbs()`, such as `burn.in.samples`, `start.value` and `thinning`, see details

## Details

We sample x ~ T(mean, Sigma, df) subject to the rectangular truncation lower <= x <= upper. Currently, two random number generation methods are implemented: rejection sampling and the Gibbs Sampler.

For rejection sampling `algorithm="rejection"`, we sample from `rmvt` and retain only samples inside the support region. The acceptance probability will be calculated with `pmvt`. `pmvt` does only accept integer degrees of freedom `df`. For non-integer `df`, `algorithm="rejection"` will throw an error, so please use `algorithm="gibbs"` instead.

The arguments to be passed along with `algorithm="gibbs"` are:

`burn.in.samples`

number of samples in Gibbs sampling to be discarded as burn-in phase, must be non-negative.

`start.value`

Start value (vector of length `length(mean)`) for the MCMC chain. If one is specified, it must lie inside the support region (lower <= start.value <= upper). If none is specified, the start value is taken componentwise as the finite lower or upper boundaries respectively, or zero if both boundaries are infinite. Defaults to NULL.

`thinning`

Thinning factor for reducing autocorrelation of random points in Gibbs sampling. Must be an integer >= 1. We create a Markov chain of length `(n*thinning)` and take only those samples `j=1:(n*thinning)` where `j %% thinning == 0` Defaults to 1 (no thinning of the chain).

## Warning

The same warnings for the Gibbs sampler apply as for the method `rtmvnorm`.

## Author(s)

Stefan Wilhelm <[email protected]>, Manjunath B G <[email protected]>

## References

Geweke, John F. (1991) Efficient Simulation from the Multivariate Normal and Student-t Distributions Subject to Linear Constraints. Computer Science and Statistics. Proceedings of the 23rd Symposium on the Interface. Seattle Washington, April 21-24, 1991, pp. 571–578 An earlier version of this paper is available at http://www.biz.uiowa.edu/faculty/jgeweke/papers/paper47/paper47.pdf

## Examples

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44``` ```########################################################### # # Example 1 # ########################################################### # Draw from multi-t distribution without truncation X1 <- rtmvt(n=10000, mean=rep(0, 2), df=2) X2 <- rtmvt(n=10000, mean=rep(0, 2), df=2, lower=c(-1,-1), upper=c(1,1)) ########################################################### # # Example 2 # ########################################################### df = 2 mu = c(1,1,1) sigma = matrix(c( 1, 0.5, 0.5, 0.5, 1, 0.5, 0.5, 0.5, 1), 3, 3) lower = c(-2,-2,-2) upper = c(2, 2, 2) # Rejection sampling X1 <- rtmvt(n=10000, mu, sigma, df, lower, upper) # Gibbs sampling without thinning X2 <- rtmvt(n=10000, mu, sigma, df, lower, upper, algorithm="gibbs") # Gibbs sampling with thinning X3 <- rtmvt(n=10000, mu, sigma, df, lower, upper, algorithm="gibbs", thinning=2) plot(density(X1[,1], from=lower[1], to=upper[1]), col="red", lwd=2, main="Gibbs vs. Rejection") lines(density(X2[,1], from=lower[1], to=upper[1]), col="blue", lwd=2) legend("topleft",legend=c("Rejection Sampling","Gibbs Sampling"), col=c("red","blue"), lwd=2) acf(X1) # no autocorrelation in Rejection sampling acf(X2) # strong autocorrelation of Gibbs samples acf(X3) # reduced autocorrelation of Gibbs samples after thinning ```

tmvtnorm documentation built on May 29, 2017, 9:36 p.m.