5
$\begingroup$

Here is a nice example of importance sampling:

%% true probability distribution
true_func = @(x) betapdf(x,1+1,1+10);

%% Do importance sampling
N = 10^6;
% uniform proposal distribution
x_samples = rand(N,1);
proposal = 1/N;
% evaluate for each sample
target = true_func(x_samples);
% calculate importance weight
w = target ./ proposal;
w = w ./ sum(w);
% resample, with replacement, according to importance weight
samples = randsample(x_samples,N,true,w);

%% plot
subplot(1,2,1)
x = linspace(0,1,1000);
plot(x, true_func(x) )
axis square

subplot(1,2,2)
hist(samples,1000)
title('importance sampling')
axis square

Result: enter image description here

I don't get it. If I already know what the target pdf looks like, then I can simply do this:

plot(true_func(linspace(0,1,N)));

I don't need to do importance sampling, I simply evaluate the target pdf at a linear space of choice.

Improve this question
$\endgroup$
9
$\begingroup$

Importance resampling is not for plotting the PDF. It is for sampling that PDF. Sometimes, even if you know the formula of the target PDF, it's hard to sample from it. For example, a typical method is inverse transform method, where you need to be able to analytically calculate the inverse of the CDF. A typical example where you cannot do so is normal distribution, in which another methods such as Box-Müller is used for sampling from it.

Improve this answer
$\endgroup$
  • $\begingroup$ I don't understand why it can be hard to sample from a PDF I know. Can you elaborate? In my example above, I essentially sampled from the pdf simply by writing true_func(linspace(0,1,N)). Plotting is one thing, but I can also calculate all possible moments, etc. $\endgroup$ – buffalo8 Dec 15 '19 at 21:12
  • 1
    $\begingroup$ How would you generate samples from $f(x)=ce^{-x^8}$ for example? Generating samples is not taking PDF values at regular intervals, you'll going to give an array of values which would produce the histogram of $f(x)$. $\endgroup$ – gunes Dec 15 '19 at 21:14
  • $\begingroup$ How about this: i=0; for x = lim_1 to lim_2; i = i + 1; samples(i) = c*exp(-x^8); end; I know something is wrong in my understanding but I can't seem to identify it. $\endgroup$ – buffalo8 Dec 15 '19 at 21:29
  • 4
    $\begingroup$ plot hist(samples), is it like $f(x)$? If not, you haven't sampled from it. $\endgroup$ – gunes Dec 15 '19 at 21:30
  • $\begingroup$ I got it. Thanks. $\endgroup$ – buffalo8 Dec 15 '19 at 21:32
9
$\begingroup$

Importance sampling is a Monte Carlo integration method that can be used to estimate the expected value of a function of a random variable. The method is useful in cases where the PDF is known, but the expected value of interest is unknown (and cannot be computed analytically from the PDF). In these cases, the method gives quite an efficient computation, so long as the generating distribution is reasonable. In your question, you have not actually finished the importance sampling, since you have not given an estimate of any function of the random variable with the distribution you are using.

It is also important to bear in mind that, for pedagogical purposes, it is useful to use importance sampling on a distribution where the expected value of interest is known (from analytic computation), to confirm that the estimate does indeed converge to the desired value.

Improve this answer
$\endgroup$
7
$\begingroup$

The code is actually producing samples by sampling importance resampling (Rubin 1987), that is by drawing points from the original iid sample $(y_1,\ldots,y_n)\sim\mathcal U(0,1)$

x_samples = rand(N,1);

according to the distribution $$\mathbb P(X=y_i)=\omega_i\big/\sum_{j=1}^n \omega_j$$ where the $\omega_j$'s are the importance weights, $\omega_j=f(y_j)/g(y_j)$ with $f$ the target true_func = @(x) betapdf(x,1+1,1+10);and $g$ the Uniform density (which should be $1$ rather than 1/N.

While this produces converging (in $n$) approximations to integrals depending on $f$ and to simulation from $f$, it is not an exact iid simulation from $f$ for a given $f$, due to (a) the bias induced by the normalisation w = w ./ sum(w); and (b) the dependence between the resampled values.

Improve this answer
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.