# Gibbs sampling deriving complete conditionals with mixture priors

My question is about the derivation of the complete conditionals for Gibbs sampling in a hierarchical model where some of the parameters are mixtures of point-masses and Normal distributions. The model is from a research paper but I've slightly simplified the notation to try and make my question easy to follow. The 4 sets of parameters that I am interested in are $\theta_{bj}$, $\pi_b$, $\mu_b$ and $\sigma^2_b$ which have the following distributions:

$$\theta_{bj} \sim \pi_b I_{[\theta_{bj} = 0]} + (1 - \pi_b) I_{[\theta_{bj} \ne 0]}N(\mu_b, \sigma^2_{b})$$

$$\pi_{b} \sim Beta(\alpha_\pi, \beta_\pi)$$

$$\mu_{b} \sim N(\mu_0, \tau^2_0)$$

$$\sigma^2_{b} \sim IG(\alpha_0, \beta_0)$$

$$b = 1,...,B; j = 1,..., k_b$$ where $I$ is the indicator function, $IG$ is the inverse-gamma distribution, $\theta_{bj}$ is a mixture of a point-mass at zero and a normal distribution.

The parameters $\mu_b$, $\pi_b$ and $\sigma^2_b$ appear in the joint posterior distribution as follows:

\begin{split} \prod_{b = 1}^B \prod_{j = 1}^{k_b} & \left[ \pi_b I_{[\theta_{bj} = 0]} + (1 - \pi_b) I_{[\theta_{bj} \ne 0]} N(\mu_b, \sigma^2_b) \right ]\\ & \times \prod_{b=1}^B \left[ \frac{1}{\sqrt{2\pi \sigma^2_b}}\exp\left(-\frac{(\mu_b - \mu_0}{2\sigma_b^2}\right)\right ]\\ & \times \prod_{b=1}^B \left(\sigma^2_{b}\right)^{-\alpha_{0} - 1} \exp \left( -\frac{\beta_{0}}{\left(\sigma^2_{b}\right)} \right) \\ & \times \prod_{b=1}^B \left[ \pi_b^{\alpha_\pi - 1}(1 - \pi_b)^{\beta_\pi - 1} \right ]\\ \end{split}

The complete conditionals given, but not derived, in the paper for $\mu_b$, $\pi_b$ and $\sigma^2_b$ are:

$$[\pi_b |...] \sim Beta\left(\alpha_\pi + \sum_{j = 1}^{k_b}I_{[\theta_{bj} = 0]}, \beta_\pi + k_b - \sum_{j = 1}^{k_b} I_{[\theta_{bj = 0]}}\right)$$

$$[\mu_{b} | ... ] \sim N\left(\frac{\mu_{ 0}\sigma^2_{b} + \tau^2_{0} \sum_{j = 1}^{k_b} \theta_{bj} }{\sigma^2_{b} + k_b\tau^2_{ 0}}, \frac{\tau^2_{ 0}\sigma^2_{ b}}{\sigma^2_{ b} + k_b\tau^2_{ 0}} \right)$$

$$[\sigma^2_{ b} | ...] \sim IG\left(\alpha_{0} + \frac{k_b}{2}, \beta_{0} + \frac{1}{2} \sum_{j = 1}^{k_b}(\theta_{bj} - \mu_{ b})^2\right)$$

The derivation for $\pi_b$ appears to be straightforward:

$$\begin{split} [\pi_b | ...] \propto & \left[ \pi_b^{\alpha_\pi - 1}(1 - \pi_b)^{\beta_\pi - 1} \right ]\prod_{j = 1}^{k_b} \left[ \pi_b I_{[\theta_{bj} = 0]} + (1 - \pi_b) I_{[\theta_{bj} \ne 0]}N(\mu_\theta b, \sigma^2_\theta b) \right ]\\ \propto & \left[ \pi_b^{\alpha_\pi - 1}(1 - \pi_b)^{\beta_\pi - 1} \right ] \left[ \pi_b^{\sum_{j=1}^{k_{bh}}I_{[\theta_{bj} = 0]}} (1 - \pi_b)^{\sum_{j = 1}^{k_b}I_{[\theta_{bj} \ne 0]}} \right ]\\ = & \pi_b^{\alpha_\pi + \sum_{j=1}^{k_{bh}}I_{[\theta_{bj} = 0]} - 1}(1 - \pi_b)^{\beta_\pi + \sum_{j = 1}^{k_b}I_{[\theta_{bj} \ne 0]} - 1} \\ = & \pi_b^{\alpha_\pi + \sum_{j=1}^{k_{bh}}I_{[\theta_{bj} = 0]} - 1}(1 - \pi_b)^{\beta_\pi + k_b - \sum_{j = 1}^{k_b}I_{[\theta_{bj} = 0]} - 1} \\ \\ [\pi_b | ...] \sim & Beta\left(\alpha_\pi + \sum_{j = 1}^{k_b}I_{[\theta_{bj} = 0]}, \beta_\pi + k_b - \sum_{j = 1}^{k_b} I_{[\theta_{bj = 0]}}\right)\\ \end{split}$$

The issue I have is that it looks to me that the complete conditionals for $\mu_b$ and $\sigma^2_b$ are derived using standard conjugate results for the normal mean and normal variance. Looking at the mean $\mu_b$ we have:

$$\begin{split} [\mu_b | ... ] \propto & \left[ \frac{1}{\sqrt{2\pi \sigma^2_b}}\exp\left(-\frac{(\mu_b - \mu_0}{2\sigma_b^2}\right)\right ] \prod_{j = 1}^{k_b} \left[ N(\mu_b, \sigma^2_b) \right ]\\ \end{split}$$ giving the answer quoted in the paper $$\begin{split} [\mu_b | ... ] \sim & N\left(\frac{\mu_{ 0}\sigma^2_{b} + \tau^2_{0} \sum_{j = 1}^{k_b} \theta_{bj} }{\sigma^2_{b} + k_b\tau^2_{ 0}}, \frac{\tau^2_{ 0}\sigma^2_{ b}}{\sigma^2_{ b} + k_b\tau^2_{ 0}} \right) \end{split}$$

It looks to me that the possibility that $\theta_{bj} = 0$ is not being handled correctly (i.e. assumed to be from a Normal distribution when in fact it is from the point-mass). I would therefore expect that the complete conditional for $\mu_b$ to be:

$$\begin{split} [\mu_b | ... ] & \propto \left[ \frac{1}{\sqrt{2\pi \sigma^2_b}}\exp\left(-\frac{(\mu_b - \mu_0}{2\sigma_b^2}\right)\right ] \prod_{j = 1}^{k_b} & \left[ \pi_b I_{[\theta_{bj} = 0]} + (1 - \pi_b)I_{[\theta_{bj} \ne 0]} N(\mu_\theta b, \sigma^2_\theta b) \right ]\\ & = \left[ \frac{1}{\sqrt{2\pi \sigma^2_b}}\exp\left(-\frac{(\mu_b - \mu_0}{2\sigma_b^2}\right)\right ] \prod_{j = 1}^{k_b} & \left[ \pi_b^{I_{[\theta_{bj} = 0]}}\left[ (1 - \pi_b) N(\mu_\theta b, \sigma^2_\theta b) \right]^{I_{[\theta_{bj} \ne 0]}} \right ] \end{split}$$

giving

$$[\mu_b|...] \sim N\left(\frac{\mu_{ 0}\sigma^2_{b} + \tau^2_{0} \sum_{j = 1}^{k_b} \theta_{bj} }{\sigma^2_{b} + K_b\tau^2_{ 0}}, \frac{\tau^2_{ 0}\sigma^2_{ b}}{\sigma^2_{ b} + K_b\tau^2_{ 0}} \right), \; K_b = \sum_{j =1}^{k_b} I_{[\theta_{bj} \ne 0]}$$

with something similar for $\sigma^2_b$.

My question is then - which complete conditional is correct for $\mu_b$, is the one quoted in the paper actually correct (and I have made an error), or does the paper contain a typo and the actual complete conditional is the one I derived?

• The likelihood for $\theta_{bj}$ does not have a $I(\theta_{bj} \ne 0)$ next to the normal. So the likelihood just adds a point mass of $0$ at $\theta_{bj} = 0$. It is not a mixture distribution. – Greenparker May 19 '16 at 12:07
• I meant point mass of $\pi_{b}$ at 0. – Greenparker May 19 '16 at 12:16
• OK - I've added in the indicator for the second part of the $\theta_{bj}$ distribution to remove any doubts. – Ray May 19 '16 at 13:08
• Can you tell us which research paper you are referring to? This will give some context. – Greenparker May 22 '16 at 16:40
• I think the question contains all the relevant information for the issue I have, for completeness though the original paper is: Berry , SM & Berry, DA. (2004).“ Accounting for Multiplicities in Assessing Drug Safety: A Three-Level Hierarchical Mixture Model” Biometrics – Ray May 23 '16 at 9:29