# Importance Sampling: What is the role of the target distribution $p(x)$ in Monte Carlo integration?

I am currently implementing a few different approaches to Monte Carlo Integration and there is a conceptual hurdle that I don’t fully understand yet.

The notation varies from one source to another, but importance sampling is commonly conveyed in a form below (or something similar), where $$p(x)$$ is the target distribution and $$q(x)$$ is the importance/proposal distribution.

$$\int f(x) \frac{ p(x) }{ q(x) } q(x) dx$$

The role of $$q(x)$$ seems pretty clear - it defines the PDF that we draw $$x$$ from, and it should (ideally) be higher in regions where additional samples would be more useful.

I’m still not completely clear on the role of $$p(x)$$, however. The fact that it’s referred to as the target, and that we want the variance in $$\frac{p(x)}{q(x)}$$ to be low, makes me wonder if $$p(x)$$ should resemble the output distribution of $$f(x)$$.

What role does $$p(x)$$ serve in this equation? And are there reliable ways to select a good distribution?

(I apologize if this has been answered already, I searched twice and wasn't able to find a specific answer to this)

I think you might be confused a bit about notation. Let me start with a brief overview. In Monte Carlo methods, we are looking to approximate integrals. For example, we are interested in approximating the mean of a continuous variable with density $$p(x)$$:

$$\int xp(x)dx.$$

We can generalize this idea to find any particular function of $$x$$, $$f(x)$$:

$$\int f(x)p(x)dx.$$

The problem is that sometimes we don't know $$p(x)$$, so we cannot compute the integral. We call $$p(x)$$ the target distribution because this is the distribution we want to take an expectation from. In Bayesian statistics, for example, this would be the posterior distribution of a parameter $$\theta$$, given some data $$\mathcal{D}$$: $$p(\theta \mid \mathcal{D})$$. Note, however, that even though the target distribution is difficult to integrate, we have to be able to evaluate a function proportional to it $$p'(x)\propto p(x)$$ [Footnote].

The way we approximate our integral of interest with importance sampling is by introducing an "auxiliary" distribution $$q(x)$$, which is easy to sample from:

$$\int f(x)p(x)\frac{q(x)}{q(x)}dx.$$

Notice that the integral is unchanged when we multiply and divide by $$q(x)$$. By defining $$w_{i}=p(x_{i})/q(x_{i})$$, and by the general principle of MC sampling, we have that:

$$\int f(x)p(x)\frac{q(x)}{q(x)}dx \approx \frac{1}{N}\sum^{N}_{i=1}w_{i}f(x_{i}),\quad \text{with}\ x_{i}\sim q(x). \tag{\ast}$$

Having explained the general idea of Importance sampling, let me go to your questions.

a) I think you might be confusing $$f(x)$$ and $$p(x)$$. As mentioned above, $$f(x)$$ is the function of $$x$$ we want to compute with respect to the density $$p(x)$$.

b) Unfortunately, there is no straightforward way to choose the sampling distribution $$q(x)$$. If you are confident of your prior, that is, if they are sensible (maybe you checked them with prior predictive checks), then you can use it as the proposal.

More recently, there has been research on choosing the proposal distribution using something called "compiled inference", where you use the proposed generative model to learn a neural network that outputs proposal distributions for a particular datum, in the spirit of amortized inference.

[Footnote] If we use an unnormalized distribution (just as in the $$p'(x)\propto p(x)$$ case), we must use the self-normalized importance sampling algorithm. The difference between this and the "vanilla" algorithm explained above, is that the weights used in $$(\ast)$$ have to be normalized. So instead of using $$w_{i}=\frac{p(x_{i})}{q(x_{i})}$$ in the sum, we must use $$\tilde{w}_{i}=\frac{w_{i}}{\sum_{i=1}^{N}w_{i}}$$. This is a technical, yet important detail.

EDIT

A worked-out example in the context of Bayesian Modeling.

Suppose you have data, $$\mathcal{D}$$, sampled from a Gaussian distribution with a known standard deviation $$\sigma$$. We write this as $$p(\mathcal{D} \mid \mu)$$. Such function is known as the likelihood function and it is a function of the parameter, $$\mu$$, and not the data.

From the data you observed and modelling assumptions, you would like to estimate the posterior distribution of the mean (of the Gaussian distribution), $$\mu$$. We write this distribution $$p(\mu \mid \mathcal{D})$$, this is our target distribution, and we would like to estimate expectations with respect to it. In particular, we will estimate the mean of such distribution, $$\mathbb{E}_{p(\mu\mid\mathcal{D})}[\mu]$$.

In Bayesian inference, we find an expression of the posterior distribution by setting a prior distribution on the $$\mu$$. We denote this fixed distribution with $$p(\mu)$$. To make it non-conjugate, we will assume the prior is an exponential distribution with known parameter $$\lambda$$.

We now have all the ingredients to find the posterior distribution we are interested in. Using Bayes rule,

$$p(\mu \mid \mathcal{D}) = \frac{p(\mathcal{D}\mid \mu)p(\mu)}{p(\mathcal{D})}.$$

Unfortunately, $$p(\mathcal{D})=\int p(\mathcal{D}\mid \mu)p(\mu) d\mu$$ has no closed-form solution, so we have to find a different way to estimate the expectations with respect to the posterior. Nevertheless, and since $$p(\mathcal{D})$$ is constant with respect to our variable of interest $$\mu$$, we can evaluate a function proportional to the posterior distribution

$$p(\mu \mid \mathcal{D}) \propto p(\mathcal{D}\mid \mu)p(\mu).$$

Now that we have introduced all this machinery let's recall our goal is to find the following expectation

$$\mathbb{E}_{p(\mu\mid\mathcal{D})}[\mu] = \int \mu p(\mu \mid \mathcal{D}) d\mu.$$

This looks very similar to the goal of importance sampling. Let's introduce two target distributions, $$q(\mu)$$, to see how it works.

a) Use the prior distribution $$p(\mu)$$ as auxiliary distribution. In this case, we have that $$p(\mu)$$ cancels out in the numerator and the denominator:

$$\mathbb{E}_{p(\mu\mid\mathcal{D})}[\mu] = \int \mu \frac{p(\mathcal{D}\mid \mu)p(\mu)p(\mu)}{p(\mu)} = \int \mu p(\mathcal{D}\mid \mu)p(\mu) d\mu,$$

and the computation is reduced to

\begin{align} \int \mu p(\mathcal{D}\mid \mu)p(\mu) d\mu \approx \frac{1}{N}\sum^{N}_{i=1}\tilde{w}_{i}\mu_{i}, \\ \text{with}\ w_{i}=p(\mathcal{D}\mid \mu_{i}),\ \tilde{w_{i}}=\frac{w_{i}}{\sum_{i}^{N}w_{i}}\ \text{and}\ \mu_{i}\sim p(\mu). \end{align}

b) Use a distribution different from the prior distribution, say another Gaussian distribution $$q(\mu)$$. In such case, the numerator and the denominator do not cancel out, and the computation is the same as in equation ($$\ast$$) above.

• Thanks Sergio! I might have been a bit unclear in my question. The definitions of $f(x)$ and $q(x)$ you provided are consistent with my understanding. The part I'm having trouble grasping is $p(x)$, or rather $p'(x)$, which I still don't have a concrete understanding of. Could you elaborate on how to approach determine $p'(x)$? Thanks again! Commented May 6, 2022 at 22:03