{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "6452580a",
   "metadata": {},
   "source": [
    "\n",
    "File: 07-lecture_solution.py\n",
    "\n",
    "\n",
    "Michel Bierlaire\n",
    "\n",
    "Wed Aug 06 2025, 12:12:22\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bc1d6ec9",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "import biogeme.biogeme_logging as blog\n",
    "from IPython.core.display_functions import display\n",
    "from biogeme.biogeme import BIOGEME\n",
    "from biogeme.expressions import Beta, log\n",
    "from biogeme.models import logit, loglogit\n",
    "from biogeme.results_processing import (\n",
    "    EstimationResults,\n",
    "    compile_estimation_results,\n",
    "    get_pandas_estimated_parameters,\n",
    ")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0beecabe",
   "metadata": {},
   "source": [
    "Variables used for the specification of the Swissmetro model are defined in the file `swissmetro_variables.py`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "834a1695",
   "metadata": {},
   "outputs": [],
   "source": [
    "from swissmetro_variables import (\n",
    "    CAR_AV_SP,\n",
    "    CAR_CO_SCALED,\n",
    "    CAR_TT_SCALED,\n",
    "    CHOICE,\n",
    "    SM_AV,\n",
    "    SM_COST_SCALED,\n",
    "    SM_HE_SCALED,\n",
    "    SM_TT_SCALED,\n",
    "    TRAIN_AV_SP,\n",
    "    TRAIN_COST_SCALED,\n",
    "    TRAIN_HE_SCALED,\n",
    "    TRAIN_TT_SCALED,\n",
    "    database,\n",
    ")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5d3b5ced",
   "metadata": {},
   "source": [
    "As the estimation time may be long, we ask Biogeme to report the details of the iterations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8434a070",
   "metadata": {},
   "outputs": [],
   "source": [
    "logger = blog.get_screen_logger(level=blog.INFO)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2592533c",
   "metadata": {},
   "source": [
    "# Parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bb8d97e0",
   "metadata": {},
   "outputs": [],
   "source": [
    "asc_car = Beta('asc_car', 0, None, None, 0)\n",
    "asc_train = Beta('asc_train', 0, None, None, 0)\n",
    "b_time = Beta('b_time', 0, None, None, 0)\n",
    "b_cost = Beta('b_cost', 0, None, None, 0)\n",
    "b_fr = Beta('b_fr', 0, None, None, 0)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "58129af4",
   "metadata": {},
   "source": [
    "# Availability conditions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "296c93bd",
   "metadata": {},
   "outputs": [],
   "source": [
    "av = {1: TRAIN_AV_SP, 2: SM_AV, 3: CAR_AV_SP}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ad8dc7b3",
   "metadata": {},
   "source": [
    "# Logit model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71fd29ae",
   "metadata": {},
   "source": [
    "## Utility functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d8b17a27",
   "metadata": {},
   "outputs": [],
   "source": [
    "v_train = (\n",
    "    asc_train\n",
    "    + b_time * TRAIN_TT_SCALED\n",
    "    + b_cost * TRAIN_COST_SCALED\n",
    "    + b_fr * TRAIN_HE_SCALED\n",
    ")\n",
    "v_swissmetro = b_time * SM_TT_SCALED + b_cost * SM_COST_SCALED + b_fr * SM_HE_SCALED\n",
    "v_car = asc_car + b_time * CAR_TT_SCALED + b_cost * CAR_CO_SCALED\n",
    "v = {1: v_train, 2: v_swissmetro, 3: v_car}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "90bec307",
   "metadata": {},
   "source": [
    "## Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "caa08a63",
   "metadata": {},
   "outputs": [],
   "source": [
    "logprob = loglogit(v, av, CHOICE)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "23ee1bae",
   "metadata": {},
   "source": [
    "## Estimation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0549e0d0",
   "metadata": {},
   "outputs": [],
   "source": [
    "biogeme = BIOGEME(database, logprob)\n",
    "biogeme.model_name = '01logit'\n",
    "results_logit: EstimationResults = biogeme.estimate(recycle=True)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7441b601",
   "metadata": {},
   "source": [
    "## Results"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "51960fc6",
   "metadata": {},
   "source": [
    "General statistics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1a939b7b",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(results_logit.short_summary())\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5bf6faee",
   "metadata": {},
   "source": [
    "Estimated parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a9bbe136",
   "metadata": {},
   "outputs": [],
   "source": [
    "param_logit = get_pandas_estimated_parameters(estimation_results=results_logit)\n",
    "display(param_logit)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a380562f",
   "metadata": {},
   "source": [
    "# Random parameter: normal distribution"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "00ca6cd8",
   "metadata": {},
   "source": [
    "Read the results from file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5a72a157",
   "metadata": {},
   "outputs": [],
   "source": [
    "results_normal = EstimationResults.from_yaml_file(filename='02normal.yaml')\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "231c4079",
   "metadata": {},
   "source": [
    "# Random parameter: lognormal distribution"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a9128e65",
   "metadata": {},
   "source": [
    "Read the results from file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a5912b6e",
   "metadata": {},
   "outputs": [],
   "source": [
    "results_lognormal = EstimationResults.from_yaml_file(filename='03lognormal.yaml')\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3124f768",
   "metadata": {},
   "source": [
    "# Latent classes"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "917275fd",
   "metadata": {},
   "source": [
    "We consider two classes in the population. The first class of individuals have considered all variables when making\n",
    "their choice. For them, the specification of the utility function is the same as for the logit model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "38f27656",
   "metadata": {},
   "outputs": [],
   "source": [
    "v_train_class_1 = (\n",
    "    asc_train\n",
    "    + b_time * TRAIN_TT_SCALED\n",
    "    + b_cost * TRAIN_COST_SCALED\n",
    "    + b_fr * TRAIN_HE_SCALED\n",
    ")\n",
    "v_swissmetro_class_1 = (\n",
    "    b_time * SM_TT_SCALED + b_cost * SM_COST_SCALED + b_fr * SM_HE_SCALED\n",
    ")\n",
    "v_car_class_1 = asc_car + b_time * CAR_TT_SCALED + b_cost * CAR_CO_SCALED\n",
    "v_class_1 = {1: v_train_class_1, 2: v_swissmetro_class_1, 3: v_car_class_1}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c5563365",
   "metadata": {},
   "source": [
    "The second class of individuals ignored the travel time variable when making the choice. Therefore, this variable\n",
    "is removed from the utility function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b4bfaf34",
   "metadata": {},
   "outputs": [],
   "source": [
    "v1_class_2 = asc_train + b_cost * TRAIN_COST_SCALED + b_fr * TRAIN_HE_SCALED\n",
    "v2_class_2 = b_cost * SM_COST_SCALED + b_fr * SM_HE_SCALED\n",
    "v3_class_2 = asc_car + b_cost * CAR_CO_SCALED\n",
    "v_class_2 = {1: v1_class_2, 2: v2_class_2, 3: v3_class_2}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a59f8f0c",
   "metadata": {},
   "source": [
    "The following parameter captures the probability to belong to class 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c95c232c",
   "metadata": {},
   "outputs": [],
   "source": [
    "omega = Beta('omega', 0.5, 0, 1, 0)\n",
    "prob_class_1 = omega\n",
    "prob_class_2 = 1 - omega\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "51ad4dad",
   "metadata": {},
   "source": [
    "## Model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6316fa5b",
   "metadata": {},
   "source": [
    "We first calculate the choice probability for each class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1dd0de88",
   "metadata": {},
   "outputs": [],
   "source": [
    "choice_prob_class_1 = logit(v_class_1, av, CHOICE)\n",
    "choice_prob_class_2 = logit(v_class_2, av, CHOICE)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3fc7d003",
   "metadata": {},
   "source": [
    "The choice probability is obtained by using the class membership model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f1b597c6",
   "metadata": {},
   "outputs": [],
   "source": [
    "choice_prob = prob_class_1 * choice_prob_class_1 + prob_class_2 * choice_prob_class_2\n",
    "logprob = log(choice_prob)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "626b2697",
   "metadata": {},
   "source": [
    "## Estimation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5a878e0f",
   "metadata": {},
   "outputs": [],
   "source": [
    "biogeme = BIOGEME(database, logprob)\n",
    "biogeme.model_name = '04latentClass'\n",
    "results_latent: EstimationResults = biogeme.estimate()\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f2f7c4f4",
   "metadata": {},
   "source": [
    "## Results"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2c390637",
   "metadata": {},
   "source": [
    "General statistics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8e0979f9",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "print(results_latent.short_summary())\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "edff4373",
   "metadata": {},
   "source": [
    "Estimated parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "79566bf0",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "param_latent = get_pandas_estimated_parameters(estimation_results=results_latent)\n",
    "display(param_latent)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6abeb497",
   "metadata": {},
   "source": [
    "# Comparison"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2e1b2794",
   "metadata": {},
   "source": [
    "We build a summary data frame."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "90f3624c",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "summary = compile_estimation_results(\n",
    "    {\n",
    "        'Logit': results_logit,\n",
    "        'Random param. (normal)': results_normal,\n",
    "        'Random param. (lognormal)': results_lognormal,\n",
    "        'Latent class': results_latent,\n",
    "    }\n",
    ")\n",
    "display(summary[0])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ef4401a",
   "metadata": {},
   "source": [
    "We note that the probability to belong to class one is..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d5bd3b8",
   "metadata": {},
   "outputs": [],
   "source": [
    "estimated_proba_class_1 = results_latent.get_parameter_value('omega')\n",
    "print(f'Probability to belong to class 1: {100 * estimated_proba_class_1:.3g}%.')"
   ]
  }
 ],
 "metadata": {},
 "nbformat": 4,
 "nbformat_minor": 5
}
