Last data update: 2014.03.03

R: Monte Carlo Standard Errors for MCMC
mcmcse-packageR Documentation

Monte Carlo Standard Errors for MCMC

Description

Provides tools for computing Monte Carlo standard errors (MCSE) in Markov chain Monte Carlo (MCMC) settings. MCSE computation for expectation and quantile estimators is supported. The package also provides functions for computing effective sample size and for plotting Monte Carlo estimates versus sample size.

Details

Package: mcmcse
Type: Package
Version: 1.1-2
Date: 2015-08-17
License: GPL (>= 2)

Author(s)

James M. Flegal <jflegal@ucr.edu>, John Hughes <hughesj@umn.edu> and Dootika Vats <vatsx007@umn.edu>

Maintainer: James M. Flegal <jflegal@ucr.edu>

References

Flegal, J. M. (2012) Applicability of subsampling bootstrap methods in Markov chain Monte Carlo. In Wozniakowski, H. and Plaskota, L., editors, Monte Carlo and Quasi-Monte Carlo Methods 2010 (to appear). Springer-Verlag.

Flegal, J. M. and Jones, G. L. (2010) Batch means and spectral variance estimators in Markov chain Monte Carlo. The Annals of Statistics, 38, 1034–1070.

Flegal, J. M. and Jones, G. L. (2011) Implementing Markov chain Monte Carlo: Estimating with confidence. In Brooks, S., Gelman, A., Jones, G. L., and Meng, X., editors, Handbook of Markov Chain Monte Carlo, pages 175–197. Chapman & Hall/CRC Press.

Flegal, J. M., Jones, G. L., and Neath, R. (2012) Markov chain Monte Carlo estimation of quantiles. University of California, Riverside, Technical Report.

Gong, L., and Flegal, J. M. A practical sequential stopping rule for high-dimensional Markov chain Monte Carlo. Journal of Computational and Graphical Statistics (to appear).

Jones, G. L., Haran, M., Caffo, B. S. and Neath, R. (2006) Fixed-width output analysis for Markov chain Monte Carlo. Journal of the American Statistical Association, 101, 1537–1547.

Vats, D., Flegal, J. M., and, Jones, G. L Multivariate Output Analysis for Markov chain Monte Carlo, arXiv preprint arXiv:1512.07713 (2015).

Examples

library(mAr)
p <- 3
n <- 1e3
omega <- 5*diag(1,p)

## Making correlation matrix var(1) model
set.seed(100)
foo <- matrix(rnorm(p^2), nrow = p)
foo <- foo %*% t(foo)
phi <- foo / (max(eigen(foo)$values) + 1)
  
out <- as.matrix(mAr.sim(rep(0,p), phi, omega, N = n))

mcse(out[,1], method = "bart")

mcse.bm <- mcse.multi(x = out)
mcse.tuk <- mcse.multi(x = out, method = "tukey")

Results