Are marginalized two-part models superior to non-marginalized two-part models for count data with excess zeroes? Estimation of marginal effects, model misspecification, and model selection

Abstract

The marginalized two-part models, including the marginalized zero-inflated Poisson and negative binomial models, have been proposed in the literature for modelling cross-sectional healthcare utilization count data with excess zeroes and overdispersion. The motivation for these proposals was to directly capture the overall marginal effects and to avoid post-modelling effect calculations that are needed for the non-marginalized conventional two-part models. However, are marginalized two-part models superior to non-marginalized two-part models because of their structural property? Is it true that the marginalized two-part models can provide direct marginal inference? This article aims to answer these questions through a comprehensive investigation. We first summarize the existing non-marginalized and marginalized two-part models and then develop marginalized hurdle Poisson and negative binomial models for cross-sectional count data with abundant zero counts. Our interest in the investigation lies particularly in the (average) marginal effect and (average) incremental effect and the comparison of these effects. The estimators of these effects are presented, and variance estimators are derived by using delta methods and Taylor series approximations. Though the marginalized models attract attention because of the alleged convenience of direct marginal inference, we provide evidence for the impact of model misspecification of the marginalized models over the conventional models, and provide evidence for the importance of goodness-of-fit evaluation and model selection in differentiating between the marginalized and non-marginalized models. An empirical analysis of the German Socioeconomic Panel data is presented.

Appendix: Gradients of marginal effects

1.1 Gradients of marginal effects in the ZIP and ZINB models

Recall that ZIP and ZINB have the identical expression of marginal expectation of response $y_i$: $\displaystyle \text{ E }(y_i|x_i,z_i) = \mu _i (1-\psi _i) = \frac{e^{x_i'\beta }}{1 + e^{{z'_i}\gamma }}$ and as a consequence share the same marginal and incremental effect formulas. The difference is the parameter $\theta$, where $\theta =(\beta ',\gamma ')'$ for ZIP models and $\theta =(\beta ',\gamma ',\alpha )'$ for ZINB models.

To simplify computation and notations, we introduce a pair of infinitely differentiable functions on the number line: $p^{\mathrm{ZIP}}(t)= e^t$ and $q(t)=\displaystyle \frac{1}{1+e^t}$ with ${\dot{p}}^{\mathrm{ZIP}}(t) = \ddot{p}^{\mathrm{ZIP}}(t) = p^{\mathrm{ZIP}}(t) = e^t$ and ${\dot{q}}(t) = -\displaystyle \frac{e^t}{(1 + e^t)^2}$, $\ddot{q}(t) = -\displaystyle \frac{e^t}{(1 + e^t)^2}\cdot \frac{1-e^t}{1 + e^t}$, for $\forall \, t\in {\mathbb {R}}$. Even “ZIP” is used in the superscript of function p and its derivatives, their expressions are exactly the same for ZINB. The $\theta$ and superscript of p will not be restated again in this subsection. The following discussions are identical for both ZIP and ZINB unless indicated otherwise.

Considering a continuous covariate $x_{ij}$ in our regression models, we adopt the simplified notations: $p^{\mathrm{ZIP}}_i, {\dot{p}}^{\mathrm{ZIP}}_i, \ddot{p}^{\mathrm{ZIP}}_i, q_i, {\dot{q}}_i,\ddot{q}_i$ which are $p^{\mathrm{ZIP}}, {\dot{p}}^{\mathrm{ZIP}}, \ddot{p}^{\mathrm{ZIP}}, q, {\dot{q}}, \ddot{q}$ evaluated at ${x'_i} \beta$ and ${z'_i}\gamma$, respectively. Then, the marginal mean of $y_i$ is $\displaystyle \text{ E }(y_i|x_i,z_i) = \mu _i (1-\psi _i) = p^{\mathrm{ZIP}}_iq_i$ and hence the marginal effect with respect to $x_{ij}$ is $\eta _j(x_i,z_i,\theta ) = \beta _j {\dot{p}}^{\mathrm{ZIP}}_iq_i + \gamma _j p^{\mathrm{ZIP}}_i{\dot{q}}_i$.

If the covariate $x_j$, or $z_j$, is categorical, to rewrite its incremental effect from level $l_1$ to $l_2$, the values of $p^{\mathrm{ZIP}}, {\dot{p}}^{\mathrm{ZIP}}, \ddot{p}^{\mathrm{ZIP}}$ at $x_{i(-j)}^\prime \beta _{(-j)} + l_2\beta _j$ and $x_{i(-j)}^\prime \beta _{(-j)} + l_1\beta _j$ will be denoted as $p_{2i}^{\mathrm{ZIP}}, {\dot{p}}_{2i}^{\mathrm{ZIP}}, \ddot{p}_{2i}^{\mathrm{ZIP}}$ and $p_{1i}^{\mathrm{ZIP}}, {\dot{p}}_{1i}^{\mathrm{ZIP}}, \ddot{p}_{1i}^{\mathrm{ZIP}}$, respectively; values of $q, {\dot{q}}, \ddot{q}$ at $z_{i(-j)}^\prime \gamma _{(-j)} + l_2\gamma _j$ and $z_{i(-j)}^\prime \gamma _{(-j)} + l_1\gamma _j$ will be represented by $q_{2i},{\dot{q}}_{2i} ,\ddot{q}_{2i}$ and $q_{1i},{\dot{q}}_{1i} ,\ddot{q}_{1i}$, respectively. Then, the incremental effect with respect to $x_{ij}$ is $\pi _j(x_{i(-j}),z_{i(-j)},\theta ) =p^{\mathrm{ZIP}}_{2i} q_{2i} -p^{\mathrm{ZIP}}_{1i} q_{1i}$.

The gradients of marginal and incremental effects are

$$\begin{aligned}&\displaystyle \nabla _\theta \eta _j(x_i,z_i,\theta )\nonumber \\&\quad = \displaystyle \left( \beta _j \ddot{p}^{\mathrm{ZIP}}_iq_i + \gamma _j{\dot{p}}^{\mathrm{ZIP}}{\dot{q}}\right) \sum \limits _{m = 0}^{J_1} x_{im}u_{(m + 1)}+\displaystyle {\dot{p}}^{\mathrm{ZIP}}_iq _i u_{(j + 1)} \nonumber \\&\displaystyle \qquad + \,\left( \beta _j {\dot{p}}^{\mathrm{ZIP}}_i{\dot{q}}_i + \gamma _j p^{\mathrm{ZIP}}_i\ddot{q}_i\right) \sum \limits _{m = 0}^{J_2}z_{im}u_{(J_1 + m + 2)} + p^{\mathrm{ZIP}}_i{\dot{q}}_i u_{(J_1 + j + 2)}, \nonumber \\&\qquad \nabla _\theta \pi _j(x_{i(-j)},z_{i(-k)},\theta ) \nonumber \\&\quad = \displaystyle \left( {\dot{p}}^{\mathrm{ZIP}}_{2i} q_{2i} -{\dot{p}}^{\mathrm{ZIP}}_{1i} q_{1i} \right) \sum \limits _{m = 0,\ne j}^{J_1} x_{im}u_{(m + 1)} + \left( l_2{\dot{p}}^{\mathrm{ZIP}}_{2i} q_{2i} -l_1{\dot{p}}^{\mathrm{ZIP}}_{1i} q_{1i} \right) \cdot u_{(j + 1)} \nonumber \\&\qquad \displaystyle \quad +\, \left( p^{\mathrm{ZIP}}_{2i} {\dot{q}}_{2i} -p^{\mathrm{ZIP}}_{1i} {\dot{q}}_{1i} \right) \sum \limits _{m = 0,\ne j}^{J_2} z_{im} u_{(J_1 + m + 2)} + \left( l_2p^{\mathrm{ZIP}}_{2i} {\dot{q}}_{2i} -l_1p^{\mathrm{ZIP}}_{1i} {\dot{q}}_{1i} \right) \cdot u_{(J_1 + j + 2)}, \end{aligned}$$

(29)

where $u_{(m)}$ is a unit vector of dimension $J_1+J_2+2$ for ZIP and dimension $J_1+J_2+3$ for ZINB with 1 in the mth component and 0 in others.

1.2 Gradients of marginal effects in the HP models

For HP models, we introduce functions: $p^{\mathrm{HP}}(t) =\displaystyle \frac{e^{t + e^t}}{e^{e^t}-1}$ and use the same q as ZIP. Then,

$$\begin{aligned}&p^{\mathrm{HP}}(t) = e^t + \sigma (t),\quad {\dot{p}}^{\mathrm{HP}}(t) = e^t + {\dot{\sigma }}(t),\quad \ddot{p}^{\mathrm{HP}}(t) = e^t + \ddot{\sigma }(t), \end{aligned}$$

where $\displaystyle \sigma (t) = \frac{e^t}{e^{e^t}-1}$, ${\dot{\sigma }}(t) =\sigma (t)\{1-e^t-\sigma (t)\}$, $\ddot{\sigma }(t) ={\dot{\sigma }}(t)\{1-e^t-2\sigma (t)\}-e^t\sigma (t).$

Using the similar notations for p, q and their derivatives as for ZIP and ZINB models in Sect. 1, the marginal mean is rewritten as $E(y_i|x_i,z_i) = p^{\mathrm{HP}}_iq_i$, the marginal effect with respect to continuous covariate $x_{ij}$ is $\eta _j(x_i,z_i,\theta ) = \beta _j {\dot{p}}^{\mathrm{HP}}_iq_i + \gamma _j p^{\mathrm{HP}}_i{\dot{q}}_i$, and the incremental effect with respect to categorical covariate $x_{ij}$ from level $l_1$ to level $l_2$ is $\pi _j(x_{i(-j}),z_{i(-j)},\theta ) =p^{\mathrm{HP}}_{2i} q_{2i} -p^{\mathrm{HP}}_{1i} q_{1i}$, where $\theta =(\beta ',\gamma ')'$. Then, the formulas of gradients of marginal and incremental effects are in the same forms as ZIP models (29) with different layouts of $p^{\mathrm{HP}}$ and its derivatives ${\dot{p}}^{\mathrm{HP}}$ and $\ddot{p}^{\mathrm{HP}}$.

1.3 Gradients of marginal effects in the HNB models

The parameter in HNB models is $\theta =(\beta ',\gamma ',\alpha )'$, and we adopt the same q function in ZIP, ZINB, and HP models but define a new function p by $p^{\mathrm{HNB}}(t,\alpha ) =\displaystyle \frac{e^t}{1-\rho (t,\alpha )}$, where $\rho (t,\alpha ) = \tau ^\alpha (t,\alpha )$, $\tau (t,\alpha ) =\displaystyle \frac{\alpha }{\alpha + e^t}$, and $\alpha >0$. We will use the same notations in terms of q and its derivatives evaluated at ${z'_i}\gamma$, $z_{i(-j)}^\prime \gamma _{(-j)} + l_2\gamma _j$ and $z_{i(-j)}^\prime \gamma _{(-j)} + l_1\gamma _j$, as introduced in previous sections.

With simple computation, we can get derivatives of $p^{\mathrm{HNB}}$ with respect to t and $\alpha >0$. In particular, ${\dot{p}}^{\mathrm{HNB}}_t(t,\alpha ) = p^{\mathrm{HNB}}\left( 1-p^{\mathrm{HNB}}\rho \tau \right)$, $\ddot{p}^{\mathrm{HNB}}_t(t,\alpha ) ={\dot{p}}^{\mathrm{HNB}}_t(t,\alpha ) \left( 1-2p^{\mathrm{HNB}}\rho \tau \right) + (\alpha + 1)\left( p^{\mathrm{HNB}}\right) ^2\rho \tau (1-\tau )$, $\displaystyle {\dot{p}}^{\mathrm{HNB}}_\alpha (t,\alpha ) = {p^{\mathrm{HNB}}\rho (\ln \tau + 1-\tau )}/{(1-\rho )}$, and $\ddot{p}^{\mathrm{HNB}}_{t\alpha }(t,\alpha ) =\big \{\displaystyle {\dot{p}}^{\mathrm{HNB}}_\alpha (t,\alpha )$ $\cdot (1-p^{\mathrm{HNB}}\rho \tau -p^{\mathrm{HNB}}\tau )\big \}-\left\{ (p^{\mathrm{HNB}})^2\rho e^t/(\alpha + e^t)^2\right\}$, where $\tau = \tau (t,\alpha )$, $\rho = \rho (t,\alpha )$ for simplicity of notations, ${\dot{\tau }}_t(t,\alpha ) = \tau (\tau -1)$, $\ddot{\tau }_{tt}(t,\alpha ) = \tau (\tau -1)(2\tau -1)$, ${\dot{\tau }}_\alpha (t,\alpha ) = {e^t}/{(\alpha + e^t)^2}$ ${\dot{\rho }}_t(t,\alpha ) = \alpha \rho (\tau -1) = -e^t\rho \tau$, $\ddot{\rho }_{tt}(t,\alpha ) = \alpha \rho (\tau -1)\{(\alpha + 1)\tau -\alpha \} = \rho \tau ^2e^t(e^t-1)$, ${\dot{\rho }}_\alpha (t,\alpha ) = \rho (\ln \tau + 1-\tau )$.

For functions $p^{\mathrm{HNB}}, {\dot{p}}_t^{\mathrm{HNB}}, \ddot{p}_{tt}^{\mathrm{HNB}}, {\dot{p}}^{\mathrm{HNB}}_{\alpha }, \ddot{p}^{\mathrm{HNB}}_{t\alpha }$ evaluated at fixed values of $({x'_i} \beta ,\alpha )$ are denoted by $p_i^{\mathrm{HNB}}, {\dot{p}}_{ti}^{\mathrm{HNB}}, \ddot{p}_{tti}^{\mathrm{HNB}}, {\dot{p}}^{\mathrm{HNB}}_{\alpha i}, \ddot{p}^{\mathrm{HNB}}_{t\alpha i}$, respectively. Values of $p^{\mathrm{HNB}}, {\dot{p}}_t^{\mathrm{HNB}}, {\dot{p}}^{\mathrm{HNB}}_{\alpha }$ at fixed values of $\big (x_{i(-j)}^\prime \beta _{(-j)} + l_2\beta _j,\alpha \big )$ and $\big (x_{i(-j)}^\prime \beta _{(-j)} + l_1\beta _j,\alpha \big )$ are denoted by $p^{\mathrm{HNB}}_{2i}, {\dot{p}}_{2ti}^{\mathrm{HNB}}, {\dot{p}}^{\mathrm{HNB}}_{2\alpha i}$, and $p^{\mathrm{HNB}}_{1i}, {\dot{p}}_{1ti}^{\mathrm{HNB}}, {\dot{p}}^{\mathrm{HNB}}_{1\alpha i}$, respectively.

By using p, q notations, the marginal mean of $y_i$ can be rewritten as $E(y_i|x_i,z_i) = p^{\mathrm{HNB}}_iq_i$, the marginal effect with respect to continuous covariate $x_{ij}$ is $\eta _j(x_i,z_i,\theta ) = \beta _j {\dot{p}}^{\mathrm{HNB}}_{ti}q_i + \gamma _j p^{\mathrm{HNB}}_i{\dot{q}}_i$, and the incremental effect with respect to categorical covariate $x_{ij}$ from level $l_1$ to level $l_2$ is $\pi _j(x_{i(-j)},z_{i(-j)},\theta ) =p^{\mathrm{HNB}}_{2i} q_{2i} -p^{\mathrm{HNB}}_{1i} q_{1i}$.

Then, the formulas of gradients of marginal and incremental effects are in the same forms as ZIP models (29) with different layouts of $p^\mathrm{HP}$ and its derivatives ${\dot{p}}^\mathrm{HP}$ and $\ddot{p}^\mathrm{HP}$. The gradients of effects with respect to parameter $\theta$ are

$$\begin{aligned}&\displaystyle \nabla _\theta \eta _j(x_i,z_i,\theta ) \\&\quad = \displaystyle + \left( \beta _j \ddot{p}^{\mathrm{HNB}}_{tti}q_i + \gamma _j{\dot{p}}^{\mathrm{HNB}}_{ti}{\dot{q}}_i\right) \sum \limits _{m = 0}^{J_1} x_{im}u_{(m + 1)} {\dot{p}}^{\mathrm{HNB}}_{ti}q _i u_{(j + 1)}\\&\quad \displaystyle \quad +\, \left( \beta _j {\dot{p}}^{\mathrm{HNB}}_{ti}{\dot{q}}_i + \gamma _j p^{\mathrm{HNB}}_i\ddot{q}_i\right) \sum \limits _{m = 0}^{J_2}z_{im}u_{(J_1 + m + 2)} + p^{\mathrm{HNB}}_i{\dot{q}}_i u_{(J_1 + j + 2)} , \\&\quad\quad \nabla _\theta \pi _j(x_{i(-j)},z_{i(-k)},\theta ) \\&\quad = \displaystyle \left( {\dot{p}}^{\mathrm{HNB}}_{2ti} q_{2i} -{\dot{p}}^{\mathrm{HNB}}_{1ti} q_{1i} \right) \sum \limits _{m = 0,\ne j}^{J_1} x_{im}u_{(m + 1)} + \left( l_2{\dot{p}}^{\mathrm{HNB}}_{2ti} q_{2i} -l_1{\dot{p}}^{\mathrm{HNB}}_{1ti} q_{1i} \right) \cdot u_{(j + 1)}\\&\quad \displaystyle \quad + \,\left( p^{\mathrm{HNB}}_{2i} {\dot{q}}_{2i} -p^{\mathrm{HNB}}_{1i} {\dot{q}}_{1i} \right) \sum \limits _{m = 0,\ne j}^{J_2} z_{im} u_{(J_1 + m + 2)} + \left( l_2p^{\mathrm{HNB}}_{2i} {\dot{q}}_{2i} -l_1p^{\mathrm{HNB}}_{1i} {\dot{q}}_{1i} \right) \cdot u_{(J_1 + j + 2)}\\&\quad\quad + \,\left( {\dot{p}}^{\mathrm{HNB}}_{2\alpha i}q_{2i}-{\dot{p}}^{\mathrm{HNB}}_{1\alpha i}q_{1i} \right) u_{(J_1 + J_2 + 3)}, \end{aligned}$$

where $u_{(m)}$ is a unit vector of dimension $J_1+J_2+3$ with 1 in the mth component and 0 in others.