= Solution
For <EM for zero-inflated negative binomial regression>, at iteration $k$ apply <Bayes' theorem> to each zero in the <zero-inflated negative binomial model>. The E step computes the <conditional expectation> of its latent indicator:
$$
\boxed{t_i^{(k)}=\mathbb E[A_i\mid y_i]=\begin{cases}\displaystyle\frac{\pi^{(k)}}{\pi^{(k)}+(1-\pi^{(k)})(r/(r+\lambda_i^{(k)}))^r},&y_i=0,\\0,&y_i>0.\end{cases}}
$$
Treating these values as fixed, the M step maximizes the expected complete-data <log-likelihood>. Its two parameter blocks separate:
$$
\boxed{\pi^{(k+1)}=\frac1n\sum_i t_i^{(k)},\qquad \beta^{(k+1)}=\arg\max_\beta\sum_i w_i^{(k)}\log f_{\rm NB}(y_i;r,e^{x_i^T\beta}),\quad w_i^{(k)}=1-t_i^{(k)}.}
$$
The second block is a weighted <negative binomial regression> with fixed size $r$. Writing $\eta_i=x_i^T\beta$ and $\lambda_i=e^{\eta_i}$, its parameter-dependent objective is
$$
Q_\beta=\sum_iw_i\{y_i\eta_i-(y_i+r)\log(r+e^{\eta_i})\}.
$$
Its <score function> and negative <Hessian matrix> are
$$
U(\beta)=\sum_iw_i\frac{r(y_i-\lambda_i)}{r+\lambda_i}x_i,\qquad -\nabla^2Q_\beta=\sum_iw_i\frac{r\lambda_i(y_i+r)}{(r+\lambda_i)^2}x_ix_i^T.
$$
The objective is <concave>, so a finite maximizer, when it exists, can be found by <Newton method> or <iteratively reweighted least squares>. The latter uses working response $\eta_i+(y_i-\lambda_i)/\lambda_i$ and weights $w_i r\lambda_i/(r+\lambda_i)$. This is not an unweighted regression on only the positive counts: zeros retain fractional weight in the count component.
Start with $0<\pi<1$ and finite coefficients, alternate the two steps, and monitor the observed <log-likelihood>. Exact maximization makes the <expectation-maximization algorithm> nondecreasing in that likelihood, although the mixture may have multiple stationary points. All-zero data can lead to boundary fits and poor identification; multiple initializations and explicit checks for boundary parameters are useful.
Back to article page
