The joint distribution for a hierarchical Bayes model

Question

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.)

Answer 1

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.