# Fitting a Bayesian Hierarchical Poisson Regression in R

I'm trying to fit a Bayesian hierarchical poisson regression. To do so, I'm using MCMChpoisson function from MCMCpack in R. Based on this package, the model is:

$$Y_i \sim Poisson(\lambda_i)$$ $$\phi(\lambda_i) = X_i\beta + W_i \beta_i + \epsilon_i$$ $$\epsilon_i \sim N(0, \sigma^2 I_{k_i})$$ $$\dots$$

In the model above, $\phi$ is the link function.

I skipped the rest of the model as only the parts of above are related to my question. My question is why they consider an measurement error ($\epsilon_i$) in the systematic component whereas in GLM we have a function of the mean; in other words, sampling from poisson will itself generate a measurement error.

Also, I think the extra $\epsilon_i$ term above causes me to get strange results. Does anyone know any other function/package in R to fit a model very similar to the model above with no measurement error in the systematic component.

• FWIW, $\lambda$ for Poisson distribution is the mean. Apr 30, 2015 at 5:59
• $\varepsilon$ is an "overdispersion term". Described at the top of page 6 in this writeup. It doesn't seem like you can "turn it off" in MCMChpoisson. Apr 30, 2015 at 6:20
• If $\sigma=0$, the random effect gets away. Apr 30, 2015 at 8:02
• @ssdecontrol, one way to "turn it off" might be putting an informative prior on $\sigma^2_{\epsilon}$ so that it always remains very close to zero.
– Sam
Apr 30, 2015 at 14:38
• As @Xi'an pointed out, if $\sigma=0$, then $e_i=0$ for all $i$, i.e. it goes away. Real data is often overdispersed and thus typically including $\epsilon_i$ provides for a better model. May 1, 2015 at 15:01