Where have all the post offices gone? Discovering neighborhood facilities with Python and OSM

Speaker 1 [00:07]

Thanks a lot for the nice introduction. What's up, PyData, PyCon DE? Great to see so many of you here, where today we're going to talk about where have all the post offices gone. We're going to do a recap of the movie The Fast and the Furious 2. Too fast, too furious. No, of course not. We're doing geospatial analysis with Python and OpenStreetMap, but Fast and Furious might as well be the subtitle of this talk because there's a lot packed in packed in, and we're going to get through it all. So I'm Katie. I'm a staff data scientist at Blue Yonder. Shout out Blue Yonder. Check out our booth to learn a little bit more about us. So I'm curious about the evolution of public spaces, how our shared spaces, much like our communities, grow and change over time. And I also enjoy writing letters. So the idea for this hobby project, It all started when I was writing a postcard and needed an international stamp, didn't have any at home, and didn't have a printer. So no worries, took a walk to my local post office, only to find out it was closed. Not for the day, for good. And this wasn't actually the first time this happened in my neighborhood. So being skeptical and suspicious, I thought, is this a systematic plot from the German postal system targeting my neighborhood? Or was this happening all over the state? When I opened up my map apps, I found only the post offices what are open around me now, not what existed six months ago or six years ago. So I need something like a database with historical information, also ideally where I can dive a little bit deeper into these places to understand what's offered, what services are offered by these facilities. So, OpenStreetMap to the rescue. I also want to do some data analysis, ideally in Python, so that's what I did. And I'm going to show you all how you can do it, too. So, first of all, what is OpenStreetMap? Has anybody here worked with OpenStreetMap before? Awesome! So many! Great to see it. Quite a lot of hands went up. So OSM is first and foremost a database. It provides open, free geographic data to everyone. OSM website also has an online map, a geodata search engine, and an editor. So community mappers like you and me can contribute to the OpenStreetMap project. It's also supported by the OpenStreetMap Foundation. This is a non-for-profit that supports the OpenStreetMap project. Okay. So what kind of data is in this database? So all the places on the map can be described in terms of three main elements. First, we have a node. A node is a single point in space, and importantly, it has longitude and latitude coordinates. Example of this would be something like a bench or a bus stop or a post office. So nodes have coordinates, and they're used to describe other elements like ways. So a way is a line or a linear feature, and it's essentially an ordered list of nodes. There's a couple of types of ways. You have an open way. This is where the first node does not equal the last node, like a road or a stream. Then you have a closed way, where the first node equals the last node. This could be something like a roundabout, a traffic circle, or a fence, a barrier around a park. And now the park itself, this is a different kind of way. This is an area. It's enclosed and it's filled. So something like a school ground could also be an area. There are some post offices that are also areas. And then finally we have relations. like their name implies, describes logical and geographical relations among a collection of other elements, so among other nodes, ways, or relations. This is where neighborhoods get classified. So your neighborhood might belong to a district, what belongs to a city, that belongs to a state, a country, a continent on the planet, you get the idea. So now we know how to figure out how do we find post offices and neighborhoods in the Database. Another important aspect in the Database are tags. So tags are key value free Format text fields that are used to describe features of Elements. So when we're talking about the Whole planet, this can be quite a lot. So currently there are over 100,000 distinct keys and over 172 million distinct tags. How do we keep all of this Straight? so conventions for tag usage Usage and their meanings these are documented in wiki pages Let's take a look at one of these so this is for the Boundary administrative tag and we can see that also for the Admin level tag where you have admin level 2 this corresponds To the country germany zooming a little bit in in this snippet Then we can see that the admin level 8 this is like a city 9 We have kind of our districts 11 even smaller we have the We have the neighborhoods. Awesome. Now we're even a huge step closer to finding neighborhoods in this data. What about post offices? How are they tagged? So for this, I'd like to introduce a few tools for diving into the data. Overpass API. Overpass API, this keeps an up-to-date copy of the main database with minute-level changes, and it provides this for search. The API is also optimized for data consumers, and Overpass Turbo is the front end for querying the Overpass API, and then it displays the results on a map. One of my favorite features of Overpass Turbo is the query wizard. So you can provide human-readable search text. Like here, I provided post office in Darmstadt, and it built this query for me. So in WR, we're looking for nodes, ways, or relations. I didn't specify that in my query. And then we have this tag, amenity equals post office. Awesome. This might be really helpful for us. Let's see what the results are. So this is the result, then, that you get in Overpass Turbo. All of the little circles represent nodes. You can click on them, dive deeper into the data, and also discover maybe other tags that it would be useful for us, like operator equals Deutsche Post. Awesome. So now we know how to get neighborhoods. We know how to get post offices. How can we get our hands on this sweet, sweet data? So full data of the planet is available from OpenStreetMap. That's a little bit more than I want. You can also download smaller data from Overpass API via Overpass Turbo as well. That's a little bit smaller than I want. I kind of want maybe my whole state. So for that, I'm going to use Geofabrik. So Geofabrik is a company that provides OSM services. Some of them are free services, like providing data extracts for subregions. Here, a subregion can be a continent, Europe, or a country, Germany. In Germany, the smallest extract or the smallest subregion, then this is broken up by state. So I'm interested in Darmstadt, so the smallest subregion I can download is Hessen. OK, cool. Also, I can get older snapshots of the database. So here, they provide historical snapshots from January 1 of the past years. So here's an example of how I would get this extract for Hessen using the latest data or from the database snapshot on January 1st, 2018. Great. But what's going on here? Dot OSM dot PBF. New format. Who dis? I don't know what this is. How do I work with these things? So if you're not familiar with this, PBF format. This is a format that's used for working with map data that's an alternative to XML. It's designed to be space efficient. In addition to that, it's a lot faster to read and write as compared to gzipped XML files. So cool, let's get started coding, right? Step back, hang on. What am I really trying to accomplish here? So I want to use this geospatial data to learn about the distribution of post offices or other resources in my neighborhood and how this distribution has changed over time. So we have really five main steps. We're going to download the data. We're going to filter it. We're going to then do some spatial operations. We've got neighborhoods. We've got post offices. Does post office, is this in neighborhood A or neighborhood B? Then we're going to enrich the data with maybe some extra information, population data if we've got it, calculate the areas, and then I want to plot it. Okay. Cool. Most importantly, I want a catchy name. So hot take, I think AI is really bad at coming up with clever names. So I used the subtitle of this talk to make an acronym and introducing Diner OSM. That's what we're going to be using. And I want to run this for a variety of regions and also a variety of versions, so different years in the database. So I created a command line utility. There's also some configurations. So Darmstadt, this points to the Hessen extract. So I put all these configurations in a TOML file. Also information about what kind of tags I'm going to be using for fetching areas or places. So all right. We got a plan. We got a structure. Let's get started. How do we work with these pbf files, right? So for this, I used pyosmium. Pyosmium, really fun to work with. Shout out pyosmium. This is a Python module that provides the Python bindings for libosmium or osmium, which is a C++ library for working with OSM data. And in just a few lines of code, we're going to be able to do all of these things. So I'm going to show you how we do this for neighborhoods. So first we instantiate a file processor to process the file that we just downloaded. Then we're going to use with areas to enable the area processing for closed ways and multipolygons. This also enables some location caching so we can build geometries. And then we're going to use different filters. And these filters are applied on objects in order. So only when an object passes through all of these filters, that's what we're going to keep in this file processor. So we're going to filter by areas, by these tags that we saw in the Wikidata page. Then we're going to apply the GeoInterface filter. This actually adds a GeoInterface attribute. GeoInterface is a Python protocol for working with GeoJSON-like spatial data. So cool. lines of code, we have gotten the neighborhoods. Post offices we'll do in much the same way, but here we're going to use with locations instead of with areas. This enables the location caching. And then we're going to filter by entity, nodes, key, name. So I don't really care about the tag for name, but the places just should have some name, ideally. Then also these tags, what we saw using Overpass Turbo, GeoInterfaceFilter, and what is this final filter, EnrichAttributes? This is a custom filter. So Pyosmium allows you to also create custom filters. So here, this one actually doesn't filter anything out. It enriches the GeoJSON properties with some attributes that I'd like to use later, specifically the original Original entity type and the osm id so that then we'll be able to construct the url In osm to then be able to see which object this points to look it up in the osm website. Okay, so Now we've got this file processor. What does it return? So it returns a python iterable Of these filtered osm objects and they all have the geo interface implemented Which means this is perfect input to another library, Geopandas. So geopandas extends panda Support to use geospatial data. And we can actually construct Data frames from this feature processor iterable using from Features. Awesome. Then we can do a spatial join. So we're going to join the areas On the nodes. We're going to do a left join. Keep all the areas, even if they don't have neighborhoods in them. And our predicate here contains means we're going to join them On areas that contain these neighborhoods. So kind of like a point and polygon Lookup is going on here. But, all right, spatial data, what Does this look like? so this is an example from Darmstadt nord. This is the area, and then these Red dots are the post offices, what are there. Then we can group by the unique areas and get the total number of post offices in areas. Oh, it seems like maybe I've lost connection here. Oh, she's back. All right. So then finally we're going to narrow the scope and calculate the area. In the Hamburg extract, for example, it also includes the Hamburg Vattenmeer. Not really interested in it. No offense, Vattenmeer. Clip this then by the bounding box just to get Hamburg City. To get Darmstadt, i'm going To use a different clip mask. So here i have a data frame of A higher admin level, so for cities. Then i'm going to filter that only for Darmstadt. Then my districts or neighborhoods data frame, then i Can use the Darmstadt clip mask to then only get the Neighborhoods that are in Darmstadt. We can also calculate areas in GeoPandas. We have to make sure that we set the appropriate coordinate reference system and then project it to one that outputs it in meters, and then we can calculate the square kilometers, the total square kilometers. So keeping in mind this also includes this isn't just land that we're doing here. All right. Awesome. Now we've got this beautiful data frame. We're finished, right? No. It's spatial data. Let's plot this. So to do this, I use one of my favorite tools for interactive visualizations, Bokeh. Bokeh provides the GeoJSON data source, so you can just create a GeoJSON by calling toJSON on our GeoDataFrame. This outputs a GeoJSON string. And if you haven't worked with Bokeh before, in its simplest form, there's four main steps. Prepare the data, call the figure, Add some type of renderers here using the patches, and then show. We can also make this a lot more complicated and interactive with visualizations, tools like tools, widgets. We can do some custom JavaScript callbacks. If you want to write JavaScript, you can, or you don't have to. Let's take a look at the final results. So here we have the post offices in Darmstadt. This is filtered just for amenity post office with some name. a few buttons. So show places. This will show the nodes on the map. A few other buttons here. So total or by area. This is going to change the data that's used for the color map for the areas. And it's normalized here. And then we have a version slider where we can change the snapshot version over time. Let's take a look at what this looks like. So we can see how more and more places are added to the database, these white nodes, how it changes over time. Also by area. Surprise, surprise. Darmstadt Mitte has the most there. Then we can zoom in and maybe see how things have changed. So here's a node. It was there in 2023, not in 2024. Click on it. Then we can go to the OSM website. Yep, it was deleted. Makes sense. We can also click through for relations. So this is why I added this custom filter to be able to use This later to take us to the website. What? Two names? That looks crazy. Data anomaly. No. Zoom in. Sure enough, we have two different places here. We can see which ones they point to. Awesome. So maybe post offices aren't your thing. So we can also use Dynor OSM to check out oak trees in Hamburg. Let's see how this has changed. How many trees have grown in the past year? No. Of course not. This is data that people are putting into OSM. And oh, my gosh, what happened in 2024? Oak trees taking over like Little Shop of Horrors. No. We had probably a big data dump then. So we can see one spec has the most total or by area or by population. That's Harburg. All right. Cool. Interesting. If you like oak trees. Or maybe you're a little bit more interested in something like surveillance in Frankfurt. So these are all places that have some value of a surveillance tag. And here I'm using a little bit smaller admin level. Let's take a look. All right. So I also haven't filtered this by any entity type, so you see a little bit different shapes, how this has changed over time. And zoom in a little bit, a straight line. That looks weird. What's going on here? Oh, it's a camera on the Aufwampf Highway. Hmm. No. Makes sense. What's that? Ah. Deutsche Bundesbank. Okay. Surveillance. Yeah. Checks out. And so, of course, then change the underlying data with these buttons. What's that? Central bank. Okay. Makes sense. Zoom in a little bit more. Street line. Oh, McDonald's. So get your McRib. You're on camera. Okay. So that's some things that you can do with exploring your neighborhood resources. But originally, the question of this talk, maybe you're wondering where did all the post Was my experience reflected in the data? Was i even asking the right question? No, not really. I wasn't. I actually needed a more complex set of filters than What i described here. If you're interested for those Details, we don't have the time, but maybe type it into Slido and we can talk about it then. So once i did apply the appropriate filters, then i do Have a little bit of a result. So these are all the districts DISTRICTS IN MY STATE AND THE NUMBER OF POST OFFICES, HOW THEY'VE CHANGED OVER TIME. THE PURPLE ONE, THAT'S MY DISTRICT, AND YOU CAN SEE, SURE ENOUGH, FROM 2024 TO 2025 AND TO NOW, THERE'S A FEW LESS TOTAL POST OFFICES, SO WE LOST TWO. NOT SO DRASTIC, RIGHT? AND ALSO, ESPECIALLY IN THIS BLUE DISTRICT, WE SEE THAT IT'S KIND OF LOST A LOT OF POSTAL facilities. If we look at this by area, wow, my community is over-served, right? What am I complaining about? It has kind of the smallest land area, so yeah, not really. But I think it's a little bit more appropriate to look at this by population. So here this is by 10,000. And so, okay, we do actually see for this blue district, it's kind of lost a lot of these access to resources. But then if we have a look in general at this trend, and especially where we are now, it seems that there's a more even distribution in the neighborhoods of our access to these facilities. So I think that's quite a positive message to leave on. So, yeah, I thank you a lot for your time and your attention, and I encourage you to get out there, go discover your neighborhood, and all the code you can find on my GitHub repo. Thanks a lot.

Speaker 2 [22:06]

Thank you so much katie for your interesting talk, and I also want to encourage you to ask your questions via slido Which is where i'm going to go to right now and read the first question from jan jan says thanks for the tutorial When analyzing places over time how do you make sure new places are actually new instead of an old place that? Was not yet there before

Speaker 1 [22:29]

So, excellent question. How can we make sure of this in open data especially? So, to really make sure that what's in OSM reflects the ground truth is kind of an impossible task, but there is a really important tag here. There are some tags that have the built year, so you can know, Okay, maybe it only showed up in the database in 2025, but it was Built in 2020. There are also tags that let you Know the last time this location was checked. So i think it's something like last checked on or last updated. I would have to check the wiki page to be sure about the Syntax there so that you can know, okay, maybe we have this Maybe we have this place in the database in 2025, but it was last checked in 2020. Maybe it's not so up to date. So there are some tags to check the reliability of the data.

Speaker 2 [23:28]

Thank you. And then the next question is, thank you for your interesting talk. Which resources or tutorials can you recommend to dive deeper into working with OpenStreetMap and Python?

Speaker 1 [23:39]

Well, of course, i can recommend my github. No, i'm just kidding. Yeah, so for diving into this, i really, for myself, what i can Recommend are the tools like overpass turbo and overpass api. You can also dive a little bit deeper into them and read their Api docs if you want to use overpass directly. But what i would really encourage you to do is just have An idea, get your hands dirty. Maybe you come from a data Science background and you want to and you like building Recommender systems, right? so then you could use osm to Design a recommender system to say, okay, these are some Restaurants that i like in berlin. I'm going to be traveling to Dresden. Let me build a recommender to recommend restaurants that I could try on my trip to Dresden. And then you're using this data, but around maybe a bigger topic like recommender system, something that you're interested in.

Speaker 2 [24:46]

Thank you. Then next question. How did you get population data and how did you map it to neighborhoods, relations or areas?

Speaker 1 [24:53]

Excellent question. I was hoping someone would ask this. So didn't have time for it in the talk. All right. So one of the tags that you can find for areas is the wikidata tag. Okay. So then that points you to the wikidata page. So what I did is then I actually used this for all of the areas what I did have this wikidata tag for. Then I queried wikidata for the population data and filled that in. Didn't look this up every time. Of course, okay, there are some Caveats to the data here. When is this population from? It's changing every day. So those are just some kind of Caveats and like a little asterisk to have in your mind How reliable is this data that you're using. But you can also do a lot more with these wiki data tags if you Want to add even enrich your data with even more qualities Qualities and features so yeah that's how i did it it's also on github if you want to check it out

Speaker 2 [25:57]

Thank you. And then the next question. Is there any indication how accurate the data on post offices or other points of interest is? Is this peer reviewed in the data?

Speaker 1 [26:06]

Okay. Yeah. So data accuracy and open data is a big question. That would actually be a talk in itself, how you can develop some tools to find out how reliable is this data. So there are some processes in OSM from community mappers and organizations to prevent having completely erroneous data or map vandalism in OSM to really cut down on that. But I think, for me, what I did here, one of the kind of tools that I use is I only use places that have this name tag, that it's not maybe just some random way or some random node that doesn't really have a post office, and also being really specific in the filters that you use, right? Because not every post office is equal. So a kiosk that sells stamps or maybe it's incorrectly labeled that it's a post box. So post boxes usually don't have names. So that was one way to cut down on it. But to have like a really systematic analysis, that's another talk in itself. And I'd be interested to go to that talk also.

Speaker 2 [27:32]

Thank you. And then the next question. Did you look into the spatial distribution of post offices, meaning how reachable is the closest post office for a random address in Darmstadt?

Speaker 1 [27:42]

Oh, that is super interesting. So short answer, no, I did not. And also when you think about routing, we're not thinking about aerial distance, too, right? So it might be really close in the aerial distance, but there's a river in between it or something. I don't know if there's rivers in Darmstadt. Sorry. I don't know. A lake, right, that you can't cross. So actually how you get there could be much longer, and this is an interesting routing problem that a lot of folks have worked on in open source projects using OpenStreetMap, but no, long answer to the short answer, no, I did not.

Speaker 2 [28:25]

Okay, thanks. And then a last question. Do you think that the old dumps have the same data quality as newer ones?

Speaker 1 [28:33]

Oh, so meaning the older places where we're imported into the database, do they have the same quality as this fresh data? So I think that's really hard to say. And I don't think it would really be possible to give kind of a yes or a no answer there. What I would tend to do is, especially for these older places, is to look at the tags that it was last checked on, and those are good indications of if the data is still fresh.

Speaker 2 [29:12]

Thank you so much, Katie.