stan.fit: Fit Stan models with cloned data

Description

Convenient functions designed to work well with cloned data arguments and Stan.

Usage

stan.fit(data, params, model, inits = NULL,
    seed = sample.int(.Machine$integer.max, 1),
    n.chains = 3,
    format = c("mcmc.list", "stanfit"),
    stan.model = TRUE, fit = NA, ...)
stan.model(object, ...)
stan.parfit(cl, data, params, model, inits = NULL,
    seed = sample.int(.Machine$integer.max, n.chains),
    n.chains = 3,
    format = c("mcmc.list", "stanfit"),
    stan.model = TRUE, fit = NA, ...)

Value

By default, an stan.fit returns an mcmc.list object. If data cloning is used via the data argument, summary returns a modified summary containing scaled data cloning standard errors (scaled by sqrt(n.clones)), and \(R_{hat}\) values (as returned by gelman.diag).

stan.model returns the stanmodel object.

stan.parfit runs chains using multiple cores when cl

is an integer. Using a cluster object leads to recompiling the model (therefore fit is ignored), and might not be very quick to run.

Arguments

data: A list (or environment) containing the data.
params: Character vector of parameters to be sampled.
model: Character string (name of the model file), a function containing the model, or a custommodel object.
inits: Optional specification of initial values in the form of a list or a function. If NULL, initial values will be generated automatically.
seed: Random seed.
n.chains: number of Markov chains.
format: Desired output format.
stan.model: Logical, if stanmodel object should be returned.
fit: Fitted Stan object.
cl: A cluster object created by makeCluster, or an integer, see parDosa and evalParallelArgument.
object: A fitted MCMC object ('mcmc.list' class for example), with "stan.model" attribute.
...: Further arguments.

Author

Peter Solymos, solymos@ualberta.ca

Examples

Run this code

if (FALSE) {
if (require(rstan)) {
    model <- custommodel("data {
          int N;
          vector[N] y;
          vector[N] x;
        }
        parameters {
          real alpha;
          real beta;
          real sigma;
        }
        model {
          alpha ~ normal(0,10);
          beta ~ normal(0,10);
          sigma ~ cauchy(0,5);
          for (n in 1:N)
            y[n] ~ normal(alpha + beta * x[n], sigma);
        }")
    N <- 100
    alpha <- 1
    beta <- -1
    sigma <- 0.5
    x <- runif(N)
    y <- rnorm(N, alpha + beta * x, sigma)
    dat <- list(N=N, y=y, x=x)
    params <- c("alpha", "beta", "sigma")
    ## compile on 1st time only
    fit0 <- stan.fit(dat, params, model)
    ## reuse compiled fit0
    fit <- stan.fit(dat, params, model, fit=fit0)
    sm <- stan.model(fit)
    summary(fit)
    sm

    ## data cloning
    dcdat <- dclone(dat, n.clones=2, multiply="N")
    dcfit <- stan.fit(dcdat, params, model, fit=fit0)
    summary(dcfit)
    nclones(dcfit)

    ## using parallel options
    cl <- makeCluster(2)
    ## cannot utilize compiled fit0
    fit2 <- stan.parfit(cl=cl, dat, params, model)
    stopCluster(cl)
    if (.Platform$OS.type != "windows") {
        ## utilize compiled fit0
        fit3 <- stan.parfit(cl=2, dat, params, model, fit=fit0)
    }
}
}