.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/swissmetro/plot_b06c_unif_mixture_integral.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_b06c_unif_mixture_integral.py: " 6c. Mixture of logit models with uniform distribution and numerical integration =============================================================================== Example of a mixture of logit models, using numerical integration. The mixing distribution is uniform. Michel Bierlaire, EPFL Fri Jun 20 2025, 10:47:24 .. GENERATED FROM PYTHON SOURCE LINES 13-32 .. code-block:: Python from IPython.core.display_functions import display import biogeme.biogeme_logging as blog from biogeme.biogeme import BIOGEME from biogeme.distributions import normalpdf from biogeme.expressions import ( Beta, IntegrateNormal, RandomVariable, exp, log, ) from biogeme.models import logit from biogeme.results_processing import ( EstimationResults, get_pandas_estimated_parameters, ) .. GENERATED FROM PYTHON SOURCE LINES 33-34 See the data processing script: :ref:`swissmetro_data`. .. GENERATED FROM PYTHON SOURCE LINES 34-51 .. code-block:: Python 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, ) logger = blog.get_screen_logger(level=blog.INFO) logger.info('Example b06unif_mixture_integral.py') .. rst-class:: sphx-glr-script-out .. code-block:: none Example b06unif_mixture_integral.py .. GENERATED FROM PYTHON SOURCE LINES 52-53 Parameters to be estimated. .. GENERATED FROM PYTHON SOURCE LINES 53-58 .. code-block:: Python 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_cost = Beta('b_cost', 0, None, None, 0) .. GENERATED FROM PYTHON SOURCE LINES 59-61 Define a random parameter, normally distributed, designed to be used for numerical integration .. GENERATED FROM PYTHON SOURCE LINES 61-65 .. code-block:: Python b_time = Beta('b_time', 0, None, None, 0) b_time_s = Beta('b_time_s', 1, None, None, 0) omega = RandomVariable('omega') .. GENERATED FROM PYTHON SOURCE LINES 66-72 .. |infinity| unicode:: U+221E :trim: As the numerical integration ranges from -|infinity| \ to + |infinity| , we need to perform a change of variable in order to integrate between -1 and 1. .. GENERATED FROM PYTHON SOURCE LINES 72-78 .. code-block:: Python LOWER_BND = -1 UPPER_BND = 1 x = LOWER_BND + (UPPER_BND - LOWER_BND) / (1 + exp(-omega)) dx = (UPPER_BND - LOWER_BND) * exp(-omega) / ((1 + exp(-omega)) ** 2) b_time_rnd = b_time + b_time_s * x .. GENERATED FROM PYTHON SOURCE LINES 79-80 Definition of the utility functions. .. GENERATED FROM PYTHON SOURCE LINES 80-84 .. code-block:: Python v_train = asc_train + b_time_rnd * TRAIN_TT_SCALED + b_cost * TRAIN_COST_SCALED v_swissmetro = asc_sm + b_time_rnd * SM_TT_SCALED + b_cost * SM_COST_SCALED v_car = asc_car + b_time_rnd * CAR_TT_SCALED + b_cost * CAR_CO_SCALED .. GENERATED FROM PYTHON SOURCE LINES 85-86 Associate utility functions with the numbering of alternatives. .. GENERATED FROM PYTHON SOURCE LINES 86-88 .. code-block:: Python v = {1: v_train, 2: v_swissmetro, 3: v_car} .. GENERATED FROM PYTHON SOURCE LINES 89-90 Associate the availability conditions with the alternatives. .. GENERATED FROM PYTHON SOURCE LINES 90-92 .. code-block:: Python av = {1: TRAIN_AV_SP, 2: SM_AV, 3: CAR_AV_SP} .. GENERATED FROM PYTHON SOURCE LINES 93-94 Conditional on omega, we have a logit model (called the kernel). .. GENERATED FROM PYTHON SOURCE LINES 94-96 .. code-block:: Python conditional_probability = logit(v, av, CHOICE) .. GENERATED FROM PYTHON SOURCE LINES 97-98 pdf of the uniform distribution .. GENERATED FROM PYTHON SOURCE LINES 98-100 .. code-block:: Python pdf_uniform = 1 / (UPPER_BND - LOWER_BND) .. GENERATED FROM PYTHON SOURCE LINES 101-103 As the `IntegrateNormal` expression is designed for a normal distribution, we need to divide by the pdf of the normal distribution, and multiply by the pdf of the uniform distribution, after applying the change of variable. .. GENERATED FROM PYTHON SOURCE LINES 103-106 .. code-block:: Python new_integrand = conditional_probability * dx * pdf_uniform / normalpdf(omega) .. GENERATED FROM PYTHON SOURCE LINES 107-109 We integrate over omega using numerical integration. To illustrate the syntax, we specific the number of quadrature points to be used. .. GENERATED FROM PYTHON SOURCE LINES 109-117 .. code-block:: Python log_probability = log( IntegrateNormal( new_integrand, 'omega', number_of_quadrature_points=60, ) ) .. GENERATED FROM PYTHON SOURCE LINES 118-119 Create the Biogeme object. .. GENERATED FROM PYTHON SOURCE LINES 119-122 .. code-block:: Python the_biogeme = BIOGEME(database, log_probability) the_biogeme.model_name = 'b06c_unif_mixture_integral' .. rst-class:: sphx-glr-script-out .. code-block:: none Biogeme parameters read from biogeme.toml. .. GENERATED FROM PYTHON SOURCE LINES 123-124 Estimate the parameters. .. GENERATED FROM PYTHON SOURCE LINES 124-131 .. code-block:: Python try: results = EstimationResults.from_yaml_file( filename=f'saved_results/{the_biogeme.model_name}.yaml' ) except FileNotFoundError: results = the_biogeme.estimate() .. rst-class:: sphx-glr-script-out .. code-block:: none *** Initial values of the parameters are obtained from the file __b06c_unif_mixture_integral.iter Cannot read file __b06c_unif_mixture_integral.iter. Statement is ignored. Starting values for the algorithm: {} As the model is rather complex, we cancel the calculation of second derivatives. If you want to control the parameters, change the algorithm from "automatic" to "simple_bounds" in the TOML file. Optimization algorithm: hybrid Newton/BFGS with simple bounds [simple_bounds] ** Optimization: BFGS with trust region for simple bounds Iter. asc_train b_time b_time_s b_cost asc_car Function Relgrad Radius Rho 0 -1 -1 2 -1 -1 5.6e+03 0.09 1 0.39 + 1 -1 -1 2 -1 -1 5.6e+03 0.09 0.5 -0.11 - 2 -1.5 -1.3 1.5 -1.5 -0.5 5.4e+03 0.074 0.5 0.32 + 3 -1.5 -1.3 1.5 -1.5 -0.5 5.4e+03 0.074 0.25 -0.18 - 4 -1.2 -1 1.8 -1.2 -0.25 5.3e+03 0.028 0.25 0.47 + 5 -1 -1.3 1.5 -1 -0.5 5.3e+03 0.04 0.25 0.16 + 6 -1 -1.3 1.5 -1 -0.5 5.3e+03 0.04 0.12 -0.059 - 7 -0.88 -1.2 1.5 -1.1 -0.38 5.3e+03 0.024 0.12 0.63 + 8 -1 -1.3 1.4 -1.2 -0.25 5.3e+03 0.017 0.12 0.33 + 9 -0.88 -1.4 1.5 -1.1 -0.24 5.3e+03 0.0093 0.12 0.88 + 10 -0.75 -1.5 1.7 -1.3 -0.11 5.2e+03 0.012 0.12 0.71 + 11 -0.68 -1.6 1.7 -1.1 -0.079 5.2e+03 0.0073 0.12 0.71 + 12 -0.56 -1.7 1.9 -1.2 -0.043 5.2e+03 0.0094 1.2 0.91 ++ 13 -0.56 -1.7 1.9 -1.2 -0.043 5.2e+03 0.0094 0.62 -2.1 - 14 -0.56 -1.7 1.9 -1.2 -0.043 5.2e+03 0.0094 0.31 -1.5 - 15 -0.56 -1.7 1.9 -1.2 -0.043 5.2e+03 0.0094 0.16 -0.15 - 16 -0.59 -1.9 2 -1.2 0.05 5.2e+03 0.01 0.16 0.38 + 17 -0.52 -1.9 2.2 -1.2 -0.011 5.2e+03 0.0047 0.16 0.81 + 18 -0.52 -1.9 2.2 -1.2 -0.011 5.2e+03 0.0047 0.078 -0.45 - 19 -0.52 -1.9 2.2 -1.2 0.067 5.2e+03 0.0079 0.078 0.31 + 20 -0.48 -2 2.3 -1.2 0.04 5.2e+03 0.0038 0.78 0.94 ++ 21 -0.48 -2 2.3 -1.2 0.04 5.2e+03 0.0038 0.39 -1.3 - 22 -0.48 -2 2.3 -1.2 0.04 5.2e+03 0.0038 0.2 -0.86 - 23 -0.48 -2 2.3 -1.2 0.04 5.2e+03 0.0038 0.098 0.032 - 24 -0.5 -2.1 2.4 -1.3 0.062 5.2e+03 0.0066 0.098 0.47 + 25 -0.45 -2.1 2.5 -1.2 0.095 5.2e+03 0.0047 0.098 0.83 + 26 -0.44 -2.2 2.6 -1.2 0.11 5.2e+03 0.0032 0.98 0.94 ++ 27 -0.44 -2.2 2.6 -1.2 0.11 5.2e+03 0.0032 0.49 -23 - 28 -0.44 -2.2 2.6 -1.2 0.11 5.2e+03 0.0032 0.24 -3.8 - 29 -0.44 -2.2 2.6 -1.2 0.11 5.2e+03 0.0032 0.12 -0.66 - 30 -0.41 -2.2 2.7 -1.3 0.1 5.2e+03 0.0028 0.12 0.28 + 31 -0.41 -2.2 2.7 -1.3 0.1 5.2e+03 0.0028 0.061 -0.31 - 32 -0.41 -2.2 2.7 -1.3 0.12 5.2e+03 0.00091 0.061 0.75 + 33 -0.38 -2.3 2.8 -1.3 0.13 5.2e+03 0.0022 0.061 0.29 + 34 -0.38 -2.3 2.8 -1.3 0.13 5.2e+03 0.0022 0.031 -0.061 - 35 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.0013 0.031 0.82 + 36 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.0013 0.015 -0.7 - 37 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.0013 0.0076 -3.2 - 38 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.0013 0.0038 -0.56 - 39 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.00038 0.0038 0.44 + 40 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.00042 0.038 0.9 ++ 41 -0.4 -2.3 2.8 -1.3 0.13 5.2e+03 0.00042 0.019 -0.73 - 42 -0.38 -2.3 2.8 -1.3 0.14 5.2e+03 0.0004 0.019 0.4 + 43 -0.39 -2.3 2.8 -1.3 0.14 5.2e+03 0.00074 0.019 0.44 + 44 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 0.00021 0.019 0.8 + 45 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 0.00011 0.019 0.81 + 46 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 0.00011 0.0065 -2.8 - 47 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 0.00011 0.0032 -0.35 - 48 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 8.5e-05 0.0032 0.56 + 49 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 2.8e-05 0.0032 0.45 + 50 -0.38 -2.3 2.9 -1.3 0.15 5.2e+03 2.8e-05 0.0012 -0.29 - 51 -0.39 -2.3 2.9 -1.3 0.14 5.2e+03 3.5e-05 0.0012 0.33 + 52 -0.39 -2.3 2.9 -1.3 0.14 5.2e+03 1.5e-06 0.0012 0.99 + Optimization algorithm has converged. Relative gradient: 1.5147804784176e-06 Cause of termination: Relative gradient = 1.5e-06 <= 6.1e-06 Number of function evaluations: 116 Number of gradient evaluations: 63 Number of hessian evaluations: 0 Algorithm: BFGS with trust region for simple bound constraints Number of iterations: 53 Proportion of Hessian calculation: 0/31 = 0.0% Optimization time: 0:00:01.077594 Calculate second derivatives and BHHH File b06c_unif_mixture_integral.html has been generated. File b06c_unif_mixture_integral.yaml has been generated. .. GENERATED FROM PYTHON SOURCE LINES 132-134 .. code-block:: Python print(results.short_summary()) .. rst-class:: sphx-glr-script-out .. code-block:: none Results for model b06c_unif_mixture_integral Nbr of parameters: 5 Sample size: 6768 Excluded data: 3960 Final log likelihood: -5215.061 Akaike Information Criterion: 10440.12 Bayesian Information Criterion: 10474.22 .. GENERATED FROM PYTHON SOURCE LINES 135-137 .. code-block:: Python pandas_results = get_pandas_estimated_parameters(estimation_results=results) display(pandas_results) .. rst-class:: sphx-glr-script-out .. code-block:: none Name Value Robust std err. Robust t-stat. Robust p-value 0 asc_train -0.385072 0.065992 -5.835159 5.373928e-09 1 b_time -2.320575 0.126118 -18.400027 0.000000e+00 2 b_time_s 2.875959 0.200170 14.367615 0.000000e+00 3 b_cost -1.277926 0.086624 -14.752624 0.000000e+00 4 asc_car 0.144969 0.053308 2.719456 6.538948e-03 .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.402 seconds) .. _sphx_glr_download_auto_examples_swissmetro_plot_b06c_unif_mixture_integral.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_b06c_unif_mixture_integral.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_b06c_unif_mixture_integral.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_b06c_unif_mixture_integral.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_