Ordinal probit model

Example of an ordinal probit 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

author:

Michel Bierlaire, EPFL

date:

Mon Apr 10 12:15:28 2023

import biogeme.biogeme_logging as blog
import biogeme.biogeme as bio
from biogeme.models import ordered_probit
from biogeme.expressions import Beta, log, Elem

See the data processing script: Data preparation for Swissmetro.

from swissmetro_data import (
    database,
    CHOICE,
    TRAIN_TT_SCALED,
    TRAIN_COST_SCALED,
)

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

Example b18ordinal_probit.py

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 probit.

\(\tau_1 \leq 0\).

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

\(\delta_2 \geq 0\).

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

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

tau2 = tau1 + delta2

Utility

U = 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\).

the_proba = ordered_probit(
    continuous_value=U,
    list_of_discrete_values=[1, 2, 3],
    tau_parameter=tau1,
)

Extract from the dict the formula associated with the observed choice.

the_chosen_proba = Elem(the_proba, CHOICE)

Definition of the model. This is the contribution of each observation to the log likelihood function.

logprob = log(the_chosen_proba)

Create the Biogeme object.

the_biogeme = bio.BIOGEME(database, logprob)
the_biogeme.modelName = 'b18ordinal_probit'

File biogeme.toml has been parsed.

Estimate the parameters.

results = the_biogeme.estimate()

*** Initial values of the parameters are obtained from the file __b18ordinal_probit.iter
Cannot read file __b18ordinal_probit.iter. Statement is ignored.
Optimization algorithm: hybrid Newton/BFGS with simple bounds [simple_bounds]
** Optimization: Newton with trust region for simple bounds
Iter.          B_COST          B_TIME            tau1     tau1_diff_2     Function    Relgrad   Radius      Rho
    0            0.56          0.0092            -0.5             1.5      5.9e+03       0.16       10      1.1   ++
    1            0.66           0.015           -0.59             1.8      5.8e+03      0.026    1e+02      1.1   ++
    2            0.69           0.018            -0.6             1.9      5.8e+03    0.00072    1e+03        1   ++
    3            0.69           0.018            -0.6             1.9      5.8e+03      6e-07    1e+03        1   ++
Results saved in file b18ordinal_probit.html
Results saved in file b18ordinal_probit.pickle

print(results.short_summary())

Results for model b18ordinal_probit
Nbr of parameters:              4
Sample size:                    6768
Excluded data:                  3960
Final log likelihood:           -5789.055
Akaike Information Criterion:   11586.11
Bayesian Information Criterion: 11613.39

pandas_results = results.getEstimatedParameters()
pandas_results

	Value	Rob. Std err	Rob. t-test	Rob. p-value
B_COST	0.687183	0.036818	18.664468	0.000000
B_TIME	0.018053	0.023389	0.771826	0.440217
tau1	-0.604796	0.038571	-15.680045	0.000000
tau1_diff_2	1.913927	0.025234	75.848255	0.000000