What is importance sampling?

Question

I'm trying to learn reinforcement learning and this topic is really confusing to me. I have taken an introduction to statistics, but I just couldn't understand this topic intuitively.

Answer 1

Importance sampling is a form of sampling from a distribution different from the distribution of interest to more easily obtain better estimates of a parameter from the distribution of interest. Typically this will provide estimates of the parameter with a lower variance than would be obtained by sampling directly from the original distribution with the same sample size.

It is applied in various contexts. In general sampling from the different distributions allows for more samples to be taken in a portion of the distribution of interest that is dictated by the application (important region).

One example might be that you want to have a sample that includes more samples from the tails of the distribution than pure random sampling from the distribution of interest would provide.

The wikipedia article that I have seen on this subject is too abstract. It is better to look at various specific examples. However, it does include links to interesting applications such as Bayesian Networks.

One example of importance sampling in the 1940s and 1950s is a variance reduction technique (a form of the Monte Carlo Method). See for example the book Monte Carlo Methods by Hammersley and Handscomb published as a Methuen Monograph/Chapman and Hall in 1964 and reprinted in 1966 and later by other publishers. Section 5.4 of the book covers Importance Sampling.

Answer 2

Importance sampling is a simulation or Monte Carlo method intended for approximating integrals. The term "sampling" is somewhat confusing in that it does not intend to provide samples from a given distribution.

The intuition behind importance sampling is that a well-defined integral, like $$\mathfrak{I}=\int_\mathfrak{X} h(x)\,\text{d}x$$ can be expressed as an expectation for a wide range of probability distributions: $$\mathfrak{I}=\mathbb{E}_f[H(X)]=\int_\mathfrak{X} H(x)f(x)\,\text{d}x$$ where $f$ denotes the density of a probability distribution and $H$ is determined by $h$ and $f$. (Note that $H(\cdot)$ is usually different from $h(\cdot)$.) Indeed, the choice $$H(x)=\dfrac{h(x)}{f(x)}$$leads to the equalities $H(x)f(x)=h(x)$ and $\mathfrak{I}=\mathbb{E}_f[H(X)]$$-$under some restrictions on the support of $f$, meaning $f(x)>0$ when $h(x)\ne 0$$-$. Hence, as pointed out by W. Huber in his comment, there is no unicity in the representation of an integral as an expectation, but on the opposite an infinite array of such representations, some of which are better than others once a criterion to compare them is adopted. For instance, Michael Chernick mentions choosing $f$ towards reducing the variance of the estimator.

Once this elementary property is understood, the implementation of the idea is to rely on the Law of Large Numbers as in other Monte Carlo methods, i.e., to simulate [via a pseudo-random generator] an iid sample $(x_1,\ldots,x_n)$ distributed from $f$ and to use the approximation $$\hat{\mathfrak{I}}=\frac{1}{n} \sum_{i=1}^n H(x_i)$$which

1. is an unbiased estimator of $\mathfrak{I}$
2. converges almost surely to $\mathfrak{I}$

Depending on the choice of the distribution $f$, the above estimator $\hat{\mathfrak{I}}$ may or may not have a finite variance. However, there always exist choices of $f$ that allow for a finite variance and even for an arbitrarily small variance (albeit those choices may be unavailable in practice). And there also exist choices of $f$ that make the importance sampling estimator $\hat{\mathfrak{I}}$ a very poor approximation of ${\mathfrak{I}}$. This includes all the choices where the variance gets infinite, even though a recent paper by Chatterjee and Diaconis studies how to compare importance samplers with infinite variance. The picture below is taken from my blog discussion of the paper and illustrates the poor convergence of infinite variance estimators.

Importance sampling with importance distribution an Exp(1) distribution target distribution an Exp(1/10) distribution, and function of interest $h(x)=x$. The true value of the integral is $10$.

[The following is reproduced from our book Monte Carlo Statistical Methods.]

The following example from Ripley (1987) shows why it may actually pay to generate from a distribution other than the (original) distribution $f$ appearing in the integral $$\int_\mathfrak{X} h(x) f(x)\,\text{d}x$$of interest or, in other words, to modify the representation of an integral as an expectation against a given density.

Example (Cauchy tail probability) Suppose that the quantity of interest is the probability, $p$, that a Cauchy ${\mathcal{C}}(0,1)$ variable is larger than $2$, that is, $$ p = \int_2^{+\infty} \; {1\over \pi(1 + x^2)} \; \text{d}x \;. $$ When $p$ is evaluated through the empirical average $$ {\hat{p}}_1 = {1\over m} \; \sum_{j=1}^m \; \mathbb{I}_{X_{j} > 2} $$ of an iid sample $X_1,\ldots,X_m$ $\sim$ $\; \mathcal{C}(0,1)$, the variance of this estimator is $p(1-p)/m$ (equal to $0.127/m$ since $p=0.15$).

This variance can be reduced by taking into account the symmetric nature of ${\mathcal{C}}(0,1)$, since the average $$ {\hat{p}}_2 = {1\over 2m} \; \sum_{j=1}^m \; \mathbb{I}_{|X_{j}| > 2} $$ has variance $p(1-2p)/2m$ equal to $0.052/m$.

The (relative) inefficiency of these methods is due to the generation of values outside the domain of interest, $[2,+\infty)$, which are, in some sense, irrelevant for the approximation of $p$. [This relates to Michael Chernick mentioning tail area estimation.] If $p$ is written as $$ p = {1\over 2} - \int_0^2 \; {1\over \pi(1 + x^2)} \; \text{d}x \;, $$ the integral above can be considered to be the expectation of $h(X) = 2/\pi(1 + X^2)$, where $X \sim {\mathcal{U}}_{[0, 2]}$. An alternative method of evaluation for $p$ is therefore $$ {\hat{p}}_3 = {1\over 2} - {1\over m} \; \sum_{j=1}^m \; h(U_j) $$ for $U_j \sim {\mathcal{U}}_{[0,2]}$. The variance of ${\hat{p}}_3$ is $(\mathbb{E}[h^2] - \mathbb{E}[h]^2)/m$ and an integration by parts shows that it is equal to $0.0285/m$. Moreover, since $p$ can be written as $$ p = \int_0^{1/2} \; {y^{-2}\over \pi(1 + y^{-2})} \; \text{d}y \;, $$ this integral can also be seen as the expectation of ${1\over 4} \; h(Y) = 1/2\pi(1 + Y^2)$ against the uniform distribution on $[0,1/2]$ and another evaluation of $p$ is $$ {\hat{p}}_4 = {1\over 4 m} \; \sum_{j=1}^m \; h(Y_j) $$ when $Y_j \sim {\mathcal{U}}_{[0,1/2]}$. The same integration by parts shows that the variance of ${\hat{p}}_{4}$ is then $0.95 \; 10^{-4}/m$.

Compared with ${\hat{p}}_1$, the reduction in variance brought by ${\hat p}_4$ is of order $10^{-3}$, which implies, in particular, that this evaluation requires $\sqrt{1000} \approx 32$ times fewer simulations than $\hat p_1$ to achieve the same precision. $\blacktriangleright$