saveJAGS: Write MCMC draws to file
In mikemeredith/saveJAGS: Run JAGS and Regularly Save Output to Files

Description Usage Arguments Details Value Author(s) Examples

View source: R/saveJags.R

Takes a list of data, initial values, a list of parameters to monitor, and the name of a text file containing the model, calls JAGS via rjags, and writes the output to a series of files.

Note: Currently, this only works with parallel processing with one chain per core.

1
2
3

saveJAGS(data, inits, params, modelFile, fileStub,
        chains=3, sample2save=1000, nSaves=3, burnin=1000, thin=1,
        modules = "glm", firstChainID="AA")

`data`	a named list or environment containing the data.
`inits`	optional specification of initial values in the form of a list or a function.
`params`	a character vector giving the names of variables to be monitored.
`modelFile`	the name of the file containing a description of the model.
`chains`	the number of chains to run; must not exceed the number of cores/threads available on the machine.
`sample2save`	the number of MCMC draws to include in each file AFTER thinning.
`nSaves`	the number of files to write per chain; the total iterations generated per chain = sample2save * thin * nSaves.
`burnin`	the number of iterations to discard at the beginning of each chain; unlike other wrappers, adaptation continues until the end of the burn-in sequence, there is no separate adaptation phase.
`thin`	the thinning interval.
`fileStub`	a character value to use as the first part of the file name; it can include the path provided the target directory exists; to this is added a letter to identify the chain, a sequence number and date-time the block was started.
`modules`	a character vector with names of the modules to load before calling JAGS.
`firstChainID`	the ID for the first chain to be saved; must be doubled upper-case letters: AA, AB, AC, ..., ZZ.

After running saveJAGS, you can extract the saved values to create an mcmc.list object with combineSaves; see the examples. recoverSaves can be used to recreate the file list if saveJAGS terminates before completion.

Note that pressing 'Esc' will interrupt the master process in the R Console, but will not stop the workers, which will continue until the task is finished, even if you close R. To stop the workers, right-click on the task in Task Manager and select 'End task' (Windows) or highlight in Activity Monitor and press the 'Force Quit' button (Mac).

The files written are *.RData files, which are compressed binary files. Each file contains four objects:

`out`	an object of class `mcmc.list` with 1 MCMC chain.
`jm`	the JAGS model; use `print(jm)` to see the model code, `jm$data()` to access the data.
`JAGSsettings`	a list with details of the modules loaded and the status of samplers.
`adaptIsAdequate`	logical, TRUE if adaptation was adequate.

An object of class 'saveJAGSfileList', a list of file names with one component per chain. There is a summary method for the class.

Mike Meredith. Based on code used to run JAGS in the wiqid package.

# An example with salamander occupancy: the number of occasions (out of 5)
#  that salamanders were detected at each of 39 sites
sal <- rep(0:4, c(21,12,1,4,1))

# JAGS code for psi(.) p(.) occupancy model
modelText <- "
model {
  for(i in 1:nSites) {
    z[i] ~ dbern(psi)
    y[i] ~ dbin(p * z[i], n)
  }
  psi ~ dbeta(1, 1)
  p ~ dbeta(1, 1)
} "
writeLines(modelText, con = "JAGSmodel.txt")

JAGSdata <- list(y = sal, n = 5, nSites = length(sal))
inits <- function(chain) list(z = rep(1, 39))
wanted <- c("p", "psi", "z")

# Create a folder for the output:
dir.create("mysaves")

res1 <- saveJAGS(JAGSdata, inits, wanted, "JAGSmodel.txt",
        chains=2, sample2save=1000, nSaves=4, burnin=1000, thin=1,
        fileStub="mysaves/testing")
str(res1)
summary(res1)

# If the saveJAGS run was terminated before completion,
#   we can recover the file list with recoverSaves:
res2 <- recoverSaves("mysaves/testing")
all.equal(res1, res2)

# Load the results into R as an mcmc.list object
mcmc <- as.mcmc.list(res1)
str(mcmc)
#... or subset by parameters (omit z), or thin:
mcmc2 <- as.mcmc.list(res1, params=c("p", "psi"), thin=4)
str(mcmc2)

# Load the same subset into R as an mcmcOutput object
( mco <- mcmcOutput(res1, params=c("p", "psi"), thin=4) )
summary(mco)
# Extract individual parameters
str(mco$psi)

# Clean up
unlink("JAGSmodel.txt")
unlink("mySaves", recursive=TRUE)