Note
Go to the end to download the full example code.
Route Choice#
In this example, we show how to perform route choice set generation using BFSLE and Link penalisation, for a city in La Serena Metropolitan Area in Chile.
References
See also
Several functions, methods, classes and modules are used in this example:
# Imports
from uuid import uuid4
from tempfile import gettempdir
from os.path import join
from aequilibrae.utils.create_example import create_example
# We create the example project inside our temp folder
fldr = join(gettempdir(), uuid4().hex)
project = create_example(fldr, "coquimbo")
import logging
import sys
# When the project opens, we can tell the logger to direct all messages to the terminal as well
logger = project.logger
stdout_handler = logging.StreamHandler(sys.stdout)
formatter = logging.Formatter("%(asctime)s;%(levelname)s ; %(message)s")
stdout_handler.setFormatter(formatter)
logger.addHandler(stdout_handler)
Model parameters#
import numpy as np
We’ll set the parameters for our route choice model. These are the parameters that will be used to calculate the utility of each path. In our example, the utility is equal to \(distance * theta\), and the path overlap factor (PSL) is equal to \(beta\).
# Distance factor
theta = 0.00011
# PSL parameter
beta = 1.1
Let’s select a set of nodes of interest
nodes_of_interest = (71645, 74089, 77011, 79385)
Let’s build all graphs
project.network.build_graphs()
# We get warnings that several fields in the project are filled with NaNs.
# This is true, but we won't use those fields.
2024-11-08 03:02:30,075;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
2024-11-08 03:02:30,155;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
2024-11-08 03:02:30,249;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
2024-11-08 03:02:30,340;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
We also see what graphs are available
project.network.graphs.keys()
dict_keys(['b', 'c', 't', 'w'])
We grab the graph for cars
graph = project.network.graphs["c"]
# Let's say that utility is just a function of distance, so we build our 'utility' field as distance * theta
graph.network = graph.network.assign(utility=graph.network.distance * theta)
# Prepare the graph with all nodes of interest as centroids
graph.prepare_graph(np.array(nodes_of_interest))
# And set the cost of the graph the as the utility field just created
graph.set_graph("utility")
# We allow flows through "centroid connectors" because our centroids are not really centroids
# If we have actual centroid connectors in the network (and more than one per centroid), then we
# should remove them from the graph
graph.set_blocked_centroid_flows(False)
2024-11-08 03:02:30,407;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
Mock demand matrix#
We’ll create a mock demand matrix with demand 1 for every zone and prepare it for computation.
from aequilibrae.matrix import AequilibraeMatrix
names_list = ["demand", "5x demand"]
mat = AequilibraeMatrix()
mat.create_empty(zones=graph.num_zones, matrix_names=names_list, memory_only=True)
mat.index = graph.centroids[:]
mat.matrices[:, :, 0] = np.full((graph.num_zones, graph.num_zones), 10.0)
mat.matrices[:, :, 1] = np.full((graph.num_zones, graph.num_zones), 50.0)
mat.computational_view()
Create plot function#
Before dive into the Route Choice class, let’s define a function to plot assignment results.
import folium
def plot_results(link_loads):
link_loads = link_loads[link_loads.tot > 0]
max_load = link_loads["tot"].max()
links = project.network.links.data
loaded_links = links.merge(link_loads, on="link_id", how="inner")
loads_lyr = folium.FeatureGroup("link_loads")
# Maximum thickness we would like is probably a 10, so let's make sure we don't go over that
factor = 10 / max_load
# Let's create the layers
for _, rec in loaded_links.iterrows():
points = rec.geometry.wkt.replace("LINESTRING ", "").replace("(", "").replace(")", "").split(", ")
points = "[[" + "],[".join([p.replace(" ", ", ") for p in points]) + "]]"
# we need to take from x/y to lat/long
points = [[x[1], x[0]] for x in eval(points)]
_ = folium.vector_layers.PolyLine(
points,
tooltip=f"link_id: {rec.link_id}, Flow: {rec.tot:.3f}",
color="red",
weight=factor * rec.tot,
).add_to(loads_lyr)
long, lat = project.conn.execute("select avg(xmin), avg(ymin) from idx_links_geometry").fetchone()
map_osm = folium.Map(location=[lat, long], tiles="Cartodb Positron", zoom_start=12)
loads_lyr.add_to(map_osm)
folium.LayerControl().add_to(map_osm)
return map_osm
Route Choice class#
Here we’ll construct and use the Route Choice class to generate our route sets
from aequilibrae.paths import RouteChoice
This object construct might take a minute depending on the size of the graph due to the construction of the compressed link to network link mapping that’s required. This is a one time operation per graph and is cached.
rc = RouteChoice(graph)
# Let's check the default parameters for the Route Choice class
print(rc.default_parameters)
{'generic': {'seed': 0, 'max_routes': 0, 'max_depth': 0, 'max_misses': 100, 'penalty': 1.01, 'cutoff_prob': 0.0, 'beta': 1.0, 'store_results': True}, 'link-penalisation': {}, 'bfsle': {'penalty': 1.0}}
Let’s add the demand. If it’s not provided, link loading cannot be preformed.
rc.add_demand(mat)
It is highly recommended to set either max_routes or max_depth to prevent runaway results.
rc.set_choice_set_generation("bfsle", max_routes=5)
We can now perform a computation for single OD pair if we’d like. Here we do one between the first and last centroid as well as an assignment.
results = rc.execute_single(77011, 74089, demand=1.0)
print(results[0])
(24222, 30332, 30333, 10435, 30068, 30069, 14198, 14199, 31161, 30928, 30929, 30930, 30931, 24172, 30878, 30879, 30880, 30881, 30882, 30883, 30884, 30885, 30886, 30887, 30888, 30889, 30890, 30891, 5179, 5180, 5181, 5182, 26463, 26462, 26461, 26460, 26459, 26458, 26457, 26456, 26480, 3341, 3342, 3339, 9509, 9510, 9511, 9512, 18487, 14972, 14973, 32692, 32693, 32694, 2300, 2301, 33715, 19978, 19979, 19977, 19976, 19975, 19974, 19973, 19972, 19971, 19970, 22082, 22080, 5351, 5352, 2280, 2281, 2282, 575, 576, 577, 578, 579, 536, 537, 538, 539, 540, 541, 15406, 15407, 15408, 553, 552, 633, 634, 635, 630, 631, 632, 623, 624, 625, 626, 471, 5363, 34169, 34170, 34171, 34785, 6466, 6465, 29938, 29939, 29940, 29941, 1446, 1447, 1448, 1449, 1450, 939, 940, 941, 9840, 9841, 26314, 26313, 26312, 26311, 26310, 26309, 26308, 26307, 26306, 26305, 26304, 26303, 26302, 26301, 26300, 34079, 34147, 29962, 26422, 26421, 26420, 765, 764, 763, 762, 761, 760, 736, 10973, 10974, 10975, 725, 10972, 727, 728, 26424, 733, 734, 29899, 20970, 20969, 20968, 20967, 20966, 20965, 20964, 20963, 20962, 9584, 9583, 20981, 21398, 20982, 34208, 35, 36, 59, 60, 61, 22363, 22364, 22365, 22366, 22367, 28958, 28959, 28960, 28961, 28962, 28805, 28806, 28807, 28808, 28809, 28810, 28827, 28828, 28829, 28830, 28874)
Because we asked it to also perform an assignment we can access the various results from that. The default return is a Pyarrow Table but Pandas is nicer for viewing.
res = rc.get_results().to_pandas()
res.head()
plot_results(rc.get_load_results()["demand"])
Batch operations#
To perform a batch operation we need to prepare the object first. We can either provide a list of tuple of the OD pairs we’d like to use, or we can provided a 1D list and the generation will be run on all permutations.
rc.prepare()
Now we can perform a batch computation with an assignment
rc.execute(perform_assignment=True)
res = rc.get_results().to_pandas()
res.head()
Since we provided a matrix initially we can also perform link loading based on our assignment results.
rc.get_load_results()
We can plot these as well
plot_results(rc.get_load_results()["demand"])
Select link analysis#
We can also enable select link analysis by providing the links and the directions that we are interested in. Here we set the select link to trigger when (7369, 1) and (20983, 1) is utilised in “sl1” and “sl2” when (7369, 1) is utilised.
rc.set_select_links({"sl1": [[(7369, 1), (20983, 1)]], "sl2": [[(7369, 1)]]})
rc.execute(perform_assignment=True)
We can get then the results in a Pandas DataFrame for both the network.
sl = rc.get_select_link_loading_results()
sl
We can also access the OD matrices for this link loading. These matrices are sparse and can be converted to SciPy sparse matrices for ease of use. They’re stored in a dictionary where the key is the matrix name concatenated with the select link set name via an underscore.
rc.get_select_link_od_matrix_results()
{'sl1': {'demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7f869e2a68c0>, '5x demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7f869e2a7400>}, 'sl2': {'demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7f869e2a4be0>, '5x demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7f869e2a4c40>}}
od_matrix = rc.get_select_link_od_matrix_results()["sl1"]["demand"]
od_matrix.to_scipy().toarray()
array([[0. , 0. , 0. , 3.04610785],
[0. , 0. , 0. , 0. ],
[0. , 0. , 0. , 0. ],
[0. , 0. , 0. , 0. ]])
project.close()
INFO:aequilibrae:Closed project on /tmp/7d56c6bae11e46c3a2538c265d416fd8
Total running time of the script: (0 minutes 4.134 seconds)