knitr::opts_chunk$set(collapse = TRUE, comment = "#>") options(digits = 4)
A demand estimate is only as reportable as its standard error. The default maximum-likelihood variance — the inverse Hessian — is correct when the model is correctly specified and choice situations are independent draws. Applied work routinely violates one or both: the same decision maker contributes several choice situations (panel and survey data), the sample was drawn by outcome (choice-based sampling), or the likelihood is at best a useful approximation. choicer treats the variance estimator as a choice you make — and can revise after estimation, without refitting.
Every MNL, MXL and NL fit exposes four estimators through one interface:
| type | Estimator | When it is the right choice |
|---|---|---|
| "hessian" | Inverse analytical Hessian | Correct specification, independent choice situations. The classical default. |
| "bhhh" | Inverse outer product of scores (BHHH/OPG) | Same asymptotic justification as the Hessian; cheap; a common companion for simulated likelihoods. |
| "robust" | Huber-White sandwich $A^{-1} B A^{-1}$ | Quasi-ML under misspecification; also the valid variance under WESML / choice-based weighting (the meat is the weight-squared score outer product). |
| "cluster" | Sandwich over within-cluster sums of scores | Dependence across choice situations from the same decision maker, market, or sampling unit. |
Here the unit $i$ is the choice situation. If the same person supplies several situations, cluster on the person. No small-sample correction is applied—not a degrees-of-freedom multiplier, leverage correction, or few-cluster reference distribution—so the number of clusters should be reasonably large. With few clusters, the asymptotic cluster sandwich can be seriously anti-conservative; choicer does not currently implement CR2, a wild-cluster bootstrap, or randomization inference. Those designs require an external small-sample procedure or a research design with a defensible higher-level sampling argument.
library(choicer) set_num_threads(2)
The cleanest way to see what clustering buys is to build data with genuine
within-person dependence. simulate_hmnl_data() draws a panel in which each
person has their own taste vector $\beta_i \sim N(b, W)$ and makes several
choices with it — a stylized version of the persistent within-person dependence
in a survey or scanner panel.
sim <- simulate_hmnl_data(N = 500, T = 8, J = 4, seed = 3)
We deliberately fit a plain multinomial logit. Because a homogeneous logit is
misspecified for this random-coefficients DGP, its point estimates target the
pseudo-true homogeneous-logit projection, not generally the population mean
$b$. Its likelihood also treats all 500 * 8 = 4000 choice situations as
independent — which they are not, because each person's taste deviation persists
across their eight tasks. Passing cluster_col
supplies the person labels at fit time and selects cluster-robust standard
errors for summary():
fit <- run_mnlogit( data = sim$data, id_col = "task", alt_col = "alt", choice_col = "choice", covariate_cols = c("x1", "x2"), cluster_col = "pid", include_outside_option = TRUE ) summary(fit)
How much did it matter? vcov(fit, type = ) recomputes any of the four
estimators from the stored per-situation scores:
se_of <- function(V) sqrt(diag(V)) tab <- cbind( hessian = se_of(vcov(fit, type = "hessian")), bhhh = se_of(vcov(fit, type = "bhhh")), robust = se_of(vcov(fit, type = "robust")), cluster = se_of(vcov(fit, type = "cluster")) ) round(tab, 4) round(tab[, "cluster"] / tab[, "hessian"], 2)
The clustered standard errors on the slope coefficients are about 25% wider
than the Hessian ones; those on the alternative-specific constants differ
little. That
pattern is the fingerprint of the data-generating process: the persistent
component is each person's taste deviation, so it is the slope scores that
are correlated within person. Note also that "robust" alone recovers almost
none of this — the Huber-White sandwich guards against misspecification of
each situation's likelihood, not against dependence across situations. If
the same decision maker appears repeatedly, robust-without-clustering is the
wrong comfort.
Nothing above required deciding at fit time. Any fit stored with
keep_data = TRUE (the default) can be re-inferenced later — useful when the
clustering level is itself a robustness question:
fit0 <- run_mnlogit( data = sim$data, id_col = "task", alt_col = "alt", choice_col = "choice", covariate_cols = c("x1", "x2"), include_outside_option = TRUE ) # one cluster label per choice situation, named by choice-situation id task_person <- unique(sim$data[, c("task", "pid")]) cl <- setNames(task_person$pid, task_person$task) se_of(vcov(fit0, type = "cluster", cluster = cl))
The cluster vector has one entry per choice situation — here, one pid per
task — and is named by the choice-situation id. choicer uses those names to
realign labels safely to its prepared order.
A caveat that matters for the mixed logit. run_mxlogit() maximizes a
cross-sectional simulated likelihood: each choice situation's probability is
integrated over the mixing distribution separately, rather than holding each
simulation draw fixed across a person's tasks before integrating their joint
probability. On panel data,
type = "cluster" makes the standard errors robust to within-person
dependence around its pseudo-true target, but the point estimates still target
that cross-sectional model —
clustering does not turn the fit into a panel mixed logit, in which one taste
draw is shared across a person's choices. When the panel structure is the
object of interest — individual-level tastes, their distribution, persistence
— estimate the panel model instead: the hierarchical Bayesian logit
(run_hmnlogit() with person_col), covered in the
hierarchical Bayes vignette.
When the sample was drawn by outcome (interviewing travellers at the terminal
of the mode they chose, oversampling rare hospitals), the weights are not
optional and the variance choice is forced: under WESML weighting the inverse
weighted Hessian and the ordinary BHHH variance are both invalid, and
type = "robust" — whose meat carries the squared weights — is the correct
estimator, conditional on the supplied population shares Q. If weighted
observations are also dependent within sampling units, use the weighted cluster
sandwich instead; if Q is estimated, its uncertainty requires an additional
resampling or sensitivity layer. The choice-based sampling vignette
develops the full workflow, from population shares to weighted fit.
Berndt, E., Hall, B., Hall, R. and Hausman, J. (1974). Estimation and inference in nonlinear structural models. Annals of Economic and Social Measurement, 3(4), 653-665.
Cameron, A. C. and Miller, D. L. (2015). A practitioner's guide to cluster-robust inference. Journal of Human Resources, 50(2), 317-372.
Manski, C. F. and Lerman, S. R. (1977). The estimation of choice probabilities from choice based samples. Econometrica, 45(8), 1977-1988.
White, H. (1982). Maximum likelihood estimation of misspecified models. Econometrica, 50(1), 1-25.
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.