How does PSCL fit zero-inflated negative binomial models?

Question

I've been using the PSCL zero-inflated negative binomial model to simply fit count data. I'm not specifying any complex formula, essentially I just provide the count data and try to find the parameters of the underlying zero-inflated NB distribution (that I think) is generating the obtained counts. Like below,

pscl::zeroinfl(formula = counts ~ 1, data = df, dist = "negbin")

I am trying to better understand how this fitting works, and I was specifically wondering how the PSCL package computes the parameters of a ZINB model? Does it use MLE or MOM? Would appreciate any pointers on this, and also even any links to the exact solutions written out for the fit. Thanks!

Answer 1

The accompanying paper from the Journal of Statistical Software (doi:10.18637/jss.v027.i08) also provided as vignette("countreg", package = "pscl") says in Section 2.3 about zero-inflation models:

Again, ML estimates of all parameters are obtained from optim(), with control options set in zeroinfl.control() and employing analytical gradients. Starting values can be user-supplied, estimated by the expectation maximization (EM) algorithm, or by glm.fit() (the default). The covariance matrix estimate is derived numerically using the Hessian matrix returned by optim(). Using EM estimation for deriving starting values is typically slower but can be numerically more stable. It already maximizes the likelihood, but a single optim() iteration is used for determining the covariance matrix estimate. See Appendix B for further technical details.

Thus, in your case this essentially corresponds to computing the mean count and the proportion of zero observations as starting values and then optimizing the zero-inflated negative binomial log-likelihood.