Package 'superdiag'

Test for Markov Chain Nonconvergence

Description

The superdiag function takes MCMC samples as input. It provides a comprehensive test suite for Markov chain nonconvergence, which integrates five standard empirical MCMC convergence diagnostics: Gelman-Rubin, Geweke, Heidelberger-Welch, Raftery-Lewis, Hellinger distance. It can also present trace plots and density histograms along with the diagnostics as options.

Usage

superdiag(
  mcmcoutput,
  burnin,
  terms = "all",
  plot = FALSE,
  confidence.gr = 0.95,
  frac1.gw = 0.1,
  frac2.gw = 0.5,
  eps.hw = 0.1,
  pvalue.hw = 0.05,
  q.rl = 0.025,
  r.rl = 0.005,
  s.rl = 0.95,
  eps.rl = 0.001,
  bins = 5,
  binwidth = NULL
)

Arguments

mcmcoutput	A mcmc or mcmc.list object.
burnin	The number of burn-in iterations. Defaults to half of the chain(s).
terms	The convergence diagnostic methods. Defaults to all five methods. Users can also specify one or more particular methods chosen from "geweke", "heidel", "raftery", "gelman", "hellinger".
plot	Logical values indicates whether a graph of superdiagPlot should be presented with a run of the function. Defaults to FALSE.
confidence.gr	(1-$\alpha$)% for the Gelman and Rubin test. The upper 95% credible interval is the default.
frac1.gw	frac1 for the Geweke test. The proportion of the early era of the chain, defaulted to 0.1.
frac2.gw	frac2 for the Geweke test. The proportion of the late era of the chain, defaulted to 0.5.
eps.hw	epsilon for the Heidelberger and Welch test. The accuracy parameter determines whether the halfwidth is passed or not, defaulted to 0.1.
pvalue.hw	p-value for the Heidelberger and Welch test. The halfwidth of the test calculates a (1-$\alpha$)% credible interval around the sample mean for each parameter dimension. The default is 0.05.
q.rl	q-parameter for the Raftery and Lewis test. The posterior tail threshold of interest, defaulted to 0.025.
r.rl	r-parameter for the Raftery and Lewis test. The tolerance for the tail threshold, defaulted to 0.0005.
s.rl	s-parameter for the Raftery and Lewis test. The desired probability of being within the tolerance, defaulted to 0.95.
eps.rl	convergence epsilon for the Raftery and Lewis test. The convergence tolerance value, which is used to determine a stopping point based on a parallel chain process, defaulted to 0.001.
bins	Number of bins for within-chain Hellinger distance test. Defaults to 5.
binwidth	Alternative specification for bins. The size of bin for each batch of the chain(s) for computing the distance between batches.

Details

If only one chain is analyzed, the default settings in boa and coda for Geweke, Raftery-Lewis, and Heidelberger-Welch diagnostics are used. If multiple chains are provided, only the first chain uses the defaults and all other chain analyses get random values as follows. For Geweke test, random non-overlapping proportions up from the start of the chain and down from the end of the chain are generated. For Heidelberger-Welch test, the value of pvalue.hw is sampled with replacement from common $\alpha$ values; the value of eps.hw is sampled uniformly in the interval [0.01:0.2]. For Raftery-Lewis test, each of these four parameters are sample from a vector (changeable by users) of values around the defaults (larger and smaller) to provide a reasonable range of alternatives. For Hellinger distance, if only one chain is analyzed, only within-chain distance will be reported.

Value

A list object including the results for all the diagnostics. A superdiagPlot including both traceplot(s) and density histogram will also be returned if plot=TRUE.

References

Tsai, Tsung-han and Gill, Jeff (2012). “superdiag: A Comprehensive Test Suite for Markov Chain Non-Convergence.” The Political Methodologist, 19 (2), 12-18.

Plummer, Martyn, Nicky Best, Kate Cowles, and Karen Vines (2006). "CODA: convergence diagnosis and output analysis for MCMC." R news, 6 (1), 7-11.

See Also

superdiagPlot, gelman.diag, heidel.diag, raftery.diag, gelman.diag

Examples

## Not run: 

data(tobit.list)
summary(tobit.list[1])

# FOR mcmc.list OBJECT
superdiag(tobit.list, burnin=0)

# FOR mcmc OBJECT
tobit.diag <- superdiag(tobit.list[[1]], burnin=0, plot=TRUE)
tobit.diag

## End(Not run)

---

Summary Plots for mcmc Objects

Description

The superdiagPlot function takes mcmc or mcmc.list as input. It provides summary plots including trace plot(s) and density histogram(s) for each variable in the MCMC chain(s).

Usage

superdiagPlot(
  mcmcoutput,
  burnin,
  params,
  title,
  col,
  trace.options,
  density.options,
  title.options
)

## S3 method for class 'mcmc'
superdiagPlot(
  mcmcoutput,
  burnin,
  params = "all",
  title,
  col = c("grey25", "grey75"),
  trace.options = list(trace.axis = TRUE, box = TRUE),
  density.options = list(),
  title.options = list()
)

## S3 method for class 'mcmc.list'
superdiagPlot(
  mcmcoutput,
  burnin,
  params = "all",
  title,
  col,
  trace.options = list(trace.axis = TRUE, box = TRUE),
  density.options = list(),
  title.options = list()
)

Arguments

mcmcoutput	A mcmc or mcmc.list object.
burnin	The number of burn-in iterations. Defaults to half of the chain(s).
params	The parameters to be summarized. Defaults to all. Users can specify a subset of variables by name or index to reduce the size of plots.
title	A character vector specifies the title of the plot.
col	A vector or list specifies the color schemes of the plot. For mcmc object, the default is a two-color grey scheme; Users can specify any numbers of color no more than the number of variables. For mcmc.list object, it requires a list or vector with multiple colors to examine the mixing. The default is a list of two sets of colors.
trace.options	Additional options for trace plots. It can take arguments as in plot and lines functions.
density.options	Additional options for density histograms. It can take arguments as in hist (e.g breaks, freq, probability) and barplot functions.
title.options	Additional options for title. It can take arguments as in mtext function.

Value

A summary plot for each variable in the MCMC chain(s).

See Also

superdiag

Examples

## Not run: 

data(tobit.list)
summary(tobit.list[1])

# FOR mcmc.list OBJECT
superdiagPlot(tobit.list, burnin=0)

# FOR mcmc OBJECT
superdiagPlot(tobit.list[[1]], burnin=0, col=c("grey25", "dodgerblue"),
  title = "Tobit Model (Chain 1)", title.options=list(cex=1.2))

## End(Not run)

---

Data: A mcmc.list from a Tobit model

Description

The data is the output from a Tobit model that was estimated with MCMC techniques. The estimation was performed with five parallel chains of 50,000 iterations each and the first 40,000 iterations were discarded.

Details

The R code for the Gibbs sampler is available at https://jeffgill.org/.

References

Gill, Jeff. (2007). Bayesian Methods for the Social and Behavior Sciences. Second Edition. New York: Chapman & Hall.

---