Rejection sampling
In numerical analysis and computational statistics, rejection sampling is a basic technique used to generate observations from a distribution. It is also commonly called the acceptance-rejection method or "accept-reject algorithm" and is a type of exact simulation method. The method works for any distribution in formula_1 with a density. Rejection sampling is based on the observation that to sample a random variable in one dimension, one can perform a uniformly random sampling of the two-dimensional Cartesian graph, and keep the samples in the region under the graph of its density function. Note that this property  can be extended to "N"-dimension functions. Description. To visualize the motivation behind rejection sampling, imagine graphing the density function of a random variable onto a large rectangular board and throwing darts at it. Assume that the darts are uniformly distributed around the board. Now remove all of the darts that are outside the area under the curve. The remaining darts will be distributed uniformly within the area under the curve, and the x-positions of these darts will be distributed according to the random variable's density. This is because there is the most room for the darts to land  where the curve is highest and thus the probability density is greatest. The visualization as just described is equivalent to a particular form of rejection sampling where the "proposal distribution" is uniform (hence its graph is a rectangle). The general form of rejection sampling assumes that the board is not necessarily rectangular but is shaped according to the density of some proposal distribution that we know how to sample from (for example, using inversion sampling), and which is at least as high at every point as the distribution we want to sample from, so that the former completely encloses the  latter. (Otherwise, there would be parts of the curved area we want to sample from that could never be reached.) Rejection sampling works as follows: This algorithm can be used to sample from the area under any curve, regardless of whether the function integrates to 1. In fact, scaling a function by a constant has no effect on the sampled x-positions. Thus, the algorithm can be used to sample from a distribution whose normalizing constant is unknown, which is common in computational statistics. Theory. The rejection sampling method generates sampling values from a target distribution formula_2 with arbitrary probability  density function formula_3 by using a proposal distribution formula_4 with probability density formula_5. The idea is that one can generate a sample value from formula_2 by instead sampling from formula_4 and accepting the sample from formula_4 with probability formula_9, repeating the draws from formula_4 until a value is accepted. formula_11 here is a constant, finite bound on the likelihood ratio formula_12, satisfying formula_13 over the support of formula_2; in other words, M must satisfy formula_15 for all values of formula_16. Note that this requires that the support of formula_4 must include the support of formula_2—in other words, formula_19 whenever  formula_20. The validation of this method is the envelope principle: when simulating the pair formula_21, one produces a uniform simulation over the subgraph of formula_22. Accepting only pairs such that formula_23 then produces pairs formula_24 uniformly distributed over the subgraph of formula_3 and thus, marginally, a simulation from formula_26 This means that, with enough replicates, the algorithm generates a sample from the desired distribution formula_3. There are a number of extensions to this algorithm, such as the Metropolis algorithm. This method relates to the general field of Monte Carlo techniques, including Markov chain Monte Carlo algorithms that also use a  proxy distribution to achieve simulation from the target distribution formula_3. It forms the basis for algorithms such as the Metropolis algorithm. The unconditional acceptance probability is the proportion of proposed samples which are accepted, which is formula_29where formula_30, and the value of formula_31 each time is generated under the density function formula_32 of the proposal distribution formula_4. The number of samples required from formula_4 to obtain an accepted value thus follows a geometric distribution with probability formula_35, which has mean formula_11. Intuitively, formula_11 is the expected number of the iterations that are needed, as a measure of the computational complexity  of the algorithm. Rewrite the above equation, formula_38 Note that formula_39, due to the above formula, where formula_40 is a probability which can only take values in the interval formula_41. When formula_11 is chosen closer to one, the unconditional acceptance probability is higher the less that ratio varies, since formula_11 is the upper bound for the likelihood ratio formula_44. In practice, a value of formula_11 closer to 1 is preferred as it implies fewer rejected samples, on average, and thus fewer iterations of the algorithm. In this sense, one prefers to have formula_11 as small as possible (while still satisfying  formula_15, which suggests that formula_5 should generally resemble formula_3 in some way. Note, however, that formula_11 cannot be equal to 1: such would imply that formula_51, i.e. that the target and proposal distributions are actually the same distribution. Rejection sampling is most often used in cases where the form of formula_3 makes sampling difficult. A single iteration of the rejection algorithm requires sampling from the proposal distribution, drawing from a uniform distribution, and evaluating the formula_9 expression. Rejection sampling is thus more efficient than some other method whenever M times the cost of these operations—which is the expected  cost of obtaining a sample with rejection sampling—is lower than the cost of obtaining a sample using the other method. Algorithm. The algorithm (used by John von Neumann and dating back to Buffon and his needle) to obtain a sample from distribution formula_2 with density formula_55 using samples from distribution formula_4 with density formula_57 is as follows: The algorithm will take an average of formula_11 iterations to obtain a sample. Advantages over sampling using naive methods. Rejection sampling can be far more efficient compared with the naive methods in some situations. For example, given a problem as sampling formula_67  conditionally on formula_2 given the set formula_69, i.e., formula_70, sometimes formula_71 can be easily simulated, using the naive methods (e.g. by inverse transform sampling): The problem is this sampling can be difficult and inefficient, if formula_75. The expected number of iterations would be formula_76, which could be close to infinity. Moreover, even when you apply the Rejection sampling method, it is always hard to optimize the bound formula_11 for the likelihood ratio. More often than not, formula_11 is large and the rejection rate is high, the algorithm can be very inefficient. The Natural Exponential Family (if it exists), also known  as exponential tilting, provides a class of proposal distributions that can lower the computation complexity, the value of formula_11 and speed up the computations (see examples: working with Natural Exponential Families). Examples: working with natural exponential families. Given a random variable formula_80, formula_81 is the target distribution. Assume for the simplicity, the density function can be explicitly written as formula_3. Choose the proposal as formula_83where formula_84 and formula_85. Clearly, formula_86, is from a natural exponential family. Moreover, the likelihood ratio is formula_87Note that formula_88 implies that it is indeed a log moment-generation function, that is, formula_89. And it is  easy to derive the log moment-generation function of the proposal and therefore the proposal's moments.formula_90 As a simple example, suppose under formula_91, formula_92, with formula_93. The goal is to sample formula_94, formula_95. The analysis goes as follows. For the above example, as the measurement of the efficiency, the expected number of the iterations the NEF-Based Rejection sampling method is of order b, that is formula_113, while under the naive method, the expected number of the iterations is formula_114, which is far more inefficient. In general, exponential tilting, a parametric class of proposal distribution, solves the optimization problems  conveniently, with its useful properties that directly characterize the distribution of the proposal. For this type of problem, to simulate formula_2 conditionally on formula_116, among the class of simple distributions, the trick is to use NEFs, which helps to gain some control over the complexity and considerably speed up the computation. Indeed, there are deep mathematical reasons for using NEFs. Drawbacks. Rejection sampling can lead to a lot of unwanted samples being taken if the function being sampled is highly concentrated in a certain region, for example a function that has a spike at some location. For many distributions, this  problem can be solved using an adaptive extension (see adaptive rejection sampling), or with an appropriate change of variables with the method of the ratio of uniforms. In addition, as the dimensions of the problem get larger, the ratio of the embedded volume to the "corners" of the embedding volume tends towards zero, thus a lot of rejections can take place before a useful sample is generated, thus making the algorithm inefficient and impractical. See curse of dimensionality. In high dimensions, it is necessary to use a different approach, typically a Markov chain Monte Carlo method such as Metropolis sampling  or Gibbs sampling. (However, Gibbs sampling, which breaks down a multi-dimensional sampling problem into a series of low-dimensional samples, may use rejection sampling as one of its steps.) Adaptive rejection sampling. For many distributions, finding a proposal distribution that includes the given distribution without a lot of wasted space is difficult. An extension of rejection sampling that can be used to overcome this difficulty and efficiently sample from a wide variety of distributions (provided that they have log-concave density functions, which is in fact the case for most of the common distributions—even those whose "density" functions  are not concave themselves) is known as adaptive rejection sampling (ARS). There are three basic ideas to this technique as ultimately introduced by Gilks in 1992: The method essentially involves successively determining an envelope of straight-line segments that approximates the logarithm better and better while still remaining above the curve, starting with a fixed number of segments (possibly just a single tangent line). Sampling from a truncated exponential random variable is straightforward. Just take the log of a uniform random variable (with appropriate interval and corresponding truncation). Unfortunately, ARS can only be applied from sampling from log-concave target  densities. For this reason, several extensions of ARS have been proposed in literature for tackling non-log-concave target distributions. Furthermore, different combinations of ARS and the Metropolis-Hastings method have been designed in order to obtain a universal sampler that builds a self-tuning proposal densities (i.e., a proposal automatically constructed and adapted to the target). This class of methods are often called as Adaptive Rejection Metropolis Sampling (ARMS) algorithms. The resulting adaptive techniques can be always applied but the generated samples are correlated in this case (although the correlation vanishes quickly to zero as the number of iterations  grows).