Q9. Gibbs sampling for a correlated normal target

Topic: Markov chains, Gibbs sampling, and Metropolis–Hastings · more questions on this topic
‹ Q5 Q10 ›

Consider the joint posterior density p(\theta_1,\theta_2\mid{\boldsymbol{x}}) \propto \exp\left[ -\frac{\theta_1^2-2\rho\theta_1\theta_2+\theta_2^2} {2(1-\rho^2)} \right], where \theta_1,\theta_2\in\mathbb{R} and -1<\rho<1 is known.

(a) Derive the two full conditional distributions, identifying their means and variances. [4 marks]

(b) Give pseudocode for a Gibbs sampler with this equilibrium density. Make clear which current or updated value is used at each step. [3 marks]

(c) Suppose the initial pair is drawn exactly from the equilibrium distribution. Are the subsequent sampled pairs necessarily independent? Explain. [2 marks]

(d) When \rho is close to 1, explain why this sampler may move slowly through the joint posterior, even though every Gibbs update is accepted. [3 marks]

Worked answer (click to show)

Definitions. The target coordinates \theta_1,\theta_2 have zero means, unit marginal variances, and correlation \rho, with -1<\rho<1. A full conditional holds the other coordinate and the data x fixed. The Normal kernel \exp\{-(u-m)^2/(2v)\} identifies mean m and variance v>0. Superscript (t) labels the Gibbs sweep.

The joint density has kernel \pi(\theta_1,\theta_2\mid x) \propto \exp\left\{ -\frac{1}{2(1-\rho^2)} \left(\theta_1^2-2\rho\theta_1\theta_2+\theta_2^2\right) \right\}.

(a) Full conditional distributions

For the first full conditional, p(\theta_1\mid\theta_2,x) = \frac{p(\theta_1,\theta_2\mid x)}{p(\theta_2\mid x)} \propto p(\theta_1,\theta_2\mid x) \quad\text{as a function of }\theta_1. The denominator is constant when \theta_2 and x are held fixed. Thus the joint kernel is sufficient: terms involving \theta_1 must be kept, while terms free of \theta_1 may be absorbed into the proportionality constant.

We aim to match the Normal kernel \exp\left\{-\frac{(\theta_1-m)^2}{2v}\right\}, where m is the mean and v is the variance. Treat the quadratic form in the joint exponent as a quadratic in \theta_1. Completing the square gives \begin{aligned} \theta_1^2-2\rho\theta_1\theta_2+\theta_2^2 &=(\theta_1-\rho\theta_2)^2 -\rho^2\theta_2^2+\theta_2^2\\ &=(\theta_1-\rho\theta_2)^2 +(1-\rho^2)\theta_2^2. \end{aligned} Substitution into the joint kernel gives \begin{aligned} p(\theta_1\mid\theta_2,x) &\propto \exp\left\{ -\frac{(\theta_1-\rho\theta_2)^2+(1-\rho^2)\theta_2^2} {2(1-\rho^2)} \right\}\\ &= \exp\left\{-\frac{(\theta_1-\rho\theta_2)^2} {2(1-\rho^2)}\right\} \exp\left\{-\frac{\theta_2^2}{2}\right\}\\ &\propto \exp\left\{-\frac{(\theta_1-\rho\theta_2)^2} {2(1-\rho^2)}\right\}. \end{aligned} The factor \exp(-\theta_2^2/2) is dropped because it is constant with respect to \theta_1. Matching the remaining expression to the Normal kernel yields \boxed{ \theta_1\mid\theta_2,x \sim N\left(\rho\theta_2,\,1-\rho^2\right). } The condition -1<\rho<1 ensures that the conditional variance 1-\rho^2 is positive.

For the other coordinate, the same calculation gives \theta_2^2-2\rho\theta_1\theta_2+\theta_1^2 = (\theta_2-\rho\theta_1)^2+(1-\rho^2)\theta_1^2. Hence, by symmetry, \boxed{ \theta_2\mid\theta_1,x \sim N\left(\rho\theta_1,\,1-\rho^2\right). }

(b) Gibbs sampler

Starting from any valid pair (\theta_1^{(0)},\theta_2^{(0)})\in\mathbb{R}^2, repeat for t=1,\ldots,T:

1. Draw \theta_1^{(t)} \sim N\left(\rho\theta_2^{(t-1)},\,1-\rho^2\right).

2. Draw \theta_2^{(t)} \sim N\left(\rho\theta_1^{(t)},\,1-\rho^2\right).

3. Store (\theta_1^{(t)},\theta_2^{(t)}).

Each update is sampled directly from a full conditional distribution, so no accept–reject step is required.

(c) Equilibrium and independence

If the chain is started in equilibrium, every stored pair has the required joint distribution. The pairs are nevertheless generally not independent, because the distribution of the next pair depends on the current pair through the Gibbs updates.

When \rho=0, the coordinates are independent and fresh conditional draws produce independent pairs across full sweeps. This exception does not apply to a correlated target.

(d) Effect of strong correlation

When \rho is close to 1, the joint density is concentrated around the diagonal \theta_1=\theta_2, while the conditional variance 1-\rho^2 is close to zero. Each Gibbs update therefore makes only a small move along the narrow high-density region. This produces strong autocorrelation and slow mixing. Thus, although every Gibbs proposal is accepted, an acceptance rate of one does not imply efficient exploration of the target distribution.