#' Matrix Gamma sampler
#'
#' Samples a matrix Gamma distribution.
#'
#' @param n sample size, a positive integer
#' @param nu shape parameter, a positive number; if \code{nu < (p-1)/2},
#' where \code{p} is the dimension (the order of \code{Sigma}), then \code{nu}
#' must be a half integer
#' @param theta scale parameter, a positive number
#' @param Sigma scale matrix, a symmetric positive definite matrix, or
#' \code{NULL} for the identity matrix of order \code{p}
#' @param p if \code{Sigma} is \code{NULL}, this sets \code{Sigma} to the
#' identity matrix of order \code{p}; ignored if \code{Sigma} is not \code{NULL}
#' @param checkSymmetry logical, whether to check that \code{Sigma} is a
#' symmetric positive definite matrix
#'
#' @return A numeric three-dimensional array;
#' simulations are stacked along the third dimension.
#' @export
#'
#' @details This is the distribution of
#' \ifelse{html}{\out{θ/2×S}}{\eqn{\frac{\theta}{2}S}} where
#' \ifelse{html}{\out{S ~ W<sub>p</sub>(2ν,Σ)}}{\eqn{S \sim \mathcal{W}_p(2\nu,\Sigma)}}.
#'
#' @references Gupta & al.
#' Properties of Matrix Variate Confluent Hypergeometric Function Distribution.
#' \emph{Journal of Probability and Statistics}
#' vol. 2016, Article ID 2374907, 12 pages, 2016.
#'
#' @examples nu <- 3; theta <- 4; Sigma <- toeplitz(2:1)
#' Gsims <- rmatrixgamma(10000, nu, theta, Sigma)
#' apply(Gsims, c(1,2), mean) # should be nu * theta * Sigma
#' nu * theta * Sigma
rmatrixgamma <- function(n, nu, theta, Sigma = NULL, p, checkSymmetry = TRUE){
stopifnot(nu > 0, theta > 0)
if(is.null(Sigma) && missing(p)){
stop("Enter `Sigma` or `p`")
}
if(!is.null(Sigma) && !(isRealScalar(Sigma) || nrow(Sigma) == ncol(Sigma))){
stop("`Sigma` is not square")
}
if(checkSymmetry && !is.null(Sigma)){
. <- isSymmetricPositive(Sigma)
}
SigmaIsId <- is.null(Sigma) || isTRUE(all.equal(Sigma, diag(nrow(Sigma))))
if(SigmaIsId){
if(!is.null(Sigma)) p <- nrow(Sigma)
theta/2 * rwishart_I(n, 2*nu, p)
}else{
theta/2 * rwishart(n, 2*nu, Sigma, checkSymmetry = FALSE)
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.