8. Box-Cox transforms

Bayesian estimation of a logit model, with a Box-Cox transform of variables.

Michel Bierlaire, EPFL Mon Nov 03 2025, 13:41:40

from pathlib import Path

from IPython.core.display_functions import display

See the data processing script: Data preparation for Swissmetro.

from swissmetro_data import (
    CAR_AV_SP,
    CAR_CO_SCALED,
    CAR_TT_SCALED,
    CHOICE,
    SM_AV,
    SM_COST_SCALED,
    SM_TT_SCALED,
    TRAIN_AV_SP,
    TRAIN_COST_SCALED,
    TRAIN_TT_SCALED,
    database,
)

import biogeme.biogeme_logging as blog
from biogeme.bayesian_estimation import (
    BayesianResults,
    BayesianResultsSummary,
    get_pandas_estimated_parameters,
)
from biogeme.biogeme import BIOGEME
from biogeme.expressions import Beta
from biogeme.models import boxcox, loglogit

logger = blog.get_screen_logger(level=blog.INFO)

Parameters to be estimated.

asc_car = Beta('asc_car', 0, None, None, 0)
asc_train = Beta('asc_train', 0, None, None, 0)
asc_sm = Beta('asc_sm', 0, None, None, 1)
b_time = Beta('b_time', 0, None, 0, 0)
b_cost = Beta('b_cost', 0, None, 0, 0)
boxcox_parameter = Beta('boxcox_parameter', 1, -2, 2, 0)

Definition of the utility functions.

v_train = (
    asc_train
    + b_time * boxcox(TRAIN_TT_SCALED, boxcox_parameter)
    + b_cost * TRAIN_COST_SCALED
)
v_swissmetro = (
    asc_sm + b_time * boxcox(SM_TT_SCALED, boxcox_parameter) + b_cost * SM_COST_SCALED
)
v_car = (
    asc_car + b_time * boxcox(CAR_TT_SCALED, boxcox_parameter) + b_cost * CAR_CO_SCALED
)

Associate utility functions with the numbering of alternatives.

v = {1: v_train, 2: v_swissmetro, 3: v_car}

Associate the availability conditions with the alternatives.

av = {1: TRAIN_AV_SP, 2: SM_AV, 3: CAR_AV_SP}

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

log_probability = loglogit(v, av, CHOICE)

Create the Biogeme object.

the_biogeme = BIOGEME(database, log_probability, bayesian_draws=10000, warmup=10000)
the_biogeme.model_name = 'b08_boxcox'

Biogeme parameters read from biogeme.toml.

Estimate the posterior distribution of the parameters, or read the results if already available.

yaml_file = Path('saved_results') / f'{the_biogeme.model_name}.yaml'
try:
    summary_results = BayesianResultsSummary.from_yaml_file(filename=yaml_file)
except FileNotFoundError:
    results: BayesianResults = the_biogeme.bayesian_estimation()
    summary_results = results.to_summary()

print(summary_results.short_summary())

Sample size                                              6768
Sampler                                                  NUTS
Number of chains                                         4
Number of draws per chain                                10000
Total number of draws                                    40000
Acceptance rate target                                   0.9
Run time                                                 0:01:28.836393
Posterior predictive log-likelihood (sum of log mean p)  -5838.32
Expected log-likelihood E[log L(Y|θ)]                    -7908.00
Best-draw log-likelihood (posterior upper bound)         -5779.02
LOO (Leave-One-Out Cross-Validation)
LOO Standard Error
Effective number of parameters (p_LOO)

Present the parameter estimates in a pandas table.

pandas_results = get_pandas_estimated_parameters(
    estimation_results=summary_results,
)
display(pandas_results)

               Name  Value (mean)  ...  ESS (bulk)  ESS (tail)
0         asc_train     -0.005509  ...    4.003203    4.003203
1           asc_car      0.074463  ...    4.003203    4.003203
2            b_time     -1.216759  ...    4.003203    4.003203
3  boxcox_parameter     -1.976255  ...    4.003203    4.003203
4            b_cost     -0.763031  ...    4.003203    4.003203

[5 rows x 12 columns]

Report the variables stored in the Bayesian estimation results.

display(summary_results.report_stored_variables())

             group           variable                dims             shape
0    constant_data          CAR_AV_SP               [obs]            [6768]
1    constant_data      CAR_CO_SCALED               [obs]            [6768]
2    constant_data      CAR_TT_SCALED               [obs]            [6768]
3    constant_data             CHOICE               [obs]            [6768]
4    constant_data              SM_AV               [obs]            [6768]
5    constant_data     SM_COST_SCALED               [obs]            [6768]
6    constant_data       SM_TT_SCALED               [obs]            [6768]
7    constant_data        TRAIN_AV_SP               [obs]            [6768]
8    constant_data  TRAIN_COST_SCALED               [obs]            [6768]
9    constant_data    TRAIN_TT_SCALED               [obs]            [6768]
10  log_likelihood            _choice  [chain, draw, obs]  [4, 10000, 6768]
11       posterior            asc_car       [chain, draw]        [4, 10000]
12       posterior          asc_train       [chain, draw]        [4, 10000]
13       posterior             b_cost       [chain, draw]        [4, 10000]
14       posterior             b_time       [chain, draw]        [4, 10000]
15       posterior   boxcox_parameter       [chain, draw]        [4, 10000]
16       posterior           log_like  [chain, draw, obs]  [4, 10000, 6768]
17           prior            asc_car       [chain, draw]        [1, 10000]
18           prior          asc_train       [chain, draw]        [1, 10000]
19           prior             b_cost       [chain, draw]        [1, 10000]
20           prior             b_time       [chain, draw]        [1, 10000]
21           prior   boxcox_parameter       [chain, draw]        [1, 10000]
22           prior           log_like  [chain, draw, obs]  [1, 10000, 6768]
23    sample_stats    acceptance_rate       [chain, draw]        [4, 10000]
24    sample_stats          diverging       [chain, draw]        [4, 10000]
25    sample_stats             energy       [chain, draw]        [4, 10000]
26    sample_stats                 lp       [chain, draw]        [4, 10000]
27    sample_stats            n_steps       [chain, draw]        [4, 10000]
28    sample_stats          step_size       [chain, draw]        [4, 10000]
29    sample_stats         tree_depth       [chain, draw]        [4, 10000]