An article by Suilou Huang for catatrophe modeler AIR-WorldWide of Boston about rejection sampling in CAT modeling got me thinking about pulling together some notes about sampling algorithms of various kinds. There are, of course, books written about this subject, including Yves Tillé’s Sampling Algorithms, 2006, Springer, which covers reservoir sampling in its section 4.4.5. Tillé does not cover rejection sampling or slice sampling, so I thought I’d put some notes and figures here, not so much to teach, but reference. And, to some degree for balance.

Note the Wikipedia article regarding slice sampling is, at least of this date and in my opinion, not up to Wikipedia‘s usual quality standards for statistical articles.

This is in large measure a technical blog, although I have a lot of material about climate activism and promotion of solar energy. Sure, I very much care about these things, but I’m at heart a quantitative problem solver, with a fondness for engineering quantitative software, and, so, this blog doesn’t reflect me fully if I leave it at activism.

The theoretical underpinnings of many of these algorithms is recounted in C. P. Robert and G. Casella’s textbook Monte Carlo Statistical Methods, 2^en edition, 2004. Slice sampling receives an extended and worthwhile treatment in their Chapter 8. Robert and Casella also treat algorithms like importance sampling, not touched here, in their section 3.3, and rejection sampling in their section 2.3.

Huang covers rejection sampling pretty well. Fewer people know about slice sampling, so the remaining references will cover that.

The idea was introduced by Professor Radford Neal of the University of Toronto in 2003 (also see accompanying discussion) and it has undergone several innovations since. He offers R code to do it. There is an R package, one MfUSampler, which offers various sampling algorithms under a supervisory framework, including slice sampling. These are intended to be used in the context of a Markov Chain Monte Carlo (MCMC) stochastic search, typically exploring a Bayesian posterior density, but the cautions of Robert and Casella regarding adaptive sampling schemes in their section 7.6.3 for this purpose, which they credit Neal for identifying, are worth a look.

There are also at least a couple of instances of code for Python, including nice tutorials by Kristiadi, and Professor Rahul Dave’s AM207 course notes from 2017.

Update, 2018-09-29, 14:44 EDT

Additional study revealed

M. M. Tibbits, M. Haran, J. C. Liechty, ``Parallel multivariate slice sampling'', Statistics and Computing, 2011, 21(13), 415-430.

which, in turn, traces the origins of this idea back to

R. M. Neal, ``Markov chain Monte Carlo methods based on 'slicing' the density function'', Technical Report, Department of Statistics, University of Toronto, 1997.

and

P. Damien, J. Wakefield, S. Walker, ``Gibbs sampling for Bayesian non-conjugate and hierarchical models by using auxiliary variables'', Journal of the Royal Statistical Society, Series B (Statistical Methods), 1999, 61(2), 331-344.

and

A. Mira, L. Tierney, ``Efficiency and convergence properties of slice samplers'', Scandinavian Journal of Statististics, 2002, 29(1), 1-12..

There is also the doctoral dissertation of M. B. Thompson, “Slice Sampling with multivariate steps”, from the Graduate Department of Statistics, University of Toronto, 2011, under the supervision of Professor Radford Neal.

Also, note the R package MCMCpack uses slice sampling at selected points for implementations.

Update, 2018-09-30, 16:54 EDT

A more careful read of Professor Rahul’s page on slice sampling reveals the claim:

One of the methods that we’ll look at is Slice Sampling. Slice sampling was invented by Radfor Neal and John Skilling. It is a very elegant and simple way to solve the problems of MH and Gibbs sampling. As Pavlos would say, its conceptually so easy you would think his (very capable) grandma invented it …

(Emphasis added.)

First — and I’m admittedly nitpicking — Professor Neal’s first name is “Radford” not “Radfor”.

Second, and more seriously, I don’t see any evidence to justify the claim that Dr John Skilling co-invented slice sampling. Skilling did invent another kind of sampling, nested sampling, which has been useful in some applications, even though it has received serious criticism from some experts, technically addressed at

N. Chopin, C. P. Robert, ``Properties of nested sampling''. Biometrika, 2010, 97(3), 741-755.

E. Higson, W. Handley, M. Hobson, A. Lasenby, ``Sampling errors in Nested Sampling parameter estimation'', Bayesian Analysis, 2018, 13(3), 873-896.

Skilling and MacKay did supply a comment on Neal’s original paper in its Discussion, proposing to use integer arithmetic for a special version of slice sampling. Neal addresses that in his Rejoinder.

Third, the discussion surrounding the Python code Professor Rahul offers is a bit weak, particularly regarding the multimodal case. This is unfortunate, because that’s the interesting one. This suggests I’m talking myself into doing a tutorial on slice sampling in R some time down the road, remedying these problems. When I do so, I should really include finding the area of the Batman shape as previously promised.