All Things Being Iso — Mapbox Isochrone API with Python

Using Python and the Mapbox Isochrone API to Generate and Visualize Isochrone and Isodistance Polygons

Folium Map of Isodistances. Image by the Author.

Background

Buffering incident locations with straight line distance (i.e. “as the crow flies”) can be useful in geospatial analyses when studying continuous phenomena like temperature, and air quality. But when studying the movement of humans, it is best to measure distances along traversable networks like roads and walking paths. Isochrones use time as the measure from an origin point along these networks, while isodistances measure distance.

In the late aughts and early 2010’s, all a GIS analyst needed to generate some “network buffers” (i.e. service areas, isochrones, etc.) was ESRI ArcMap GIS software (yes, we used to capitalize ‘Esri’) and their Network Analyst extension. This commercial software extension came with access to ESRI Data and Maps (via DVDs), and contained a national, ready-built street network dataset that one could use to power the extension’s tools. Over the years, access to this street network went away, along with the company’s acronymic branding. The current generation version of the Network Analyst extension for Esri’s ArcGIS Pro now consumes credits (pay as you go currency) every time you use their cloud-based street network in the software’s tools. This is in addition, of course, to the costs associated with procurement and annual maintenance of the product.

Don’t get me wrong. I’m still an Esri dude. If I’m not directly using ArcGIS Pro or developing ArcGIS Online content, I’m certainly accessing their arcpy Python package in my day-to-day workflows. But as the industry continues to evolve, so has my desire to diversify and seek alternative methods to capturing, analyzing, and visualizing spatial data.

After exploring routingpy and OSMnx Python packages, which I couldn’t get to operate and execute the way that I had hoped, I decided to write code that accesses the Mapbox Isochrone API directly to perform the functions that I was interested in executing. To use the Mapbox Isochrone API, you must sign up for an account and agree to their pay-as-you go model. The good news is that the free tier includes up to 100,000 API request a month without incurring any expenses.

Code

As you may infer from my code, I am not a developer by training, but rather by necessity/desire. I’m certain my code includes flaws and I have not implemented any error handling, however, I was able to pull together something that operated the way that I wanted.

I created a Python class called ‘Isopleth’ which contains methods that generate isochrones and isodistances from input latitude and longitude coordinates, generate Folium maps of the GeoPandas results, and dissolve shapes when isopleths of the same value overlap.

The following code snippet, when executed in a Jupyter Notebook, generates a map for isodistances of 1, 3, and 10 miles from the input latitude and longitude coordinates for my employer’s main campus in Rockville, MD.

#import python packages
import geopandas as gpd
import pandas as pd
from py_code import isopleth as iso # isopleth is a .py file containing the Isopleth class
from py_code import cred # cred is a .py file that contains my Mapbox API credentials

#Create point dataframe
df = pd.DataFrame({"Address": ["1600 Research Blvd, Rockville, MD"],
                   "Latitude": [39.096892],
                   "Longitude": [-77.183063],})
gdf = gpd.GeoDataFrame(df, geometry=gpd.points_from_xy(df.Longitude, df.Latitude), crs="EPSG:4326")

# distance in miles, up to four values
distance = [1, 3, 10]
# colors used in Folium, must be same len as distance/minutes
colors = ['#984ea3', '#e41a1c', '#ff7f00']

# Creates an instance of the Isopleth class
iso_map = iso.Isopleth(cred.MY_API_KEY)

# Uses the add_isodist and get_isodist methods to return a gpd of isodistances
iso_map.add_isodist(gdf, id_col="Address", distance_miles=distance, colors=colors)
isodists = iso_map.get_isodists()

# Map isolines with the map_iso method
iso_map.map_iso(isodists)

Folium Map of Isodistances. Image by the Author.

While the above example shows isodistances around a single location, multiple locations can be processed with the Isopleth class. The Isochrone API is limited to 300 requests per minute, so prepare accordingly.

Another method in the Isopleth class is a function that uses the Folium DualMap plugin to show side by side results, comparing isodistances to isochrones. Below I generated an isochrome map, covering 5, 15, and 20 minute drivetimes to display next to the isodistance map generated earlier.

# Minutes for isochrone
minutes = [5, 15, 20]

#Isochrone colors
c_colors = ['#ffff33','#4daf4a', '#8dd3c7']

# Use same Isopleth class instance to generate and return isochrones gpd
iso_map.add_isochrone(gdf, id_col="Address", minute_list=minutes, colors=c_colors, mode="driving")
isochrones = iso_map.get_isochrones()

# Map isodistance and isochrones side by side with the map_dual_iso method
iso_map.map_dual_iso(isodists, isochrones)

Isodistance on Left and Isochrones on Right Using Folium DualMap Plugin. Image by the Author.

I find this to be very useful, especially in the exploratory stages of analysis.

The ability to dissolve isopleths when multiple origin points are processed was also built into the class. For example, when two point locations are input and they are geographically close to each other, the process will result in overlapping service areas. There may be cases where you want to keep these shapes separate (if performing spatial analysis on the individual locations), but for visualizing them, it may make more sense to dissolve the shapes based on isopleth values. In the below code, two point locations are processed.

#Create point dataframe with two points
df2 = pd.DataFrame({"Address": ["1600 Research Blvd, Rockville, MD", "50 Southlawn Ct, Rockville, MD 20850" ],
                   "Latitude": [39.096892, 39.106171],
                   "Longitude": [-77.183063, -77.12947],})
gdf = gpd.GeoDataFrame(df, geometry=gpd.points_from_xy(df.Longitude, df.Latitude), crs="EPSG:4326")

# distance in miles, up to four values
distance2 = [1, 3, 5]
# colors used in Folium, must be same len as distance/minutes
colors2 = ['#000000', '#6a51a3', '#feb24c']

# Creates a new instance of the Isopleth class
iso_map2 = iso.Isopleth(cred.MY_API_KEY)

# Uses the add_isodist and get_isodist methods to return a gpd of isodistances
iso_map2.add_isodist(gdf2, id_col="Address", distance_miles=distance2, colors=colors2)
isodists2 = iso_map2.get_isodists()

# Dissolves the isodistance poly by distance values
dissolved_isodist2 = iso_map2.get_dissolve_isodists()

# Map isodistance and idissolved isodistance side by side with the map_dual_iso method
iso_map2.map_dual_iso(isodists2, dissolved_isodist2)

The map produced shows the non-dissolved, overlapping isodistances on the left, and the dissolved version on the right.

Isodistance on Left and the Dissolved Isodistance Shapes on Right. Image by the Author.

Summary

Many different methods, programming languages, and data sources can be used to generate and visualize isochrones and isodistances. If you would like to give my approach a whirl, access to the Jupyter Notebook used in this post and the supporting .py code containing the Isopleth class can be found as a directory in my public GitHub repository (https://github.com/mikegwiz/public_repo/tree/main/mapbox_iso). You will, however, need your own Mapbox API key for the code to run.

If you’d like more information about isochrones and isodistances, including their uses and limitations, check out the Mapbox post on the subject.

What are Isochrones? - Mapbox

Get mapping and thanks for reading!