[R-sig-ME] glmmTMB: how to calculate posterior prob. of structural zero?

Ben Bolker bbolker @ending from gm@il@com
Wed Oct 3 21:01:09 CEST 2018


  If you're fitting a zero-inflated Poisson model, it would be
especially easy because the zero probability is simply exp(-lambda), so
I believe the "posterior"(ish) probability that the zero is due to the
structural rather than conditional part of the model would be

 P(zero|structural)/P(zero) =
P(zero|structural)/(P(zero_structural)+P(zero_conditional))

or

prob_struc <- function(model) {
   strucprob <- predict(model, type="zprob")
   condprob <- exp(-predict(model, type="conditional"))
   return(strucprob/(strucprob+condprob))
}

For a negative binomial ("nbinom2", or quadratic parameterization) the
zero probability is (k/(k+mu))^k (I think ... e.g.
dnbinom(0,mu=0.5,size=0.5) is sqrt(0.5)).  sigma(model) returns the
overdispersion parameter k, so a function as above but with

  k <- sigma(model)
  condprob <- (k/(k+predict(model,type="conditional")))^k

 inserted should do the trick.  (Please test these yourself and make
sure they are sensible before proceeding!)

On 2018-10-03 02:18 PM, Aaron Mackey wrote:
> I'm happily using glmmTMB to fit zero-inflated count models with my data,
> but I'd like to also know which zeroes in my data are more likely (or not)
> to be structural vs. expected from the conditional distribution. I know how
> to use Bayes formula to calculate the posterior, and predict(zinb,
> type="zprob") gives me the prior probabilites for each data point being
> structural or not (respecting the zero inflation part of the model), and
> the likelihoods for the structural components are 1 (if the data point is a
> zero) or 0 (if the data point is not a zero) -- but is there a way to
> extract the likelihood for each zero data point with respect to the
> conditional part of the model?
> 
> thanks,
> -Aaron
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>



More information about the R-sig-mixed-models mailing list