I fitted a Bayesian hurdle model assuming a Poisson distribution:

$ P(Y_i=y_i) = \begin{cases} 1-\pi_{i} & \mbox{if }y_i=0\\[0.5em] \pi_{i} \frac{\lambda_{i}^{y_i} e^{-\lambda_i}}{y_i!(1-e^{-\lambda_i})} &\mbox{if }y_i > 0 \end{cases} $

Because the log-likelihood is separable with respect to $\pi$ and $\lambda$, I fitted two separate regressions (i.e., logistic and zero-truncated Poisson regression) using a set of predictors. In other words:

$ \lambda = \log(\beta_{1}^T X) \\ \pi = \mbox{logit}(\beta_{2}^T X) $

My issue arises in trying to predict new values. In a Bayesian context, I have posterior distributions for the coefficients from the logistic and zero-truncated Poisson regressions. For a new predicted value, I will have a posterior distribution for $\pi$ and $\lambda$. To get the $E(\pi \cdot \lambda)$ is easy, I can just take the product of the posterior means. But what about getting the 2.5 and 97.5% percentiles for the credible interval around for the new predicted value?


1 Answer 1


Let's call the posterior predictive distribution $\tilde{Y}$.

$$\begin{align}\tilde{Y}_{\pi} &\sim\text{Bernoulli}(\pi)\\ \tilde{Y}_{\lambda} &\sim\text{Poisson}(\lambda)\\ \tilde{Y} &\sim\tilde{Y}_{\pi}*\tilde{Y}_{\lambda}\end{align}$$

You would need to truncate $\tilde{Y}_{\lambda}$ to be non-zero. Since you would then have a full predictive distribution, it would be easy to calculate intervals, expectations, etc.

Computation wise, this could be implemented in a number of ways. One possibility (in pseudo-code):

y_new = 0
y_pi = bernoulli_rng(pi)
    y_new = poisson_rng(lambda)
return y_new
  • $\begingroup$ Nice response @c-r-peterson! I guess I should have clarified this, but my real issue is computational. Both $\pi$ and $\lambda$ will have distributions, so you would need to draw from both to construct the posterior predictive distribution of $\tilde{Y}$. That can be a lot, perhaps just random sampling would work? $\endgroup$
    – user13317
    Sep 9, 2016 at 13:18
  • $\begingroup$ Assuming you're using an MCMC method to calculate your posterior, then you could just use the samples from that. This would be fairly simple to implement in Stan, for example. $\endgroup$ Sep 9, 2016 at 15:06

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.