mcmcr

Lifecycle: maturing Travis-CI Build Status AppVeyor build status Coverage Status License: MIT status CRAN_Status_Badge CRAN Downloads

Introduction

mcmcr is an R package to manipulate Monte Carlo Markov Chain (MCMC) samples (Brooks et al. 2011).

For the purposes of this discussion, an MCMC sample represents the value of a term from a single iteration of a single chain. While a simple parameter such as an intercept corresponds to a single term, more complex parameters such as an interaction between two factors consists of multiple terms with their own inherent dimensionality - in this case a matrix. A set of MCMC samples can be stored in different ways.

Existing Classes

The three most common S3 classes store MCMC samples as follows:

In the first two cases the terms/parameters are represented by a single dimension which means that the dimensionality inherent in the parameters is stored in the labelling of the variables, ie, "bIntercept", "bInteraction[1,2]", "bInteraction[2,1]", .... The structure of the mcmc and mcmc.list objects emphasizes the time-series nature of MCMC samples and is optimized for thining. In contrast mcarray objects preserve the dimensionality of the parameters.

New Classes

The mcmcr package defines three related S3 classes which also preserve the dimensionality of the parameters:

Why mcmcr?

mcmcarray objects were developed to facilitate manipulation of the MCMC samples (although they are just one aperm away from mcarray objects they are more intuitive to program with - at least for this programmer!). mcmcr objects were developed to allow a set of dimensionality preserving parameters from a single analysis to be manipulated as a whole. mcmcrs objects were developed to allow the results of multiple analyses using the same model to be manipulated together.

In addition the mcmcr package defines the term vector to store and manipulate the term labels, ie, "bIntercept", "bInteraction[1,2]", "bInteraction[2,1]", when the MCMC samples are summarised in tabular form.

The mcmcr package also introduces a variety of (often) generic functions to manipulate and query mcmcarray, mcmcr and mcmcrs objects. In particular it provides functions to

The code is opinionated which has the advantage of providing a small set of stream-lined functions. For example the only ‘convergence’ metric is the uncorrected, untransformed, univariate split R-hat (potential scale reduction factor). If you can convince me that additional features are important I will add them or accept a pull request (see below). Alternatively you might want to use the mcmcr package to manipulate your samples before coercing them to an mcmc.list to take advantage of all the summary functions in packages such as coda.

Demonstration

library(mcmcr)

mcmcr_example
#> $alpha
#> [1] 3.718025 4.718025
#> 
#> nchains:  2 
#> niters:  400 
#> 
#> $beta
#>           [,1]     [,2]
#> [1,] 0.9716535 1.971654
#> [2,] 1.9716535 2.971654
#> 
#> nchains:  2 
#> niters:  400 
#> 
#> $sigma
#> [1] 0.7911975
#> 
#> nchains:  2 
#> niters:  400

coef(mcmcr_example)
#>        term  estimate        sd   zscore     lower    upper pvalue
#> 1  alpha[1] 3.7180250 0.9007167 4.149545 2.2120540 5.232403 0.0012
#> 2  alpha[2] 4.7180250 0.9007167 5.259772 3.2120540 6.232403 0.0012
#> 3 beta[1,1] 0.9716535 0.3747971 2.572555 0.2514796 1.713996 0.0225
#> 4 beta[2,1] 1.9716535 0.3747971 5.240666 1.2514796 2.713996 0.0050
#> 5 beta[1,2] 1.9716535 0.3747971 5.240666 1.2514796 2.713996 0.0050
#> 6 beta[2,2] 2.9716535 0.3747971 7.908776 2.2514796 3.713996 0.0012
#> 7     sigma 0.7911975 0.7408373 1.306700 0.4249618 2.559520 0.0012
rhat(mcmcr_example, by = "term")
#> $alpha
#> [1] 2.002 2.002
#> 
#> $beta
#>       [,1]  [,2]
#> [1,] 1.147 1.147
#> [2,] 1.147 1.147
#> 
#> $sigma
#> [1] 1
plot(mcmcr_example[["alpha"]])

Installation

To install the latest release version from CRAN

install.packages("mcmcr")

To install the latest development version from GitHub

if(!"remotes" %in% installed.packages()[,1]) 
  install.packages("remotes")
remotes::install_github("poissonconsulting/mcmcr")

To install the latest development version from the Poisson drat repository

if(!"drat" %in% installed.packages()[,1]) 
  install.packages("drat")
drat::addRepo("poissonconsulting")
install.packages("mcmcr")

Contribution

Please report any issues.

Pull requests are always welcome.

Please note that this project is released with a Contributor Code of Conduct. By participating in this project you agree to abide by its terms.

Inspiration

coda and rjags

References

Brooks, S., Gelman, A., Jones, G.L., and Meng, X.-L. (Editors). 2011. Handbook for Markov Chain Monte Carlo. Taylor & Francis, Boca Raton. ISBN: 978-1-4200-7941-8.