# Sampling from marginal using integrated conditional

I would like to sample from: $$p(\theta_2|x)=\int p(\theta_2|\theta_1,x) . p(\theta_1|x) . d\theta_1$$ knowing that I can easily sample from $p(\theta_1|x)$ and (less easily) from $p(\theta_2|\theta_1,x)$. I have some “naive-brute-force-computationally-intractable” ideas but maybe there exists some more elegant ones?

EDIT:

1) Several answers suggest me to use the following procedure:

for i=1:N
draw $\theta_1^i$ from $p(\theta_1|x)$
draw $\theta_2^i$ from $p(\theta_2|\theta_1^i,x)$


Is this provide me i.i.d $\theta_2^i$ from $p(\theta_2|x)$ ?

2) In practice, I would like to compute an HPD interval over $\theta_2$ (e.g. using boa.hpd). Thus I do not really need i.i.d RVs. Is the following procedure valid :

for i=1:P
draw $\theta_1^{i}$ from $p(\theta_1|x)$
for j=1:N/P
draw $\theta_2^{(i-1)*P+j}$ from $p(\theta_2|\theta_1^i,x)$


i.e. drawing many $\theta_2^{(i-1)*P+j}$ from a single $\theta_1^{i}$, $P$ times ($P$ being a divisor of $N$ and being fixed and $N$ the number of samples that I want to generate) and use these $N$ samples to compute hpd. In practice such a solution would be largely more efficient than the first one in my application.

• This is kind of like Rao-Blackwellisation. Note that by sampling in the way you suggest, you already get a sample from $p(\theta_2)$. Just ignore $\theta_1$. Commented Jul 5, 2013 at 10:25
• Also, you will get much better answers if you write down what the distributions are. Commented Jul 5, 2013 at 10:28
• I’m not 100% sure but maybe the answer is, as you stated, sample from $p(θ_{1}|x)$ and then from $p(θ_{2}|θ_{1},x)$. Ignore the integral over $θ_{1}$. Commented Jul 5, 2013 at 13:58
• @probabilityislogic: from which argument can I validate to forget $\theta_1$ ? Commented Jul 5, 2013 at 14:08

As one of the comments stated (guessed), if you sample from the marginal distribution $p(\theta_1\vert x)$ and then use this value to sample from the conditional distribution $p(\theta_2\vert \theta_1,x)$, you will end up with a sample from the joint distribution $p(\theta_1,\theta_2\vert x)$. Take the projection on the second column (simply take the column of samples associated to $\theta_2$), and you get a sample from $\theta_2\vert x$. This obviously provides an answer to 1). This is a common technique in Bayesian statistics.
The answer to 2) is "maybe", it depends on the specific target distributions, $P$ and $N$. You will need to conduct an sensitivity analysis to check the robustness of the method which will probably be more complicated than using the first approach.