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.
2025-06-17 01:06:41,104;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
/home/runner/work/aequilibrae/aequilibrae/aequilibrae/paths/graph.py:248: UserWarning: Found centroids not present in the graph!
[ 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54
55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90
91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108
109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126
127 128 129 130 131 132 133]
warnings.warn("Found centroids not present in the graph!\n" + str(centroids[~present_centroids]))
2025-06-17 01:06:41,177;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
2025-06-17 01:06:41,265;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
2025-06-17 01:06:41,352;WARNING ; Field(s) speed, travel_time, capacity, osm_id, lanes has(ve) at least one NaN value. Check your computations
/home/runner/work/aequilibrae/aequilibrae/aequilibrae/paths/graph.py:248: UserWarning: Found centroids not present in the graph!
[ 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54
55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90
91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108
109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126
127 128 129 130 131 132 133]
warnings.warn("Found centroids not present in the graph!\n" + str(centroids[~present_centroids]))
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")
2025-06-17 01:06:41,416;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["demand_tot"] > 0]
max_load = link_loads["demand_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
return loaded_links.explore(
color="red",
style_kwds={
"style_function": lambda x: {
"weight": x["properties"]["demand_tot"] * factor,
}
},
)
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.
res = rc.get_results()
res.head()
plot_results(rc.get_load_results())
DEBUG:aequilibrae:Called __enter__ with self.__manual_transaction=False, self.__depth=0
DEBUG:aequilibrae:Called __exit__ with self.__manual_transaction=False, self.__depth=0
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()
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())
DEBUG:aequilibrae:Called __enter__ with self.__manual_transaction=False, self.__depth=0
DEBUG:aequilibrae:Called __exit__ with self.__manual_transaction=False, self.__depth=0
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 0x7fb2ed4ade40>, '5x demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7fb2ed4acf40>}, 'sl2': {'demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7fb2ed4af220>, '5x demand': <aequilibrae.matrix.sparse_matrix.COO object at 0x7fb2ed4afe20>}}
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()
WARNING:aequilibrae:This project at /tmp/bcb4624f1cef45f096ecbfd74b4d8983 is already closed
Total running time of the script: (0 minutes 3.288 seconds)