[Stata] Comparing Coefficients Across Nested Logit Models: khb
A common step in applied work is to fit a model with the exposure of interest, add a set of variables that might account for its association with the outcome, and compare the two coefficients. In linear regression the change in the coefficient answers that question directly. In logit and probit models it does not, because the coefficients are scaled by the residual variance of the latent outcome, and that variance changes whenever new variables enter the model. Part of the drop in a logit coefficient after adding controls can come from rescaling rather than from the controls, as Karlson, Holm, and Breen showed.
The KHB method puts the two models on the same scale. The reduced model includes the residuals of the added variables, after regressing them on the exposure, in place of the variables themselves, so it has the same fit and the same scale as the full model while leaving the exposure’s coefficient free of the added variables. The Stata command khb by Ulrich Kohler, Kristian Karlson, and Anders Holm implements it, along with a breakdown of how much each added variable contributes.
The example asks how much of the gap between adults with limited English proficiency and other adults in having no usual source of care is accounted for by health insurance and family income.
Installation
ssc install khbThe current version on SSC is 2.16 (August 2024), and its ado-file requires Stata 17 or newer.
Data
The examples use the 2024 California Health Interview Survey (CHIS) adult public use file, which has 24,810 adult respondents. The CHIS public use files are free and available through a brief application on the UCLA Center for Health Policy Research site. Every post in this series builds its variables with the same short do-file, 00_chis_prep.do:
use "chis_adult_2024.dta", clear
* outcomes
gen k6 = distress if distress >= 0 // K6 score, past month (0-24)
gen nousual = usual == 2 // no usual source of care
gen delay = ah22 == 1 // delayed or did not get medical care, past 12 months
gen erv = er == 1 // any ER visit, past 12 months
gen uninsured = ins == 2
gen needhelp = af81 == 1 if af81 > 0 // needed help for emotional/mental or substance problem
gen fairpoor = ab1 >= 4 // fair or poor self-rated health
gen spd = dstrs12 == 1 // likely serious psychological distress, past year
* groups
gen lep = spk_eng == 3 // speaks English not well or not at all
gen foreign = ah33new == 2
gen female = srsex == 2
* covariates
gen age = srage_p1
recode ahedc_p1 (1 2 = 1 "Less than high school") (3 = 2 "High school") ///
(4 5 6 = 3 "Some college") (7 8 9 = 4 "Bachelor's or more"), gen(educ)
gen pov = povll
label define povllf 4 "300%+ FPL", modify
label values pov povllf
recode ombsrr_p1 (1 = 1 "Latino") (2 = 2 "White") (3 = 3 "Black") (5 = 4 "Asian") ///
(4 6 = 5 "Other or multiracial"), gen(race)
* CHIS design: final weight plus 80 jackknife replicate weights
svyset [pw=rakedw0], jkrweight(rakedw1-rakedw80, multiplier(1)) vce(jackknife) mseLimited English proficiency (lep) means the respondent reported speaking English not well or not at all; 1,520 respondents (6.1 percent) meet that definition. Commands that accept the svy prefix use the jackknife replicate weights; commands that do not accept it use the final weight rakedw0 as a probability weight, which gives the right point estimates and robust standard errors that ignore the replicate design.
The naive comparison
First, the usual approach: a reduced model with the exposure and two background variables, and a full model that adds insurance and income.
tab pov, gen(pv)
logit nousual lep age female [pw=rakedw0]
estimates store reduced
logit nousual lep age female uninsured pv2 pv3 pv4 [pw=rakedw0]
estimates store full
estimates table reduced full, keep(lep) b(%9.3f) se--------------------------------------
Variable | reduced full
-------------+------------------------
lep | 0.921 0.435
| 0.108 0.123
--------------------------------------
Legend: b/seThe coefficient for limited English proficiency falls from 0.921 to 0.435, a drop of 53 percent. Read naively, insurance and income account for about half of the gap. That figure mixes two things: the part of the association that runs through insurance and income, and the change of scale between the two models.
KHB decomposition
In khb, the outcome and the exposure come before the ||, the variables whose contribution you want go after it, and concomitant() holds the controls that appear in both models.
khb logit nousual lep || uninsured pv2 pv3 pv4 [pw=rakedw0], ///
concomitant(age female) summary disentangleDecomposition using the KHB-Method
Model-Type: logit Number of obs = 24810
Variables of Interest: lep Pseudo R2 = 0.12
Z-variable(s): uninsured pv2 pv3 pv4
Concomitant: age female
------------------------------------------------------------------------------
| Robust
nousual | Coefficient std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
lep |
Reduced | .8905438 .1160857 7.67 0.000 .6630199 1.118068
Full | .4351783 .1231229 3.53 0.000 .1938618 .6764949
Diff | .4553654 .0439641 10.36 0.000 .3691973 .5415335
------------------------------------------------------------------------------
Summary of confounding
Variable | Conf_ratio Conf_Pct Resc_Fact
-------------+-------------------------------------
lep | 2.0463881 51.13 .96703317
---------------------------------------------------
Components of Difference
Z-Variable | Coef Std_Err P_Diff P_Reduced
-------------+---------------------------------------------
lep |
uninsured | .1884706 .0281099 41.39 21.16
pv2 | -.0416943 .0170353 -9.16 -4.68
pv3 | -.0078804 .0065694 -1.73 -0.88
pv4 | .3164696 .0430024 69.50 35.54
-----------------------------------------------------------The first table has three rows. Reduced is the exposure coefficient from the reduced model on the full model’s scale (0.891, slightly smaller than the 0.921 from the naive reduced model). Full is the coefficient from the full model (0.435). Diff is their difference (0.455), with a standard error, so you can test whether insurance and income account for a statistically significant share of the association (z = 10.36).
The summary table gives the same result as ratios. Conf_ratio is Reduced divided by Full (2.05), Conf_Pct is the share of the reduced coefficient accounted for by the added variables (51.1 percent, compared with 53 percent from the naive comparison), and Resc_Fact is the rescaling factor. A value of 0.967 means the change of scale was small here, so the naive comparison was off by only a couple of percentage points. With other outcomes and stronger controls the factor can be much further from 1.
The disentangle option splits the difference by variable. P_Diff is each variable’s share of the difference: being uninsured accounts for 41.4 percent, and the income indicators together account for the rest (69.5 − 9.2 − 1.7 = 58.6 percent). The share for the top income category is large because adults with limited English proficiency are much less likely to have family incomes at 300 percent of the poverty line or higher (21 percent, weighted, compared with 62 percent of other adults), and that category has the lowest probability of having no usual source of care. The negative shares for the two middle categories mean those categories push slightly the other way. P_Reduced expresses the same contributions as a share of the reduced coefficient.
Average partial effects
Logit coefficients are hard to explain to most readers. The ape option reports the decomposition as average partial effects on the probability scale.
khb logit nousual lep || uninsured pv2 pv3 pv4 [pw=rakedw0], ///
concomitant(age female) ape summaryDecomposition using the APE-Method
------------------------------------------------------------------------------
| Robust
nousual | Coefficient std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
lep |
Reduced | .102212 .0136914 7.47 0.000 .0753774 .1290467
Full | .0499475 .0143145 3.49 0.000 .0218917 .0780034
Diff | .0522645 . . . . .
------------------------------------------------------------------------------
Note: Standard errors of difference not known for APE methodAdjusting for age and gender, adults with limited English proficiency are 10.2 percentage points more likely to have no usual source of care. After accounting for insurance and income the gap is 5.0 percentage points, so about 5.2 points of the gap go with differences in insurance and income. The APE version does not give a standard error for the difference, so report the significance test from the coefficient version.
Notes
khb accepts probability weights but not the svy prefix, so the standard errors above are robust standard errors that do not use the CHIS replicate weights. The decomposition describes how much of an association is statistically accounted for by other variables. With cross-sectional survey data it does not establish that insurance and income are causal mediators, and the write-up should say “accounted for” rather than “explained by” or “mediated by” unless the design supports it.
References
Karlson, K. B., Holm, A., & Breen, R. (2012). Comparing regression coefficients between same-sample nested models using logit and probit: A new method. Sociological Methodology, 42(1), 286–313. https://doi.org/10.1177/0081175012444861
Kohler, U., Karlson, K. B., & Holm, A. (2011). Comparing coefficients of nested nonlinear probability models. The Stata Journal, 11(3), 420–438. https://doi.org/10.1177/1536867X1101100306
