18a. Ordinal logit model

Bayesian estimation of an ordinal logit model. This is just to illustrate the syntax, as the data are not ordered. But the example assume, for the sake of it, that the alternatives are ordered as 1->2->3

Michel Bierlaire, EPFL Mon Nov 17 2025, 16:38:41

from IPython.core.display_functions import display

import biogeme.biogeme_logging as blog
from biogeme.bayesian_estimation import BayesianResults, get_pandas_estimated_parameters
from biogeme.biogeme import BIOGEME
from biogeme.expressions import Beta, OrderedLogLogit

See the data processing script: Data preparation for Swissmetro.

from swissmetro_data import CHOICE, TRAIN_COST_SCALED, TRAIN_TT_SCALED, database

logger = blog.get_screen_logger(level=blog.INFO)
logger.info('Example b18a_ordinal_logit.py')
Example b18a_ordinal_logit.py

We define a small but positive lower bound

POSITIVE_LOWER_BOUND = 1.0e-5

Parameters to be estimated

b_time = Beta('b_time', 0, None, None, 0)
b_cost = Beta('b_cost', 0, None, None, 0)

Threshold parameters for the ordered logit.

\(\tau_1 \leq 0\).

tau1 = Beta('tau1', -1, None, 0, 0)

\(\delta_2 \geq 0\).

delta2 = Beta('delta2', 2, POSITIVE_LOWER_BOUND, None, 0)

\(\tau_2 = \tau_1 + \delta_2\)

tau2 = tau1 + delta2

Utility.

utility = b_time * TRAIN_TT_SCALED + b_cost * TRAIN_COST_SCALED

Associate each discrete indicator with an interval.

  1. \(-\infty \to \tau_1\),

  2. \(\tau_1 \to \tau_2\),

  3. \(\tau_2 \to +\infty\).

log_probability = OrderedLogLogit(
    eta=utility,
    cutpoints=[tau1, tau2],
    y=CHOICE,
    categories=[1, 2, 3],
    neutral_labels=[],
)

Create the Biogeme object.

the_biogeme = BIOGEME(database, log_probability)
the_biogeme.model_name = 'b18a_ordinal_logit'
Biogeme parameters read from biogeme.toml.

Estimate the parameters.

try:
    results = BayesianResults.from_netcdf(
        filename=f'saved_results/{the_biogeme.model_name}.nc'
    )
except FileNotFoundError:
    results = the_biogeme.bayesian_estimation()
Loaded NetCDF file size: 728.3 MB
load finished in 4212 ms (4.21 s)
print(results.short_summary())
posterior_predictive_loglike finished in 239 ms
expected_log_likelihood finished in 11 ms
best_draw_log_likelihood finished in 11 ms
waic_res finished in 632 ms
waic finished in 632 ms
loo_res finished in 7954 ms (7.95 s)
loo finished in 7955 ms (7.95 s)
Sample size                                              6768
Sampler                                                  NUTS
Number of chains                                         4
Number of draws per chain                                2000
Total number of draws                                    8000
Acceptance rate target                                   0.9
Run time                                                 0:00:35.703617
Posterior predictive log-likelihood (sum of log mean p)  -5789.12
Expected log-likelihood E[log L(Y|θ)]                    -5791.31
Best-draw log-likelihood (posterior upper bound)         -5789.32
WAIC (Widely Applicable Information Criterion)           -5793.51
WAIC Standard Error                                      48.80
Effective number of parameters (p_WAIC)                  4.39
LOO (Leave-One-Out Cross-Validation)                     -5793.52
LOO Standard Error                                       48.80
Effective number of parameters (p_LOO)                   4.40

Get the results in a pandas table

pandas_results = get_pandas_estimated_parameters(
    estimation_results=results,
)
display(pandas_results)
Diagnostics computation took 24.1 seconds (cached).
     Name  Value (mean)  Value (median)  ...     R hat   ESS (bulk)   ESS (tail)
0  b_time     -0.022481       -0.022118  ...  1.000724  3744.636441  4360.599389
1  b_cost      1.264078        1.263626  ...  1.000943  4576.045578  4726.793230
2    tau1     -1.030518       -1.031250  ...  1.001045  3877.677075  4491.405430
3  delta2      3.193749        3.193473  ...  1.001145  4844.804684  4321.703036

[4 rows x 12 columns]

Total running time of the script: (0 minutes 37.233 seconds)

Gallery generated by Sphinx-Gallery