.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/swissmetro/plot_b13_panel_simul.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_swissmetro_plot_b13_panel_simul.py: 13. Simulation of panel model ============================= Calculates each contribution to the log likelihood function using simulation. We also calculate the individual parameters. Michel Bierlaire, EPFL Sat Jun 21 2025, 17:06:31 .. GENERATED FROM PYTHON SOURCE LINES 13-28 .. code-block:: Python import sys from IPython.core.display_functions import display import biogeme.biogeme_logging as blog from biogeme.biogeme import BIOGEME from biogeme.expressions import Beta, Draws, MonteCarlo, PanelLikelihoodTrajectory, log from biogeme.jax_calculator.single_formula import ( calculate_single_formula_from_expression, ) from biogeme.models import logit from biogeme.results_processing import EstimationResults from biogeme.second_derivatives import SecondDerivativesMode .. GENERATED FROM PYTHON SOURCE LINES 29-30 See the data processing script: :ref:`swissmetro_panel`. .. GENERATED FROM PYTHON SOURCE LINES 30-44 .. code-block:: Python from swissmetro_panel 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, ) .. GENERATED FROM PYTHON SOURCE LINES 45-46 Define the number of draws to be used for Monte-Carlo integration. .. GENERATED FROM PYTHON SOURCE LINES 46-51 .. code-block:: Python NUMBER_OF_DRAWS = 100_000 logger = blog.get_screen_logger(level=blog.INFO) logger.info('Example b13_panel_simul.py') .. rst-class:: sphx-glr-script-out .. code-block:: none Example b13_panel_simul.py .. GENERATED FROM PYTHON SOURCE LINES 52-53 Parameters to be estimated. .. GENERATED FROM PYTHON SOURCE LINES 53-55 .. code-block:: Python b_cost = Beta('b_cost', 0, None, None, 0) .. GENERATED FROM PYTHON SOURCE LINES 56-58 Define a random parameter, normally distributed across individuals, designed to be used for Monte-Carlo simulation. .. GENERATED FROM PYTHON SOURCE LINES 58-60 .. code-block:: Python b_time = Beta('b_time', 0, None, None, 0) .. GENERATED FROM PYTHON SOURCE LINES 61-62 It is advised not to use 0 as starting value for the following parameter. .. GENERATED FROM PYTHON SOURCE LINES 62-65 .. code-block:: Python b_time_s = Beta('b_time_s', 1, None, None, 0) b_time_rnd = b_time + b_time_s * Draws('b_time_rnd', 'NORMAL_ANTI') .. GENERATED FROM PYTHON SOURCE LINES 66-67 We do the same for the constants, to address serial correlation. .. GENERATED FROM PYTHON SOURCE LINES 67-79 .. code-block:: Python asc_car = Beta('asc_car', 0, None, None, 0) asc_car_s = Beta('asc_car_s', 1, None, None, 0) asc_car_rnd = asc_car + asc_car_s * Draws('asc_car_rnd', 'NORMAL_ANTI') asc_train = Beta('asc_train', 0, None, None, 0) asc_train_s = Beta('asc_train_s', 1, None, None, 0) asc_train_rnd = asc_train + asc_train_s * Draws('asc_train_rnd', 'NORMAL_ANTI') asc_sm = Beta('asc_sm', 0, None, None, 0) asc_sm_s = Beta('asc_sm_s', 1, None, None, 0) asc_sm_rnd = asc_sm + asc_sm_s * Draws('asc_sm_rnd', 'NORMAL_ANTI') .. GENERATED FROM PYTHON SOURCE LINES 80-81 Definition of the utility functions. .. GENERATED FROM PYTHON SOURCE LINES 81-85 .. code-block:: Python v_train = asc_train_rnd + b_time_rnd * TRAIN_TT_SCALED + b_cost * TRAIN_COST_SCALED v_swissmetro = asc_sm_rnd + b_time_rnd * SM_TT_SCALED + b_cost * SM_COST_SCALED v_car = asc_car_rnd + b_time_rnd * CAR_TT_SCALED + b_cost * CAR_CO_SCALED .. GENERATED FROM PYTHON SOURCE LINES 86-87 Associate utility functions with the numbering of alternatives. .. GENERATED FROM PYTHON SOURCE LINES 87-89 .. code-block:: Python v = {1: v_train, 2: v_swissmetro, 3: v_car} .. GENERATED FROM PYTHON SOURCE LINES 90-91 Associate the availability conditions with the alternatives. .. GENERATED FROM PYTHON SOURCE LINES 91-93 .. code-block:: Python av = {1: TRAIN_AV_SP, 2: SM_AV, 3: CAR_AV_SP} .. GENERATED FROM PYTHON SOURCE LINES 94-96 Conditional on the random parameters, the likelihood of one observation is given by the logit model (called the kernel). .. GENERATED FROM PYTHON SOURCE LINES 96-98 .. code-block:: Python choice_probability_one_observation = logit(v, av, CHOICE) .. GENERATED FROM PYTHON SOURCE LINES 99-102 Conditional on the random parameters, the likelihood of all observations for one individual (the trajectory) is the product of the likelihood of each observation. .. GENERATED FROM PYTHON SOURCE LINES 102-106 .. code-block:: Python conditional_trajectory_probability = PanelLikelihoodTrajectory( choice_probability_one_observation ) .. GENERATED FROM PYTHON SOURCE LINES 107-108 We integrate over the random parameters using Monte-Carlo .. GENERATED FROM PYTHON SOURCE LINES 108-110 .. code-block:: Python log_probability = log(MonteCarlo(conditional_trajectory_probability)) .. GENERATED FROM PYTHON SOURCE LINES 111-112 We retrieve the parameters estimates. .. GENERATED FROM PYTHON SOURCE LINES 112-121 .. code-block:: Python try: results = EstimationResults.from_yaml_file(filename='saved_results/b12_panel.yaml') except FileNotFoundError: sys.exit( 'Run first the script plot_b12_panel.py ' 'in order to generate the ' 'file b12_panel.yaml and move it to the directory saved_results.' ) .. GENERATED FROM PYTHON SOURCE LINES 122-124 Simulate to recalculate the log likelihood directly from the formula, without the Biogeme object .. GENERATED FROM PYTHON SOURCE LINES 124-134 .. code-block:: Python simulated_loglike = calculate_single_formula_from_expression( expression=log_probability, database=database, number_of_draws=NUMBER_OF_DRAWS, the_betas=results.get_beta_values(), numerically_safe=False, second_derivatives_mode=SecondDerivativesMode.NEVER, use_jit=True, ) .. rst-class:: sphx-glr-script-out .. code-block:: none Flattening database [(6768, 38)]. Database flattened [(752, 362)] .. GENERATED FROM PYTHON SOURCE LINES 135-137 .. code-block:: Python print(f'Simulated log likelihood: {simulated_loglike}') .. rst-class:: sphx-glr-script-out .. code-block:: none Simulated log likelihood: -3577.5783861418536 .. GENERATED FROM PYTHON SOURCE LINES 138-139 We also calculate the individual parameters for the time coefficient. .. GENERATED FROM PYTHON SOURCE LINES 139-147 .. code-block:: Python numerator = MonteCarlo(b_time_rnd * conditional_trajectory_probability) denominator = MonteCarlo(conditional_trajectory_probability) simulate = { 'Numerator': numerator, 'Denominator': denominator, } .. GENERATED FROM PYTHON SOURCE LINES 148-149 Creation of the Biogeme object. .. GENERATED FROM PYTHON SOURCE LINES 149-151 .. code-block:: Python biosim = BIOGEME(database, simulate, number_of_draws=NUMBER_OF_DRAWS, seed=1223) .. rst-class:: sphx-glr-script-out .. code-block:: none Biogeme parameters read from biogeme.toml. .. GENERATED FROM PYTHON SOURCE LINES 152-153 Simulation. .. GENERATED FROM PYTHON SOURCE LINES 153-155 .. code-block:: Python sim = biosim.simulate(results.get_beta_values()) .. GENERATED FROM PYTHON SOURCE LINES 156-158 .. code-block:: Python sim['Individual-level parameters'] = sim['Numerator'] / sim['Denominator'] .. GENERATED FROM PYTHON SOURCE LINES 159-161 .. code-block:: Python print(f'{sim.shape=}') .. rst-class:: sphx-glr-script-out .. code-block:: none sim.shape=(6768, 3) .. GENERATED FROM PYTHON SOURCE LINES 162-163 .. code-block:: Python display(sim) .. rst-class:: sphx-glr-script-out .. code-block:: none Numerator Denominator Individual-level parameters 0 NaN NaN NaN 1 NaN NaN NaN 2 NaN NaN NaN 3 NaN NaN NaN 4 NaN NaN NaN ... ... ... ... 6763 NaN NaN NaN 6764 NaN NaN NaN 6765 NaN NaN NaN 6766 NaN NaN NaN 6767 NaN NaN NaN [6768 rows x 3 columns] .. rst-class:: sphx-glr-timing **Total running time of the script:** (6 minutes 56.200 seconds) .. _sphx_glr_download_auto_examples_swissmetro_plot_b13_panel_simul.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_b13_panel_simul.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_b13_panel_simul.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_b13_panel_simul.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_