5
$\begingroup$

I was looking at the Gibbs Sampler when I stumbled upon the following example:

Suppose $y = (y_{1}, y_{2}, \ldots, y_{n})$ are iid observations from an $N(\mu, \tau^{-1})$

Furthermore, suppose there exists prior information on $\mu$ and $\tau$ indicating that they are independent and that $\mu$ follows an $N(\mu_{0}, \sigma_{0}^{2})$ distribution while $\tau$ follows a $\text{Ga}(\alpha_{0}, \beta_{0})$ distribution.

For this example,

$$ p(\mu, \tau|y_{1}, y_{2}, \ldots, y_{n}) \propto f(y|\mu, \tau)p(\mu, \tau)$$

Also,

$$ p(\mu|\tau, y) \sim N(\frac{n\bar{y}\tau+\mu_{0}\tau_{0}}{n\tau+\tau_{0}}, (n\tau+\tau_{0})^{-1}) $$

where $\tau_{0} = (\sigma_{0}^{2})^{-1}$

and

$$ p(\tau|\mu, y) \sim \text{Ga}(\alpha_{0}+\frac{n}{2}, \beta_{0}+\frac{S_{n}}{2}) $$

where $S_{n} = \sum_{i=1}^{n}(y_{i}-\mu)^{2}$

Because it is not easy to compute directly from this distribution, the Gibbs sampler can be used provided the conditional distributions are known.

Could anybody demonstrate how (provide the derivations for) the conditional distributions ($p(\mu|\tau, y)$ and $p(\tau|\mu, y)$) given above?

EDIT:

For example, to achieve $p(\mu|\tau, y)$, is the following valid?

$$ p(\mu|\tau, y) \propto p(\tau, y|\mu)p(\mu)$$

If so, what form will $p(\tau, y|\mu)$ have?

$\endgroup$
  • $\begingroup$ Hello, is your $S_n$ supposed to be $\sum_{i=1}^n (y_i - \mu)^2$? $\endgroup$ – Sam Livingstone Nov 9 '12 at 10:59
  • $\begingroup$ @SamLivingstone: Yes, it is. Thanks for informing me. I've fixed it now. $\endgroup$ – user9171 Nov 9 '12 at 11:05
10
$\begingroup$

Ok, what you need to do is compute the joint posterior up to a constant, i.e. $f(y_1,...,y_n|\mu,\tau)p(\mu,\tau)$. Then to compute the conditional posterior $\pi(\mu|\mathbf{y},\tau)$ you just treat the $\tau$ terms as fixed and known, so that some of them can be cancelled out. Then you do the same thing with $\mu$ to get $\pi(\tau|\mathbf{y},\mu)$.

So what I mean is, the joint posterior is given by: \begin{equation} \pi(\mu,\tau|\mathbf{y}) \propto \tau^{\frac{n}{2}+\alpha_0-1} \exp\{-\frac{\tau}{2} \sum_{i=1}^{n}(y_i-\mu)^2-\frac{\tau_0}{2}(\mu-\mu_0)^2-\beta_0\tau\}. \end{equation}

Now, to get the conditional posterior $\pi(\tau|\mathbf{y},\mu)$ you just remove the $\mu$ terms that just multiply the expression. So we could re-write the joint posterior as: \begin{equation} \pi(\mu,\tau|\mathbf{y}) \propto \tau^{\frac{n}{2}+\alpha_0-1} \exp\{-\frac{\tau_0}{2}(\mu-\mu_0)^2\}\exp\{-\frac{\tau}{2} \sum_{i=1}^{n}(y_i-\mu)^2-\beta_0\tau\}. \end{equation}

Now, since we are treating $\mu$ as fixed and known, for the conditonal posterior $\pi(\tau|\mu,\mathbf{y})$ we can remove the first exponential term from the equation and it will still be true (owing to the $\propto$ sign rather than the equals sign). It's important that you get this: the conditional posterior says given that we know $\mu$ what is $\tau$, so $\mu$ is known (and hence fixed). So the kernel for the conditional posterior for $\tau$ becomes (after some re-arranging): \begin{equation} \pi(\tau|\mu,\mathbf{y}) \propto \tau^{\frac{n}{2}+\alpha_0-1}\exp\{-\tau(\frac{1}{2} \sum_{i=1}^{n}(y_i-\mu)^2+\beta_0)\}, \end{equation} which is the kernel for the Gamma distribution with the parameters you state.

Now, for $\pi(\mu|\tau,\mathbf{y})$ you do the same thing, but it's a bit trickier because you have to complete the square. After cancelling the relevant terms not involving $\mu$ you get: \begin{equation} \pi(\mu|\tau,\mathbf{y}) \propto \exp\{ -\frac{\tau}{2} \sum_{i=1}^{n} (y_i-\mu)^2 - \frac{\tau_0}{2}(\mu-\mu_0)^2\}. \end{equation} If you multiply out the squared brackets and take the relevant summations, then remove all terms not involving $\mu$ (as they are all just multiplying the kernel by a constant) this becomes: \begin{equation} \pi(\mu|\tau,\mathbf{y}) \propto \exp\{ -\frac{1}{2}[\mu^2(n\tau+\tau_0) - 2\mu(\tau n\bar{y} - \tau_0\mu_0)]\}. \end{equation} If you complete the square here you get the kernel for a Gaussian with the mean and variance you state.

Hope that helps.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Thank you so much! That's a superb answer. I really appreciate it. $\endgroup$ – user9171 Nov 9 '12 at 12:02
  • $\begingroup$ :-) no problem. $\endgroup$ – Sam Livingstone Nov 9 '12 at 12:17
  • 1
    $\begingroup$ Yeah, so what you have is $a\mu^2 + b\mu$. But the $c$ term is constant with respect to $\mu$, so you can add or subtract it as you please without losing the validity of the expression. If you have $\exp\{ a\mu^2 + b\mu + c\} = \exp\{c\}\exp\{ a\mu^2 + b\mu\}$ then you can remove and add in the $\exp\{c\}$ as you like, since it is fixed. Does that make sense/help? $\endgroup$ – Sam Livingstone Nov 9 '12 at 13:56
  • 1
    $\begingroup$ So to clarify, if you write it in the form $a(\mu + h)^2$, you will end up with an extra constant term at the end, but that's ok because we just need something proportional to what we had. $\endgroup$ – Sam Livingstone Nov 9 '12 at 13:59
  • 1
    $\begingroup$ Hello. For the uninformative priors you want $\alpha_0 = 1$, not $0.00001$. You can't determine what the informative priors should be the way you said (or in any way like that). In practice you would choose $\alpha_0$ and $\beta_0$ so that the shape of the prior reflected your beliefs. For simulation purposes, I think they just mean try some different values of $\alpha_0$ and $\beta_0$ that will give you a less flat prior and see how your posterior changes. $\endgroup$ – Sam Livingstone Nov 22 '12 at 17:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.