forked from SandaD/MCMCEnsembleSampler
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
161 lines (118 loc) · 5.03 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
---
output: github_document
bibliography: readme.bib
link-citations: true
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/README-",
out.width = "100%"
)
set.seed(18052020)
```
# mcmcensemble
<!-- badges: start -->
[](https://CRAN.R-project.org/package=mcmcensemble)
[](https://github.com/Bisaloo/mcmcensemble/actions)
[](https://app.codecov.io/gh/Bisaloo/mcmcensemble?branch=main)
[](https://lifecycle.r-lib.org/articles/stages.html)
This R package provides ensemble samplers for affine-invariant Monte Carlo
Markov Chain, which allow a faster convergence for badly scaled estimation
problems. Two samplers are proposed: the 'differential.evolution' sampler from
@terBraak2008 and the 'stretch' sampler from @Goodman2010.
For theoretical background about Ensemble MCMC (what are the benefits over
simple MCMC? How do they work? What are the pitfalls?), please refer for example to
[this lecture](https://doi.org/10.26207/46za-m573) from Eric B. Ford (Penn State).
## Installation
You can install the stable version of this package from
[CRAN](https://cran.r-project.org/package=mcmcensemble):
```{r, eval = FALSE}
install.packages("mcmcensemble")
```
or the development version from [GitHub](https://github.com/bisaloo), via my
[r-universe](https://bisaloo.r-universe.dev/packages):
```{r, eval = FALSE}
install.packages("mcmcensemble", repos = "https://bisaloo.r-universe.dev")
```
## Usage
```{r, include = FALSE}
set.seed(20210220)
```
```{r}
library(mcmcensemble)
## a log-pdf to sample from
p.log <- function(x) {
B <- 0.03 # controls 'bananacity'
-x[1]^2 / 200 - 1/2 * (x[2] + B * x[1]^2 - 100 * B)^2
}
## set options and starting point
n_walkers <- 10
unif_inits <- data.frame(
"a" = runif(n_walkers, 0, 1),
"b" = runif(n_walkers, 0, 1)
)
## use stretch move
res1 <- MCMCEnsemble(p.log, inits = unif_inits,
max.iter = 5000, n.walkers = n_walkers,
method = "stretch")
attr(res1, "ensemble.sampler")
str(res1)
```
If the [coda](https://cran.r-project.org/package=coda) package is installed,
you can then use the `coda = TRUE` argument to get objects of class `mcmc.list`.
The coda package then allows you to call `summary()` and `plot()` to get
informative and nicely formatted results and plots:
```{r example-stretch, dev = 'svglite'}
## use stretch move, return samples as 'coda' object
res2 <- MCMCEnsemble(p.log, inits = unif_inits,
max.iter = 5000, n.walkers = n_walkers,
method = "stretch", coda = TRUE)
attr(res2, "ensemble.sampler")
summary(res2$samples)
plot(res2$samples)
```
```{r example-de, dev = 'svglite'}
## use different evolution move, return samples as 'coda' object
res3 <- MCMCEnsemble(p.log, inits = unif_inits,
max.iter = 5000, n.walkers = n_walkers,
method = "differential.evolution", coda = TRUE)
attr(res3, "ensemble.sampler")
summary(res3$samples)
plot(res3$samples)
```
To see more plotting and MCMC diagnostic options, please refer to the relevant
vignette:
[`vignette("diagnostic-pkgs", package = "mcmcensemble")`](https://hugogruson.fr/mcmcensemble/articles/diagnostic-pkgs.html)
## Progress bar
You can choose to enable a progress bar thanks to the
[progressr](https://cran.r-project.org/package=progressr) package. This can be
done by adding the following line to your script before running
`MCMCEnsemble()`:
```{r, eval = FALSE}
progressr::handlers(global = TRUE) # requires R >= 4.0
progressr::handlers("progress")
MCMCEnsemble(p.log, inits = unif_inits,
max.iter = 5000, n.walkers = n_walkers,
method = "differential.evolution", coda = TRUE)
```
## Parallel processing
This package is set up to allow transparent parallel processing when requested
by the user thanks to the framework provided by the
[future](https://cran.r-project.org/package=future) package. To enable parallel
processing, you must run:
```{r, eval = FALSE}
future::plan("multiprocess")
```
at the start of your session.
## Similar projects
The Goodman-Weare 'stretch' sampler is also available in the [tonic R package](https://github.com/SimonVaughanDataAndCode/tonic).
The methods used in this package also have (independent) implementations in
other languages:
- [emcee v3: A Python ensemble sampling toolkit for affine-invariant MCMC](https://doi.org/10.21105/joss.01864)
- [GWMCMC which implements the Goodman-Weare 'stretch' sampler in Matlab](https://github.com/grinsted/gwmcmc)
## Who is talking about this package?
- [R View from October 2020](https://rviews.rstudio.com/2020/11/19/october-2020-top-40-new-cran-packages/)
## References