dcmle-package: Hierarchical Models Made Easy with Data Cloning

Description

S4 classes around infrastructure provided by the dclone package to make package development with data cloning for hierarchical models easy as a breeze.

Arguments

Author

Peter Solymos

Maintainer: Peter Solymos <psolymos@gmail.com>

Details

The package defines S4 object classes for plain BUGS models ("gsFit", after BU*GS*/JA*GS*), and BUGS models made ready for data cloning ("dcFit"). It also defines virtual classes for S3 object classes defined in the dclone and coda packages.

The S4 class "dcmle" is a fitted model object containing MCMC results as returned by the dcmle function. These object classes are easily extensible to allow inclusion into functions fitting specific models to the data (see Examples).

References

Forum: https://groups.google.com/forum/#!forum/dclone-users

Issues: https://github.com/datacloning/dcmle/issues

Data cloning website: https://datacloning.org

Examples

Run this code

## Data and model taken from Ponciano et al. 2009
## Ecology 90, 356-362.

## Function to create template object for the Beverton-Holt model
## R CMD check will not choke on character representation of model
## the convenient makeDcFit creator function is used here
bevholtFit <- 
function(y) {
makeDcFit(
  data = list(ncl=1, n=length(y), Y=dcdim(data.matrix(y))),
  model = structure(
    c("model {",
      "  for (k in 1:ncl) {",
      "    for(i in 2:(n+1)) {",
      "      Y[(i-1), k] ~ dpois(exp(X[i, k]))",
      "      X[i, k] ~ dnorm(mu[i, k], 1 / sigma^2)",
      "      mu[i,k] <- X[(i-1),k]+log(lambda)-log(1+beta*exp(X[(i-1),k]))",
      "    }",
      "    X[1, k] ~ dnorm(mu0, 1 / sigma^2)",
      "  }",
      "  beta ~ dlnorm(-1, 1)",
      "  sigma ~ dlnorm(0, 1)",
      "  tmp ~ dlnorm(0, 1)",
      "  lambda <- tmp + 1",
      "  mu0 <- log(2)  + log(lambda) - log(1 + beta * 2)",
      "}"),
      class = "custommodel"),
  multiply = "ncl",
  unchanged = "n",
  params <- c("lambda","beta","sigma"))
}
## S4 class 'bevholtMle' extends the 'dcmle' class
## it can have additional slots
setClass("bevholtMle", 
    representation(y="numeric", title="character"), 
    contains = "dcmle")
## Function to fit the Beverton-Holt model to data
bevholt <- function(y, n.clones, ...) {
    new("bevholtMle",
        dcmle(bevholtFit(y), n.clones=n.clones, ...),
        y = y,
        title = "Beverton-Holt Model")
}
## Show method with appropriate heading
setMethod("show", "bevholtMle", function(object)
    show(summary(as(object, "dcmle"), object@title)))
paurelia <- c(17,29,39,63,185,258,267,392,510,
    570,650,560,575,650,550,480,520,500)
if (FALSE) {
(m <- bevholt(paurelia, n.clones=2, n.iter=1000))
vcov(m)
m@y
}