Run a Julia HMC Sampler for a BUGS-like Probabilistic Model

Executes a Hamiltonian Monte Carlo (HMC) sampler in Julia from R, using a model specified in Julia or in BUGS syntax. It compiles the model, converts data, sets sampler parameters, and returns posterior samples in various formats. The setup for the HMC sampler uses Not-U-Turn Sampler (NUTS) with the target acceptance probability δ=0.8) for step size adaptation.

Usage

juliaBUGS(
  data,
  model_def,
  params_to_save,
  initializations = NULL,
  name = "sampler_juliaBUGS",
  n_iter = 2000,
  n_warmup = floor(n_iter/2),
  n_discard = n_warmup,
  n_thin = 1,
  n_chain = 1,
  use_parallel = TRUE,
  posterior_type = "array",
  force_setup_juliaBUGS = FALSE,
  control = NULL,
  progress = TRUE,
  ...
)

Arguments

data

A named list of numeric values (integer or double). All elements must be named.

model_def

A character string with the model definition, either in Julia-compatible format or BUGS syntax.

params_to_save

Character vector with the names of model parameters to extract from the sampler output.

initializations

A named list of parameter names for which you may wish to set corresponding initial values for the sampler. The default is NULL, which means no default initial values are used.

name

Character. Name for the sampler object created in Julia (must be a valid Julia variable name).

n_iter

Integer. Total number of MCMC iterations. Default is 2000.

n_warmup

Integer. Number of iterations used warm-up or tuning (e.g., adaption steps in NUTS). Default is floor(n_iter / 2).

n_discard

Integer. Number of initial samples to be completely discarded. Default is n_warmup, i.e: discard all the iterations used as adaptation steps.

n_thin

Integer. Thinning interval. Default is 1 (no thinning).

n_chain

Integer. Number of MCMC chains. Default is 1.

use_parallel

Logical. Whether to use AbstractMCMC.MCMCThreads() for parallel sampling. Default is TRUE.

posterior_type

Character. Format of the posterior samples. One of "array", "rvar", "mcmc", or "draws". Default is "array".

force_setup_juliaBUGS

Logical. If TRUE, forces reinitialization of the Julia environment via setup_juliaBUGS(). Default is FALSE.

control

Optional list of control parameters. Supported entries:

data_convert_int: Logical. If TRUE, coerces numeric values to integers when possible. Default is TRUE.
convert_var_name: Logical. If TRUE, automatically renames variables in the BUGS model. Default is FALSE.
julia_model: Logical. If TRUE, assumes the model is already in Julia format used by the models under the Turing.jl appraoch, not a BUGS model. Default is FALSE.

progress

Logical. If TRUE, a progress bar will be displayed; if FALSE, it will not. The default is TRUE.

...

Additional arguments passed to setup_juliaBUGS().

Value

An object of class "rjuliabugs" containing:

params: Posterior samples, in the format specified by posterior_type.
name: Character string identifying the Julia sampler object.
sampler: Sampler object returned by AbstractMCMC.sample.
n_threads: Number of Julia threads detected.
mcmc: List of MCMC configuration parameters.
control: Control options used in the sampler setup.

Details

This function relies on Julia packages LogDensityProblems, AdvancedHMC, and AbstractMCMC. Gradients are computed via ReverseDiff. The model is compiled before sampling.

The posterior_type argument determines the return format:

"array": 3D numeric array (iterations × chains × parameters).
"rvar": posterior::rvar object.
"mcmc": coda::mcmc (single chain) or mcmc.list (multiple chains).
"draws": posterior::draws_array.

Note

You must call setup_juliaBUGS() at least once before using this function. If parallel sampling is requested but only one Julia thread is available, a warning is issued and sampling will run serially.

Examples

if (FALSE) { # \dontrun{
model_def <- "model = @model ... end"
data <- list(N = 10, x = rnorm(10))
result <- juliaBUGS(
  data = data,
  model_def = model_def,
  params_to_save = c("mu"),
  name = "my_sampler"
)
} # }