IDCeMPy Package¶
Description¶
IDCeMPy is a Python package that provides functions to fit and assess the performance of the following distinct sets of “inflated” discrete choice models.
Fit the Zero-Inflated Ordered Probit (ZIOP) model without and with correlated errors (ZIOPC model) to evaluate zero-inflated ordered choice outcomes that results from a dual data generating process (d.g.p.).
Fit the Middle-Inflated Ordered Probit (MIOP) model without and with correlated errors (MIOPC) to account for the inflated middle-category in ordered choice measures that relates to a dual d.g.p.
Fit Generalized Inflated Multinomial Logit (GIMNL) models that account for the preponderant and heterogeneous share of observations in the baseline or any lower category in unordered polytomous choice outcomes.
Compute AIC and Log-likelihood statistics and the Vuong Test statistic to assess the performance of each inflated discrete choice model in the package.
IDCeMPy uses Newton numerical optimization methods to estimate the models listed above via Maximum Likelihood Estimation (MLE).
When should you use IDCeMPy?¶
An excessive (“inflated”) share of observations—stemming from two distinct d.g.p’s—fall into a single choice category in many ordered and unordered polytomous outcome variables. Standard Ordered Probit and Multinomial Logit models cannot account for such category inflation which leads to biased inferences. Examples include,
The inflated zero-category of “no smoking” in ordered measures of self-reported smoking behavior is generated from nonsmokers who never smoke cigarettes and those who smoked previously but temporarily stopped smoking because of high cigarette prices.
The inflated “indifference” middle-category in ordered measures of immigration attitudes includes respondents truly indifferent to immigration and those that choose indifference for social desirability reasons.
The inflated baseline or other lower outcome categories of unordered polytomous outcome measures of vote choice include nonvoters who temporarily abstain from voting and routine nonvoters who always abstain.
IDCeMPy includes the ZIOP(C) models for evaluating zero-inflated ordered choice outcomes that results from a dual d.g.p, the MIOP(C) models that address inflated middle-category ordered outcome measures arising from distinct d.g.p’s, and GIMNL models that account for inflated baseline or other categories for unordered polytomous outcomes.
Each inflated discrete choice model in this package addresses category inflation in one’s discrete outcome—unordered or unordered polytomous—of interest by jointly estimating a binary split-stage equation and an ordered or multinomial discrete choice outcome equation.
Installation¶
The package can be installed in two different ways.
From PyPi:
# Import the package
pip install idcempy
From its GitHub Repository:
# Import the package
git clone https://github.com/hknd23/idcempy.git
cd idcempy
python setup.py install
Examples¶
The Standard Ordered Probit (OP) model¶
The package also includes a function that estimates a standard Ordered Probit (OP) model. The OP model does not account for the “zero inflation”, so it does not have a split-probit stage.
We first import the required libraries, set up the package and import the dataset:
# Import the necessary libraries and package
import numpy as np
import pandas as pd
import urllib
from idcempy import zmiopc
# Import the "Youth Tobacco Consumption" dataset.
url='https://github.com/hknd23/zmiopc/blob/main/data/tobacco_cons.csv'
# Read the dataset
data=pd.read_csv(url)
# Define a list of variable names X,Y:
# X = Column names of covariates (from `Data.Frame`) in the OP equation
# Y = Column name of outcome variable (from `Data.Frame`).
X = ['age', 'grade', 'gender_dum']
Y = ['cig_count']
Your data is not ready for estimation.
# Starting parameters for optimization:
pstartop = np.array([.01, .01, .01, .01, .01, .01, .01])
# Model estimation:
op_tob = zmiopc.opmod(pstartop, data, X, Y, method='bfgs', weights=1, offsetx=0)
# See estimates:
print(ziop_tob.coefs)
Results from the model:
The following message will appear when the model has converged:
Warning: Desired error not necessarily achieved due to precision loss.
Current function value: 4411.710049
Iterations: 10
Function evaluations: 976
Gradient evaluations: 121
zmiopc.OpModel stores results from model estimation and other information in its attributes.
The following line of code to see the estimates of coefficients:
print(op_tob.coefs)
Coef SE tscore p 2.5% 97.5%
cut1 1.696175 0.047320 35.844532 0.000000 1.603427 1.788922
cut2 -0.705037 0.031650 -22.276182 0.000000 -0.767071 -0.643004
cut3 -2.304405 0.121410 -18.980329 0.000000 -2.542369 -2.066441
cut4 2.197381 0.235338 9.337141 0.000000 1.736119 2.658643
age -0.070615 0.007581 -9.314701 0.000000 -0.085474 -0.055756
grade 0.233741 0.010336 22.614440 0.000000 0.213483 0.254000
gender_dum 0.020245 0.032263 0.627501 0.530331 -0.042991 0.083482
Log-likelihood, AIC, and Variance-Covariance matrix can be extracted with:
print(op_tob.llik)
print(op_tob.AIC)
print(op_tob.vcov)
The Vuong Test¶
Harris and Zhao (2007) suggest that a variant of the Vuong (1989) Test (with a v statistic) can be used to compare the performance of the ZiOP versus the standard Ordered Probit (OP) model using zmiopc.vuong_opiop().
The formula to estimate the Vuong test is:
where v < -1.96 favors the more general (ZiOP/ZiOPC) model, -1.96 < v < 1.96 lends no support to either model, and v > 1.96 supports the simpler (OP) model.
The OP and ZiOP models must have the same number of observations, and the OP must have the same number of covariates as ZiOP’s OP stage. The statistic reveals that the OP model is preferred over the ZiOP model.
zmiopc.vuong_opiop(op_tob, ziop_tob)
6.624742132792222
The Vuong test can also be implemented to compare the ZiOPC, MiOP and MiOPC models and the OP model.
Generalized Inflated Multinomial Logit (GiMNL) Model¶
The IDCeMPy package also includes a function that estimates General “inflated” Multinomial Logit models (GiMNL). GiMNL models minimize issues present when unordered polytomous outcome variables have an excessive share and heterogeneous pool of observations in the lower category. Failing to account for such inflation could lead to inaccurate inferences.
To estimate the GiMNL model, we first import the library and the dataset introduced above.
from idcempy import gimnl
url= 'https://github.com/hknd23/zmiopc/raw/main/data/replicationdata.dta'
data= pd.read_stata(url)
We the define the list of covariates in the split-stage (z), the second-stage (x) and the outcome variable (y).
# Define the variable list for x, y and z.
# x = Column names of covariates (from `Data.Frame`) in the outcome-stage.
# z = Column names of covariates (from `Data.Frame`) in the split-stage.
# y = Column names of outcome variable (from `Data.Frame`).
x = ['educ', 'party7', 'agegroup2']
z = ['educ', 'agegroup2']
y = ['vote_turn']
Users can employ the argument inflatecat to specify any unordered category as the inflated category (dictated by the distribution) in their unordered-polytomous outcome measure. If a higher category (say 1) is inflated in a 0,1,2 unordered outcome measure.
We first need to specify the order of the outcome variable. Then, you need to define which category is “inflated.”
order = [0, 1, 2]
inflatecat = "baseline"
Further, employing the argument reference, users can select which category of the unordered outcome variable is the baseline (“reference”) category by placing it first. Since the baseline (“0”) category in the Presidential vote choice outcome measure is inflated, the following code fits the BIMNL Model.
gimnl_2004vote = gimnl.gimnlmod(data, x, y, z, order, inflatecat)
The following line of code prints the coefficients of the covariates.
print(gimnl_2004vote.coefs)
Coef SE tscore p 2.5% 97.5%
Inflation: int -4.935 2.777 -1.777 0.076 -10.379 0.508
Inflation: educ 1.886 0.293 6.441 0.000 1.312 2.460
Inflation: agegroup2 1.295 0.768 1.685 0.092 -0.211 2.800
1: int -4.180 1.636 -2.556 0.011 -7.387 -0.974
1: educ 0.334 0.185 1.803 0.071 -0.029 0.697
1: party7 0.454 0.057 7.994 0.000 0.343 0.566
1: agegroup2 0.954 0.248 3.842 0.000 0.467 1.441
2: int 0.900 1.564 0.576 0.565 -2.166 3.966
2: educ 0.157 0.203 0.772 0.440 -0.241 0.554
2: party7 -0.577 0.058 -9.928 0.000 -0.691 -0.463
2: agegroup2 0.916 0.235 3.905 0.000 0.456 1.376
The results from the model are stored in a gimnlModel with the following attributes:
coefs: Model coefficients and standard errors
llik: Log-likelihood
AIC: Akaike information criterion
vcov: Variance-covariance matrix
You can, for example, print the AIC as follows.
print(gimnl_2004vote.AIC)
1656.8324085039708
Using the function gimnl.mnlmod(), users can fit a standard Multinomial Logit Model (MNL) by specifying the list of X, Y, and baseline (using reference).
mnl_2004vote = gimnl.mnlmod(data, x, y, z, order)
print(mnl_2004vote.coefs)
Coef SE tscore p 2.5% 97.5%
1: int -4.914 0.164 -29.980 0.000 -5.235 -4.593
1: educ 0.455 0.043 10.542 0.000 0.371 0.540
1: party7 0.462 0.083 5.571 0.000 0.300 0.625
1: agegroup2 0.951 0.029 32.769 0.000 0.894 1.008
2: int 0.172 0.082 2.092 0.036 0.011 0.334
2: educ 0.282 0.031 9.011 0.000 0.221 0.343
2: party7 -0.567 0.085 -6.641 0.000 -0.734 -0.399
2: agegroup2 0.899 0.138 6.514 0.000 0.629 1.170
Similar to the GiMNL model, the AIC for the MNL model can also be given by:
print(mnl_2004vote.AIC)
1657.192925769978