Acceptance-Rejection

Assume PDF $f:\mathbb{R}\to [0,\infty )$ is our target PDF and that

• $f=0$ outside $[a,b]$
• $f(x)\le c$, $\forall x\in [a,b]$ for some $c\in \mathbb{R}$¹ The goal is to sample the RV $Z\sim f$.
• Here, $c$ is green and $[a,b]$ are the red vertical bars.²

Naive Algorithm

1. Generate $X\sim U[a,b]$. For example, let $U\sim U[0,1]$ and $X=a+(b-a)U$.
2. Generate $Y\sim U[0,1]$ independent of $X$.
3. If $Y\le f(X)$, accept $X$ and set $Z\leftarrow X$.²

Why does this work?

• $(X,Y)$ is uniformly distributed in the box $[a,b]\times [0,1]$. This means that the accepted pair $(X,Y)$ is uniformly distributed on

$$B=\left\{(x,y)|a\le x\le b,0\le y\le f(x)\right\}$$

• The area of $B$ is $1$, since the area of $f$ is $1$.
• This implies that the marginal PDF for accepted $X$ is then

$${\int }_{-\infty }^{\infty }{f}_{X,Y}(x,y)dy={\int }_0^{f(x)}1dy=f(x)$$

• This tells us the distribution of accepted $X$ is $f(x)$.
• This integral works because we are integrating over $y$, and the PDF is only defined between $0$ and $f(x)$. The joint PDF collapses to $1$, because it is the value of the accepted points during sampling.

General Algorithm

We can write a more general A-R method. The goal now is to sample from a target PDF

$$f:{\mathbb{R}}^n\to [0,\infty )$$

For intuition sake, we can imagine $n=1$.

The top part is $\alpha g$ and the purple part is $f$.

We want to find a PDF $g:{\mathbb{R}}^n\to [0,\infty )$ and $\alpha \ge 1$ such that

• We know how to sample from $g$. $g$ is called the proposed PDF. It is close to $f$ and easy to sample from.
• $\alpha g(x)\ge f(x)$ for all $x\in {\mathbb{R}}^n$. $\alpha g(x)$ is called a majorizing function of $f$.
  • Note that $\alpha \ge 1$ is required because $g$ is a PDF which integrates to $1$.

$$\alpha ={\int }_{{\mathbb{R}}^n}\alpha g(x)d\underline{x}\ge {\int }_{{\mathbb{R}}^n}f(x)d\underline{x}=1$$

• We want $\alpha$ to be as small as possible.

Acceptance-Rejection Algorithm

Given target PDF $f:{\mathbb{R}}^n\to [0,\infty )$, find proposal PDF $g$ and constant $\alpha$ such that $\alpha g(x)\ge f(x)$ for all $x\in {\mathbb{R}}^n$.

1. Generate $X\sim g$ where $X\in {\mathbb{R}}^n$.
2. Generate $Y\sim U[0,\alpha g(x)]$ where $Y\in \mathbb{R}$.
3. If $y\le f(x)$, then accept $x$ and set $z\leftarrow x$. Otherwise, we reject and go step 1.

We want to measure the following rates:

Efficiency

The probability of acceptance is: $P\left((X,Y)\text{ is accepted}\right)=\frac{v(B)}{v(A)}$ Where $v(B)$ and $v(A)$ correspond to the volumes under the curves:

$$\begin{array}{rl}v(B)& ={\int }_{{\mathbb{R}}^n}f(x)dx=1\\ v(A)& ={\int }_{{\mathbb{R}}^n}\alpha g(x)dx=\alpha \end{array}$$

Thus, the expected acceptance rate is $\frac{1}{\alpha }$.

Theorem (Uniformity on Region $A$)

Let the relevant geometric regions be defined as:

$$\begin{array}{rl}A& =\{(x,y)\mid 0\le y\le \alpha g(x)\}\\ B& =\{(x,y)\mid 0\le y\le f(x)\}\subseteq A\end{array}$$

The pair $(X,Y)$ obtained from Steps 1 & 2 is uniformly distributed on $A$.

Proof:

• Let $f_{X,Y}:{\mathbb{R}}^n\times \mathbb{R}\to [0,\infty )$ be the joint PDF of $(X,Y)$.
• Let $f_{Y|X}(y|x)$ be the conditional PDF of $Y$ given $X=x$.

Recall from theorem that the joint PDF is equal to the conditional times marginal. Since $X\sim g$ (Step 1), we have:

$$f_{X,Y}(x,y)=\{\begin{array}{ll}g(x)f_{Y|X}(y|x)& (x,y)\in A\\ 0& \text{otherwise}\end{array}$$

From Step 2, we generate $Y\sim U[0,\alpha g(x)]$. The PDF of a uniform distribution on $[0,k]$ is $1/k$, so:

$$f_{Y|X}(y|x)=\{\begin{array}{ll}\frac{1}{\alpha g(x)}& y\in [0,\alpha g(x)]\\ 0& \text{otherwise}\end{array}$$

Substituting this back into the joint PDF:

$$(\star )f_{X,Y}(x,y)=g(x)\cdot \frac{1}{\alpha g(x)}=\frac{1}{\alpha }$$

Since the PDF is constant on the support $A$, $(X,Y)$ is uniformly distributed on $A$.

Theorem (From Uniformity on $A$ to Target $f(x)$)

We established that the candidate points $(X,Y)$ are uniform on $A$.

$$⟹f_{X,Y}(x,y)=\frac{1}{\alpha }\forall (x,y)\in A$$

How does $X\sim f$?

Proof:

Let $(X^{\ast },Y^{\ast })$ be the first accepted point. Because $(X,Y)\sim U(A)$ and we strictly accept points within $B\subseteq A$, the accepted points are uniformly distributed on $B$:

$$(X,Y)\sim U(A)⟹(X^{\ast },Y^{\ast })\sim U(B)$$

$B$ is the area under the PDF $f$, so its total volume is $\int f(x)dx=1$. The joint PDF of the accepted point $(X^{\ast },Y^{\ast })$ is therefore:

$$f_{X^{\ast },Y^{\ast }}(x,y)=\{\begin{array}{ll}1& (x,y)\in B\\ 0& \text{otherwise}\end{array}$$

To find the distribution of just the $x$-values, we take the marginal of $f_{X^{\ast },Y^{\ast }}$:

$$\begin{array}{rl}{f}_{X^{\ast }}(x)& ={\int }_{-\infty }^{\infty }{f}_{X^{\ast },Y^{\ast }}(x,y)dy\\ & ={\int }_0^{f(x)}1dy\\ & =f(x)\end{array}$$

Thus the accepted samples $X^{\ast }$ follow the target distribution $f(x)$.

Theorem (A-R Generates $Z\sim f$)

This theorem explains why this works via geometric intuition.

Consider the graph of $\alpha g(x)$ enveloping $f(x)$.

• Step 1 ($X\sim g$): $X$ is more likely to fall where $g$ is big. This leads to a “mushing” (clustering) of points along the x-axis at the peaks.
• Step 2 ($Y\sim U$): Choosing $Y\sim U[0,\alpha g(x)]$ gives “more space” for $Y$ to be uniform where $g$ is large.

This vertical expansion compensates for the horizontal “mushing” of $X$. The high horizontal density (large $g$) $\times$ low vertical density (large range $\alpha g$) gives us a near 2D constant density in $A$. These effects cancel out perfectly, see $(\star )$ from Theorem (Uniformity on Region A). This allows us to cut out the shape of $f(x)$ or region $B$. The remaining points (accepted points) form $f$.

In general, this is easy to visualize in ${\mathbb{R}}^2$, but this method provides justification for ${\mathbb{R}}^n$.

Modified A-R

Let $f$ be our target PDF and $g$ be our proposal PDF. Let $\alpha \in \mathbb{R}$ where $f\ge \alpha g$ for all $x\in \mathbb{R}$.

1. Generate $X\sim g$
2. Generate $U\sim U[0,1]$.
3. If $$U\le \frac{f(X)}{\alpha g(X)}$$ then accept $X$. Set $Z\leftarrow X$. Else, reject $X$ and go to Step 1.

Example 1 (Sampling Standard Normal via A-R)

Let suppose we wanted to sample $N(0,1)$, the standard normal. So,

$$f(x)=\frac{1}{\sqrt{2\pi }}\exp \left(-\frac{1}{2}{x}^2\right)x\in \mathbb{R}$$

and let

$$g(x)=\frac{1}{\alpha }\exp (-|x|)\text{on }\mathbb{R}$$

Visually,

Note that $g$ is the exponential distribution $\text{Exp}(1)$ which has PDF $e^{-x}$ for $x\ge 0$ multiplied randomly $\pm 1$.

To do Step 1 (generate $X\sim g$) we need to:

1. Generate $U_1,U_2\sim [0,1]$ independently. Let $U_1$ be transformed by Inversion and $U_2$ be $\pm 1$.
2. Set $$X=\{\begin{array}{ll}-\ln U_2& U_1\le \frac{1}{2}\\ +\ln U_2& U_1\ge \frac{1}{2}\end{array}$$
3. Now we find $\alpha \ge 1$ such that $\alpha g(x)\ge f(x)$ for all $x\in \mathbb{R}$.

For our case, we can set

$$\alpha =\underset{x\in \mathbb{R}}{\sup }\frac{f(x)}{g(x)}=\underset{x\ge 0}{\sup }\sqrt{\frac{2}{\pi }}\exp \left(x-\frac{x^2}{2}\right)=\sqrt{\frac{2e}{\pi }}\approx 1.315$$

This allows us to sample $N(0,1)$ computationally.

Footnotes

1. In practice, you want $c$ as small as possible so that sampling can converge faster. ↩

2. Ideally, $c$ is as tight as possible. The textbook denotes $c=\sup \{f(x)|x\in [a,b]\}$. ↩ ↩²