Markov chain Monte Carlo methods are a class of algorithms for sampling from probability distributions based on constructing a Markov chain that has the desired distribution as its equilibrium distribution. The state of the chain after a large number of steps is then used as a sample of the...
In mathematics and computer science, an algorithm is an effective method expressed as a finite list of well-defined instructions for calculating a function. Algorithms are used for calculation, data processing, and automated reasoning...
Pseudo-random number sampling or non-uniform pseudo-random variate generation is the numerical practice of generating pseudo-random numbers that are distributed according to a given probability distribution....
, i.e. for drawing random samples from a statistical distribution. The method is based on the observation that to sample a random variable
Random variable
In probability and statistics, a random variable or stochastic variable is, roughly speaking, a variable whose value results from a measurement on some type of random process. Formally, it is a function from a probability space, typically to the real numbers, which is measurable functionmeasurable...
one can sample uniformly from the region under the graph of its density function.
To visualize this motivation, imagine printing out a simple bell curve and throwing darts at it. Assume that the darts are uniformly distributed around the board. Now take off all of the darts that are outside the curve. The x-positions of the remaining darts will be distributed according to the bell curve. This is because there is the most room for the darts to land where curve is highest and thus the probability density is greatest. Note that this motivating darts example tells us why it is reasonable to sample in this way, but it doesn't describe how slice sampling is actually performed (the darts are actually just performing rejection sampling
Rejection sampling
In mathematics, rejection sampling is a basic pseudo-random number sampling technique used to generate observations from a distribution. It is also commonly called the acceptance-rejection method or "accept-reject algorithm"....
In probability theory, a probability density function , or density of a continuous random variable is a function that describes the relative likelihood for this random variable to occur at a given point. The probability for the random variable to fall within a particular region is given by the...
Computational statistics, or statistical computing, is the interface between statistics and computer science. It is the area of computational science specific to the mathematical science of statistics....
.
Implementation
Slice sampling gets its name from the first step: defining a slice by sampling from an auxiliary variable . This variable is sampled from , where is either the probability density function (pdf) of X or is at least proportional to its pdf. This defines a slice of X where . In other words, we are now looking at a region of X where the probability density is at least . Then the next value of X is sampled uniformly from this slice. A new value of is sampled, then X, and so on. This can be visualized as alternatively sampling the y-position and then the x-position of points under pdf, thus the Xs are from the desired distribution. The values have no particular consequences or interpretations outside of their usefulness for the procedure.
If both the pdf and its inverse are available, and the distribution is unimodal, then finding the slice and sampling from it are simple. If not, a stepping-out procedure can be used to find a region whose endpoints fall outside the slice. Then, a sample can be drawn from the slice using rejection sampling
Rejection sampling
In mathematics, rejection sampling is a basic pseudo-random number sampling technique used to generate observations from a distribution. It is also commonly called the acceptance-rejection method or "accept-reject algorithm"....
. Various procedures for this are described in detail by Neal.
Note that, in contrast to many available methods for generating random numbers from non-uniform distributions, random variates generated directly by this approach will exhibit serial statistical dependence. This is because to draw the next sample, we define the slice based on the value of f(x) for the current sample.
Motivation
Suppose you want to sample some random variable X with distribution f(x). Suppose that the following is the graph of f(x). The height of f(x) corresponds to the likelihood at that point.
If you were to uniformly sample X, each value would have the same likelihood of being sampled, and your distribution would be of the form f(x)=y for some y value instead of some non-uniform function f(x). Instead of the original black line, your new distribution would look more like the blue line.
In order to sample X in a manner which will retain the distribution f(x), some sampling technique must be used which takes into account the varied likelihoods for each range of f(x).
Compared to Other Methods
Slice sampling is a Markov chain method and as such serves the same purpose as Gibbs sampling and Metropolis. Unlike Metropolis, there is no need to manually tune the candidate function or candidate standard deviation.
Recall that Metropolis is sensitive to step size. If the step size is too small random walk causes slow decorrelation. If the step size is too large there is great inefficiency due to a high rejection rate.
In contrast to Metropolis, slice sampling automatically adjusts the step size to match the local shape of the density function. Implementation is arguably easier and more efficient than Gibbs sampling or simple Metropolis updates.
Note that, in contrast to many available methods for generating random numbers from non-uniform distributions, random variates generated directly by this approach will exhibit serial statistical dependence. In other words, not all points have the same independent likelihood of selection.
Slice Sampling requires that the distribution to be sampled be evaluable. One way to relax this requirement is to substitute an evaluable distribution which is proportional to the true unevaluable distribution.
Univariate Case
To sample a random variable X with density f(x) we introduce an auxiliary variable Y and iterate as follows:
Given a sample x we choose y uniformly at random from the interval [0, f(x)];
given y we choose x uniformly at random from the set .
The sample of x is obtained by ignoring the y values.
Our auxiliary variable Y represents a horizontal "slice" of the distribution. The rest of each iteration is dedicated to sampling an x value from the slice which is representative of the density of the region being considered.
In practice, sampling from a horizontal slice of a multimodal distribution is difficult. There is a tension between obtaining a large sampling region and thereby making possible large moves in the distribution space, and obtaining a simpler sampling region to increase efficiency. One option for simplifying this process is regional expansion and contraction.
First, a width parameter w is used to define the area containing the given x value. Each endpoint of this area is tested to see if it lies outside the given slice. If not, the region is extended in the appropriate direction(s) by w until the end both endpoints lie outside the slice.
A candidate sample is selected uniformly from within this region. If the candidate sample lies inside of the slice, then it is accepted as the new sample. If it lies outside of the slice, the candidate point becomes the new boundary for the region. A new candidate sample is taken uniformly. The process repeats until the candidate sample is within the slice. (See diagram for a visual example).
Treating each variable independently
Single variable slice sampling can be used in the multivariate case by sampling each variable in turn repeatedly, as in Gibbs sampling. To do so requires that we can compute, for each component a function that is proportional to .
To prevent random walk behavior, overrelaxation methods can be used to update each variable in turn. Overrelaxation chooses a new value on the opposite side of the mode from the current value, as opposed to choosing a new independent value from the distribution as done in Gibbs.
Hyperrectangle slice sampling
This method adapts the univariate algorithm to the multivariate case by substituting a hyperrectangle for the one-dimensional w region used in the original. The hyperrectangle H is initialized to a random position over the slice. H is then shrunken as points from it are rejected.
Reflective slice sampling
Reflective slice sampling is a technique to suppress random walk behavior in which the successive candidate samples of distribution f(x) are kept within the bounds of the slice by "reflecting" the direction of sampling inward toward the slice once the boundary has been hit.
In this graphical representation of reflective sampling, the shape indicates the bounds of a sampling slice. The dots indicate start and stopping points of a sampling walk. When the samples hit the bounds of the slice, the direction of sampling is "reflected" back into the slice.
Example
Let us consider a single variable example. Suppose our true distribution . So:
We first draw a uniform random value y from the range of f(x) in order to define our slice(es). Suppose y=0.01.
Next, we set our width parameter w which we will use to expand our region of consideration. Suppose w=2.
Next, we need an initial value for x. We draw x from the uniform distribution within the domain of f(x) at our current y. Suppose x=2.
Because x=2 and w=2, our current region of interest is bounded by (1,3).
Now, each endpoint of this area is tested to see if it lies outside the given slice. Our right bound lies outside our slice, but the left value does not. We expand the left bound by adding w to it until it extends past the limit of the slice. After this process, the new bounds of our region of interest are (-4,3).
Next, we take a uniform sample within (-4,3). Suppose this sample yields x=-3.9. Though this sample is within our region of interest, it does not lie within our slice, so we modify the left bound of our region of interest to this point. Now we take a uniform sample from (-3.9,3). This time our sample yields x=1, which is within our slice, and thus is our accepted sample. Had our new x not been within our slice, we would continue the shrinking/resampling process until a valid x within bounds is found.
Another Example
To sample from the normal distribution we first choose an initial x -- say 0. After each sample of x we choose y uniformly at random from , which is bounded the pdf of . After each y sample we choose x uniformly at random from where . This is the slice where .
. This variable is sampled from 
, where 
is either the probability density function (pdf) of X or is at least proportional to its pdf. This defines a slice of X where 
. In other words, we are now looking at a region of X where the probability density is at least 
. Then the next value of X is sampled uniformly from this slice. A new value of 
is sampled, then X, and so on. This can be visualized as alternatively sampling the y-position and then the x-position of points under pdf, thus the Xs are from the desired distribution. The 
values have no particular consequences or interpretations outside of their usefulness for the procedure.
a function that is proportional to 
.
. So:
we first choose an initial x -- say 0. After each sample of x we choose y uniformly at random from 
, which is bounded the pdf of 
. After each y sample we choose x uniformly at random from 
where 
. This is the slice where 
.
.