Week 3: Accept-Reject Sampling

STAT 536 — Statistical Computing

Soumojit Das

Fall 2026

The Problem

Motivation for Accept-Reject

The Challenge. What if:

  • No inverse transform available?
  • No simple transformation exists?
  • Distribution known only up to normalizing constant?

Examples.

  • Gamma with non-integer shape
  • Truncated distributions
  • Posterior distributions in Bayesian inference
  • Any \(f(x) \propto h(x)\) where \(\int h(x)dx\) is unknown

Solution. Accept-Reject: a general method using an easier “instrumental” distribution.

The Algorithm and Why It Works

Accept-Reject Algorithm

Theorem (Accept-Reject; von Neumann, 1951). Let \(f\) be the target density and \(g\) be an instrumental density with \[f(x) \leq M g(x) \quad \text{for all } x \text{ in the support of } f.\] Then the following algorithm generates \(X \sim f\).

Algorithm.

  1. Generate \(Y \sim g\) (instrumental)
  2. Generate \(U \sim \Uniform(0,1)\) (independent of \(Y\))
  3. If \(U \leq \dfrac{f(Y)}{M g(Y)}\), accept \(X = Y\)
  4. Otherwise, reject and return to step 1

Acceptance probability \(= 1/M\).

The Idea in One Sentence

Essentially, if we can bound the target \(f\) by a constant multiple of an easier density \(g\), proposing from \(g\) and randomly thinning the proposals produces exact draws from \(f\).

Recall Marsaglia’s polar method from last week? Propose uniformly on the square, keep only the disc. That was accept-reject for the uniform distribution on the disc, with \(g\) uniform on \([-1,1]^2\) and acceptance region \(\{V_1^2 + V_2^2 \leq 1\}\), giving \(M = 4/\pi\).

Accept-Reject: What We Need to Prove

The Claim. The Accept-Reject algorithm produces samples from the target distribution \(f\).

If \(X\) is the output of the algorithm, we need to show \[P[X \in A] = \int_A f(x)dx \quad \text{for any set } A.\] Equivalently, show the CDF: \(P[X \leq x] = \int_{-\infty}^x f(t)dt\).

\(X\) is a \(Y \sim g\) that was accepted. So we need \[P[X \leq x] = P[Y \leq x \mid Y \text{ was accepted}].\]

Accept-Reject: Understanding the Acceptance Event

The Acceptance Criterion. We accept \(Y\) when \(U \leq \dfrac{f(Y)}{Mg(Y)}\). This is a random event that depends on both \(Y\) and \(U\).

For a fixed \(y\)

\[ \begin{aligned} P[\text{accept} \mid Y = y] &= P\left[U \leq \frac{f(y)}{Mg(y)}\right]\\ &= \frac{f(y)}{Mg(y)} \end{aligned} \]

Since \(U \sim \Uniform(0,1)\).

Important.

  • High \(f(y)\): likely to accept
  • Low \(f(y)\): likely to reject
  • This biases \(Y\) toward \(f\).

The Complete Proof: The Numerator

Using rules of conditional probability: \[P[X \leq x] = P[Y \leq x \mid \text{accept}] = \frac{P[Y \leq x \text{ and accept}]}{P[\text{accept}]}\]

The numerator.

\[ \begin{aligned} P[Y \leq x \text{ and accept}] &= \int_{-\infty}^x P[\text{accept} \mid Y=y] \cdot g(y)dy\\ &= \int_{-\infty}^x \frac{f(y)}{Mg(y)} \cdot g(y)dy\\ &= \frac{1}{M}\int_{-\infty}^x f(y)dy \end{aligned} \]

The Complete Proof: Denominator and Conclusion

The denominator. Same computation over the whole line:

\[ P[\text{accept}] = \int_{-\infty}^{\infty} \frac{f(y)}{Mg(y)} \cdot g(y)dy = \frac{1}{M}\int_{-\infty}^{\infty} f(y)dy = \frac{1}{M} \]

Combining the two,

\[ \begin{aligned} P[X \leq x] &= \frac{P[Y \leq x \text{ and accept}]}{P[\text{accept}]}\\ &= \frac{\frac{1}{M}\int_{-\infty}^x f(y)dy}{\frac{1}{M}}\\ &= \int_{-\infty}^x f(y)dy = F(x) \end{aligned} \]

Correctness of Accept-Reject

Therefore, \(X\) has CDF \(F\) \(\implies\) \(X \sim f\). \(\square\)

Accept-Reject: Alternative View via Joint Distribution

Consider the joint distribution of \((Y, U)\) conditional on acceptance. The accepted points satisfy \(0 \leq U \leq \dfrac{f(Y)}{Mg(Y)}\), so \[f_{Y,U \mid \text{accept}}(y,u) = \frac{g(y) \cdot 1 \cdot \mathbf{1}\left[0 \leq u \leq \frac{f(y)}{Mg(y)}\right]}{P[\text{accept}]}.\]

Marginal density of \(Y\) given acceptance:

\[ f_{Y \mid \text{accept}}(y) = \int_0^{f(y)/(Mg(y))} \frac{g(y)}{1/M} \, du = \frac{g(y) \cdot \frac{f(y)}{Mg(y)}}{1/M} = f(y) \]

Accepted \(Y\) values are distributed exactly according to \(f\).

Accept-Reject: Why the Algorithm Makes Sense

The Fundamental Principle. To sample from \(f\), we need to see values \(y\) with frequency proportional to \(f(y)\).

What A-R Does

  1. Propose \(y\) with frequency \(g(y)\)
  2. Accept with probability \(\dfrac{f(y)}{Mg(y)}\)
  3. Net frequency: \[g(y) \times \frac{f(y)}{Mg(y)} \propto f(y)\]

Example. If \(g(y_0) = 2f(y_0)\):

  • See \(y_0\) twice as often as needed
  • Accept at half the rate of a point where \(f = g\): probability \(1/(2M)\) against \(1/M\)
  • Net effect: correct relative frequency.

The Role of \(M\). Ensures \(\dfrac{f(y)}{Mg(y)} \leq 1\) so we have valid probabilities.

Proposals Under \(Mg\), Accepted When They Fall Under \(f\)

Geometric Interpretation

Accept-Reject uniformly samples points under the curve \(Mg(x)\) and keeps those under \(f(x)\).

Theorem (Fundamental Theorem of Simulation). To generate \(X \sim f\), it suffices to:

  1. Generate \((X,V)\) uniformly on \(\{(x,v): 0 \leq v \leq f(x)\}\)
  2. Return the \(X\)-coordinate

Essentially, this says that scattering points uniformly under a density curve and reading off their \(x\)-coordinates produces draws from the density.

Accept-Reject implements this by:

  1. Generating \((Y,V)\) uniformly on \(\{(y,v): 0 \leq v \leq Mg(y)\}\)
  2. Keeping only points where \(V \leq f(Y)\)

Uniform Points Under the Curve

Choosing \(g\) and \(M\)

The proof asks nothing of \(g\) beyond \(f \leq Mg\), so what is left is a design problem: pick a \(g\) you can sample from, then make \(M\) as small as it will go.

Choosing the Instrumental Distribution

Requirements.

  • Easy to sample from \(g\)
  • Can evaluate \(g(x)\)
  • \(\text{supp}(f) \subseteq \text{supp}(g)\)
  • Can find \(M < \infty\)

Critical. If \(g\) has thinner tails than \(f\): \[\lim_{x \to \infty} \frac{f(x)}{g(x)} = \infty.\] Then \(M = \infty\).

If the instrumental density drops off faster than the target in the tails, the ratio \(f/g\) grows without bound, no finite \(M\) exists, and the algorithm cannot run.

Making \(M\) Small

Optimization. Want small \(M\) for efficiency:

  • \(M \geq 1\) always
  • \(M = 1 \iff f = g\)
  • Shape similarity helps

Common Choices.

  • Normal \(\leftarrow\) \(t\)-distribution
  • Gamma \(\leftarrow\) Exponential
  • Beta \(\leftarrow\) Uniform (for \(\alpha, \beta \geq 1\); otherwise the density is unbounded)

Finding the Optimal \(M\)

General Approach.

  1. Define \(h(x) = \dfrac{f(x)}{g(x)}\)
  2. Find \(x^* = \argmax_x h(x)\) (often by calculus)
  3. Set \(M = h(x^*) = \sup_x \dfrac{f(x)}{g(x)}\)

Example: Normal from Double Exponential

Target: \(f(x) = \dfrac{1}{\sqrt{2\pi}} e^{-x^2/2}\). Instrumental: \(g(x) = \dfrac{1}{2} e^{-|x|}\).

  • For \(x > 0\): \(h(x) = \sqrt{\dfrac{2}{\pi}} e^{-x^2/2 + x}\)
  • Maximize: \(h'(x) = 0 \Rightarrow x^* = 1\)
  • \(M = h(1) = \sqrt{\dfrac{2e}{\pi}} \approx 1.315\)

Example: Gamma Distribution (Non-Integer Shape)

Problem. Generate \(X \sim \Gam(\alpha, \beta)\) where \(\alpha \notin \mathbb{N}\) and \(\alpha > 1\).

Instrumental: an integer-shape Gamma. Let \(a = \lfloor\alpha\rfloor \geq 1\) and \(d = \alpha - a \in (0,1)\), and use \(g \sim \Gam(a, b)\), which is a sum of \(a\) exponentials.

  • Ratio: \(\dfrac{f(x)}{g(x)} \propto x^d e^{-(\beta - b)x}\)
  • The instrumental needs the smaller rate, \(b < \beta\). At \(b = \beta\) the ratio is \(x^d\), which is unbounded, so \(M = \infty\)
  • With \(b < \beta\) the ratio is maximized at \(x^* = \dfrac{d}{\beta - b}\)

Gamma: The Optimal Rate

Minimizing \(M\) over \(b\) gives \(b^* = \dfrac{a\beta}{\alpha}\), and then

\[M = \frac{\Gamma(a)}{\Gamma(\alpha)} \cdot \frac{\alpha^\alpha}{a^a} \cdot e^{-(\alpha - a)}.\]

Efficiency. \(M\) is free of the rate \(\beta\): only the shape matters. For \(\alpha = 2.5\), \(M \approx 1.13\) and the sampler accepts 89% of proposals.

Example: Truncated Normal Distribution

Problem. Generate \(X \sim f(x) \propto e^{-x^2/2} \mathbf{1}[x > a]\) (truncated normal).

Method 1: Naive

  1. Generate \(Z \sim \Normal(0,1)\)
  2. If \(Z > a\), accept
  3. Otherwise, reject

Efficiency: \(P[Z > a] = 1 - \Phi(a)\). Very inefficient for large \(a\).

Method 2: Shifted Exponential

Use \(g(x) = \lambda e^{-\lambda(x-a)}\mathbf{1}[x > a]\).

Optimal \(\lambda = \dfrac{a + \sqrt{a^2 + 4}}{2}\).

As \(a \to \infty\): \(M \to 1\), taking \(M\) as the bound for the normalized target. At \(a = 4\), \(M = 1.026\) and the sampler accepts 97% of proposals; the naive method accepts 0.003%.

The Shifted Exponential Envelope for a Truncated Normal

Unnormalized Targets

Working Without Normalizing Constants

The normalizing constant never enters the algorithm

Accept-Reject works even when \(f(x) = c \cdot f_1(x)\) with unknown \(c\).

Modified Algorithm. If \(f(x) = c f_1(x)\) and \(f_1(x) \leq M' g(x)\):

  1. Generate \(Y \sim g\), \(U \sim \Uniform(0,1)\)
  2. Accept if \(U \leq \dfrac{f_1(Y)}{M' g(Y)}\)

Multiply numerator and denominator by \(c\): \[\frac{f_1(y)}{M' g(y)} = \frac{c\,f_1(y)}{c M' g(y)} = \frac{f(y)}{M g(y)}, \qquad M = cM'.\]

This is the original algorithm with bound \(M = cM'\) for \(f\), so the proof already covers it and \(c\) never appears in the code.

Unnormalized Targets: The Bayesian Case

Bayesian Application. Posterior: \(\pi(\theta \mid x) \propto f(x \mid \theta)\pi(\theta) = f_1(\theta)\).

Don’t need \(\int f(x \mid \theta)\pi(\theta)d\theta\).

Note. The acceptance rate is \(\dfrac{1}{cM'}\): still constant, but not computable, since \(c\) is the constant you avoided. The observed acceptance rate estimates it, and so estimates \(c\).

This is the property the Bayesian half of the course runs on: the Metropolis-Hastings acceptance ratio in week 10 touches the posterior only through \(f_1\), for exactly the reason above.

What You Actually Do in Practice

The theory says: if \(f \leq Mg\), accepted proposals follow \(f\) exactly. In practice you execute a recipe. Each step draws on a specific result from this lecture.

  1. Pick \(g\) with cheap draws and tails at least as heavy as \(f\). [Thin-tails condition: \(g\) lighter than \(f\) means \(M = \infty\).]
  2. Compute \(M = \sup f/g\), usually by calculus on \(\log(f/g)\). [Optimal \(M\) procedure.]
  3. Code the loop: draw \(Y \sim g\), \(U \sim \Uniform(0,1)\); accept \(Y\) when \(U \leq f(Y)/(Mg(Y))\). [Algorithm + correctness theorem.]
  4. Check: compare the observed acceptance rate to \(1/M\). A large discrepancy is a bug. For an unnormalized target the rate is \(1/(cM')\), so the observed rate gives an estimate of \(c\). [Acceptance probability + unnormalized-target result.]

Efficiency

Envelope Accept-Reject

Motivation. Evaluating \(f(x)\) may be computationally expensive.

Theorem (Envelope Accept-Reject). If \(g_l(x) \leq f(x) \leq g_u(x) \leq Mg(x)\), then:

  1. Generate \(Y \sim g\), \(U \sim \Uniform(0,1)\)
  2. If \(U \leq \dfrac{g_l(Y)}{Mg(Y)}\), accept immediately
  3. Else if \(U > \dfrac{g_u(Y)}{Mg(Y)}\), reject immediately
  4. Else evaluate \(f(Y)\) and accept if \(U \leq \dfrac{f(Y)}{Mg(Y)}\)

Essentially, cheap bounds \(g_l\) and \(g_u\) bracket \(f\) and decide most proposals without evaluating \(f\). Only proposals falling between the bounds require the expensive evaluation.

Efficiency Analysis

Expected Attempts.

  • Number of attempts (including the accepted one) \(N \sim \text{Geom}(1/M)\) on \(\{1, 2, \ldots\}\)
  • \(\E[N] = M\)
  • \(\Var(N) = M(M-1)\)

On average, you draw \(M\) proposals from \(g\) for every accepted sample from \(f\).

Computational Cost. Per iteration: generate from \(g\), evaluate \(f(Y)/g(Y)\), generate a uniform.

Comparison.

Method Cost
Inverse \(O(1)\)
Box-Muller \(O(1)\)
Accept-Reject \(O(M)\)

Trade-off. Flexibility against efficiency.

Accept-Reject: Strengths and Limitations

Strengths.

  • Universal: works for any bounded \(f/g\)
  • Flexible: no need for \(F^{-1}\) or special structure
  • Unnormalized: works with \(f(x) \propto h(x)\)
  • Exact: produces exact samples

Limitations.

  • Efficiency depends critically on \(M\)
  • Finding good \(g\) can be challenging
  • May be slow for high-dimensional problems

Three Results, Three Roles

Result Statement Role
Correctness Accepted draws \(\sim f\) The method produces exact samples
Acceptance probability \(P[\text{accept}] = 1/M\) The cost: \(M\) proposals per accepted draw
Unnormalized targets Only \(f_1/g\) enters The normalizing constant never appears

What to Carry Out of the Room

  • Accept-Reject converts any instrumental with \(f \leq Mg\) for finite \(M\) into exact draws from \(f\).
  • The price is \(M\): on average \(M\) attempts per accepted draw, whatever \(f\) and \(g\) are.
  • The design problem is \(g\): tails at least as heavy as \(f\), cheap draws, and \(M = \sup f/g\) as small as it will go.

Because only the ratio \(f_1/g\) is ever evaluated, the target needs no normalizing constant, which is what makes the method reach posteriors at all.

Next week takes the two obvious complaints in turn. Adaptive rejection sampling builds the envelope for you when \(\log f\) is concave, ratio of uniforms sidesteps the envelope by moving to a bounded region in the plane, and importance sampling drops rejection altogether: keep every draw, weight it by \(f/g\), and spend the rest of the week making the variance of that estimator small.