{ "cells": [ { "cell_type": "code", "metadata": {}, "source": [ "#%%\n", "\"\"\"File 03nestedElasticities.py\n", "\n", ":author: Michel Bierlaire, EPFL\n", ":date: Wed Sep 11 10:28:11 2019\n", "\n", " We use a previously estimated nested logit model.\n", " Three alternatives: public transporation, car and slow modes.\n", " RP data.\n", " We calculate disaggregate and aggregate direct point elasticities.\n", "\"\"\"\n", "\n", "import sys\n", "import pandas as pd\n", "import biogeme.database as db\n", "import biogeme.biogeme as bio\n", "from biogeme import models\n", "import biogeme.results as res\n", "from biogeme.expressions import Beta, Derive\n", "\n", "\n", "# Read the data\n", "df = pd.read_csv('optima.dat', sep='\\t')\n", "database = db.Database('optima', df)\n", "\n", "# The following statement allows you to use the names of the variable\n", "# as Python variable.\n", "globals().update(database.variables)\n", "\n", "# Exclude observations such that the chosen alternative is -1\n", "database.remove(Choice == -1.0)\n", "\n", "# Normalize the weights\n", "sumWeight = database.data['Weight'].sum()\n", "numberOfRows = database.data.shape[0]\n", "normalizedWeight = Weight * numberOfRows / sumWeight\n", "\n", "# Calculate the number of accurences of a value in the database\n", "numberOfMales = database.count('Gender', 1)\n", "print(f'Number of males: {numberOfMales}')\n", "\n", "numberOfFemales = database.count('Gender', 2)\n", "print(f'Number of females: {numberOfFemales}')\n", "\n", "# For more complex conditions, we use Pandas\n", "unreportedGender = database.data[\n", " (database.data['Gender'] != 1) & (database.data['Gender'] != 2)\n", "].count()['Gender']\n", "print(f'Unreported gender: {unreportedGender}')\n", "\n", "# List of parameters. Their value will be set later.\n", "ASC_CAR = Beta('ASC_CAR', 0, None, None, 0)\n", "ASC_PT = Beta('ASC_PT', 0, None, None, 1)\n", "ASC_SM = Beta('ASC_SM', 0, None, None, 0)\n", "BETA_TIME_FULLTIME = Beta('BETA_TIME_FULLTIME', 0, None, None, 0)\n", "BETA_TIME_OTHER = Beta('BETA_TIME_OTHER', 0, None, None, 0)\n", "BETA_DIST_MALE = Beta('BETA_DIST_MALE', 0, None, None, 0)\n", "BETA_DIST_FEMALE = Beta('BETA_DIST_FEMALE', 0, None, None, 0)\n", "BETA_DIST_UNREPORTED = Beta('BETA_DIST_UNREPORTED', 0, None, None, 0)\n", "BETA_COST = Beta('BETA_COST', 0, None, None, 0)\n", "\n", "# Define new variables. Must be consistent with estimation results.\n", "TimePT_scaled = TimePT / 200\n", "TimeCar_scaled = TimeCar / 200\n", "MarginalCostPT_scaled = MarginalCostPT / 10\n", "CostCarCHF_scaled = CostCarCHF / 10\n", "distance_km_scaled = distance_km / 5\n", "male = Gender == 1\n", "female = Gender == 2\n", "unreportedGender = Gender == -1\n", "fulltime = OccupStat == 1\n", "notfulltime = OccupStat != 1\n", "\n", "# Definition of utility functions:\n", "V_PT = (\n", " ASC_PT\n", " + BETA_TIME_FULLTIME * TimePT_scaled * fulltime\n", " + BETA_TIME_OTHER * TimePT_scaled * notfulltime\n", " + BETA_COST * MarginalCostPT_scaled\n", ")\n", "V_CAR = (\n", " ASC_CAR\n", " + BETA_TIME_FULLTIME * TimeCar_scaled * fulltime\n", " + BETA_TIME_OTHER * TimeCar_scaled * notfulltime\n", " + BETA_COST * CostCarCHF_scaled\n", ")\n", "V_SM = (\n", " ASC_SM\n", " + BETA_DIST_MALE * distance_km_scaled * male\n", " + BETA_DIST_FEMALE * distance_km_scaled * female\n", " + BETA_DIST_UNREPORTED * distance_km_scaled * unreportedGender\n", ")\n", "\n", "# Associate utility functions with the numbering of alternatives\n", "V = {0: V_PT, 1: V_CAR, 2: V_SM}\n", "\n", "# Definition of the nests:\n", "# 1: nests parameter\n", "# 2: list of alternatives\n", "\n", "MU_NOCAR = Beta('MU_NOCAR', 1.0, 1.0, None, 0)\n", "\n", "CAR_NEST = 1.0, [1]\n", "NO_CAR_NEST = MU_NOCAR, [0, 2]\n", "nests = CAR_NEST, NO_CAR_NEST\n", "\n", "# The choice model is a nested logit\n", "prob_pt = models.nested(V, None, nests, 0)\n", "prob_car = models.nested(V, None, nests, 1)\n", "prob_sm = models.nested(V, None, nests, 2)\n", "\n", "# Calculation of the direct elasticities.\n", "# We use the 'Derive' operator to calculate the derivatives.\n", "direct_elas_pt_time = Derive(prob_pt, 'TimePT') * TimePT / prob_pt\n", "direct_elas_pt_cost = (\n", " Derive(prob_pt, 'MarginalCostPT') * MarginalCostPT / prob_pt\n", ")\n", "direct_elas_car_time = Derive(prob_car, 'TimeCar') * TimeCar / prob_car\n", "direct_elas_car_cost = Derive(prob_car, 'CostCarCHF') * CostCarCHF / prob_car\n", "direct_elas_sm_dist = Derive(prob_sm, 'distance_km') * distance_km / prob_sm\n", "\n", "simulate = {\n", " 'weight': normalizedWeight,\n", " 'Prob. car': prob_car,\n", " 'Prob. public transportation': prob_pt,\n", " 'Prob. slow modes': prob_sm,\n", " 'direct_elas_pt_time': direct_elas_pt_time,\n", " 'direct_elas_pt_cost': direct_elas_pt_cost,\n", " 'direct_elas_car_time': direct_elas_car_time,\n", " 'direct_elas_car_cost': direct_elas_car_cost,\n", " 'direct_elas_sm_dist': direct_elas_sm_dist,\n", "}\n", "\n", "biogeme = bio.BIOGEME(database, simulate)\n", "biogeme.modelName = '03nestedElasticties'\n", "\n", "# Read the estimation results from the file\n", "try:\n", " results = res.bioResults(pickleFile='01nestedEstimation.pickle')\n", "except FileNotFoundError:\n", " sys.exit(\n", " 'Run first the script 01nestedEstimation.py in order to generate '\n", " 'the file 01nestedEstimation.pickle.'\n", " )\n", "\n", "# simulatedValues is a Panda dataframe with the same number of rows as\n", "# the database, and as many columns as formulas to simulate.\n", "simulatedValues = biogeme.simulate(results.getBetaValues())\n", "\n", "# We calculate the elasticities\n", "simulatedValues['Weighted prob. car'] = (\n", " simulatedValues['weight'] * simulatedValues['Prob. car']\n", ")\n", "simulatedValues['Weighted prob. PT'] = (\n", " simulatedValues['weight'] * simulatedValues['Prob. public transportation']\n", ")\n", "simulatedValues['Weighted prob. SM'] = (\n", " simulatedValues['weight'] * simulatedValues['Prob. slow modes']\n", ")\n", "\n", "denominator_car = simulatedValues['Weighted prob. car'].sum()\n", "denominator_pt = simulatedValues['Weighted prob. PT'].sum()\n", "denominator_sm = simulatedValues['Weighted prob. SM'].sum()\n", "\n", "direct_elas_term_car_time = (\n", " simulatedValues['Weighted prob. car']\n", " * simulatedValues['direct_elas_car_time']\n", " / denominator_car\n", ").sum()\n", "print(\n", " f'Aggregate direct point elasticity of car wrt time: '\n", " f'{direct_elas_term_car_time:.3g}'\n", ")\n", "\n", "direct_elas_term_car_cost = (\n", " simulatedValues['Weighted prob. car']\n", " * simulatedValues['direct_elas_car_cost']\n", " / denominator_car\n", ").sum()\n", "print(\n", " f'Aggregate direct point elasticity of car wrt cost: '\n", " f'{direct_elas_term_car_cost:.3g}'\n", ")\n", "\n", "direct_elas_term_pt_time = (\n", " simulatedValues['Weighted prob. PT']\n", " * simulatedValues['direct_elas_pt_time']\n", " / denominator_pt\n", ").sum()\n", "print(\n", " f'Aggregate direct point elasticity of PT wrt time: '\n", " f'{direct_elas_term_pt_time:.3g}'\n", ")\n", "\n", "direct_elas_term_pt_cost = (\n", " simulatedValues['Weighted prob. PT']\n", " * simulatedValues['direct_elas_pt_cost']\n", " / denominator_pt\n", ").sum()\n", "print(\n", " f'Aggregate direct point elasticity of PT wrt cost: '\n", " f'{direct_elas_term_pt_cost:.3g}'\n", ")\n", "\n", "direct_elas_term_sm_dist = (\n", " simulatedValues['Weighted prob. SM']\n", " * simulatedValues['direct_elas_sm_dist']\n", " / denominator_sm\n", ").sum()\n", "print(\n", " f'Aggregate direct point elasticity of SM wrt distance: '\n", " f'{direct_elas_term_sm_dist:.3g}'\n", ")\n" ], "outputs": [], "execution_count": null } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.1" } }, "nbformat": 4, "nbformat_minor": 4 }