Travel-time isochrones with Mapbox, Python, and GeoPandas

Travel-time isochrones are powerful analytical tools that represent the reachable area from a location for a given time and travel mode. In R, my package mapboxapi seamlessly integrates with R’s GIS infrastructure to allow for the use of Mapbox’s isochrones in spatial analysis workflows. In Python, there isn’t a package that directly connects Mapbox’s navigation toolkit to GeoPandas for spatial data analysis. However, these services are accessible via the routingpy Python package.

In this blog post, I’ll present a workflow to help you connect GeoPandas with Mapbox’s isochrone services via routingpy. You’ll be able to use GeoPandas POINT geometries as inputs and get back isochrone polygons as GeoDataFrames. The goal is to replicate some of the functionality of R’s mb_isochrone() function in Python.

To get started, let’s import a few libraries we’ll need, and make a connection to Mapbox’s navigation services which are named MapboxOSRM in routingpy. You’ll need a Mapbox account and a Mapbox access token to get this to work; you get 100,000 isochrones for free each month, so you shouldn’t have to worry about getting charged.

import geopandas as gp
import pandas as pd
from shapely.geometry import Polygon
from routingpy.routers import MapboxOSRM
import numpy as np

mb = MapboxOSRM(api_key = "YOUR KEY GOES HERE")

Let’s use a public libraries dataset in the city of Dallas, Texas as an example. Mapbox’s routing services will run for any location in the world covered by OpenStreetMap, so you can try this out for other datasets of interest as well.

dallas_libraries = gp.read_file("https://egis.dallascityhall.com/resources/Downloads/ShpZip/Library/Libraries.zip")

dallas_libraries.explore()

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

Next comes the mb_isochrone() function. I’ve written this to work in a similar way to mb_ischrone() in R, though it is much more limited. Read through the comments in the code to get a sense of how it works.

def mb_isochrone(gdf, time = [5, 10, 15], profile = "driving"):

    # Grab X and Y values in 4326
    gdf['LON_VALUE'] = gdf.to_crs(4326).geometry.x
    gdf['LAT_VALUE'] = gdf.to_crs(4326).geometry.y

    coordinates = gdf[['LON_VALUE', 'LAT_VALUE']].values.tolist()

    # Build a list of shapes
    isochrone_shapes = []

    if type(time) is not list:
        time = [time]

    # Use minutes as input, but the API requires seconds
    time_seconds = [60 * x for x in time]

    # Given the way that routingpy works, we need to iterate through the list of 
    # coordinate pairs, then iterate through the object returned and extract the 
    # isochrone geometries.  
    for c in coordinates:
        iso_request = mb.isochrones(locations = c, profile = profile,
                                    intervals = time_seconds, polygons = "true")

        for i in iso_request:
            iso_geom = Polygon(i.geometry[0])
            isochrone_shapes.append(iso_geom)

    # Here, we re-build the dataset but with isochrone geometries
    df_values = gdf.drop(columns = ['geometry', 'LON_VALUE', 'LAT_VALUE'])

    time_col = time * len(df_values)

    # We'll need to repeat the dataframe to account for multiple time intervals
    df_values_rep = pd.DataFrame(np.repeat(df_values.values, len(time_seconds), axis = 0))
    df_values_rep.columns = df_values.columns

    isochrone_gdf = gp.GeoDataFrame(
        data = df_values_rep,
        geometry = isochrone_shapes,
        crs = 4326
    )

    isochrone_gdf['time'] = time_col

    # We are sorting the dataframe in descending order of time to improve visualization
    # (the smallest isochrones should go on top, which means they are plotted last)
    isochrone_gdf = isochrone_gdf.sort_values('time', ascending = False)

    return(isochrone_gdf)

Let’s try it out! The function runs seamlessly over a dataset of 29 input points, returning 29 5-minute isochrones that we can visualize with .explore().

library_isos = mb_isochrone(dallas_libraries, time = 5, 
                            profile = "driving-traffic")

library_isos.explore()

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

One feature the R version of mb_isochrone() includes is the ability to get isochrones directly from an address. Here, we’ll need to geocode the address first and pass it to mb_isochrone(). Our result is an interactive map of multiple travel-times around our input address.

dickies = gp.tools.geocode("1911 Montgomery St, Fort Worth, TX 76107")

dickies_isos = mb_isochrone(dickies, time = [5, 10, 15], profile = "driving-traffic")

dickies_isos.explore(column = "time")

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

A great benefit of using isochrones with GeoPandas is that the isochrone shapes can be integrated into all sorts of spatial analysis workflows. If you are interested in integrating isochrones into your application or business workflow, or if you’d like a custom workshop to help you get up and running with these tools, please don’t hesistate to reach out to kyle@walker-data.com!