Exploring NYC Taxi Geospatial Data With DuckDB, H3, and KeplerGL
Success and Failures Working in New Spaces
The Initial Spark
I started exploring the NYC Taxi and Limo Commission (TLC) data earlier this month in a broader effort to familiarize myself with available public datasets for an interesting data story to tell. However, as I was poking around the possibilities, I was initially disappointed to find that the geospatial information was rounded up to rather large regions, limiting some of my initial plans.
Just as I was about to turn to other datasets, I stumbled upon Stephen Kent’s article about his journey exploring map creation with Python and SQL. What caught my eye was his building-level map of taxi pickups from the TLC data. This led me to discover that the pre-2010 TLC data contained actual latitude and longitude coordinates. While this is relatively outdated data, I saw it as an opportunity to explore some data I have not extensively worked with and learn some new tools along the way. Plus the original map was rather pretty and I wanted to replicate it with my own spin.
Goals of My Exploration
As I set off in my exploratory work I had a few goals in mind. I knew that in the short time I had (a quiet weekend), I was unlikely to find any deep new insight but I saw this as a great opportunity for learning.
- I aimed to improve my experience working with DuckDB and Geospatial data
- Practice working with a new GIS visualization tool. I explored lonboard, pydeck, and checked some findings in QGIS but ultimately settled on KeplerGL as it was a tool I had not previously worked with and seemed like it could handle the data set and provide the possibility of providing an interactive dashboard.
- Host the interactive dashboard that other users could use to explore the dataset.
What I Ended Up Using
With a short project, not everything goes as planned. In the end, the implementation combined several tools:
- DuckDB for efficient data processing
- H3 for geospatial indexing
- KeplerGL for interactive visualization
- Shapely for geometry handling
- GeoPandas to ease a pass into KeplerGL
Goal 1: Setup DuckDB as the Foundation
Data Processing Pipeline
At first, I thought I might look at the variation in the number of taxi pickups per building by day. I set up my aggregation like so.
SELECT
cell_id,
boundary,
MAX(CASE WHEN day_of_week = 'Monday' THEN trip_count ELSE 0 END) as cnt_monday,
MAX(CASE WHEN day_of_week = 'Tuesday' THEN trip_count ELSE 0 END) as cnt_tuesday,
MAX(CASE WHEN day_of_week = 'Wednesday' THEN trip_count ELSE 0 END) as cnt_wednesday,
MAX(CASE WHEN day_of_week = 'Thursday' THEN trip_count ELSE 0 END) as cnt_thursday,
MAX(CASE WHEN day_of_week = 'Friday' THEN trip_count ELSE 0 END) as cnt_friday,
MAX(CASE WHEN day_of_week = 'Saturday' THEN trip_count ELSE 0 END) as cnt_saturday,
MAX(CASE WHEN day_of_week = 'Sunday' THEN trip_count ELSE 0 END) as cnt_sunday
FROM (
SELECT
cell_id,
boundary,
CASE EXTRACT(DOW FROM CAST(pickup_datetime AS TIMESTAMP))
WHEN 0 THEN 'Sunday'
WHEN 1 THEN 'Monday'
WHEN 2 THEN 'Tuesday'
WHEN 3 THEN 'Wednesday'
WHEN 4 THEN 'Thursday'
WHEN 5 THEN 'Friday'
WHEN 6 THEN 'Saturday'
END as day_of_week,
COUNT(*) as trip_count
FROM taxi_data
GROUP BY 1,2,3
HAVING COUNT(*) >= {min_cnt}
)
GROUP BY 1,2
ORDER BY cell_idThen I joined the geometries with the Overture building dataset.
SELECT
b.geometry as geometry,
c.cell_id,
c.cnt_monday,
c.cnt_tuesday,
c.cnt_wednesday,
c.cnt_thursday,
c.cnt_friday,
c.cnt_saturday,
c.cnt_sunday
FROM buildings AS b
JOIN {temporal_table_name} AS c
ON ST_Intersects(ST_GeomFromText(b.geometry), ST_GeomFromText(c.boundary))I stored everything in a local DuckDB database for easy access for future use. I had also experimented with geopandas but ran into hardware constraints; at least on my hardware DuckDB was able to handle the parquet files much more efficiently.
Goal 2: Visualization with KeplerGL
My goal was to create an interactive visualization showing how pickup locations shifted over time. I was drawn to KeplerGL as it was touted as being able to handle a large amount of data well and I wanted to get hands-on with a new geospatial tool as previously my GIS work was mostly with ArcGIS. To implement the visualization layers I had to edit the basic configs to be able to handle both the spatial and temporal aspects of the data. My end config looked like this:
{
'version': 'v1',
'config': {
'visState': {
'filters': [],
'layers': [
{
'id': f'buildings-layer-{day}',
'type': 'geojson',
'config': {
'dataId': 'NYC Taxi Pickup by Day',
'label': f'Pickups - {day.capitalize()}',
'columns': {'geojson': 'geometry'},
'isVisible': day == 'monday', # Only Monday visible by default
'visConfig': {
'opacity': 0.8,
'filled': True,
'enable3d': False,
'colorRange': {
'name': 'Global Warming',
'type': 'sequential',
'category': 'Uber',
'colors': ['#FFC300', '#F1920E', '#E3611C', '#C70039', '#900C3F', '#5A1846']
}
}
},
'visualChannels': {
'colorField': {'name': f'cnt_{day}', 'type': 'integer'},
'colorScale': 'quantile'
}
} for day in ['monday', 'tuesday', 'wednesday', 'thursday', 'friday', 'saturday', 'sunday']
],
'interactionConfig': {
'tooltip': {
'fieldsToShow': {
'NYC Taxi Pickup by Day': [
{'name': 'cnt_monday', 'format': None},
{'name': 'cnt_tuesday', 'format': None},
{'name': 'cnt_wednesday', 'format': None},
{'name': 'cnt_thursday', 'format': None},
{'name': 'cnt_friday', 'format': None},
{'name': 'cnt_saturday', 'format': None},
{'name': 'cnt_sunday', 'format': None}
]
},
'enabled': True
}
}
},
'mapState': {
'latitude': 40.7128, # NYC center latitude
'longitude': -74.0060, # NYC center longitude
'zoom': 11,
'pitch': 0,
'bearing': 0
}
}
}Disappointing Initial Results
However, after examining the output it was evident that the only story this was telling was that when looking at count the map only showed that popular buildings stayed popular. We see some differences in the total number of pickups in workplace locations vs residential, but not in a significant way.
If I binned data slightly differently or used a metric such as difference from average I may be able to glean something useful, however first I wanted to examine the data by the hour because intuitively I felt at a minimum I should be able to see a difference in pickup location across hours in a day.
Slightly More Interesting Hour Aggregation
So I rewrote my aggregation to sum up the counts by hour over a month. Hoping to at minimum see some sort of flow in the taxi pickups as people commute to and from work.
SELECT
cell_id,
boundary,
MAX(CASE WHEN hour_of_day = 0 THEN trip_count ELSE 0 END) as cnt_0,
MAX(CASE WHEN hour_of_day = 1 THEN trip_count ELSE 0 END) as cnt_1,
MAX(CASE WHEN hour_of_day = 2 THEN trip_count ELSE 0 END) as cnt_2,
MAX(CASE WHEN hour_of_day = 3 THEN trip_count ELSE 0 END) as cnt_3,
MAX(CASE WHEN hour_of_day = 4 THEN trip_count ELSE 0 END) as cnt_4,
MAX(CASE WHEN hour_of_day = 5 THEN trip_count ELSE 0 END) as cnt_5,
MAX(CASE WHEN hour_of_day = 6 THEN trip_count ELSE 0 END) as cnt_6,
MAX(CASE WHEN hour_of_day = 7 THEN trip_count ELSE 0 END) as cnt_7,
MAX(CASE WHEN hour_of_day = 8 THEN trip_count ELSE 0 END) as cnt_8,
MAX(CASE WHEN hour_of_day = 9 THEN trip_count ELSE 0 END) as cnt_9,
MAX(CASE WHEN hour_of_day = 10 THEN trip_count ELSE 0 END) as cnt_10,
MAX(CASE WHEN hour_of_day = 11 THEN trip_count ELSE 0 END) as cnt_11,
MAX(CASE WHEN hour_of_day = 12 THEN trip_count ELSE 0 END) as cnt_12,
MAX(CASE WHEN hour_of_day = 13 THEN trip_count ELSE 0 END) as cnt_13,
MAX(CASE WHEN hour_of_day = 14 THEN trip_count ELSE 0 END) as cnt_14,
MAX(CASE WHEN hour_of_day = 15 THEN trip_count ELSE 0 END) as cnt_15,
MAX(CASE WHEN hour_of_day = 16 THEN trip_count ELSE 0 END) as cnt_16,
MAX(CASE WHEN hour_of_day = 17 THEN trip_count ELSE 0 END) as cnt_17,
MAX(CASE WHEN hour_of_day = 18 THEN trip_count ELSE 0 END) as cnt_18,
MAX(CASE WHEN hour_of_day = 19 THEN trip_count ELSE 0 END) as cnt_19,
MAX(CASE WHEN hour_of_day = 20 THEN trip_count ELSE 0 END) as cnt_20,
MAX(CASE WHEN hour_of_day = 21 THEN trip_count ELSE 0 END) as cnt_21,
MAX(CASE WHEN hour_of_day = 22 THEN trip_count ELSE 0 END) as cnt_22,
MAX(CASE WHEN hour_of_day = 23 THEN trip_count ELSE 0 END) as cnt_23
FROM (
SELECT
cell_id,
boundary,
EXTRACT(HOUR FROM CAST(pickup_datetime AS TIMESTAMP)) as hour_of_day,
COUNT(*) as trip_count
FROM taxi_data
GROUP BY 1,2,3
HAVING COUNT(*) >= {min_cnt}
)
GROUP BY 1,2
ORDER BY cell_idThis created a slightly more informative map allowing me to clearly identify the more residential areas of Brooklyn Heights compared to the corporate Financial District by the change in the density of taxi pickups. In the morning (Hour 7) note the deeper red and purple shades in Brooklyn Heights west shore as everyone leaves for work vs the lighter yellow shades in the afternoon (Hour 18) as everyone is departing from work.
Goal 3: Interactive Dashboard — It All Goes South
KeplerGL can easily export the data to a html for hosting on a dashboard or other site however it seems to have to store a great deal of information to effectively work, ending up with the user facing long load times. While the implementation is usable deeper explorations are limited compared to tools like QGIS.
Get Tibeggs’s stories in your inbox
Join Medium for free to get updates from this writer.
If you would like to look at a KeplerGL page with the hourly data populated I have a version currently up on my github pages here https://tibeggs.github.io/TLC_map_exploration/ (I make no promises of the efficiency of that page… still learning how to optimize with KeplerGL)
Things might be smoother if I imposed more limits on the data such as a tighter geographic region or not displaying as detailed building data. Another tool may be the answer or perhaps providing a server to pipe only the currently needed data to a front end, but I leave that open to future improvements.
That is All for Now
While my findings are certainly not groundbreaking, the work highlighted some important considerations for future analysis. I would have liked to do a deeper dive with a broader time scale rather than a month and spend more time analyzing the different neighborhoods within the dataset. Pure count is probably a poor choice for this analysis as highlighted by the daily aggregation which tends to just show what is a popular building to get picked up from rather than the individual impact of a a time frame. Ideally one could re-approach this with a different statistical measure, perhaps a difference from the daily or hourly average for a building.
Despite these limitations, this project proved to be a learning experience in several ways. It provided hands-on experience with modern geospatial tools like KeplerGL and H3, while also demonstrating the importance of iterative exploration in data analysis. The process of moving from day-of-week to hourly analysis showed how different temporal aggregations can reveal distinct patterns in the same underlying data. Furthermore, working with DuckDB for large-scale data processing highlighted the importance of efficient query optimization when dealing with geospatial datasets.
If you are interested in replicating or expanding on this you can see my full code at my GitHub repo linked below.
Technical Resources
For those interested in replicating or building upon this work:
- My Github Repo: https://github.com/tibeggs/TLC_map_exploration
- DuckDB Documentation: https://duckdb.org/docs/
- KeplerGL GitHub: https://github.com/keplergl/kepler.gl
- NYC TLC Data: https://www.nyc.gov/site/tlc/about/tlc-trip-record-data.page
- Kent’s Article: https://medium.com/@stephen.kent.data/how-i-got-started-making-maps-with-python-and-sql-5aaefcfb2b27
