Agent-based Transportaion Model¶

Let’s create a simple agent-based transportation model. First, start with a transportation network. Our transportation network can be created as a Networkx Graph or a (Geo)Pandas (Geo)DataFrame. There are some sample transportation networks available in dpd.mapping.samples. Real-life networks can be imported from OpenStreetMap via pyrosm.

[1]:

from networkx import draw

from dpd.mapping.samples import cross_box

graph = cross_box

pos = {
    node: (graph.nodes()[node]["geometry"].x, graph.nodes()[node]["geometry"].y)
    for node in graph.nodes
}

node_color = [
    (
        "red"
        if graph.nodes()[node].get("type") == "stop_sign"
        else (
            "yellow"
            if graph.nodes()[node].get("type") == "yield_sign"
            else (
                "blue"
                if graph.nodes()[node].get("type") == "stop_light"
                else "orange" if graph.nodes()[node].get("type") == "stop" else "green"
            )
        )
    )
    for node in graph.nodes()
]

draw(graph, pos=pos, node_color=node_color, with_labels=True)

../_images/notebooks_agent-based_transportation_model_2_0.png

[2]:

from geopandas import GeoDataFrame
from matplotlib import pyplot as plt

edges_df = GeoDataFrame(
    [graph.edges[edge] for edge in graph.edges], index=list(graph.edges)
)
nodes_df = GeoDataFrame([graph.nodes[node] for node in graph.nodes])

fig = plt.figure(figsize=(9, 8))
ax = fig.add_subplot(111)

nodes_df["geometry"].plot(ax=ax, color=node_color, markersize=1000)
nodes_df.apply(
    lambda x: ax.annotate(text=x.name, xy=x.geometry.coords[0], ha="center", size=20),
    axis=1,
)
edges_df["geometry"].plot(ax=ax)

[2]:

<Axes: >

../_images/notebooks_agent-based_transportation_model_3_1.png

networkx provides the ability to compute a path from any node to another node. When using OpenStreetMap, the same can be accomplished via the Open Source Routing Machine.

[3]:

from networkx import shortest_path

node_ids = shortest_path(graph, 0, 1)
print("Node IDs:", node_ids)

Node IDs: [0, 1]

To create a couple agents, we can create some transportation zones. For simplicity, we will create one zone per node with three Production and three Attraction so each zone has one person that goes to each other zone. Below are some other DataFrames we can generate from the Zones DataFrame. Also, we will create a path for each person using networkx.

[4]:

from dpd.modeling import Zones

zones = Zones(
    data=[
        {"Name": "Zone 0", "Production": 3, "Attraction": 3},
        {"Name": "Zone 1", "Production": 3, "Attraction": 3},
        {"Name": "Zone 2", "Production": 3, "Attraction": 3},
        {"Name": "Zone 3", "Production": 3, "Attraction": 3},
    ],
    index=nodes_df.index,
)
zones["geometry"] = nodes_df["geometry"]
zones

/tmp/ipykernel_2317/4093024280.py:12: FutureWarning: You are adding a column named 'geometry' to a GeoDataFrame constructed without an active geometry column. Currently, this automatically sets the active geometry column to 'geometry' but in the future that will no longer happen. Instead, either provide geometry to the GeoDataFrame constructor (GeoDataFrame(... geometry=GeoSeries()) or use `set_geometry('geometry')` to explicitly set the active geometry column.
  zones["geometry"] = nodes_df["geometry"]

[4]:

	Name	Production	Attraction	geometry
0	Zone 0	3	3	POINT (0.00000 0.00000)
1	Zone 1	3	3	POINT (0.00000 1.00000)
2	Zone 2	3	3	POINT (1.00000 0.00000)
3	Zone 3	3	3	POINT (1.00000 1.00000)

[5]:

distance_dataframe = zones.calculate_distance_dataframe()
distance_dataframe

[5]:

	0	1	2	3
0	0.000000	111195.080234	111195.080234	157249.598474
1	111195.080234	0.000000	157249.598474	111195.080234
2	111195.080234	157249.598474	0.000000	111178.144254
3	157249.598474	111195.080234	111178.144254	0.000000

[6]:

import numpy

from dpd.modeling import TripDataFrame

trip_dataframe = TripDataFrame(
    data=numpy.ones([4, 4]), index=zones.index, columns=zones.index
).map(int)

trip_dataframe

[6]:

	0	1	2	3
0	1	1	1	1
1	1	1	1	1
2	1	1	1	1
3	1	1	1	1

[7]:

from dpd.modeling import Population

population = Population.from_trip_dataframe(trip_dataframe)
population = population[
    population.origin != population.destination
]  # remove rows where the origin and destination zone are equal
population

[7]:

	origin	destination
1	0	1
2	0	2
3	0	3
4	1	0
6	1	2
7	1	3
8	2	0
9	2	1
11	2	3
12	3	0
13	3	1
14	3	2

[8]:

population["node_ids"] = population.apply(
    lambda x: shortest_path(graph, x["origin"], x["destination"]), axis=1
)
population

[8]:

	origin	destination	node_ids
1	0	1	[0, 1]
2	0	2	[0, 2]
3	0	3	[0, 3]
4	1	0	[1, 0]
6	1	2	[1, 2]
7	1	3	[1, 3]
8	2	0	[2, 0]
9	2	1	[2, 1]
11	2	3	[2, 3]
12	3	0	[3, 0]
13	3	1	[3, 1]
14	3	2	[3, 2]

Next, we setup and run our agent-based model. We need to tranform our transportation network in to Python objects for Edges and Nodes. Again, this can be done with either a Graph or a (Geo)DataFrame.

[9]:

from uuid import uuid4

from dpd.driving import EdgesLanesNodesDriver
from dpd.mapping import add_object_to_edges_and_nodes
from dpd.mapping.nodes import NodeModel
from dpd.mechanics import KinematicBodyWithAcceleration
from dpd.mechanics.datacollection import BODY_AGENT_REPORTERS
from dpd.modeling import TransportationModel

body_model = TransportationModel(
    agent_reporters=BODY_AGENT_REPORTERS | {"geometry": "geometry"}
)

node_model = NodeModel()

graph = add_object_to_edges_and_nodes(graph, node_model)

for index, row in population.iterrows():
    kbwas = KinematicBodyWithAcceleration(
        initial_acceleration=0.1,
        initial_velocity=0.1,
        initial_position=0,
        max_deceleration=0.1,
        min_velocity=0,
        unique_id=uuid4(),
        model=body_model,
    )
    elnd = EdgesLanesNodesDriver.from_node_ids(
        nodes_dict=graph.nodes,
        edges_dict=graph.edges,
        node_ids=row["node_ids"],
        body=kbwas,
        driver_final_velocity=0,
        unique_id=uuid4(),
        model=body_model,
    )
    body_model.schedule.add(elnd)

while body_model.running:
    body_model.step()
    node_model.step()

df = GeoDataFrame(body_model.get_dataframe())
df

[9]:

		position	geometry	time
Step	AgentID
0	05f18d17-4cc7-4689-bcc4-3a8a4010f09a	0.000000	POINT (0.00000 0.00000)	0
	f45a6d4d-74fc-4a13-b9ee-0b6569d5edfb	0.000000	POINT (0.00000 0.00000)	0
	5278ae06-7b75-4777-ac2d-71a8aa968713	0.000000	POINT (0.00000 0.00000)	0
	ac073fd6-6f11-4100-9dda-91e91dd83aeb	0.000000	POINT (0.00000 1.00000)	0
	c64d8a6b-8649-4779-8e74-44f04905db0a	0.000000	POINT (0.00000 1.00000)	0
	545df317-81b4-4bfd-b909-52df90c0e6a8	0.000000	POINT (0.00000 1.00000)	0
	31377bb5-0242-434f-8d3a-50a6e71dc217	0.000000	POINT (1.00000 0.00000)	0
	2300db3a-1ec9-49a4-963f-f2c2623273ad	0.000000	POINT (1.00000 0.00000)	0
	ac0c4b34-ae4d-4056-83da-fb936aa40487	0.000000	POINT (1.00000 0.00000)	0
	cfd50172-03ac-46c6-9a09-d913ea49932b	0.000000	POINT (1.00000 1.00000)	0
	a93f8bc1-4504-4e3f-8f9d-0c3b0de4f111	0.000000	POINT (1.00000 1.00000)	0
	b7304cff-d53b-455b-b068-8c8bbf83bcdc	0.000000	POINT (1.00000 1.00000)	0
1	05f18d17-4cc7-4689-bcc4-3a8a4010f09a	0.200000	POINT (0.00000 0.20000)	1
	f45a6d4d-74fc-4a13-b9ee-0b6569d5edfb	0.200000	POINT (0.20000 0.00000)	1
	5278ae06-7b75-4777-ac2d-71a8aa968713	0.200000	POINT (0.14142 0.14142)	1
	ac073fd6-6f11-4100-9dda-91e91dd83aeb	0.200000	POINT (0.00000 0.80000)	1
	c64d8a6b-8649-4779-8e74-44f04905db0a	0.200000	POINT (0.14142 0.85858)	1
	545df317-81b4-4bfd-b909-52df90c0e6a8	0.200000	POINT (0.20000 1.00000)	1
	31377bb5-0242-434f-8d3a-50a6e71dc217	0.200000	POINT (0.80000 0.00000)	1
	2300db3a-1ec9-49a4-963f-f2c2623273ad	0.200000	POINT (0.85858 0.14142)	1
	ac0c4b34-ae4d-4056-83da-fb936aa40487	0.200000	POINT (1.00000 0.20000)	1
	cfd50172-03ac-46c6-9a09-d913ea49932b	0.200000	POINT (0.85858 0.85858)	1
	a93f8bc1-4504-4e3f-8f9d-0c3b0de4f111	0.200000	POINT (0.80000 1.00000)	1
	b7304cff-d53b-455b-b068-8c8bbf83bcdc	0.200000	POINT (1.00000 0.80000)	1
2	05f18d17-4cc7-4689-bcc4-3a8a4010f09a	0.500000	POINT (0.00000 0.50000)	2
	f45a6d4d-74fc-4a13-b9ee-0b6569d5edfb	0.500000	POINT (0.50000 0.00000)	2
	5278ae06-7b75-4777-ac2d-71a8aa968713	0.500000	POINT (0.35355 0.35355)	2
	ac073fd6-6f11-4100-9dda-91e91dd83aeb	0.500000	POINT (0.00000 0.50000)	2
	c64d8a6b-8649-4779-8e74-44f04905db0a	0.500000	POINT (0.35355 0.64645)	2
	545df317-81b4-4bfd-b909-52df90c0e6a8	0.500000	POINT (0.50000 1.00000)	2
	31377bb5-0242-434f-8d3a-50a6e71dc217	0.500000	POINT (0.50000 0.00000)	2
	2300db3a-1ec9-49a4-963f-f2c2623273ad	0.500000	POINT (0.64645 0.35355)	2
	ac0c4b34-ae4d-4056-83da-fb936aa40487	0.500000	POINT (1.00000 0.50000)	2
	cfd50172-03ac-46c6-9a09-d913ea49932b	0.500000	POINT (0.64645 0.64645)	2
	a93f8bc1-4504-4e3f-8f9d-0c3b0de4f111	0.500000	POINT (0.50000 1.00000)	2
	b7304cff-d53b-455b-b068-8c8bbf83bcdc	0.500000	POINT (1.00000 0.50000)	2
3	05f18d17-4cc7-4689-bcc4-3a8a4010f09a	0.816228	POINT (0.00000 0.81623)	3
	f45a6d4d-74fc-4a13-b9ee-0b6569d5edfb	0.816228	POINT (0.81623 0.00000)	3
	5278ae06-7b75-4777-ac2d-71a8aa968713	0.900000	POINT (0.63640 0.63640)	3
	ac073fd6-6f11-4100-9dda-91e91dd83aeb	0.816228	POINT (0.00000 0.18377)	3
	c64d8a6b-8649-4779-8e74-44f04905db0a	0.900000	POINT (0.63640 0.36360)	3
	545df317-81b4-4bfd-b909-52df90c0e6a8	0.816228	POINT (0.81623 1.00000)	3
	31377bb5-0242-434f-8d3a-50a6e71dc217	0.816228	POINT (0.18377 0.00000)	3
	2300db3a-1ec9-49a4-963f-f2c2623273ad	0.900000	POINT (0.36360 0.63640)	3
	ac0c4b34-ae4d-4056-83da-fb936aa40487	0.816228	POINT (1.00000 0.81623)	3
	cfd50172-03ac-46c6-9a09-d913ea49932b	0.900000	POINT (0.36360 0.36360)	3
	a93f8bc1-4504-4e3f-8f9d-0c3b0de4f111	0.816228	POINT (0.18377 1.00000)	3
	b7304cff-d53b-455b-b068-8c8bbf83bcdc	0.816228	POINT (1.00000 0.18377)	3
4	5278ae06-7b75-4777-ac2d-71a8aa968713	1.220691	POINT (0.86316 0.86316)	4
	c64d8a6b-8649-4779-8e74-44f04905db0a	1.220691	POINT (0.86316 0.13684)	4
	2300db3a-1ec9-49a4-963f-f2c2623273ad	1.220691	POINT (0.13684 0.86316)	4
	cfd50172-03ac-46c6-9a09-d913ea49932b	1.220691	POINT (0.13684 0.13684)	4

[10]:

from datetime import datetime

from geopandas import GeoDataFrame
from movingpandas import TrajectoryCollection
from pandas import to_timedelta

START_TIME = datetime(1970, 1, 1, 0, 0, 0)
timedelta = to_timedelta(df.index.levels[0], unit="s")
index = START_TIME + timedelta
gdf = df
gdf.index = gdf.index.set_levels(index, level=0)
gdf.reset_index(level="AgentID", inplace=True)
tc = TrajectoryCollection(gdf, "AgentID")
tc.add_speed()
tc.hvplot(line_width=10)

/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/dpd/envs/latest/lib/python3.9/site-packages/movingpandas/trajectory.py:166: MissingCRSWarning: Trajectory generated without CRS. Computations will use Euclidean distances.
warnings.warn(

[10]:

[11]:

import random


def random_color():
    return (
        "#"
        + hex(random.randint(0, 0xFF))[2:]
        + hex(random.randint(0, 0xFF))[2:]
        + hex(random.randint(0, 0xFF))[2:]
    )


features = []
for trajectory in tc.trajectories:
    color = random_color()
    df = trajectory.df.copy()
    df["previous_geometry"] = df["geometry"].shift()
    df["time"] = df.index
    df["previous_time"] = df["time"].shift()
    for _, row in df.iloc[1:].iterrows():
        coordinates = [
            [row["previous_geometry"].xy[0][0], row["previous_geometry"].xy[1][0]],
            [row["geometry"].xy[0][0], row["geometry"].xy[1][0]],
        ]
        times = [row["previous_time"].isoformat(), row["time"].isoformat()]
        features.append(
            {
                "type": "Feature",
                "geometry": {
                    "type": "LineString",
                    "coordinates": coordinates,
                },
                "properties": {
                    "times": times,
                    "style": {
                        "color": color,
                        "weight": 5,
                    },
                },
            }
        )

import folium
from folium.plugins import TimestampedGeoJson

m = folium.Map(location=[0, 0], zoom_start=16)

TimestampedGeoJson(
    {
        "type": "FeatureCollection",
        "features": features,
    },
    period="PT1S",
    add_last_point=True,
    transition_time=1000,
).add_to(m)

m

[11]:

Make this Notebook Trusted to load map: File -> Trust Notebook

[ ]:

Agent-based Transportaion Model¶

dpd

Navigation

Related Topics