# The joint distribution for a hierarchical Bayes model

Consider the following hierarchical Bayes model. We observe random variables $$X_1, X_2, ..., X_n$$ conditionally independent and having Poisson distributions $$X_i \sim Poiss(m_i\theta_i)$$, where $$m_i$$ are known multipliers and $$\theta_i$$ are unknown parameters. Assume that $$\theta_1, ..., \theta_n$$ are independet random variables with the same distribution $$Gamma(\alpha, \lambda)$$. Parameter $$\lambda$$ is a random variable with $$Gamma(\beta, \sigma)$$ distribution. Hyperparameters $$\alpha, \beta, \sigma$$ are known.

I would like to compute the the joint distribution.

Is the following computation correct?

$$P(X_1, X_2, ..., X_n|\theta_1, \theta_2, ..., \theta_n) P(\theta_1, \theta_2, ..., \theta_n| \lambda) P(\lambda) =$$ $$= \prod \limits_{i=1}^n P(X_i|\theta_1, \theta_2, ..., \theta_n) \prod \limits_{i=1}^n P(\theta_i|\lambda) P(\lambda) =$$ $$\prod \limits_{i=1}^n \frac{(m_i\theta_i)^{x_i}}{x_i!} e^{-m_i\theta_i} \prod \limits_{i=1}^n \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \theta_i^{\alpha-1}e^{-\lambda\theta_i} \frac{\sigma^{\beta}}{\Gamma(\beta)} \lambda^{\beta-1} e^{-\sigma \lambda}$$

(Next step is to propose Gibbs sampler.)

• How do you define the likelihood? I ask because the likelihood is typically expressed as $p(\bf x | \Theta)$ where $\bf x$ represents the observations and $\Theta$ represents the model parameters. Jan 18, 2020 at 19:16
• I was hoping to follow the definition you mentioned. Next step is to propose Gibbs sampler. Jan 18, 2020 at 19:24
• I see. To implement the Gibbs sampler (in my experience) we often use the joint distributions of the observations and the parameters, $p( \bf x , \Theta)$, which looks like what you've found. Jan 18, 2020 at 20:03
• The next step would be to obtain the full conditional distributions for $\theta_1, \ldots, \theta_n$, and $\lambda$. Jan 18, 2020 at 20:07
• I see. You are of course right :) Thanks :) I will edit the post. Jan 18, 2020 at 20:45

You are on the right way. Observe that

$$p(\theta_{1:n}, \lambda \vert x_{1:n}) = \prod \limits_{i=1}^n \frac{(m_i\theta_i)^{x_i}}{x_i!} e^{-m_i\theta_i} \prod \limits_{i=1}^n \frac{\lambda^{\alpha}}{\Gamma(\alpha)} \theta_i^{\alpha-1}e^{-\lambda\theta_i} \frac{\sigma^{\beta}}{\Gamma(\beta)} \lambda^{\beta-1} e^{-\sigma \lambda} \\ \propto \prod \limits_{i=1}^n (m_i\theta_i)^{x_i} e^{-m_i\theta_i} \prod \limits_{i=1}^n \lambda^{\alpha} \theta_i^{\alpha-1}e^{-\lambda\theta_i} \lambda^{\beta-1} e^{-\sigma \lambda}.$$

To sample $$\lambda | \theta_{1:n}$$, it is easy to see that $$p(\lambda \vert x_{1:n}, \theta_{1:n}) \propto \prod \limits_{i=1}^n \lambda^{\alpha} e^{-\lambda\theta_i} \lambda^{\beta-1} e^{-\sigma \lambda} = \lambda^{\alpha+\beta-1} e^{-\lambda(\sigma + \sum_{i=1}^n \theta_i)},$$ which is a Gamma distribution with parameters $$(\alpha + \beta, \sigma + \sum_{i=1}^n \theta_i)$$. Note that the data does not inform us about $$\lambda$$.

Conversely, to sample $$\theta_{1:n} | \lambda, x_{1:n}$$ we have

$$p(\theta_{1:n} | \lambda, x_{1:n}) \propto \prod \limits_{i=1}^n (m_i\theta_i)^{x_i} e^{-m_i\theta_i} \theta_i^{\alpha-1}e^{-\lambda\theta_i} = \frac{1}{\prod_{i=1}^n m_i^{\alpha-1}} \prod \limits_{i=1}^n (m_i\theta_i)^{x_i} e^{-m_i\theta_i} (m_i\theta_i)^{\alpha-1}e^{-\lambda\theta_i} \\ \propto \prod \limits_{i=1}^n (m_i\theta_i)^{x_i + \alpha - 1} e^{-(m_i + \lambda)\theta_i} \\ \propto \prod \limits_{i=1}^n \theta_i^{x_i + \alpha - 1} e^{-(m_i + \lambda)\theta_i},$$ such that each $$\theta_i | \lambda, x_i$$ is proportional to a Gamma distribution with parameters $$(x_i + \alpha, m_i + \lambda)$$. These can be sampled independently across multiple processors.

• Thanks so much! Unfortunately I do not have high enough reputation to give you an UpVote :) Jan 18, 2020 at 21:26
• No worries, just make sure that you verify I've gotten it anywhere near correct :-) Jan 18, 2020 at 21:34
• It's very typical to confuse the product notation: $$\prod_{i=1}^n c = c^n$$ Using some parentheses will reveal your mistakes. Jan 18, 2020 at 23:16