New feature: Retrieving OSM POIs with osmnx?

[HTenkanen]

Hi,

first of all, awesome work with the package! Highly useful and needed.

I was wondering if there would be interest to include separate functionalities to retrieve OpenStreetMap POIs with osmnx? I am currently starting to write such functionalities for my own project, and I think I will be using osmnx as a basis for fetching Points of interest from OpenStreetMap (such as restaurants etc.). I know that it is possible to work with other infrastructures than roads, but if I have understood correctly that function works only with graphs, not with point data such as POIs.

If you think this would be interesting, I would be happy to help implementing these functions.

Cheers,
Henrikki

P.S. FYI: I introduce osmnx also with the fully open GIS course (Lesson 7) that can be found here http://autogis.github.io.

[gboeing]

@HTenkanen your course looks really cool. Nice work. Yes I think it might be interesting to work with POIs directly in OSMnx... perhaps with a separate new POI module.

[HTenkanen]

@gboeing thanks! 👍

Nice, I think having POIs included in OSMnx would be a nice addition! And yes, I think a separate module would be good, and it could be written e.g. into pois.py in a similar manner as you have now the buildings.py. Or what would you suggest?

Basically the idea could be that you can pass a list of amenities to a function pois_from_poly() that you want to retrieve from Overpass API, and then it would return a GeoDataFrame out of those amenities. This would follow the style that have, and there would be separate functions for pois_from_bbox(), pois_from_place(), pois_from_address() etc.

I actually already wrote a first version of module with this style (only pois_from_poly() now) for my own purposes, and it works quite nicely. There are, however, still some things that needs a bit of thinking:

Should the module provide information about what amenities are available?
- e.g. ox.pois_available could return a list of possible amenity values
- This would mean that the amenity list should be updated every now and then...

The function searches now amenities from nodes, ways and relations:

query_template = ('[out:json][timeout:{timeout}]{maxsize};((node["amenity"~"{amenities}"]({south:.8f},'
                  '{west:.8f},{north:.8f},{east:.8f});(._;>;););(way["amenity"~"{amenities}"]({south:.8f},'
                  '{west:.8f},{north:.8f},{east:.8f});(._;>;););(relation["amenity"~"{amenities}"]'
                  '({south:.8f},{west:.8f},{north:.8f},{east:.8f});(._;>;);));out;')

How should we deal with the amenities that are of type way? E.g.

{'type': 'way', 'id': 531711920, 'nodes': [5161070021, 5161088574, 5161088573, 5161088572, 5161070022, 5161070023, 5161070024, 5161070021], 'tags': {'amenity': 'kindergarten', 'operator': 'Boställs Daghem'}}

Should the functions return two GeoDataFrames where first one are the points and the second one contains edges? There are for example a lot of parking area features and schools that seem to be ways. I think that might work.

And then the last thing, how should the functions deal with relations?

{'type': 'relation', 'id': 7158483, 'members': [{'type': 'way', 'ref': 29875899, 'role': 'outer'}, {'type': 'way', 'ref': 486843935, 'role': 'inner'}], 'tags': {'addr:housenumber': '4', 'addr:street': 'Talvelantie', 'amenity': 'hospital', 'building': 'public', 'name': 'Malmin terveysasema ja hammashoitola', 'type': 'multipolygon'}} 

These can be fairly tricky to handle. Any ideas?

Okey, so these are some of the thoughts that I have now after testing a little bit what kind of amenities are returned by the Overpass. I would really appreciate all the thoughts and ideas from you about this, as you know the best how you want to develop the package.

I need to be cleaning my code a little bit and implement some of the missing features, but after that I will push my changes to a fork I made where you can take a look https://github.com/HTenkanen/osmnx.

Cheers,
Henrikki

[HTenkanen]

@gboeing

I now wrote a pois.py which is a first version of the POI module that works in a similar manner as the other modules. At the moment it works only for node tags and returns a GeoDataFrame of POIs with points of interests. I wrote it into this branch:
https://github.com/HTenkanen/osmnx/tree/1-osm-poi-dev

Please take a look, and try it out. E.g. following should work:

>>> import osmnx as ox
>>> place = "Kamppi, Helsinki, Finland"
>>> ox.available_amenities
>>> selected_amenities = ['restaurant', 'bar', 'cafe']
>>> poi_gdf = ox.pois_from_place(place=place, amenities=selected_amenities)
>>> poi_gdf['name'].head()
60062502                Kabuki
60068035             Cafe Java
60133792          Ateljé Finne
62965963          Empire Plaza
62967659    Ravintola Pääposti
Name: name, dtype: object

Cheers,
Henrikki

[gboeing]

@HTenkanen this is looking like a great start! A few ideas regarding your questions:

I'm not sure if the module needs to provide information about which amenities are available, as this would create a maintenance task on the OSMnx side... is it possible to just provide a link in the docstring to an authoritative documentation source from OSM?

Regarding amenities of type "way"... we could return a GeoDataFrame with both node-based and way-based amenities as individual rows. The former would have point geometries, and the latter would have polygon or multipolygon geometries. Would that work? If not, I'd be inclined to take the centroid of the way-based amenity to convert it to a point.

Regarding relations, do you have any examples? What kinds of amenities are these, and must they be handled as relations?

[HTenkanen]

@gboeing Thanks for your thoughts!

I now implemented your ideas and now the poi-functions return a single GeoDataFrame consisting of node, way and relation elements. I added a column called element_type that tells which kind of element a specified row is. Possible values for that column are: 'node', 'way', and 'relation'. In this way, it is possible to easily select only such POIs that are points (i.e. nodes), or polygons (ways), or multipolygons (relations).

Regarding your question about relations, there are quite many cases where for example a school consists of many Polygons and hence those are dealt with a relation. And it is typically the case that the individual 'way' elements do not contain any attribute information about the school (like name, address, etc) as those attributes are assigned to the 'relation' element. Hence, in these cases it is necessary to handle the relation and create a MultiPolygon object having all the attributes. I have now wrote a function called parse_osm_relations() which does exactly that.

I also modularized the code a bit, and I think there are few functions that could be used directly e.g. in buildings.py module so that there wouldn't be repetitive code in different modules. But this would maybe need another module for those functions that are shared commonly between different modules (I was thinking to put them in utils.py but I'm not sure if that would be the place to put them or to make a separate module for such functions that are related to handling of geometries and parsing osm responses).

Overall, I think the poi-module is now almost there, but it still needs some tests, documentation and a look from you.
You can check it out from https://github.com/HTenkanen/osmnx/tree/1-osm-poi-dev.
If you have any suggestions or ideas, I'm all ears. 👂 😊

Now e.g. following should work:

>>> import osmnx as ox
>>> place = "Manhattan Island, Manhattan, New York City, New York, United States of America"

>>> all_pois = ox.pois_from_place(place=place)
>>> len(all_pois)
13182
>>> all_pois['amenity'].value_counts()
bicycle_parking          3443
restaurant               2119
parking                   768
cafe                      762
school                    611
place_of_worship          595
fast_food                 529
bar                       385
bank                      345
bicycle_rental            342
drinking_water            334
bench                     308
post_box                  258
embassy                   212
pharmacy                  202
pub                       147
theatre                   132
...
# Get restaurants from Manhattan, New York
>>> restaurants = ox.pois_from_place(place=place, amenities=['restaurant'])
>>> restaurants['element_type'].value_counts()
node        1941
way          179
relation       2
Name: element_type, dtype: int64

# Get schools from Manhattan, New York
>>> schools = os.pois_from_place(place=place, amenities=['school']
>>> schools['name'].head()
42473959                                                   NaN
42883022                                                   NaN
357544551                                   All Hallows School
357544577                               All Saints High School
357548361    High School for Health Profesion and Human Ser...
Name: name, dtype: object

>>> schools.plot(markersize=2)

# Get all Universities from Boston area (with 2 km buffer to cover also Cambridge)
>>> boston_q = "Boston, Massachusetts, United States of America"
>>> boston_poly = ox.gdf_from_place(boston_q, buffer_dist=2000)
>>> universities = ox.pois_from_polygon(boston_poly.geometry.values[0], amenities=['university'])
>>> list(universities.name.unique())
[nan, 'Picower Institute for Learning and Memory (MIT)', 'McGovern Institute for Brain Research (MIT)', 'Boston ollege', 'Boston College (Newton Campus)', 'Harvard Business School', 'John F Kennedy School of overnment', 'Radcliffe Institute for Advanced Studies', 'Soldiers Field Park Apartments', 'MIT Museum', 'Winthrop House', '250 Mt Vernon St @ Geiger Gibson Health Center', 'Boston University School of edicine', 'Citibank', 'Cambridge Innovation Center', 'Harvard Law School', 'Rubinstein Hall', 'Boston niversity', 'Suffolk University', 'Conant Hall', 'Longfellow Building', 'Putnam Building', 'Quincy College', 'Engineering Science Lab', 'Lowell Lecture Hall', 'Littauer Center', 'Memorial Hall', 'University Hall', 'Sever all', 'Boylston Hall', 'Myles Standish Hall', 'Stone Hearth', 'Swiss Bakers', 'Benjamin Franklin Institute of echnology', 'Randolph Hall', 'Lehman Dudley House', '7/11', 'Byerly Hall', 'Divinity Hall', 'Smith Campus enter', 'CGIS South', 'Tufts University Hirsh Health Sciences Library', 'Mather House', '38 Oxford', 'Farlow', 'Radcliffe University', 'CGIS Knafel', 'Phillips Brooks House', 'Holden Chapel', 'Vanserg Building', 'Biological aboratory', 'Larsen Building', 'Buckingham House', 'Gutman Building', 'Claverly Hall', 'Knowles Center', 'Tufts University Dowling Hall', 'Blackstone South', 'Harvard University Information Systems', 'Pforzheimer ouse', 'Emerson Hall', 'Currier House', 'Tufts Dental School', 'Mount Ida College', 'Teele Hall', 'Harvard Ed ortal', 'Hastings Building', 'Robinson Hall', 'University Herbaria', 'Harvard Hall', 'White Hall', 'Eliot House', 'Harvard University', 'Radcliffe Institute at Harvard University', 'Radcliffe Gym', 'Agassiz House', 'Boston ollege (Brighton Campus)', 'Boston College - Newton Campus', 'Harvard University Northwest Building', 'University of Massachusetts Boston', 'Harvard Medical School', 'Harvard School of Public Health', Harvard chool of Dental Medicine', '25 Travis Street', 'Northeastern University', 'Boston University Medical Campus', 'Leverett House', 'Quincy House', 'Cabot House', 'Dunster House', 'Lowell House', 'Harvard Science Center', 'Massachusetts Institute of Technology', 'Westmorely Hall', 'Kirkland House', 'Tufts University']
>>> universities.name.value_counts().head()
Suffolk University                    8
Harvard University                    7
Boston University                     7
Boston University Medical Campus      3
University of Massachusetts Boston    3
Name: name, dtype: int64
>>> universities.plot(alpha=0.7)

[gboeing]

Looks good! Yes I think converting relations like schools -> multipolygons is a good solution. For now I think putting shared code in utils.py is a good idea.

[gboeing]

@HTenkanen just checking back in on this

[HTenkanen]

@gboeing Hi! Sorry for silence, have been busy with funding proposals and teaching last couple of weeks. I try to finish this up and make a pull request in a week or so. I leave the main readme.md documentation for you as I don't know how you would like to present this stuff. How does that sound?

By the way having a generic filter that you are discussing in #90 sounds really good and useful, and it was actually wished from few persons who I taught these stuff today. So that indeed would be a nice feature. E.g. related to buildings, being able to fetch buildings from an area that were built in 1954, or between years 2000-2015 are examples of such filtering that might be highly useful, and should be possible to do with OverpassAPI.

[gboeing]

Sounds good. I can edit the readme file as needed, but please add documentation via function docstrings so that it appears automatically at https://osmnx.readthedocs.io/en/latest/

[gboeing]

@HTenkanen just checking back if this is still an open item.

[BlackArbsCEO]

Great project here, has the poi functionality been added yet or that still in testing?

[HTenkanen]

@gboeing @BlackArbsCEO:

Sorry for leaving this hanging for quite a long time. Just made a pull request for the added features.

[TRomijn]

Great work guys! I've tried using the code, but it doesn't work for me. I'm not sure whether I've made a mistake or maybe something on osm's side has updated.

When I try to run this code:

import osmnx as ox 
gdf = ox.pois_from_place(place='Kamppi, Helsinki, Finland')

I'm receiving these errors:

JSONDecodeError: Expecting value: line 1 column 1 (char 0)

and

Exception: Server returned no JSON data.

Am I doing something wrong, or is something wrong in the code?

[gboeing]

@TRomijn the server started enforcing changed syntax. Should be fixed in 1413217.

[claytongroth]

Hello guys,

Amazing stuff! I've been loving these modules. I am getting stuck with a particular part of the OSM POIs module. I was wondering if you had any suggestions.

ox.pois_from_point(point1, distance = 500)

For some reason whenever I try to use a distance parameter higher than 500 I get the following error.

raise ValueError("No Shapely geometry can be created from null value") ValueError: No Shapely geometry can be created from null value

The error occurs at this line.
POIshapes = all_pois1['geometry'].representative_point()

For my current project, I really would like to run it at 1000.

I can post more code if necessary.

Thanks!

[gboeing]

@claytongroth can you provide a complete minimal working example so I can reproduce it?

[claytongroth]

@gboeing My apologies for the long response. I was away for the weekend.
Here is the minimal working example required. It throws the error only for some (lat,lon)s and for some only when the radius for consideration is greater than 500.

Also, I apologize if this is a bit more than minimal...

import osmnx as ox
import networkx as nx
import sys
import pandas as pd
import matplotlib.cm as cm
def getXY(pt):
    return (pt.y, pt.x)

lat1 = (43.801357)
lon1 =  (-91.239578)
point1 = (lat1,lon1)


G = ox.core.graph_from_point(point1, distance = 1000, network_type='walk')
all_pois1 = ox.pois_from_point(point1, distance = 1000)
stats = ox.basic_stats(G, area=500)


POIshapes =  all_pois1['geometry'].representative_point()
POICoordslist = map(getXY, POIshapes)
print POICoordslist
POIdistances = []
for coordset in POICoordslist:
    orig_node = ox.get_nearest_node(G, (point1))
    dest_node = ox.get_nearest_node(G, (coordset))
    Dpoi = nx.shortest_path_length(G, orig_node, dest_node, weight='length')
    POIdistances.append(Dpoi)
print POIdistances

• works for 40.758896,-73.985130 and 43.801357, -91.239578 and 35.2000506, -105.1422244
• throws error for 43.073051, -89.401230
• throws errors for SOME (lat ,lon) at 1000m but not at 500m

[HTenkanen]

Hi @claytongroth! Can you still share us the error you receive?

Could you also check what version of Pandas you are using? I at least cannot reproduce the error.
Newest version of pandas is required (0.23.1).

import pandas as pd
pd.__version__

[gboeing]

@claytongroth I don't see any OSMnx error. I had to change your code snippet to turn it into a working example that reproduces the error you mentioned (this took some guesswork as your snippet works fine out of the box).

MWE:

import osmnx as ox
point = (43.073051, -89.401230)
dist = 1000
G = ox.graph_from_point(point, distance=dist, network_type='walk')
pois = ox.pois_from_point(point, distance=dist)
POIshapes = pois['geometry'].representative_point()

This throws a geopandas ValueError as you're calling geopandas.GeoSeries.representative_point() on a DataFrame of mixed types (Point, Polygon, MultiPolygon). The "within" logic of representative_point doesn't make sense on a Point. If you call .centroid instead, this snippet works fine.

[claytongroth]

@gboeing, @HTenkanen
Thank you guys for your help!
My pandas version is 0.23.1

When I run @gboeing's above code with .centroid I get :

Traceback (most recent call last): File "C:\Users\clayt\Desktop\PYProject\testShapely.py", line 8, in <module> POIshapes = pois['geometry'].centroid() TypeError: 'GeoSeries' object is not callable

When I run the above code with .representative_point I get:

No handlers could be found for logger "shapely.geos" Traceback (most recent call last):

File "C:\Users\clayt\Desktop\PYProject\testShapely.py", line 8, in <module> POIshapes = pois['geometry'].representative_point()

File "C:\Users\clayt\Anaconda2\lib\site-packages\geopandas\base.py", line 180, in representative_point for geom in self.geometry],

File "C:\Users\clayt\Anaconda2\lib\site-packages\shapely\impl.py", line 37, in wrapper return func(*args, **kwargs)

File "C:\Users\clayt\Anaconda2\lib\site-packages\shapely\geometry\base.py", line 476, in representative_point return geom_factory(self.impl['representative_point'](self))

File "C:\Users\clayt\Anaconda2\lib\site-packages\shapely\geometry\base.py", line 76, in geom_factory raise ValueError("No Shapely geometry can be created from null value")

ValueError: No Shapely geometry can be created from null value

[gboeing]

@claytongroth please see the geopandas documentation. centroid is an attribute, not a method (it should not have parentheses), and representative_point doesn't make sense to call on a GDF that contains Point type geometries.

[claytongroth]

@gboeing I just read it and was about to delete my post haha...

Thank you for your help and patience! I'm all set!

[chesterharvey]

@gboeing What do you think about generalizing the amenities parameter so you can retrieve features with other tags, too? I wrote a generalization in this branch that switches out amenities for tags, which takes a dictionary to query specific key:value pairs.

https://github.com/chesterharvey/osmnx/blob/master/osmnx/pois.py

For example, let's say you want to get all features with an amenity tag (specified with a 'None' value in the dictionary), but only park within the leisure tag, a couple different types of landuse, and bus_stop within the highway tag.

tags = {
    'amenity':None,
    'leisure':'park',
    'landuse':['retail','commercial'],
    'highway':'bus_stop'
}

gdf = ox.pois_from_place(place='Santa Clara, Santa Clara County, California, USA', tags=tags)

Since there are a lot of POIs without consistent amenity tags (e.g., parks), this generalization was necessary for me to retrieve a broad array of POI types and seemed like low hanging fruit.

A **kwargs parameter in the pois_from* functions continues to accept the amenities parameter and automatically converts it to a tag dictionary to handle any existing uses.

Happy Friday,
Chester

[HTenkanen]

@chesterharvey There is actually an open issue for this feature, see #190. So definitely this is something that would be good and important to include. Great that you worked on it! 👍 And the way you solved this issue follows more or less the discussion that we had with @gboeing.

[gboeing]

@chesterharvey thanks! Following up on @HTenkanen's comment, would you take a quick look at #190 to see how your fork corresponds to the discussion there, and I'd be happy to review a pull request if you'll open one.

[chesterharvey]

Thanks @HTenkanen and @gboeing for pointing out #190! Yes, my implementation is surprisingly similar, dictionaries and all. The biggest difference seems to be @HTenkanen's use of booleans to specify all or no values for a given tag key. I really like this; True is more intuitive than my use of None to specify all values, and False also intuitivly facilitates a NOT query.

Switching None to True is easy, but adding the False option will take a little more work to translate into an appropriate Overpass QL query string. I can work on that, hopefully within the next week. Will open a pull request once this is done.

[chesterharvey]

@gboeing Status update: I'm having trouble implementing a reliable NOT query, I think related to this known issue with Overpass that keeps popping up: drolbr/Overpass-API#192

What this means in practice is that a tag dictionary like the one below, which is supposed to exclude any element with a leisure tag, ends up excluding most of them, but not all.

tags = {
    'amenity':True,
    'leisure':False,
    'landuse':['retail','commercial'],
    'highway':'bus_stop',
}

For example, if I use this dictionary to query Piedmont, CA, I end up with a couple of stray elements with leisure tags:

gdf = ox.pois_from_place(place='Piedmont, Alameda County, California, 94611, USA', tags=tags)
gdf['leisure'].value_counts()

park              1
fitness_centre    1
Name: leisure, dtype: int64

One the one hand, the API is supposed to facilitate this sort of query, so maybe it makes sense to assume it will someday and include the functionality. On the other hand, I don't want to lead people into thinking it will work perfectly.

Thoughts?

I still think TRUE is the most intuitive way to ask for all values.

[gboeing]

I'm having trouble implementing a reliable NOT query, I think related to this known issue with Overpass that keeps popping up: drolbr/Overpass-API#192 What this means in practice is that a tag dictionary like the one below, which is supposed to exclude any element with a leisure tag, ends up excluding most of them, but not all.

@mmd-osm do you have any insights on this?

[mmd-osm]

drolbr/Overpass-API#192 was resolved in release 0.7.55.5, so this shouldn't happen anymore, if you're using overpass-api.de as a backend. In any case, we'd need a query in Overpass QL format, if this still needs further research.

[gboeing]

drolbr/Overpass-API#192 was resolved in release 0.7.55.5, so this shouldn't happen anymore, if you're using overpass-api.de as a backend. In any case, we'd need a query in Overpass QL format, if this still needs further research.

@chesterharvey if you're still experiencing this issue, could you produce an MWE in Overpass QL?

[gboeing]

@chesterharvey also, feel free to open a PR in the meantime. We can always add more commits there, and it'll make discussion a little easier.

[gboeing]

@chesterharvey just checking to see if this is still active and if you're interested in opening a PR?

[chesterharvey]

Sorry @gboeing, I let this slide for a while. I'm opening a PR now. My modified osm_poi_download and pois_from_place functions have optional parameters to return the Overpass queries constructed from the input tag dictionaries so they can the troubleshooting, if necessary.

[corneliailin]

Hi!

I am trying to retrieve name and geometry data for all the schools in Nairobi, Kenya. I am using the pois_from_place() function available in the osmnx package for this. When I plot the data on a map I see schools outside the boundaries of Nairobi and I don't understand why (please see graph attached). Doing the same thing with QGIS tells me that there are also more schools inside the boundaries of Nairobi that I don't see in my osmnx data.

Am I missing something here? Any help appreciated.

[elharchaoui]

Hi all,

I'm really impressed with the work done since the begining of this module project that was started from a needs ! I would like to thank you for this amazing job done and that collaborative work. it's really good to see that!

[xxC3llxx]

Hi,
I tried the below code but got an error.

Code :
import osmnx as ox
place = "Manhattan Island, Manhattan, New York City, New York, United States of America"

all_pois = ox.pois_from_place(place=place)

Error:
AttributeError: module 'osmnx' has no attribute 'pois_from_place'

Why this could be?

[gboeing]

@xxC3llxx per the contributing guidelines, the issue tracker is for tracking bugs or proposing pull requests: please read the package documentation for usage details, review the examples for usage demonstrations, and ask any remaining "how-to" usage questions on StackOverflow.