Applying a location and mode choice model to populations¶

This notebook applies a simple location and mode choice model to a PAM population.

The pam.planner.choice.ChoiceMNL class allows the user to apply an MNL specification for selecting the location of activities and the mode for accessing them, given person characteristics, network conditions and/or zone attraction data.

The typical workflow goes as follows:

choice_model = ChoiceMNL(population, od, zones) # initialize the model and point to the data objects 
choice_model.configure(u, scope) # configure the model by specifying a utility function and the scope of application.
choice_model.apply() # apply the model and update the population with the results.

In [1]:

import logging
import os
import random

import numpy as np
import pandas as pd
from pam.operations.cropping import link_population
from pam.planner import choice_location as choice
from pam.planner.od import ODFactory, ODMatrix
from pam.read import read_matsim
from prettytable import PrettyTable

logging.basicConfig(level=logging.DEBUG)
random.seed(0)

import logging import os import random import numpy as np import pandas as pd from pam.operations.cropping import link_population from pam.planner import choice_location as choice from pam.planner.od import ODFactory, ODMatrix from pam.read import read_matsim from prettytable import PrettyTable logging.basicConfig(level=logging.DEBUG) random.seed(0)

/var/folders/6n/0h9tynqn581fxsytcc863h94tm217b/T/ipykernel_95065/106397807.py:6: DeprecationWarning: 
Pyarrow will become a required dependency of pandas in the next major release of pandas (pandas 3.0),
(to allow more performant data types, such as the Arrow string type, and better interoperability with other libraries)
but was not found to be installed on your system.
If this would cause problems for you,
please provide us feedback at https://github.com/pandas-dev/pandas/issues/54466
        
  import pandas as pd

Data¶

We read an example population, and set the location of all activities to zone a:

In [2]:

population = read_matsim(os.path.join("..", "tests", "test_data", "test_matsim_plansv12.xml"))
link_population(population)
for _hid, _pid, person in population.people():
    for act in person.activities:
        act.location.area = "a"

population = read_matsim(os.path.join("..", "tests", "test_data", "test_matsim_plansv12.xml")) link_population(population) for _hid, _pid, person in population.people(): for act in person.activities: act.location.area = "a"

In [3]:

def print_activity_locs(population, act_scope="work"):
    summary = PrettyTable(["pid", "seq", "location", "mode"])
    for _hid, pid, person in population.people():
        for seq, act in enumerate(person.plan.activities):
            if (act.act == act_scope) or (act_scope == "all"):
                trmode = act.previous.mode if act.previous is not None else "NA"
                summary.add_row([pid, seq, act.location.area, trmode])

    print(summary)


print("Work locations and travel modes:")
print_activity_locs(population, act_scope="work")

def print_activity_locs(population, act_scope="work"): summary = PrettyTable(["pid", "seq", "location", "mode"]) for _hid, pid, person in population.people(): for seq, act in enumerate(person.plan.activities): if (act.act == act_scope) or (act_scope == "all"): trmode = act.previous.mode if act.previous is not None else "NA" summary.add_row([pid, seq, act.location.area, trmode]) print(summary) print("Work locations and travel modes:") print_activity_locs(population, act_scope="work")

Work locations and travel modes:
+--------+-----+----------+------+
|  pid   | seq | location | mode |
+--------+-----+----------+------+
| chris  |  1  |    a     | car  |
| fatema |  1  |    a     | bike |
|  fred  |  3  |    a     | walk |
| gerry  |  3  |    a     | walk |
|  nick  |  1  |    a     | car  |
+--------+-----+----------+------+

Our zones dataset includes destination attraction data, for example the number of jobs or schools in each likely destination zone:

In [4]:

data_zones = pd.DataFrame({"zone": ["a", "b"], "jobs": [100, 200], "schools": [3, 1]}).set_index(
    "zone"
)
data_zones

data_zones = pd.DataFrame({"zone": ["a", "b"], "jobs": [100, 200], "schools": [3, 1]}).set_index( "zone" ) data_zones

Out[4]:

	jobs	schools
zone
a	100	3
b	200	1

The od object holds origin-destination data, for example travel time and travel distance between each origin and destination, for each travel mode:

In [5]:

zone_labels = ("a", "b")
od = ODFactory.from_matrices(
    [
        ODMatrix("time", "car", zone_labels, zone_labels, np.array([[20, 40], [40, 20]])),
        ODMatrix("time", "bus", zone_labels, zone_labels, np.array([[30, 45], [45, 30]])),
        ODMatrix("distance", "car", zone_labels, zone_labels, np.array([[5, 8], [8, 5]])),
        ODMatrix("distance", "bus", zone_labels, zone_labels, np.array([[5, 9], [9, 5]])),
    ]
)
od

zone_labels = ("a", "b") od = ODFactory.from_matrices( [ ODMatrix("time", "car", zone_labels, zone_labels, np.array([[20, 40], [40, 20]])), ODMatrix("time", "bus", zone_labels, zone_labels, np.array([[30, 45], [45, 30]])), ODMatrix("distance", "car", zone_labels, zone_labels, np.array([[5, 8], [8, 5]])), ODMatrix("distance", "bus", zone_labels, zone_labels, np.array([[5, 9], [9, 5]])), ] ) od

Out[5]:

Origin-destination dataset 
--------------------------------------------------
Labels(vars=['time', 'distance'], origin_zones=('a', 'b'), destination_zones=('a', 'b'), mode=['car', 'bus'])
--------------------------------------------------
time - car:
[[20. 40.]
 [40. 20.]]
--------------------------------------------------
time - bus:
[[30. 45.]
 [45. 30.]]
--------------------------------------------------
distance - car:
[[5. 8.]
 [8. 5.]]
--------------------------------------------------
distance - bus:
[[5. 9.]
 [9. 5.]]
--------------------------------------------------

The dimensions of the od object are always (in order): variables, origin zone, destination zone, and mode. It can be sliced using the respective labels under od.labels, for example:

In [6]:

od["time", "a", "b", "bus"]

od["time", "a", "b", "bus"]

Out[6]:

45.0

Choice model¶

In [7]:

planner = choice.ChoiceMNL(population, od, data_zones)

planner = choice.ChoiceMNL(population, od, data_zones)

INFO:pam.planner.choice_location:Updated model configuration

INFO:pam.planner.choice_location:ChoiceConfiguration(u=None, scope=None, func_probabilities=<function calculate_mnl_probabilities at 0x12a98dea0>, func_sampling=<function sample_weighted at 0x12a98de10>)

We configure the model by specifying:

• the scope of the model. For example, work activities.
• the utility formulation of each alternative.

Both settings are defined as strings. The stings may comprise mathematical operators, coefficients, planner data objects (od / zones), and/or PAM population objects (person/ act).

Coefficients can be passed either as a number, or as a list, with each element in the list corresponding to one of the modes in the od object.

In [8]:

scope = "act.act=='work'"
asc = [0, -1]  # one value for each mode, 0->car, -1->
asc_shift_poor = [0, 2]  # one value for each mode
beta_time = [-0.05, -0.07]  # one value for each mode
beta_zones = 0.4
u = f""" \
    {asc} + \
    (np.array({asc_shift_poor}) * (person.attributes['subpopulation']=='poor')) + \
    ({beta_time} * od['time', person.home.area]) + \
    ({beta_zones} * np.log(zones['jobs']))
"""

planner.configure(u=u, scope=scope)

scope = "act.act=='work'" asc = [0, -1] # one value for each mode, 0->car, -1-> asc_shift_poor = [0, 2] # one value for each mode beta_time = [-0.05, -0.07] # one value for each mode beta_zones = 0.4 u = f""" \ {asc} + \ (np.array({asc_shift_poor}) * (person.attributes['subpopulation']=='poor')) + \ ({beta_time} * od['time', person.home.area]) + \ ({beta_zones} * np.log(zones['jobs'])) """ planner.configure(u=u, scope=scope)

INFO:pam.planner.choice_location:Updated model configuration

INFO:pam.planner.choice_location:ChoiceConfiguration(u="[0,-1]+(np.array([0,2])*(person.attributes['subpopulation']=='poor'))+([-0.05,-0.07]*od['time',person.home.area])+(0.4*np.log(zones['jobs']))\n", scope="act.act=='work'", func_probabilities=<function calculate_mnl_probabilities at 0x12a98dea0>, func_sampling=<function sample_weighted at 0x12a98de10>)

The .get_choice_set() provides with with the utilities of each alternative, as perceived by each agent.

In [9]:

choice_set = planner.get_choice_set()
print("Activities in scope: \n", choice_set.idxs)
print("\nAlternatives: \n", choice_set.choice_labels)
print("\nChoice set utilities: \n", choice_set.u_choices)

choice_set = planner.get_choice_set() print("Activities in scope: \n", choice_set.idxs) print("\nAlternatives: \n", choice_set.choice_labels) print("\nChoice set utilities: \n", choice_set.u_choices)

Activities in scope: 
 [ChoiceIdx(pid='chris', hid='chris', seq=1, act=<pam.activity.Activity object at 0x12a9af700>), ChoiceIdx(pid='fatema', hid='fatema', seq=1, act=<pam.activity.Activity object at 0x12a9afaf0>), ChoiceIdx(pid='fred', hid='fred', seq=3, act=<pam.activity.Activity object at 0x12a9925f0>), ChoiceIdx(pid='gerry', hid='gerry', seq=3, act=<pam.activity.Activity object at 0x12a9902e0>), ChoiceIdx(pid='nick', hid='nick', seq=1, act=<pam.activity.Activity object at 0x12a992560>)]

Alternatives: 
 [ChoiceLabel(destination='a', mode='car'), ChoiceLabel(destination='a', mode='bus'), ChoiceLabel(destination='b', mode='car'), ChoiceLabel(destination='b', mode='bus')]

Choice set utilities: 
 [[ 0.84206807 -1.25793193  0.11932695 -2.03067305]
 [ 0.84206807  0.74206807  0.11932695 -0.03067305]
 [ 0.84206807  0.74206807  0.11932695 -0.03067305]
 [ 0.84206807  0.74206807  0.11932695 -0.03067305]
 [ 0.84206807 -1.25793193  0.11932695 -2.03067305]]

The .apply() method samples from the alternatives, and updates the location and mode of each activity accordingly:

In [10]:

planner.apply()
print("Sampled choices: \n", planner._selections.selections)

planner.apply() print("Sampled choices: \n", planner._selections.selections)

INFO:pam.planner.choice_location:Applying choice model...

INFO:pam.planner.choice_location:Configuration: 
ChoiceConfiguration(u="[0,-1]+(np.array([0,2])*(person.attributes['subpopulation']=='poor'))+([-0.05,-0.07]*od['time',person.home.area])+(0.4*np.log(zones['jobs']))\n", scope="act.act=='work'", func_probabilities=<function calculate_mnl_probabilities at 0x12a98dea0>, func_sampling=<function sample_weighted at 0x12a98de10>)

INFO:pam.planner.choice_location:Choice model application complete.

Sampled choices: 
 [ChoiceLabel(destination='b', mode='car'), ChoiceLabel(destination='b', mode='car'), ChoiceLabel(destination='a', mode='bus'), ChoiceLabel(destination='a', mode='car'), ChoiceLabel(destination='a', mode='car')]

The population's activity locations and travel modes have now been updated accordingly:

In [11]:

print_activity_locs(planner.population)

print_activity_locs(planner.population)

+--------+-----+----------+------+
|  pid   | seq | location | mode |
+--------+-----+----------+------+
| chris  |  1  |    b     | car  |
| fatema |  1  |    b     | car  |
|  fred  |  3  |    a     | bus  |
| gerry  |  3  |    a     | car  |
|  nick  |  1  |    a     | car  |
+--------+-----+----------+------+