Inflated Ordered Probit and Logit Models

Estimation, inference, and quantities of interest for ordered probit and ordered logit models whose outcome contains an inflated category: a single ordered category (bottom, middle, top, or any other) that mixes observations generated by the ordered process with observations generated by a distinct split-population process. Fits the zero-inflated ordered probit of Harris and Zhao (2007) and its middle- and top-inflated extensions (Bagozzi and Mukherjee 2012 ; Bagozzi, Hill, Moore and Mukherjee 2015 ; Bagozzi, Joo and Mukherjee 2024 ), generalized to an arbitrary inflated category and to the logit link, with optional correlated errors for the probit form, plus the standard ordered probit and logit and their partial proportional-odds (non-parallel) variants on the same footing. Provides analytic, robust, and cluster-robust standard errors, survey weights and offsets, model comparison (Vuong, likelihood-ratio, information criteria), regime-specific predicted probabilities and first differences, simulation for residual diagnostics, and tidy/table-package integration. The likelihood, its gradient, and the bivariate-normal probabilities are implemented in C++.

iop: Inflated Ordered Probit and Logit Models

Documentation site: https://bagozzib.github.io/iop/.

iop estimates ordered probit and ordered logit models whose outcome contains an inflated category -- a single ordered category (bottom, middle, top, or any other) that mixes observations generated by the ordered process with observations generated by a distinct split-population process that places them in that category regardless of the ordered mechanism. It fits the zero-inflated ordered probit of Harris and Zhao (2007) and its middle- and top-inflated extensions (Bagozzi and Mukherjee 2012; Bagozzi, Hill, Moore and Mukherjee 2015; Bagozzi, Joo and Mukherjee 2024), generalized to any inflated category and to the logit link, alongside the plain ordered models on the same engine. It provides:

• oprobit(), ologit() -- the standard ordered models, with survey weights, offsets, partial proportional-odds (non-parallel) effects (parallel =, parallel_test(), parallel = "auto"), and analytic, robust, or cluster-robust standard errors;
• iop() -- the inflated ordered probit for any single inflated category (inflate = "bottom" | "middle" | "top" or a category label), with optional correlated errors (correlated = TRUE: ZiOPC / MiOPC / TiOPC) and optional category-specific split equations (split = "category": the generalised GZiOP / GMiOP of Brown, Harris and Spencer 2020, tested against the common split by split_test());
• iol() -- the inflated ordered logit counterpart;
• panel tools on every model: unit random intercepts (re =, adaptive Gauss--Hermite quadrature, optionally in both equations), unit fixed effects (fe =, with the split-panel jackknife fe_correction = "jackknife"), the mundlak() correlated-random-effects device, and cluster-robust standard errors;
• one set of methods for all four: summary(), predict() by regime and by type of zero (with delta-method se.fit), first_difference() and ame() (average marginal effects, by stage -- outcome equation, inflation equation, or both -- and with the two-types-of-zeros decomposition) and their plot() methods, the diord()/piord()/qiord()/riord() distribution family of the response, vuong(), lr_test(), inflation_test() (Vuong and a parametric-bootstrap likelihood-ratio test), split_test(), classification() (confusion table, Brier and ranked probability scores), compare_models(), simulate() (DHARMa-ready), ranef(), analytic / robust / cluster-robust / bootstrap standard errors, and broom / modelsummary / texreg integration.

The likelihood, its analytic gradient, and the bivariate-normal probabilities are implemented in C++.

Installation

# development version
# install.packages("remotes")
remotes::install_github("bagozzib/iop", build_vignettes = TRUE)

Documentation

• vignette("iop") -- getting started: the four estimators on the bundled political-violence (zero-inflated), trade-agreement (top-inflated), and simulated survey (middle-inflated) data.
• vignette("quantities") -- predicted probabilities with standard errors, first differences, average marginal effects, plots, regression tables, and simulated-residual diagnostics.
• vignette("panels") -- random intercepts, fixed effects and the jackknife, the Mundlak device, cluster-robust standard errors, and the package's Monte Carlo on short panels.
• vignette("model") -- the likelihood and sign conventions, identification, multi-start estimation and boundary cases, which test for which comparison, and validation.
• ?iop for the estimators' full argument list; ?"iop-package" for an index.

Quick start

library(iop)

## political violence (none / repression / civil war): zero-inflated ordered probit
data(bp)
m_op   <- oprobit(violence ~ loggdppc + parliament + disaster + major_oil + major_primary, data = bp)
m_ziop <- iop(violence ~ loggdppc + parliament + disaster + major_oil + major_primary |
                loggdppc + parliament + disaster + major_oil + major_primary,
              data = bp, inflate = "bottom")
summary(m_ziop)
inflation_test(m_ziop)                     # Vuong (raw / AIC / BIC) vs the plain ordered probit
compare_models(op = m_op, ziop = m_ziop)

## regime-specific quantities
head(predict(m_ziop, type = "inflated"))   # P(structurally peaceful)
predict(m_ziop, newdata = data.frame(loggdppc = 8, parliament = 1, disaster = 0,
                                     major_oil = 0, major_primary = 0), se.fit = TRUE)
first_difference(m_ziop, "loggdppc", from = 7, to = 9)           # total effect
first_difference(m_ziop, "loggdppc", from = 7, to = 9, stage = "inflation")
plot(ame(m_ziop, vars = c("loggdppc", "disaster")))              # average marginal effects

## correlated errors (ZiOPC)
m_ziopc <- iop(violence ~ loggdppc + parliament + disaster + major_oil + major_primary |
                 loggdppc + parliament + disaster + major_oil + major_primary,
               data = bp, inflate = "bottom", correlated = TRUE)
lr_test(m_ziop, m_ziopc)                   # rho = 0?

## top-inflated: escape-flexibility provisions in trade agreements
data(pta)
m_tiop <- iop(flexibility ~ depth * democracy + gdp + gdppc + trade + gattwto + members +
                democratization | gdp + gdppc + democracy + democratization,
              data = pta, inflate = "top")

## partial proportional odds for the plain ordered logit
m_pl <- ologit(flexibility ~ depth + democracy + gdp, data = pta, parallel = "auto")
m_pl$autofit$relaxed

## panels: cluster-robust SEs, random intercepts, the Mundlak device
m_cl <- iop(violence ~ loggdppc + parliament + disaster | loggdppc + parliament + disaster,
            data = bp, inflate = "bottom", cluster = "country")
m_re <- oprobit(violence ~ loggdppc + parliament + disaster, data = bp, re = "country")
head(ranef(m_re))
md <- mundlak(violence ~ loggdppc + disaster, data = bp, unit = "country")
oprobit(md$formula, md$data, cluster = "country")

## tables
texreg::screenreg(list(m_op, m_ziop))
modelsummary::modelsummary(list(OP = m_op, ZiOP = m_ziop))

Model

Outcome stage (ordered probit/logit with cutpoints tau): P(y <= j | s = 1) = F(tau_j - x'beta). Inflation stage: P(s = 1) = F(z'gamma) is the probability of the ordered (non-inflated) regime, so positive inflation coefficients raise the probability of the ordered regime; units in the inflated regime (s = 0) are observed in the inflated category k with certainty. Hence

P(y = j) = F(z'gamma) * pi_j(x) + 1{j = k} * [1 - F(z'gamma)].

With correlated = TRUE (probit only) the two latent errors are bivariate normal with correlation rho, and the regime-1 cell probabilities are rectangle probabilities of the bivariate normal. See vignette("model").

Validation

The package's tests check the bivariate-normal routine against pbivnorm and mvtnorm, every analytic gradient against numerical differentiation, oprobit()/ologit() against MASS::polr, ordinal::clm, and VGAM::vglm, ordinal::clmm for the random-intercept models, parameter recovery on simulated inflated data, and -- on the bundled bp data -- Table 1 of Bagozzi, Hill, Moore and Mukherjee (2015): OP, ZiOP, ZiOPC, and ZiOPC2 log-likelihoods and the error correlations. Beyond the test suite, iop() reproduces Table 1 and the marginal-effects figures of Bagozzi and Mukherjee (2012), the TiOP/TiOPC log-likelihoods of Bagozzi, Joo and Mukherjee (2024), and the Python idcempy implementation where it converges; the original likelihood code of each paper evaluated at iop's estimates returns iop's log-likelihood to machine precision. See data-raw/oracles/README.md in the source repository, and inst/mc/ for the Monte Carlo behind the panel recommendations.

Citation

citation("iop") gives the package citation and the articles that introduced the models.

References

Harris, M.N. and Zhao, X. (2007). A zero-inflated ordered probit model, with an application to modelling tobacco consumption. Journal of Econometrics, 141, 1073-1099.

Bagozzi, B.E. and Mukherjee, B. (2012). A mixture model for middle category inflation in ordered survey responses. Political Analysis, 20, 369-386.

Bagozzi, B.E., Hill, D.W., Moore, W.H. and Mukherjee, B. (2015). Modeling two types of peace: The zero-inflated ordered probit (ZiOP) model in conflict research. Journal of Conflict Resolution, 59, 728-752.

Bagozzi, B.E., Joo, M.M. and Mukherjee, B. (2024). Top-category inflation in ordered international relations outcomes. Foreign Policy Analysis, 20, orae006.