openbugs: Wrapper to run OpenBUGS

Description

The openbugs function takes data and starting values as input. It automatically calls the package BRugs and runs something similar to BRugsFit. Not available in S-PLUS.

Usage

openbugs(data, inits, parameters.to.save,
    model.file = "model.txt", n.chains = 3, n.iter = 2000,
    n.burnin = floor(n.iter/2),
    n.thin = max(1, floor(n.chains * (n.iter - n.burnin) / n.sims)),
    n.sims = 1000,  DIC = TRUE, 
    bugs.directory = "c:/Program Files/OpenBUGS/",
    working.directory = NULL, digits = 5, over.relax = FALSE, seed=NULL)

Value

A bugs object.

Arguments

data: either a named list (names corresponding to variable names in the model.file) of the data for the OpenBUGS model, or a vector or list of the names of the data objects used by the model. If data is a one element character vector (such as "data.txt"), it is assumed that data have already been written to the working directory into that file, e.g. by the function bugs.data.
inits: a list with n.chains elements; each element of the list is itself a list of starting values for the OpenBUGS model, or a function creating (possibly random) initial values. Alternatively, if inits are missing or inits = NULL, initial values are generated by OpenBUGS.
parameters.to.save: character vector of the names of the parameters to save which should be monitored
model.file: file containing the model written in OpenBUGS code. The extension can be either .bug or .txt. If .bug, a copy of the file with extension .txt will be created in the bugs() call and removed afterwards. Note that similarly named .txt files will be overwritten.
n.chains: number of Markov chains (default: 3)
n.iter: number of total iterations per chain (including burn in; default: 2000)
n.burnin: length of burn in, i.e. number of iterations to discard at the beginning. Default is n.iter/2, that is, discarding the first half of the simulations.
n.thin: thinning rate. Must be a positive integer. Set n.thin > 1 to save memory and computation time if n.iter is large. Default is max(1, floor(n.chains * (n.iter-n.burnin) / 1000)) which will only thin if there are at least 2000 simulations.
n.sims: The approximate number of simulations to keep after thinning.
DIC: logical; if TRUE (default), compute deviance, pD, and DIC. This is done in BRugs directly.
digits: number of significant digits used for OpenBUGS input, see formatC
bugs.directory: directory that contains the OpenBUGS executable - currently unused
working.directory: sets working directory during execution of this function; WinBUGS in- and output will be stored in this directory; if NULL, a temporary working directory via tempdir is used.
over.relax: If TRUE, over-relaxed form of MCMC is used if available from OpenBUGS.
seed: random seed (default is no seed)

Author

Andrew Gelman, gelman@stat.columbia.edu; modifications and packaged by Sibylle Sturtz, sturtz@statistik.tu-dortmund.de, and Uwe Ligges.

Examples

Run this code

# An example model file is given in:
model.file <- system.file(package = "R2WinBUGS", "model", "schools.txt")
# Let's take a look:
file.show(model.file)

# Some example data (see ?schools for details):
data(schools)
schools

J <- nrow(schools)
y <- schools$estimate
sigma.y <- schools$sd
data <- list ("J", "y", "sigma.y")
inits <- function(){
    list(theta = rnorm(J, 0, 100), mu.theta = rnorm(1, 0, 100),
         sigma.theta = runif(1, 0, 100))
}
## or alternatively something like:
# inits <- list(
#   list(theta = rnorm(J, 0, 90), mu.theta = rnorm(1, 0, 90),
#        sigma.theta = runif(1, 0, 90)),
#   list(theta = rnorm(J, 0, 100), mu.theta = rnorm(1, 0, 100),
#        sigma.theta = runif(1, 0, 100))
#   list(theta = rnorm(J, 0, 110), mu.theta = rnorm(1, 0, 110),
#        sigma.theta = runif(1, 0, 110)))

parameters <- c("theta", "mu.theta", "sigma.theta")

if (FALSE) {
## both write access in the working directory and package BRugs required:
schools.sim <- bugs(data, inits, parameters, model.file,
    n.chains = 3, n.iter = 5000,
    program = "openbugs", working.directory = NULL)
print(schools.sim)
plot(schools.sim)
}