Importance Sampling

Importance Sampling is a way to compute an expectation/integral more efficiently by sampling from a different distribution that puts more probability mass on the parts of the domain that matter the most, and then reweighing them to keep the estimator correct.

Motivation: Monte Carlo Integration

For example, consider MC Integration:

$$I={\int }_{D}g(x)f(x)d\underline{x}$$

where

• $D\subseteq {\mathbb{R}}^{\ell }$ where $D$ is bounded, i.e. $D=[0,1{]}^{\ell }$.
• $g:D\to \mathbb{R}$ is bounded and integrable (continuous)
• $f:D\to \mathbb{R}$ is a PDF where $f\ge 0$ with ${\int }_{D}f(x)d\underline{x}=1$.

Note that $I={\mathbb{E}}_f[g(x)]=\mathbb{E}\left\{g(x)\right\}$ and

$$\begin{array}{rl}{\sigma }^2& :=\text{Var}(g(x))\\ & =\mathbb{E}\left\{g(x{)}^2\right\}-{\left[\mathbb{E}\left\{g(x)\right\}\right]}^2\\ & ={\int }_{D}g(x{)}^{2}f(x)dx-I^2\\ & ={\int }_{D}g(x{)}^{2}f(x)dx-2I^2+I^2\\ & ={\int }_{D}g(x{)}^{2}f(x)dx-{\int }_{D}2Ig(x)f(x)dx+{\int }_D{I}^{2}f(x)dx\\ & ={\int }_D\left(g(x{)}^2-2Ig(x)+I^2\right)f(x)dx\\ & ={\int }_D(g(x)-I{)}^{2}f(x)dx\end{array}$$

Special Case: CMC via a Uniform Choice

More generally, any integral on $D$ can be rewritten as an expectation by introducing a PDF $p$ on $D$:

$${\int }_{D}h(x)dx={\int }_D\frac{h(x)}{p(x)}p(x)dx={\mathbb{E}}_p\left[\frac{h(X)}{p(X)}\right].$$

One simple choice is the Uniform density on $D$:

$$p(x)=\frac{1}{\text{vol}(D)}$$

If we plug this is in,

$$\frac{h(x)}{p(x)}=\frac{h(x)}{1/\text{vol}(D)}=\text{vol}(D)\cdot h(x),$$

such that

$${\int }_{D}h(x)dx=\mathbb{E}\left\{\text{vol}(D)\cdot h(X)\right\}=\text{vol}(D)\cdot \mathbb{E}\left\{h(X)\right\}.$$

as desired.

Importance Sampling (Variance Reduction)

Instead of sampling from the “default” density, choose a proposal density $\varphi$ on $D$ (with $\varphi (x)>0$ wherever $g(x)f(x)\ne 0$) and reweigh:

$$I={\int }_{D}g(x)f(x)dx={\int }_D\frac{g(x)f(x)}{\varphi (x)}\varphi (x)dx$$

This gives

$$I={\mathbb{E}}_f[g(X)]={\mathbb{E}}_{\varphi }\left[\frac{g(Y)f(Y)}{\varphi (Y)}\right]$$

for $X\sim f$ and $Y\sim \varphi$. Some remarks:

• $\varphi$ is called the importance function
• $\frac{f(x)}{\varphi (x)}$ is the likelihood ratio

Algorithm

1. Generate $Y_1,\dots ,Y_N\sim \varphi$ iid.
2. Set $$J_N=\frac{1}{N}\sum _{i=1}^N\frac{g(Y_i)f(Y_i)}{\varphi (Y_i)}.$$

Choosing $\varphi$

Pick $\varphi$ to satisfy:

• Support: $\varphi (x)>0$ anywhere $g(x)f(x)\ne 0$ (so the weight $f(x)/\varphi (x)$ is well-defined).
• Computable: easy to sample from.
• Low-variance weights: avoid making $\varphi (x)$ tiny where $|g(x)f(x)|$ is large (that creates huge weights and high variance).

Variance Reduction means choosing $\varphi$ such that

$${\int }_D\frac{g(x{)}^{2}f(x{)}^2}{\varphi (x)}dx$$

is small. The idea is to make $\varphi (x)$ larger when $(g(x)f(x){)}^2$ is large. Why do we want that integral small? This comes from the variance formula. Let

$$W:=\frac{g(Y)f(Y)}{\varphi (Y)},Y\sim \varphi .$$

Then $J_N=\frac{1}{N}{\sum }_{i=1}^N{W}_i$ with $W_i$ i.i.d., so

$$\text{Var}(J_N)=\frac{1}{N}\text{Var}(W).$$

Using $\text{Var}(W)=\mathbb{E}[W^2]-\mathbb{E}[W{]}^2$:

$${\mathbb{E}}_{\varphi }[W]={\int }_D\frac{g(x)f(x)}{\varphi (x)}\varphi (x)dx={\int }_{D}g(x)f(x)dx=I,$$

and

$${\mathbb{E}}_{\varphi }[W^2]={\int }_D{\left(\frac{g(x)f(x)}{\varphi (x)}\right)}^2\varphi (x)dx={\int }_D\frac{g(x{)}^{2}f(x{)}^2}{\varphi (x)}dx.$$

Therefore

$$\text{Var}(J_N)=\frac{1}{N}\left({\int }_D\frac{g(x{)}^{2}f(x{)}^2}{\varphi (x)}dx-I^2\right).$$

Since $I^2$ does not depend on $\varphi$, making ${\int }_D\frac{g^2{f}^2}{\varphi }$ small makes the estimator variance small.

Example 1

Consider $D=[0,1]\subseteq \mathbb{R},f(x)=1$ on $D$ where $(f(x)=0,x\notin D)$ and $g(x)=4\sqrt{1-x^2}$ on $D$.

By Crude MC,

$$I_N=\frac{1}{N}\sum _{k=1}^{N}g(U_k)$$

the variance is

$$\frac{1}{N}\left[{\int }_0^{1}16(1-x^2)dx-I^2\right]\approx \frac{0.797\cdots }{N}.$$

Note that

$$I={\int }_0^{1}4\sqrt{1-x^2}dx=\pi .$$

Example 2 (Using Variance Reduction)

Let’s try

$$\varphi (x)=\frac{4}{3}-\frac{2}{3}xx\in [0,1].$$

This gives the estimator

$$J_N=\frac{1}{N}\sum _{k=1}^N\frac{g(Y_k)f(Y_k)}{\varphi (Y_k)}$$

for $Y_1,\dots ,Y_N\sim \varphi$ iid. So,

$$J_N=\frac{1}{N}\sum _{k=1}^N\frac{6\sqrt{1-Y_k^2}}{2-Y_k}.$$

Then

$$\begin{array}{rl}{\text{Var}}_{\varphi }(J_N)& =\frac{1}{N}\left[{\int }_0^1\frac{g(x{)}^2}{\varphi (x)}dx-I^2\right]\\ & =\frac{1}{N}\left[{\int }_0^1\frac{24(1-x^2)}{2-x}dx-{\pi }^2\right]\\ & \approx \frac{0.224\cdots }{N}.\end{array}$$

by applying the variance formula for estimators above. What if we overdo it with a large $\varphi$?

Let’s try $\psi (x)=2-2x$ for $x\in [0,1]$. $\psi$ being $0$ when $x=1$ is not ideal. What does that look like?

The estimator

$$K_N=\frac{1}{N}\sum _{k=1}^N\frac{g(Z_k)f(Z_k)}{\psi (Z_k)}$$

with $Z_1,\dots ,Z_N\sim \psi$. Then

$${\text{Var}}_{\psi }(K_N):=\frac{1}{N}\left[{\int }_0^1\frac{8(1-x^2)}{1-x}dx-{\pi }^2\right]\approx \frac{2.1\cdots }{N},$$

which shows it’s a bad estimator.