Q51. Three-variable Gibbs sampling

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

[10 marks]

Consider the joint probability mass/density function \begin{aligned} f(x,y,z)&\propto\binom nx\binom my z^{x+y+\alpha-1}(1-z)^{n+m-x-y+\beta-1}\\ &\qquad\times\gamma^{m-y}(1-\gamma)^{n-x} \delta^x(1-\delta)^y, \end{aligned} where x\in\{0,1,\ldots,n\},\qquad y\in\{0,1,\ldots,m\},\qquad 0<z<1. Here n,m are known positive integers, and \alpha,\beta>0 and 0<\gamma,\delta<1 are known constants.

(a) Derive the following full conditional distributions, identifying each by name and specifying its parameters and support: X\mid Y=y,Z=z,\qquad Y\mid X=x,Z=z,\qquad Z\mid X=x,Y=y. [6 marks]

(b) Using your answers to part (a), give pseudocode for a Gibbs sampler targeting f(x,y,z). Specify a valid initial state and make clear which iteration’s values are used in each update. [3 marks]

(c) Suppose the chain has reached equilibrium. Are successive sampled triples independent observations from f(x,y,z)? Briefly justify your answer. [1 mark]

No implementation or numerical simulation is required.

Worked answer (click to show)

Definitions. A weighted Binomial kernel \binom nx a^xb^{n-x} normalises to success probability a/(a+b) for a,b>0. A \operatorname{Beta}(A,B) kernel is z^{A-1}(1-z)^{B-1} on (0,1). Here X\in\{0,\ldots,n\}, Y\in\{0,\ldots,m\}, and 0<Z<1; n,m are positive integers, \alpha,\beta>0, and 0<\gamma,\delta<1. Superscript (t) indexes complete Gibbs sweeps.

(a) Full conditional distributions.

Conditional distribution of X. Holding y and z fixed, \begin{aligned} \mathbb{P}(X=x\mid Y=y,Z=z) &\propto\binom nx z^x(1-z)^{n-x}(1-\gamma)^{n-x}\delta^x\\ &=\binom nx(z\delta)^x\{(1-z)(1-\gamma)\}^{n-x}. \end{aligned} For positive weights a,b, the binomial theorem gives \sum_{j=0}^n\binom nj a^j b^{n-j}=(a+b)^n. Dividing the kernel by this sum yields a Binomial mass function with success probability a/(a+b). Here define p_X(z)=\frac{z\delta}{z\delta+(1-z)(1-\gamma)}. Consequently, \boxed{X\mid Y=y,Z=z\sim\operatorname{Binomial}\{n,p_X(z)\}, \qquad x\in\{0,\ldots,n\}.} The conditional does not depend on y.

Conditional distribution of Y. Holding x and z fixed, the y-dependent kernel is \begin{aligned} \mathbb{P}(Y=y\mid X=x,Z=z) &\propto\binom my z^y(1-z)^{m-y}\gamma^{m-y}(1-\delta)^y\\ &=\binom my\{z(1-\delta)\}^y\{(1-z)\gamma\}^{m-y}. \end{aligned} The normalising sum is \{z(1-\delta)+(1-z)\gamma\}^m. Thus, with p_Y(z)=\frac{z(1-\delta)}{z(1-\delta)+(1-z)\gamma}, we obtain \boxed{Y\mid X=x,Z=z\sim\operatorname{Binomial}\{m,p_Y(z)\}, \qquad y\in\{0,\ldots,m\}.} This conditional does not depend on x. In fact, the joint kernel for (X,Y) at fixed z factors into these two Binomial kernels, so X and Y are conditionally independent given Z.

Conditional distribution of Z. Holding x and y fixed leaves f(z\mid x,y)\propto z^{\alpha+x+y-1}(1-z)^{\beta+n+m-x-y-1},\qquad0<z<1. Matching the Beta kernel gives \boxed{Z\mid X=x,Y=y\sim \operatorname{Beta}(\alpha+x+y,\ \beta+n+m-x-y),\qquad0<z<1.} Both Beta parameters are positive. Also, 0<p_X(z),p_Y(z)<1 under the stated assumptions, so all three conditional distributions are valid.

(b) Gibbs sampler. Use the functions p_X and p_Y from part (a). A valid initial state is \boxed{(x^{(0)},y^{(0)},z^{(0)})=(0,0,1/2).} Choose a run length T. For t=1,\ldots,T, perform one complete sweep:

1. Update X using the previous Z: x^{(t)}\sim\operatorname{Binomial}\!\left( n,\frac{z^{(t-1)}\delta} {z^{(t-1)}\delta+(1-z^{(t-1)})(1-\gamma)}\right). This is the full conditional given y^{(t-1)},z^{(t-1)}; its formula does not depend on y^{(t-1)}.

2. Update Y using the previous Z: y^{(t)}\sim\operatorname{Binomial}\!\left( m,\frac{z^{(t-1)}(1-\delta)} {z^{(t-1)}(1-\delta)+(1-z^{(t-1)})\gamma}\right). The latest available conditioning values are x^{(t)},z^{(t-1)}; the formula does not depend on x^{(t)}.

3. Update Z using both newly sampled counts: z^{(t)}\sim\operatorname{Beta}\!\left( \alpha+x^{(t)}+y^{(t)},\ \beta+n+m-x^{(t)}-y^{(t)}\right).

4. Store the completed triple (x^{(t)},y^{(t)},z^{(t)}).

Use fresh random numbers for each conditional draw. Each full-conditional update leaves the target distribution invariant, and therefore so does the complete sweep. No accept/reject step is needed because the draws are directly from the full conditionals.

An arbitrary valid starting state is not necessarily an equilibrium draw; allow for convergence before treating the output as representative of the target distribution.

(c) Independence at equilibrium. No. At equilibrium each completed triple has the target distribution, but successive triples are dependent: the new counts are sampled using the preceding value of Z, and the next Z is then sampled using those counts. Reaching equilibrium does not remove this dependence.

Key result. Conditional independence of X and Y given Z is not independence between successive Gibbs-sampler triples. Stationary output has the correct distribution, but it is not an iid sample.