Gibbs sampling in R
I have the following model:
$y_{it}=\alpha + x'{it}\beta{i} + \epsilon_{it}, \text{ } i=1,2,...,N, \text{ } t=1,2,...,T$ (1)
$\beta_{i}= z'{i}\gamma+\eta{i}$ (2)
with $\epsilon_{it} \sim N(0,\sigma_{\epsilon_{i}}^{2})$ and $\eta_{i} \sim N(0,\sigma_{\eta}^{2})$
Also i consider the following prior specifications:
$p(\beta,\gamma) \propto 1 \ p(\sigma_{\epsilon_{i}}^{2}) \propto \sigma_{\epsilon_{i}}^{-2} \ p(\sigma_{eta}^{2}) \propto \sigma_{\eta}^{-2}$
which has the following likelihood function:
$p(Y|\theta)=\prod_{i=1}^{N}\int_{\beta_{i}} \left( \prod_{t=1}^{T}\frac{1}{\sigma_{\epsilon_{i}}\sqrt{2\pi}}exp \left( -\frac{1}{2\sigma_{\epsilon_{i}}^{2}}(y_{it}-\alpha-x'{it}\beta{i})^{2}\right) \right) \times \frac{1}{\sigma_{\eta}\sqrt{2\pi}}exp \left( -\frac{1}{2{\sigma_{\eta}^{2}}}(\beta_{i}-z'{i}\gamma)^{2} \right) d\beta{i}$
I would like to estimate the parameters $\theta=(\alpha,\gamma,\sigma_{\epsilon_{i}}^{2},\sigma_{\eta}^{2})$. So $\theta$ is the parameter vector. Only $\beta_{i}$ is random. So I consider $\beta_{i}$ as the latent variables.
I would like to estimate these by using the following Gibbs sampling schema:
sample $\alpha$ given ${\beta_{i}}{i=1}^{N},\gamma,\sigma{\epsilon}^{2},\sigma_{\eta}^{2},Y$
sample $\gamma$ given ${\beta_{i}}{i=1}^{N},\alpha,\sigma{\epsilon}^{2},\sigma_{\eta}^{2},Y$
sample $\sigma_{\epsilon}^{2}$ given ${\beta_{i}}{i=1}^{N},\gamma,\alpha,\sigma{\eta}^{2},Y$
sample $\sigma_{\eta}^{2}$ given ${\beta_{i}}{i=1}^{N},\gamma,\alpha,\sigma{\epsilon}^{2},Y$
sample ${\beta_{i}}{i=1}^{N}$ given $\gamma,\alpha,\sigma{\epsilon}^{2},\sigma_{\eta}^{2},Y$
I know that:
To sample $\alpha$ I re-write (1) as $y_{it} - x'{it}\beta{i}=\alpha + \epsilon_{it}$, so the full conditional distribution of $\alpha$ is normal with mean $\hat{\alpha}=\frac{1}{T\times N} \sum_{i=1}^{N}\sum_{t=1}^{T} (y_{it} - x'{it}\beta{i})$ and variance $\frac{\sigma_{\epsilon}^{2}}{T \times N}$.
To sample $\sigma_{\epsilon}^{2}$ I consider the regression model in (1) , hence it can be sampled from an inverted Gamma-2 distribution with parameter $ \sum_{i=1}^{N}\sum_{t=1}^{T} (y_{it} - \alpha - x'{it}\beta{i})^{2}$ with $N \times T$ degrees of freedom. 402