Mixed logit and preference heterogeneity • choicer

The mixed (random-coefficients) logit lets tastes vary across decision units. Instead of a single coefficient on each random attribute, choicer estimates a distribution of coefficients. Substitution then reflects both observed covariates and the estimated distribution of tastes. Estimation is by simulated maximum likelihood using Halton draws, with the likelihood, gradient and Hessian evaluated in parallel C++. By default run_mxlogit() estimates a cross-sectional simulated likelihood: each choice situation’s probability is integrated over the taste distribution separately. When the same decision maker is observed in several choice situations, person_col switches to the panel likelihood instead (Revelt and Train 1998): one taste draw is held fixed across that person’s situations and their conditional probabilities are integrated jointly, so repeated choices discipline the taste distribution directly — see Panel data: one taste draw per decision maker below.

A useful way to see the mechanism is to write the mixed-logit probability as a logit kernel averaged over the taste distribution, $P_{ij} = \int L_{ij}(\beta)\, f(\beta)\, d\beta$ , where $L_{ij}(\beta) = \exp(X_{ij}\beta) / \sum_k \exp(X_{ik}\beta)$ . Conditional on a draw $\beta$ , the kernel is an ordinary logit. The empirical content of the mixed logit is in the averaging: people with different $\beta$ ’s place different values on the same attributes, so demand leaving an alternative need not go to the same destinations for all consumers. For counterfactual work, the distinction between observed heterogeneity and the estimated mixing distribution is central: flexibility that is not disciplined by the available variation is supplied by the maintained distribution $f(\beta)$ .

library(choicer)
set_num_threads(2)

Simulate correlated random coefficients

simulate_mxl_data() draws choices in which the coefficients on w1 and w2 are themselves random and correlated across the population.

sim <- simulate_mxl_data(N = 2000, J = 4, seed = 1)
sim
#> <choicer_sim: mxl>
#>   settings:
#>     N = 2000
#>     T = 1
#>     J = 4
#>     K_x = 2
#>     K_w = 2
#>     outside_option = TRUE
#>     vary_choice_set = TRUE
#>   rows in $data: 7982
#>   true_params: beta, delta, Sigma, L_params, mu, rc_dist, rc_correlation, gamma_i

Fit

A robust recipe for mixed logit: for a large model, warm-start from a plain MNL (as inst/simulations/mxl_simulation.R does); let the optimizer work in coordinates scaled by the curvature at the start values (the default, scale_vars = "bhhh"); and use enough Halton draws. Here we estimate a full (correlated) covariance of the random coefficients.

fit <- run_mxlogit(
  data            = sim$data,
  id_col          = "id",
  alt_col         = "alt",
  choice_col      = "choice",
  covariate_cols  = c("x1", "x2"),  # fixed coefficients
  random_var_cols = c("w1", "w2"),  # random coefficients
  rc_correlation  = TRUE,           # estimate their full covariance
  S               = 100L,           # Halton draws per choice situation (no person_col here)
  draws           = "generate",     # generate draws on the fly (low memory)
  seed            = 7L,
  scale_vars      = "bhhh",         # the optimizer's coordinates (the default)
  se_method       = "bhhh"
)
#> Optimization run time 0h:0m:0.26s
summary(fit)
#> Mixed Logit (MXL) model
#> 
#> Parameter    Estimate  Std.Error  z-value  Pr(>|z|)  
#> x1           0.776408   0.070138  11.0697  0.00e+00  ***
#> x2          -0.675643   0.068318  -9.8897  0.00e+00  ***
#> Sigma_11     0.905183   0.451785   2.0036  4.51e-02  *
#> Sigma_21     0.894700   0.310679   2.8798  3.98e-03  **
#> Sigma_22     1.884033   0.575835   3.2718  1.07e-03  **
#> ASC_1        0.536159   0.078972   6.7893  1.13e-11  ***
#> ASC_2       -0.477738   0.105075  -4.5466  5.45e-06  ***
#> ASC_3        0.591089   0.077545   7.6225  2.49e-14  ***
#> ASC_4       -0.513880   0.105365  -4.8771  1.08e-06  ***
#> ---
#> Signif. codes:  '***' 0.001 '**' 0.01 '*' 0.05
#> 
#> Random coefficient covariance (Sigma):
#>        w1     w2
#> w1 0.9052 0.8947
#> w2 0.8947 1.8840
#> 
#> Std. Errors: BHHH (OPG) 
#> Log-likelihood: -2435.84 
#> AIC: 4889.67  | BIC: 4940.08 
#> McFadden R2: 0.106 (adj: 0.102) | Hit rate: 0.452 
#> N: 2000  | Parameters: 9 
#> Optimization time: 0.26 s
#> Convergence: 1 ( NLOPT_SUCCESS: Generic success return value. )

This vignette uses se_method = "bhhh" because the outer-product calculation is fast and keeps package-build time short. For final empirical work, compare it with the default analytical-Hessian standard errors; when the data come from a choice-based or otherwise weighted sample, use se_method = "sandwich" so the reported covariance is the robust WESML sandwich rather than the inverse weighted Hessian. If WESML observations are also dependent within people or markets, use the weighted cluster sandwich (se_method = "cluster", with cluster_col) and have enough independent clusters for its asymptotics.

Tip. For real applications increase the number of draws (S) until your estimates are stable, and keep the default scale_vars = "bhhh". Stability in S is a statistical requirement, not just numerical hygiene: simulated ML with a fixed number of draws is biased — the log of an unbiased probability simulator is not unbiased — and the classical asymptotics require S to grow with the sample (see the math companion for the conditions). If the solver struggles, pass an explicit theta_init: place the MNL slope and ASC estimates in the corresponding MXL parameter blocks and initialize the additional mean and Cholesky coordinates deliberately. Bounds can keep Cholesky diagonals away from numerically pathological regions. See inst/simulations/mxl_simulation.R for a fully hardened example.

Parameter recovery

recovery_table(fit, sim$true_params)
#> <choicer_recovery> model=choicer_mxl level=0.95
#>    parameter  group    true estimate     se    bias rel_bias_pct z_vs_true
#>       <char> <char>   <num>    <num>  <num>   <num>        <num>     <num>
#> 1:        x1   beta  0.8000   0.7764 0.0701 -0.0236       -2.949   -0.3364
#> 2:        x2   beta -0.6000  -0.6756 0.0683 -0.0756       12.607   -1.1072
#> 3:      L_11  sigma  0.0000  -0.0498 0.2496 -0.0498           NA   -0.1996
#> 4:      L_21  sigma  0.5000   0.9404 0.3057  0.4404       88.078    1.4406
#> 5:      L_22  sigma  0.1116  -0.0002 0.3566 -0.1117     -100.136   -0.3133
#> 6:     ASC_1    asc  0.5000   0.5362 0.0790  0.0362        7.232    0.4579
#> 7:     ASC_2    asc -0.5000  -0.4777 0.1051  0.0223       -4.452    0.2119
#> 8:     ASC_3    asc  0.5000   0.5911 0.0775  0.0911       18.218    1.1747
#> 9:     ASC_4    asc -0.5000  -0.5139 0.1054 -0.0139        2.776   -0.1317
#>    lower_ci upper_ci covers
#>       <num>    <num> <lgcl>
#> 1:   0.6389   0.9139   TRUE
#> 2:  -0.8095  -0.5417   TRUE
#> 3:  -0.5389   0.4393   TRUE
#> 4:   0.3412   1.5396   TRUE
#> 5:  -0.6991   0.6988   TRUE
#> 6:   0.3814   0.6909   TRUE
#> 7:  -0.6837  -0.2718   TRUE
#> 8:   0.4391   0.7431   TRUE
#> 9:  -0.7204  -0.3074   TRUE

The beta rows are the fixed coefficients, the sigma rows describe the covariance of the random coefficients (its Cholesky elements), and the asc rows are the alternative-specific constants.

Substitution through taste heterogeneity

With random coefficients, diversion is mediated by the distribution of tastes. If the estimated mixing distribution captures economically meaningful heterogeneity, people who value one alternative tend to value nearby substitutes on that latent margin. Diversion can therefore depend on which attribute is changing, so diversion_ratios() takes a wrt_var:

elasticities(fit, elast_var = "x2")
#>   0         1         2         3         4
#> 0 0 -0.015656 -0.009793 -0.025914 -0.009414
#> 1 0 -0.023929 -0.005506 -0.017254 -0.005718
#> 2 0 -0.010723 -0.016649 -0.017014 -0.007051
#> 3 0 -0.008498 -0.006243 -0.006971 -0.006645
#> 4 0 -0.011956 -0.006111 -0.021578 -0.020221
diversion_ratios(fit, wrt_var = "x2")
#>        0      1      2      3      4
#> 0 0.0000 0.3435 0.3183 0.3533 0.3177
#> 1 0.3112 0.0000 0.2616 0.3195 0.2629
#> 2 0.1789 0.1623 0.0000 0.1658 0.1434
#> 3 0.3342 0.3336 0.2790 0.0000 0.2760
#> 4 0.1757 0.1605 0.1412 0.1614 0.0000

# For a random-coefficient attribute the perturbation coordinate matters.
elasticities(fit, elast_var = "w2", is_random_coef = TRUE)
#>   0        1        2        3        4
#> 0 0 -0.01809 -0.02890 -0.01330 -0.02617
#> 1 0  0.05861 -0.02036 -0.02040 -0.01905
#> 2 0 -0.02383  0.13609 -0.02394 -0.02000
#> 3 0 -0.02050 -0.01863  0.04831 -0.01820
#> 4 0 -0.02493 -0.02164 -0.02110  0.13556
diversion_ratios(fit, wrt_var = "w2", is_random_coef = TRUE)
#>         0       1      2        3       4
#> 0  0.0000  0.5513 -6.311  0.53192 -0.4087
#> 1 -2.0961  0.0000  8.744  0.62491  1.8795
#> 2  1.3819  0.5036  0.000 -0.08724 -0.7105
#> 3  1.4015 -0.4331  1.050  0.00000  0.2397
#> 4  0.3127  0.3782 -2.482 -0.06959  0.0000

The rest of the toolkit — predict(), wtp(), consumer_surplus(), blp() — uses the same fitted object as in the getting-started vignette. Prediction and share inversion integrate over the fitted taste distribution; wtp() respects choicer’s supported coefficient parameterizations. In v0.2.0, consumer_surplus() requires a fixed price coefficient because inverse moments of a random price coefficient need not exist.

Identification and tails

The mixed logit is a genuine generalization of the MNL — if tastes are in fact homogeneous, the estimator simply returns a near-zero variance and you are back to a logit. The issue is not that random coefficients are intrinsically fragile. The issue is identification: the additional substitution structure is carried by the mixing distribution $f(\beta)$ , and $f(\beta)$ is often hardest to pin down where it matters most for welfare and diversion — in the tails.

Two consequences are worth keeping in front of you:

The tails drive the economics you report. A lognormal price coefficient puts a slice of the population at near-zero price sensitivity, which can give a willingness-to-pay distribution with no finite mean and explosive welfare numbers. An unbounded normal coefficient implies a fraction of consumers with the wrong sign (who prefer paying more). These artifacts come from the assumed shape of $f$ , not from the data, and the estimator will happily contort a tail to match an aggregate moment.
$f(\beta)$ is hard to identify and estimate in practice. A single cross-section of choices — one decision per person, a fixed menu — carries little information about the spread of tastes. In mixed-logit research, reliable estimation typically uses repeated choices from the same individual (panel data), substantial variation in choice sets or attributes across markets (the BLP setting), or the rich, designed attribute variation of a stated-preference experiment. Without one of these, the random-coefficient variances are weakly identified and the estimates can be fragile. run_mxlogit() exploits the first source when person_col is supplied: one taste draw is then held fixed across a person’s choice situations, so it is their joint probability that gets integrated over $f(\beta)$ (Revelt and Train 1998) — see Panel data: one taste draw per decision maker below. Without person_col, the default cross-sectional likelihood integrates each situation separately and does not use that within-person link; cluster-robust standard errors (cluster_col=, or vcov(fit, type = "cluster"); see Which standard errors, and when) then repair the inference for within-person dependence without changing what the point estimates target. The hierarchical Bayesian logit, run_hmnlogit(), is the Bayesian alternative to the panel likelihood: the same taste-sharing idea, plus partial pooling of alternative-level effects and posterior rather than asymptotic inference.

Panel data: one taste draw per decision maker

Estimand. For decision maker $n$ observed in choice situations $t = 1, \dots, T_n$ , the panel likelihood holds one taste draw $\beta_n \sim f(\beta \mid \theta)$ fixed across all of $n$ ’s situations and integrates their joint probability (Revelt and Train 1998): $L_n(\theta) = \int \prod_{t=1}^{T_n} P_{nt}(j_{nt} \mid \beta)\, f(\beta \mid \theta)\, d\beta.$ The target is the same taste distribution $f(\beta \mid \theta)$ ; what differs is the maintained model of a person’s repeated choices. The default likelihood integrates every situation on its own — the full likelihood when tastes are redrawn each situation or each person is observed once.

Why repeated choices identify $\Sigma$ . A person who keeps picking the same kind of alternative, situation after situation, is revealing that their tastes sit away from the population average and stay there — that persistence is exactly the variation that identifies $\Sigma$ most sharply (the marginals identify it too, but only through the shape of the mixture). The cross-sectional likelihood cannot see it: integrating each situation on its own treats a person’s five consistent choices the same as five different people who each chose once. It still targets the same population mixing distribution — situation by situation, each probability is correctly specified, so it stays consistent (a composite, or pseudo-, likelihood) — but it is statistically inefficient, and the loss falls squarely on $\Sigma$ .

Fitting the panel likelihood. Point person_col at the respondent identifier; id_col still identifies choice situations, and must stay unique across the whole data set (situation ids do not restart within person). S becomes draws per decision maker rather than per choice situation, since one draw block is now shared across all of a person’s situations. A product of $T_n$ per-situation probabilities is more sharply peaked in $\beta$ than a single choice probability and needs more draws to simulate accurately — the recovery script in inst/simulations/mxl_panel_simulation.R uses S = 1000 at T = 10; this vignette keeps S = 100 for speed.

panel_sim <- simulate_mxl_data(N = 300, T = 6, J = 4, seed = 1)

T = 6 gives each of the 300 decision makers six choice situations, sharing one realized taste vector (panel_sim$true_params$gamma_i[, n]) across them.

panel_fit <- run_mxlogit(
  data            = panel_sim$data,
  id_col          = "id",
  alt_col         = "alt",
  choice_col      = "choice",
  covariate_cols  = c("x1", "x2"),
  random_var_cols = c("w1", "w2"),
  rc_correlation  = TRUE,
  person_col      = "pid",          # one taste draw per pid, shared across its T tasks
  S               = 100L,           # Halton draws per decision maker
  draws           = "generate",
  seed            = 7L
)
#> Optimization run time 0h:0m:0.11s
summary(panel_fit)
#> Mixed Logit (MXL) model
#> 
#> Parameter    Estimate  Std.Error  z-value  Pr(>|z|)  
#> x1           0.867807   0.067285  12.8974  0.00e+00  ***
#> x2          -0.723011   0.069041 -10.4721  0.00e+00  ***
#> Sigma_11     1.241486   0.251382   4.9386  7.87e-07  ***
#> Sigma_21     0.532926   0.187874   2.8366  4.56e-03  **
#> Sigma_22     1.693020   0.301606   5.6134  1.98e-08  ***
#> ASC_1        0.537943   0.080112   6.7149  1.88e-11  ***
#> ASC_2       -0.414377   0.094427  -4.3883  1.14e-05  ***
#> ASC_3        0.485624   0.081056   5.9912  2.08e-09  ***
#> ASC_4       -0.590391   0.098788  -5.9763  2.28e-09  ***
#> ---
#> Signif. codes:  '***' 0.001 '**' 0.01 '*' 0.05
#> 
#> Random coefficient covariance (Sigma):
#>        w1     w2
#> w1 1.2415 0.5329
#> w2 0.5329 1.6930
#> 
#> Std. Errors: Analytical Hessian 
#> Log-likelihood: -2126.78 
#> AIC: 4271.57  | BIC: 4321.03 
#> McFadden R2: 0.132 (adj: 0.128) | Hit rate: 0.436 
#> N: 1800  | Respondents: 300  | Parameters: 9 
#> Optimization time: 0.11 s
#> Convergence: 3 ( NLOPT_FTOL_REACHED: Optimization stopped because ftol_rel or ftol_abs (above) was reached. )

recovery_table(panel_fit, panel_sim$true_params)
#> <choicer_recovery> model=choicer_mxl level=0.95
#>    parameter  group    true estimate     se    bias rel_bias_pct z_vs_true
#>       <char> <char>   <num>    <num>  <num>   <num>        <num>     <num>
#> 1:        x1   beta  0.8000   0.8678 0.0673  0.0678        8.476    1.0078
#> 2:        x2   beta -0.6000  -0.7230 0.0690 -0.1230       20.502   -1.7817
#> 3:      L_11  sigma  0.0000   0.1082 0.1012  0.1082           NA    1.0683
#> 4:      L_21  sigma  0.5000   0.4783 0.1600 -0.0217       -4.341   -0.1356
#> 5:      L_22  sigma  0.1116   0.1907 0.1010  0.0791       70.897    0.7829
#> 6:     ASC_1    asc  0.5000   0.5379 0.0801  0.0379        7.589    0.4736
#> 7:     ASC_2    asc -0.5000  -0.4144 0.0944  0.0856      -17.125    0.9068
#> 8:     ASC_3    asc  0.5000   0.4856 0.0811 -0.0144       -2.875   -0.1774
#> 9:     ASC_4    asc -0.5000  -0.5904 0.0988 -0.0904       18.078   -0.9150
#>    lower_ci upper_ci covers
#>       <num>    <num> <lgcl>
#> 1:   0.7359   0.9997   TRUE
#> 2:  -0.8583  -0.5877   TRUE
#> 3:  -0.0903   0.3066   TRUE
#> 4:   0.1646   0.7920   TRUE
#> 5:  -0.0073   0.3887   TRUE
#> 6:   0.3809   0.6950   TRUE
#> 7:  -0.5995  -0.2293   TRUE
#> 8:   0.3268   0.6445   TRUE
#> 9:  -0.7840  -0.3968   TRUE

What changes in inference. With person_col, scores are already per decision maker, so se_method = "sandwich" (equivalently vcov(panel_fit, type = "robust")) is robust to within-person dependence out of the box — there is nothing left within a person to cluster on. cluster_col is for coarser groups that nest decision makers (households, markets); it must be constant within each person. Weights, if supplied, are decision-maker weights, constant within person, since the objective is $\sum_n w_n \log L_n$ ; WESML (choice-based) weights are not supported with person_col, because a choice-based weight varies with the alternative a person chose, which can differ across their situations. nobs() still counts choice situations; a panel fit additionally reports n_persons, printed as “Respondents” above.

As a check on the efficiency argument above: refitting the same data cross-sectionally, clustering on pid for valid inference, repairs the standard errors’ validity but not their size.

cs_fit <- run_mxlogit(
  data            = panel_sim$data,
  id_col          = "id",
  alt_col         = "alt",
  choice_col      = "choice",
  covariate_cols  = c("x1", "x2"),
  random_var_cols = c("w1", "w2"),
  rc_correlation  = TRUE,
  cluster_col     = "pid",          # inference fix only: still one draw per situation
  S               = 100L,
  draws           = "generate",
  seed            = 7L
)
#> Optimization run time 0h:0m:0.19s
se_ratio <- sqrt(diag(vcov(cs_fit, type = "cluster"))) / sqrt(diag(vcov(panel_fit)))
round(se_ratio, 2)
#>    x1    x2  L_11  L_21  L_22 ASC_1 ASC_2 ASC_3 ASC_4 
#>  1.08  1.11  1.88  1.68  1.64  1.03  1.13  1.05  1.18

The ratio sits close to one on the mean-utility parameters (beta, ASCs — about 1.0-1.2 here) but visibly above one on the Cholesky block — the cost of never sharing a taste draw within a person, exactly where a composite likelihood’s efficiency loss concentrates. (cs_fit’s numerator is the cluster-robust sandwich; panel_fit’s denominator is the default analytical-Hessian SE, so the ratio also reflects that choice of estimator, not only the likelihood.) Because it targets a different (composite) likelihood, the cross-sectional fit’s log-likelihood, AIC and BIC are not comparable to the panel fit’s.

Individual-level tastes. conditional_tastes() computes, for each decision maker, the mean and SD of $\beta_n$ conditional on the choices they made (Revelt and Train 2000; Train 2009, ch. 11) — the frequentist analogue of the hierarchical Bayesian beta_i summaries of run_hmnlogit():

ct <- conditional_tastes(panel_fit)
ct
#> Conditional tastes E[beta | observed choices] for 300 decision makers (S = 100 draws each)
#> Unit-level means and SDs in $mean and $sd.
#> Population distribution vs conditional means (revealed + residual should be close to 1):
#>    model_mean model_sd mean_cond sd_cond revealed residual
#> w1          0    1.114     0.064   0.767    0.473    0.497
#> w2          0    1.301     0.143   0.952    0.535    0.438

cs_ct <- conditional_tastes(cs_fit)

Read this carefully: ct$mean[, n] is not decision maker $n$ ’s own $\beta_n$ . It is the average taste among people who would make $n$ ’s choices, so it shrinks toward the population mean — more so for people with fewer situations — and its spread across people understates the population spread in $\Sigma$ . The printed revealed and residual columns decompose population variance by the law of total variance, $\mathrm{Var}(\beta) = \mathrm{Var}(E[\beta \mid y]) + E[\mathrm{Var}(\beta \mid y)]$ ; they should sum to about one. This is a numerical check, not a specification test: the random-coefficient covariance $\Sigma$ is always freely estimated, so this identity is the likelihood’s first-order condition for $\Sigma$ and holds at a converged fit whatever the true taste distribution — a departure points to too few draws S or an unconverged fit. (Here rc_mean = FALSE fixes the random-coefficient means at zero rather than estimating them freely; a mean_cond far from zero below would flag that restriction.) In this example the conditional means correlate 0.67 (w1) and 0.74 (w2) with the realized tastes in panel_sim$true_params$gamma_i. At the true $\theta$ , shrinkage sets this correlation exactly: $\mathrm{cor}(E[\beta \mid y], \beta) = sd(E[\beta \mid y]) / sd(\beta) = \sqrt{\texttt{revealed}}$ , so the same shrinkage that keeps revealed below one keeps the correlation below one too. A single choice reveals little about a person’s tastes, so fitting without person_col leaves revealed much smaller — 0.17 vs. 0.47 for w1, 0.2 vs. 0.54 for w2 — than with repeated choices; conditional means only earn their keep with repeated choices.

Using conditional means as regressors in a second-stage analysis is worse than merely ignoring their sampling and simulation error: the draw weights behind each conditional mean are not exported, so regressing on ct$mean directly measures tastes with error and biases the second-stage coefficient toward zero (attenuation), on top of understating its standard error. Prefer entering the observables in utility instead — interactions with the random-coefficient columns, or rc_mean shifts — or, if a second-stage regression is unavoidable, propagate the uncertainty with a decision-maker-level bootstrap rather than treating conditional means as if they were data.

A research checklist

For a reportable run_mxlogit() application:

Fit a transparent MNL benchmark and use its estimates to inform explicit theta_init values. The default start moves with the covariates’ units (a normal random coefficient’s standard deviation starts at one over its covariate’s typical step: 1 for a dummy, the standard deviation for a dense covariate), so the optimizer’s path barely depends on them (exactly so, in exact arithmetic, under scale_vars = "sd"); a shifted log-normal coefficient’s random part, though, has median one in its covariate’s units: put such a covariate in units where the coefficient is of order one.
Run several manually supplied starting values and retain their objective, convergence code and substantive outputs. choicer does not yet automate multistart estimation.
Increase S until coefficients, covariance elements, standard errors, WTP, elasticities and diversion ratios are stable—not merely until the optimizer reports convergence.
Compare choicer’s supported normal and shifted-lognormal restrictions where economically meaningful, and show what their tails imply for signs and WTP.
Compare the analytical-Hessian and BHHH results; use robust or WESML inference for weighting/misspecification and cluster inference for repeated sampling units as appropriate.
On panel data, decide what you can defend about tastes within a person before choosing the standard error: fit person_col when tastes can be defended as constant within a decision maker’s situations, and keep the clustered cross-sectional fit when they cannot — it relies only on each situation’s marginal and pays in precision, mostly on $\Sigma$ . Add run_hmnlogit() when partially-pooled alternative effects or posterior uncertainty are wanted too.
Compare substitution and welfare with the MNL or NL benchmark and state which data variation, rather than the mixing distribution alone, earns the change.

Bounded/censored mixing distributions, WTP-space estimation, and latent-class logit can be useful alternative modeling strategies, but they are not implemented in choicer v0.2.0. They require another implementation today and are roadmap candidates here; they should not be presented as options to run_mxlogit(). The broader tradeoff is laid out in Choosing among choice models.

References

McFadden, D. and Train, K. (2000). Mixed MNL models for discrete response. Journal of Applied Econometrics, 15(5), 447-470.

Revelt, D. and Train, K. (1998). Mixed logit with repeated choices: households’ choices of appliance efficiency level. Review of Economics and Statistics, 80(4), 647-657.

Revelt, D. and Train, K. (2000). Customer-specific taste parameters and mixed logit: households’ choice of electricity supplier. Working Paper E00-274, Department of Economics, University of California, Berkeley.

Train, K. E. (2009). Discrete Choice Methods with Simulation (2nd ed.). Cambridge University Press.

Train, K. and Weeks, M. (2005). Discrete choice models in preference space and willingness-to-pay space. In R. Scarpa and A. Alberini (Eds.), Applications of Simulation Methods in Environmental and Resource Economics. Springer.