Monte Carlo or other general sampling approaches for conditional distributions? Suppose we have a sampler – eg a Monte Carlo sampler – for the posterior probability distribution $p(x \vert D)$ of a quantity $x=(y,z)$ consisting of a pair of continuous and multidimensional quantities. $D$ is a set of previously observed values of $x$.
Are there any general methods to obtain samples from the conditional distribution $p(y \vert z^*, D)$, for a given $z^*$, from the previous sampler or from samples of $x$ obtained therefrom?
(To be clearer, what I'm interested in is not just the conditional mean $\mathrm{E}(y \vert z^*, D)$ or similar statistics, but the full conditional distribution.)
 A: In the event the joint density $p(y,z|D)$ is available in closed form (possibly up to a multiplicative constant), a sample $x_1,\ldots,x_n$ can be exploited to produce a sample from $p(y|z^*,D)$ by

*

*acceptance-rejection, assuming $p(y,z^*|D)/p(y,z_i|D)$ is upper bounded [exact sampling, random sample size]

*sampling importance resampling, selecting with replacement $y_i$'s with probability proportional to $p(y,z^*|D)/p(y,z_i|D)$ [biased resampling, due to the renormalisation of the importance weights]

*Markov chain Monte Carlo, with, at iteration $t$, selecting one index i at random and replacing the current value $y^{(t)}$ by $y_i$ with Metropolis-Hastings probability$$1 \wedge \dfrac{p(y_i,z^*|D)}{p(y_i,z_i|D)} \dfrac{p(y^{(t)},z_i|D)}{p(y^{(t)},z^*|D)}$$[exact sampling assuming enough iterations for the Markov chain to reach stationarity]
Here is a basic R illustration
#generate conditional from joint
#x ~ Exp(1) & y|x ~ N(0,x)
#joint sample
N=1e5
smpl=cbind(z<-rexp(N),y<-z*rnorm(N))
#conditional on 
zo=2
we=dnorm(y,sd=zo)/dnorm(y,sd=z)
# 1. Accept-reject
yo=smpl[z>zo,2]
ya=yo*(runif(length(yo))<sqrt(zo/z[z>zo])*we[z>zo])
ya=ya[ya!=0]
# 2. sampling importance resampling
ess=round(sum(we)^2/sum(we^2)) #effective sample size
yi=sample(y,ess,prob=we,rep=TRUE)
# 3. MCMC
ym=rep(yi[1],N/10)
for (t in 2:length(ym)){
  i=sample(1:N,1)
  ym[t]=ifelse(runif(1)<we[i]*dnorm(ym[t-1],sd=z[i])/dnorm(ym[t-1],sd=zo),
               y[i],ym[t-1])}

Without such items of information, i.e. in a generative model when the sampler is a black box, there is no generic solution and approximate methods such as ABC or synthetic likelihood need be used.
