Title:	No-U-Turn MCMC Sampling for 'ADMB' Models
Version:	1.1.2
Description:	Bayesian inference using the no-U-turn (NUTS) algorithm by Hoffman and Gelman (2014) https://www.jmlr.org/papers/v15/hoffman14a.html. Designed for 'AD Model Builder' ('ADMB') models, or when R functions for log-density and log-density gradient are available, such as 'Template Model Builder' models and other special cases. Functionality is similar to 'Stan', and the 'rstan' and 'shinystan' packages are used for diagnostics and inference.
Depends:	R (≥ 3.6.0), snowfall (≥ 1.84.6.1)
URL:	https://github.com/Cole-Monnahan-NOAA/adnuts
BugReports:	https://github.com/Cole-Monnahan-NOAA/adnuts/issues
License:	GPL-3 | file LICENSE
Encoding:	UTF-8
LazyData:	true
RoxygenNote:	7.1.1
ByteCompile:	true
Suggests:	shinystan (≥ 2.5.0), matrixcalc (≥ 1.0.3), stats, knitr, TMB, rmarkdown, withr, testthat (≥ 2.1.0)
Imports:	ellipse, rstan, R2admb, ggplot2, rlang
VignetteBuilder:	knitr
NeedsCompilation:	no
Packaged:	2021-03-02 19:04:08 UTC; cole.monnahan
Author:	Cole Monnahan [aut, cre]
Maintainer:	Cole Monnahan <monnahc@uw.edu>
Repository:	CRAN
Date/Publication:	2021-03-02 19:50:09 UTC

Check that the model is compiled with the right version of ADMB which is 12.0 or later

Description

Check that the model is compiled with the right version of ADMB which is 12.0 or later

Usage

.check_ADMB_version(model, path = getwd(), min.version = 12, warn = TRUE)

Arguments

Details

Some functionality of packages adnuts is imbedded in the ADMB source code so that when a model is compiled it is contained in the model executable. If this code does not exist adnuts will fail. The solution is to update ADMB and recompile the model.

Value

Nothing, errors out if either model could not be run or the version is incompatible. If compatible nothing happens.

Check if the session is interactive or Rstudio which has implications for parallel output

Description

Check if the session is interactive or Rstudio which has implications for parallel output

Usage

.check_console_printing(parallel)

Arguments

Details

When using RStudio and RGui, the parallel output does not show on the console. As a workaround it is captured in each cluster into a file and then read in and printed.

Value

Boolean whether output should be printed to console progressively, or saved to file and printed at the end.

Check that the file can be found

Description

Check that the file can be found

Usage

.check_model_path(model, path)

Arguments

Read in admodel.hes file

Description

Read in admodel.hes file

Usage

.getADMBHessian(path)

Arguments

Value

The Hessian matrix

Hidden wrapper function for sampling from ADMB models

Description

Hidden wrapper function for sampling from ADMB models

Usage

.sample_admb(
  model,
  path = getwd(),
  iter = 2000,
  init = NULL,
  chains = 3,
  warmup = NULL,
  seeds = NULL,
  thin = 1,
  mceval = FALSE,
  duration = NULL,
  cores = NULL,
  control = NULL,
  verbose = TRUE,
  algorithm = "NUTS",
  skip_optimization = TRUE,
  skip_monitor = FALSE,
  skip_unbounded = TRUE,
  admb_args = NULL
)

Arguments

Convert model name depending on system

Description

Convert model name depending on system

Usage

.update_model(model)

Arguments

Value

Updated model name to use with system call

Constructor for the "adfit" (A-D fit) class

Description

Constructor for the "adfit" (A-D fit) class

Usage

adfit(x)

Arguments

Value

An object of class "adfit"

adnuts: No-U-turn sampling for AD Model Builder (ADMB)

Description

Draw Bayesian posterior samples from an ADMB model using the no-U-turn MCMC sampler. Adaptation schemes are used so specifying tuning parameters is not necessary, and parallel execution reduces overall run time.

Details

The software package Stan pioneered the use of no-U-turn (NUTS) sampling for Bayesian models (Hoffman and Gelman 2014, Carpenter et al. 2017). This algorithm provides fast, efficient sampling across a wide range of models, including hierarchical ones, and thus can be used as a generic modeling tool (Monnahan et al. 2017). The functionality provided by adnuts is based loosely off Stan and R package rstan

The adnuts R package provides an R workflow for NUTS sampling for ADMB models (Fournier et al. 2011), including adaptation of step size and metric (mass matrix), parallel execution, and links to diagnostic and inference tools provided by rstan and shinystan. The ADMB implementation of NUTS code is bundled into the ADMB source itself (as of version 12.0). Thus, when a user builds an ADMB model the NUTS code is incorporated into the model executable. Thus, adnuts simply provides a convenient set of wrappers to more easily execute, diagnose, and make inference on a model. More details can be found in the package vignette.

Note that previous versions of adnuts included functionality for TMB models, but this has been replaced by tmbstan (Kristensen et al. 2016, Monnahan and Kristensen 2018).

References

Carpenter, B., Gelman, A., Hoffman, M.D., Lee, D., Goodrich, B., Betancourt, M., Riddell, A., Guo, J.Q., Li, P., Riddell, A., 2017. Stan: A Probabilistic Programming Language. J Stat Softw. 76:1-29.

Fournier, D.A., Skaug, H.J., Ancheta, J., Ianelli, J., Magnusson, A., Maunder, M.N., Nielsen, A., Sibert, J., 2012. AD Model Builder: using automatic differentiation for statistical inference of highly parameterized complex nonlinear models. Optim Method Softw. 27:233-249.

Hoffman, M.D., Gelman, A., 2014. The no-U-turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo. J Mach Learn Res. 15:1593-1623.

Kristensen, K., Nielsen, A., Berg, C.W., Skaug, H., Bell, B.M., 2016. TMB: Automatic differentiation and Laplace approximation. J Stat Softw. 70:21.

Kristensen, K., 2017. TMB: General random effect model builder tool inspired by ADMB. R package version 1.7.11.

Monnahan, C.C., Thorson, J.T., Branch, T.A., 2017. Faster estimation of Bayesian models in ecology using Hamiltonian Monte Carlo. Methods in Ecology and Evolution. 8:339-348.

Monnahan C.C., Kristensen K. (2018). No-U-turn sampling for fast Bayesian inference in ADMB and TMB: Introducing the adnuts and tmbstan R packages PLoS ONE 13(5): e0197954. https://doi.org/10.1371/journal.pone.0197954

Stan Development Team, 2016. Stan modeling language users guide and reference manual, version 2.11.0.

Stan Development Team, 2016. RStan: The R interface to Stan. R package version 2.14.1. http://mc-stan.org.

Convert object of class adfit to data.frame. Calls extract_samples

Description

Convert object of class adfit to data.frame. Calls extract_samples

Usage

## S3 method for class 'adfit'
as.data.frame(x, row.names = NULL, optional = FALSE, ...)

Arguments

Details

This calls the default settings of extract_samples, no warmup samples and no column for the log-posterior (lp__). Use this function directly for finer control.

Value

A data frame with parameters as columns and samples as rows.

Check identifiability from model Hessian

Description

Check identifiability from model Hessian

Usage

check_identifiable(model, path = getwd())

Arguments

Details

Read in the admodel.hes file and check the eigenvalues to determine which parameters are not identifiable and thus cause the Hessian to be non-invertible. Use this to identify which parameters are problematic. This function was converted from a version in the FishStatsUtils package.

Value

Prints output of bad parameters and invisibly returns it.

Extract sampler parameters from a fit.

Description

Extract information about NUTS trajectories, such as acceptance ratio and treedepth, from a fitted object.

Usage

extract_sampler_params(fit, inc_warmup = FALSE)

Arguments

Details

Each trajectory (iteration) in NUTS has associated information about the trajectory: stepsize, acceptance ratio, treedepth, and number of leapfrog steps. This function extracts these into a data.frame, which may be useful for diagnosing issues in certain cases. In general, the user should not need to examine them, or preferably should via plot_sampler_params or launch_shinyadmb.

Value

An invisible data.frame containing samples (rows) of each parameter (columns). If multiple chains exist they will be rbinded together.

See Also

launch_shinyadmb.

Examples

fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
sp <- extract_sampler_params(fit, inc_warmup=TRUE)
str(sp)

Extract posterior samples from a model fit.

Description

A helper function to extract posterior samples across multiple chains into a single data.frame.

Usage

extract_samples(
  fit,
  inc_warmup = FALSE,
  inc_lp = FALSE,
  as.list = FALSE,
  unbounded = FALSE
)

Arguments

Details

This function is loosely based on the rstan function extract. Merging samples across chains should only be used for inference after appropriate diagnostic checks. Do not calculate diagnostics like Rhat or effective sample size after using this function, instead, use monitor. Likewise, warmup samples are not valid and should never be used for inference, but may be useful in some cases for diagnosing issues.

Value

If as.list is FALSE, an invisible data.frame containing samples (rows) of each parameter (columns). If multiple chains exist they will be rbinded together, maintaining order within each chain. If as.list is TRUE, samples are returned as a list of matrices.

Examples

## A previously run fitted ADMB model
fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
post <- extract_samples(fit)
tail(apply(post, 2, median))

Check object of class adfit

Description

Check object of class adfit

Usage

is.adfit(x)

Arguments

Launch shinystan for an ADMB fit.

Description

Launch shinystan for an ADMB fit.

Usage

launch_shinyadmb(fit)

Arguments

See Also

launch_shinytmb

Launch shinystan for a TMB fit.

Description

Launch shinystan for a TMB fit.

Usage

launch_shinytmb(fit)

Arguments

See Also

launch_shinyadmb

Plot pairwise parameter posteriors and optionally the MLE points and confidence ellipses.

Description

Plot pairwise parameter posteriors and optionally the MLE points and confidence ellipses.

Usage

pairs_admb(
  fit,
  order = NULL,
  diag = c("trace", "acf", "hist"),
  acf.ylim = c(-1, 1),
  ymult = NULL,
  axis.col = gray(0.5),
  pars = NULL,
  label.cex = 0.8,
  limits = NULL,
  add.mle = TRUE,
  add.monitor = TRUE,
  unbounded = FALSE,
  ...
)

Arguments

Details

This function is modified from the base pairs code to work specifically with fits from the 'adnuts' package using either the NUTS or RWM MCMC algorithms. If an invertible Hessian was found (in fit$mle) then estimated covariances are available to compare and added automatically (red ellipses). Likewise, a "monitor" object from rstan::monitor is attached as fit$monitor and provides effective sample sizes (ESS) and Rhat values. The ESS are used to potentially order the parameters via argument order, but also printed on the diagonal.

Value

Produces a plot, and returns nothing.

Author(s)

Cole Monnahan

Examples

fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
pairs_admb(fit)
pairs_admb(fit, pars=1:2)
pairs_admb(fit, pars=c('b', 'a'))
pairs_admb(fit, pars=1:2, order='slow')
pairs_admb(fit, pars=1:2, order='fast')

Plot object of class adfit

Description

Plot object of class adfit

Usage

## S3 method for class 'adfit'
plot(x, y, ...)

Arguments

Value

Plot created

Plot marginal distributions for a fitted model

Description

Plot marginal distributions for a fitted model

Usage

plot_marginals(
  fit,
  pars = NULL,
  mfrow = NULL,
  add.mle = TRUE,
  add.monitor = TRUE,
  breaks = 30
)

Arguments

Details

This function plots grid cells of all parameters in a model, comparing the marginal posterior histogram vs the asymptotic normal (red lines) from the inverse Hessian. Its intended use is to quickly gauge differences between frequentist and Bayesian inference on the same model.

If fit$monitor exists the effective sample size (ESS) and R-hat estimates are printed in the top right corner. See https://mc-stan.org/rstan/reference/Rhat.html for more information. Generally Rhat>1.05 or ESS<100 (per chain) suggest inference may be unreliable.

This function is customized to work with multipage PDFs, specifically: pdf('marginals.pdf', onefile=TRUE, width=7,height=5) produces a nice readable file.

Examples

fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
plot_marginals(fit, pars=1:2)

Plot adaptation metrics for a fitted model.

Description

Plot adaptation metrics for a fitted model.

Usage

plot_sampler_params(fit, plot = TRUE)

Arguments

Details

This utility function quickly plots the adaptation output of NUTS chains.

Value

Prints and invisibly returns a ggplot object

Examples

fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
plot_sampler_params(fit)

Plot MLE vs MCMC marginal standard deviations for each parameter

Description

Plot MLE vs MCMC marginal standard deviations for each parameter

Usage

plot_uncertainties(fit, log = TRUE, plot = TRUE)

Arguments

Details

It can be helpful to compare uncertainty estimates between the two paradigms. This plots the marginal posterior standard deviation vs the frequentist standard error estimated from the .cor file. Large differences often indicate issues with one estimation method.

Value

Invisibly returns data.frame with parameter name and estimated uncertainties.

Examples

fit <- readRDS(system.file('examples', 'fit.RDS', package='adnuts'))
x <- plot_uncertainties(fit, plot=FALSE)
head(x)

Print summary of adfit object

Description

Print summary of adfit object

Usage

## S3 method for class 'adfit'
print(x, ...)

Arguments

Value

Summary printed to console

Deprecated version of wrapper function. Use sample_nuts or sample_rwm instead.

Description

Deprecated version of wrapper function. Use sample_nuts or sample_rwm instead.

Usage

sample_admb(
  model,
  path = getwd(),
  iter = 2000,
  init = NULL,
  chains = 3,
  warmup = NULL,
  seeds = NULL,
  thin = 1,
  mceval = FALSE,
  duration = NULL,
  parallel = FALSE,
  cores = NULL,
  control = NULL,
  skip_optimization = TRUE,
  algorithm = "NUTS",
  skip_monitor = FALSE,
  skip_unbounded = TRUE,
  admb_args = NULL
)

Arguments

Warning

This is deprecated and will cease to exist in future releases

Function to generate random initial values from a previous fit using adnuts

Description

Function to generate random initial values from a previous fit using adnuts

Usage

sample_inits(fit, chains)

Arguments

Value

A list of lists which can be passed back into sample_admb.

Bayesian inference of an ADMB model using the no-U-turn sampler (NUTS) or random walk Metropolis (RWM) algorithms.

Description

Draw Bayesian posterior samples from an AD Model Builder (ADMB) model using an MCMC algorithm. 'sample_nuts' and 'sample_rwm' generates posterior samples from which inference can be made.

Usage

sample_nuts(
  model,
  path = getwd(),
  iter = 2000,
  init = NULL,
  chains = 3,
  warmup = NULL,
  seeds = NULL,
  thin = 1,
  mceval = FALSE,
  duration = NULL,
  parallel = FALSE,
  cores = NULL,
  control = NULL,
  skip_optimization = TRUE,
  verbose = TRUE,
  skip_monitor = FALSE,
  skip_unbounded = TRUE,
  admb_args = NULL,
  extra.args = NULL
)

sample_rwm(
  model,
  path = getwd(),
  iter = 2000,
  init = NULL,
  chains = 3,
  warmup = NULL,
  seeds = NULL,
  thin = 1,
  mceval = FALSE,
  duration = NULL,
  parallel = FALSE,
  cores = NULL,
  control = NULL,
  skip_optimization = TRUE,
  verbose = TRUE,
  skip_monitor = FALSE,
  skip_unbounded = TRUE,
  admb_args = NULL,
  extra.args = NULL
)

Arguments

Details

Adaptation schemes are used with NUTS so specifying tuning parameters is not necessary. See vignette for options for adaptation of step size and mass matrix. The RWM algorithm provides no new functionality not available from previous versions of ADMB. However, 'sample_rwm' has an improved console output, is setup for parallel execution, and a smooth workflow for diagnostics.

Parallel chains will be run if argument 'cores' is greater than one. This entails copying the folder, and starting a new R session to run that chain, which are then merged back together. Note that console output is inconsistent when using parallel, and may not show. On Windows the R terminal shows output live, but the GUI does not. RStudio is a special case and will not show live, and instead is captured and returned at the end. It is strongly recommended to start with serial execution as debugging parallel chains is very difficult.

Note that the algorithm code is in the ADMB source code, and 'adnuts' provides a wrapper for it. The command line arguments are returned and can be examined by the user. See vignette for more information.

This function implements algorithm 6 of Hoffman and Gelman (2014), and loosely follows package rstan. The step size can be adapted or specified manually. The metric (i.e., mass matrix) can be unit diagonal, adapted diagonal (default and recommended), a dense matrix specified by the user, or an adapted dense matrix. Further control of algorithms can be specified with the control argument. Elements are:

adapt_delta

The target acceptance rate. D

metric

The mass metric to use. Options are: "unit" for a unit diagonal matrix; NULL to estimate a diagonal matrix during warmup; a matrix to be used directly (in untransformed space).

adapt_delta

Whether adaptation of step size is turned on.

adapt_mass

Whether adaptation of mass matrix is turned on. Currently only allowed for diagonal metric.

adapt_mass_dense

Whether dense adaptation of mass matrix is turned on.

max_treedepth

Maximum treedepth for the NUTS algorithm.

stepsize

The stepsize for the NUTS algorithm. If NULL it will be adapted during warmup.

adapt_init_buffer

The initial buffer size during mass matrix adaptation where sample information is not used (default 50)

adapt_term_buffer

The terminal buffer size (default 75) during mass matrix adaptation (final fast phase)

adapt_window

The initial size of the mass matrix adaptation window, which gets doubled each time thereafter.

refresh

The rate at which to refresh progress to the console. Defaults to even 10 progress updates.

The adaptation scheme (step size and mass matrix) is based heavily on those by the software Stan, and more details can be found in that documentation and this vignette.

Warning

The user is responsible for specifying the model properly (priors, starting values, desired parameters fixed, etc.), as well as assessing the convergence and validity of the resulting samples (e.g., through the coda package), or with function launch_shinytmb before making inference. Specifically, priors must be specified in the template file for each parameter. Unspecified priors will be implicitly uniform.

Author(s)

Cole Monnahan

Examples

## Not run: 
## This is the packaged simple regression model
path.simple <- system.file('examples', 'simple', package='adnuts')
## It is best to have your ADMB files in a separate folder and provide that
## path, so make a copy of the model folder locally.
path <- 'simple'
dir.create(path)
trash <- file.copy(from=list.files(path.simple, full.names=TRUE), to=path)
## Compile and run model
oldwd <- getwd()
setwd(path)
system('admb simple.tpl')
system('simple')
setwd('..')
init <- function() rnorm(2)
## Run NUTS with defaults
fit <- sample_nuts(model='simple', init=init, path=path)
unlink(path, TRUE) # cleanup folder
setwd(oldwd)

## End(Not run)

Bayesian inference of a TMB model using the no-U-turn sampler.

Description

Draw Bayesian posterior samples from a Template Model Builder (TMB) model using an MCMC algorithm. This function generates posterior samples from which inference can be made. Adaptation schemes are used so specification tuning parameters are not necessary, and parallel execution reduces overall run time.

Usage

sample_tmb(
  obj,
  iter = 2000,
  init,
  chains = 3,
  seeds = NULL,
  warmup = floor(iter/2),
  lower = NULL,
  upper = NULL,
  thin = 1,
  parallel = FALSE,
  cores = NULL,
  path = NULL,
  algorithm = "NUTS",
  laplace = FALSE,
  control = NULL,
  ...
)

Arguments

Details

This function implements algorithm 6 of Hoffman and Gelman (2014), and loosely follows package rstan. The step size can be adapted or specified manually. The metric (i.e., mass matrix) can be unit diagonal, adapted diagonal (default and recommended), or a dense matrix specified by the user. Further control of algorithms can be specified with the control argument. Elements are:

adapt_delta

The target acceptance rate.

metric

The mass metric to use. Options are: "unit" for a unit diagonal matrix; "diag" to estimate a diagonal matrix during warmup; a matrix to be used directly (in untransformed space).

adapt_engaged

Whether adaptation of step size and metric is turned on.

max_treedepth

Maximum treedepth for the NUTS algorithm.

stepsize

The stepsize for the NUTS algorithm. If NULL it will be adapted during warmup.

Value

A list containing the samples, and properties of the sampler useful for diagnosing behavior and efficiency.

Warning

This is deprecated and will cease to exist in future releases

Author(s)

Cole Monnahan

See Also

extract_samples to extract samples and launch_shinytmb to explore the results graphically which is a wrapper for the launch_shinystan function.

Examples

## Build a fake TMB object with objective & gradient functions and some
## other flags
## Not run: 
f <- function(x, order=0){
  if(order != 1) # negative log density
    -sum(dnorm(x=x, mean=0, sd=1, log=TRUE))
  else x # gradient of negative log density
}
init <- function() rnorm(2)
obj <- list(env=list(DLL='demo', last.par.best=c(x=init()), f=f,
  beSilent=function() NULL))
## Run NUTS for this object
fit <- sample_tmb(obj, iter=1000, chains=3, init=init)
## Check basic diagnostics
mon <- rstan::monitor(fit$samples, print=FALSE)
Rhat <- mon[,"Rhat"]
max(Rhat)
ess <- mon[, 'n_eff']
min(ess)
## Or do it interactively with ShinyStan
launch_shinytmb(fit)

## End(Not run)

Draw MCMC samples from a model posterior using a static HMC sampler.

Description

Draw MCMC samples from a model posterior using a static HMC sampler.

Usage

sample_tmb_hmc(
  iter,
  fn,
  gr,
  init,
  L,
  eps,
  warmup = floor(iter/2),
  seed = NULL,
  chain = 1,
  thin = 1,
  control = NULL
)

Arguments

Details

This function implements algorithm 5 of Hoffman and Gelman (2014), which includes adaptive step sizes (eps) via an algorithm called dual averaging.

Value

A list containing samples ('par') and algorithm details such as step size adaptation and acceptance probabilities per iteration ('sampler_params').

References

• Neal, R. M. (2011). MCMC using Hamiltonian dynamics. Handbook of Markov Chain Monte Carlo.

• Hoffman and Gelman (2014). The No-U-Turn sampler: Adaptively setting path lengths in Hamiltonian Monte Carlo. J. Mach. Learn. Res. 15:1593-1623.

Hoffman and Gelman (2014). The No-U-Turn sampler: Adaptively setting path lengths in Hamiltonian Monte Carlo. J. Mach. Learn. Res. 15:1593-1623.

See Also

sample_tmb

sample_tmb

Draw MCMC samples from a model posterior using the No-U-Turn (NUTS) sampler with dual averaging.

Description

Draw MCMC samples from a model posterior using the No-U-Turn (NUTS) sampler with dual averaging.

Usage

sample_tmb_nuts(
  iter,
  fn,
  gr,
  init,
  warmup = floor(iter/2),
  chain = 1,
  thin = 1,
  seed = NULL,
  control = NULL
)

Arguments

Details

This function implements algorithm 6 of Hoffman and Gelman (2014), which includes adaptive step sizes (eps) via an algorithm called dual averaging. It also includes an adaptation scheme to tune a diagonal mass matrix (metric) during warmup.

These fn and gr functions must have Jacobians already applied if there are transformations used.

References

Hoffman and Gelman (2014). The No-U-Turn sampler: Adaptively setting path lengths in Hamiltonian Monte Carlo. J. Mach. Learn. Res. 15:1593-1623.

See Also

sample_tmb

[Deprecated] Draw MCMC samples from a model posterior using a Random Walk Metropolis (RWM) sampler.

Description

[Deprecated] Draw MCMC samples from a model posterior using a Random Walk Metropolis (RWM) sampler.

Usage

sample_tmb_rwm(
  iter,
  fn,
  init,
  alpha = 1,
  chain = 1,
  warmup = floor(iter/2),
  thin = 1,
  seed = NULL,
  control = NULL
)

Arguments

Details

This algorithm does not yet contain adaptation of alpha so some trial and error may be required for efficient sampling.

Value

A list containing samples and other metadata.

References

Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H., Teller, E., 1953. Equation of state calculations by fast computing machines. J Chem Phys. 21:1087-1092.

See Also

sample_tmb

Print summary of object of class adfit

Description

Print summary of object of class adfit

Usage

## S3 method for class 'adfit'
summary(object, ...)

Arguments

Value

Summary printed to screen