R/expo.gibbs.R

In BaM: Functions and Datasets for "Bayesian Methods: A Social and Behavioral Sciences Approach"

#' expo.gibbs
#'
#' Simple Gibbs sampler demonstration on conditional exponentials from Chapter 1 (pages 25-27).
#' 
#' @usage expo.gibbs(B,k,m)
#' 
#' @param B an upper bound
#' @param k length of the subchains
#' @param m number of iterations
#'
#' @author Jeff Gill
#' @importFrom stats runif rexp
#' @export
expo.gibbs <- function(B=5, k=15, m=5000)  {
  x <- y <- NULL
  while (length(x) < m)  {
    x.val <- c(runif(1,0,B),rep((B+1),length=k))
    y.val <- c(runif(1,0,B),rep((B+1),length=k))
    for (j in 2:(k+1))  {
      while(x.val[j] > B) x.val[j] <- rexp(1,y.val[j-1])
      while(y.val[j] > B) y.val[j] <- rexp(1,x.val[j])
    }
    x <- c(x,x.val[(k+1)])
    y <- c(y,y.val[(k+1)])
  }
  return(cbind(x,y))
}