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 prettytable import PrettyTable

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

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

import logging import os import random import numpy as np import pandas as pd from prettytable import PrettyTable 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 logging.basicConfig(level=logging.DEBUG) random.seed(0)

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 0x16565d300>, func_sampling=<function sample_weighted at 0x16565d3a0>)

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 0x16565d300>, func_sampling=<function sample_weighted at 0x16565d3a0>)

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 0x165694b90>), ChoiceIdx(pid='fatema', hid='fatema', seq=1, act=<pam.activity.Activity object at 0x165695690>), ChoiceIdx(pid='fred', hid='fred', seq=3, act=<pam.activity.Activity object at 0x165696ad0>), ChoiceIdx(pid='gerry', hid='gerry', seq=3, act=<pam.activity.Activity object at 0x1656a0250>), ChoiceIdx(pid='nick', hid='nick', seq=1, act=<pam.activity.Activity object at 0x1656a1250>)]

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 0x16565d300>, func_sampling=<function sample_weighted at 0x16565d3a0>)
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  |
+--------+-----+----------+------+