Skip to contents

Computes common convergence diagnostics for selected parameters from a JAGS/BUGS fit and returns a compact, report-ready table. The diagnostics include Gelman–Rubin \(\hat{R}\), Geweke z-scores, Heidelberger-Welch stationarity p-values, and autocorrelation at lag 50.

Usage

mcmcDiag(bml.out, parameters)

Arguments

bml.out

A model fit object containing JAGS output, typically as returned by R2jags::jags(), with component $jags.out$BUGSoutput.

parameters

Character vector of parameter names (or patterns) to extract, passed to crMCMC. Depending on your crMCMC() implementation, these may be exact names or patterns (e.g., a prefix like "b.l2" that matches "b.l2[1]", "b.l2[2]", …).

Value

A data.frame with one row per diagnostic and one column per parameter; cell entries are the average diagnostic values across chains. Row names include: "Gelman/Rubin convergence statistic", "Geweke z-score", "Heidelberger/Welch p-value", "Autocorrelation (lag 50)".

Details

Internally, the function converts the BUGS/JAGS output to a coda::mcmc.list via crMCMC, then computes per-chain diagnostics and averages them across chains for each parameter:

  • Gelman–Rubin (\(\hat{R}\)): coda::gelman.diag(). Values close to 1 indicate convergence; a common heuristic is \(\hat{R} \le 1.1\).

  • Geweke z-score: coda::geweke.diag(). Large absolute values (e.g., \(|z|>2\)) suggest lack of convergence.

  • Heidelberger–Welch p-value: coda::heidel.diag() tests the null of stationarity in the chain segment.

  • Autocorrelation (lag 50): coda::autocorr() at lag 50, averaged across chains.

All statistics are rounded to three decimals. The returned table is transposed so that rows are diagnostics and columns are parameters.

Author

Benjamin Rosche benrosche@nyu.edu

Examples

if (FALSE) { # \dontrun{
mcmcdiag(fit, parameters = c("beta", "sigma"))
} # }