Fragment: Components for Rolling Your Own GIS Inside Jupyter Notebooks

I’ve been tinkering with maps again… I think we really do need to start paying more attention to them. In particular, how we can help students to create them in concert with data that has geospatial relevance.

Last night, preparing for a Skyped 9am meeting this morning, I was looking for a recipe for allowing students to draw a shape on a map, get polygon data describing the shape back from the map, and then use that polygon data to do a lookup for crime figures or OSM streets within that boundary.

My recipe was, as ever, a bit hacky, but it was a first step and it was gone past 1am…

The recipe took the following form:

#gdf is a (geo)pandas dataframe
#Get the lat/lon of the the first point (i.e. an arbitrary point)
example_coords = gdf[['lat','lon']].iloc[0].values
#Generate a URL to open a geojson.net at that location and zoom into it
print('https://geojson.net/#12/{}/{}'.format(example_coords[0],example_coords[1]))

The printed out link is actually clickable…

Draw the shape and grab the geoJSON:

#Paste the JSON data below so we can assign it to the jsondata variable
jsondata='''
{
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "properties": {},
      "geometry": {
        "type": "Polygon",
        "coordinates": [
          [
            [-0.960474, 51.165091],
            [-0.988426, 51.154801],
            [-1.007309, 51.141879],
            [-0.971955, 51.135133],
            [-0.945854, 51.151786],
            [-0.960474, 51.165091]
          ]
        ]
      }
    }
  ]
}
'''

We can then parse that into a Shapely object:

#If the json data is a FeatureCollection with a single features list item
# we can parse it as follows...
import json

#Load the json string as a python dict and extract the geometty
gj = json.loads(jsondata)['features'][0]['geometry']

from shapely.geometry import shape
#Generate and preview the shape
shape( gj )

#Alternatively, we can let geopandas handle the conversion
gpd = geopandas.GeoDataFrame.from_features( gj )

#Preview the first geometry item
gpd['geometry'].iloc[0]

So that… works… but it’s a bit clunky….

Let’s recap what we’re doing: go to geojson.net, draw a shape, get the JSON, turn the JSON into a Shapely shape, derive boundary data from the shape that we can use for searches within the shape.

But can we do better?

The ipyleaflet Jupyter notebook extension has been coming on leaps and bounds recently, and now includes the ability to embed ipywidgets within the map.

It also includes a draw control that lets us draw shapes on the map and get the corresponding geojson back into the notebook (code) environment.

Via StackOverflow, I found a handy recipe for getting the data back into the notebook from the shape drawn on the map. I also used the new ipyleaflet ipywidget integration to add a button to the map that lets us preview the shape:

#via https://gis.stackexchange.com/a/312462/119781
from ipyleaflet import Map, basemaps, basemap_to_tiles, DrawControl, WidgetControl
from ipywidgets import Button

from shapely.geometry import shape

watercolor = basemap_to_tiles(basemaps.Stamen.Watercolor)
m = Map(layers=(watercolor, ), center=(50, -1), zoom=5)
draw_control = DrawControl()

draw_control.circle = {
    "shapeOptions": {
        "fillColor": "#efed69",
        "color": "#efed69",
        "fillOpacity": 1.0
    }
}

feature_collection = {
    'type': 'FeatureCollection',
    'features': []
}

def handle_draw(self, action, geo_json):
    """Do something with the GeoJSON when it's drawn on the map"""    
    feature_collection['features'].append(geo_json)

draw_control.on_draw(handle_draw)

m.add_control(draw_control)

button = Button(description="Shape co-ords")

def on_button_clicked(b):
    #Generate and preview the shape
    display(shape(feature_collection['features'][0]['geometry']))
    #print(feature_collection)

button.on_click(on_button_clicked)

widget_control = WidgetControl(widget=button, position='bottomright')
m.add_control(widget_control)
m

#Draw a shape on the map below... then press the button...
#The data is available inside the notebook via the variable: feature_collection

This is the result:

The green shape is the rendering of a shapely object generated from the shapefile extracted from the drawn shape on the map. It looks slightly different to the one on the map because “projections, innit?”. But it is the same. Just displayed in a differently projected co-ordinate grid.

So that’s maybe an hour last night working on my own hacky recipe and finding the SO question last night, half an hour today trying out the SO recipe and another half hour on this blog post.

And the result is a recipe that would allow us to easily draw a shape on a map and grab some data back from it….

Like this for example…

Start off by getting the boundary into the form required by the Police API:

drawn_area = shape(feature_collection['features'][0]['geometry'])

#Get the co-ords
lon, lat = drawn_area.exterior.coords.xy
drawn_boundary = list(zip(lat, lon))

Then we can look-up crimes in that area and pop the data into a geopandas dataframe:

import pandas as pd
import geopandas
from shapely.geometry import Point

#get crime in area
from police_api import PoliceAPI
api = PoliceAPI()


def setCrimesAsGeoDataFrame(crimes, df=None):
    ''' Convert crimes result to geodataframe. '''
    if df is None:
        df=pd.DataFrame(columns = ['cid','type', 'month','name','lid','location','lat','lon'])
        #[int, object, object, int, object, float, float]
    for c in crimes:
        df = df.append({'cid':c.id,'type':c.category.id, 'month':c.month, 'name':c.category.name,
                            'lat':c.location.latitude,'lon':c.location.longitude, 
                        'lid':c.location.id,'location':c.location.name }, ignore_index=True)

    df['lat']=df['lat'].astype(float)
    df['lon']=df['lon'].astype(float)
    
    df['Coordinates'] = list(zip(df['lon'], df['lat']))
    df['Coordinates'] = df['Coordinates'].apply(Point)
    
    return geopandas.GeoDataFrame(df, geometry='Coordinates')

crimes_df = setCrimesAsGeoDataFrame( api.get_crimes_area(drawn_boundary, date='2019-01') )
crimes_df.head()

Then we can plot the data on a map:

from ipyleaflet import Marker, MarkerCluster
import geopandas

m = Map( zoom=2)
m.add_layer(MarkerCluster(
    markers=[Marker(location=geolocation.coords[0][::-1]) for geolocation in crimes_df.geometry])
    )
m

(Cue another 90 mins that should only have been 20 mins because I made a stupid mistake of centering the map on (51, 359), as cut and pasted from code-on-the-web, rather than a sensible (51, -1) and then spent ages trying to debug the wrong thing…)

We can then do another iteration and create a map with a draw control that lets us draw a shape on the map and button that will take the shape, look up crimes within it, and plot corresponding markers on the same map…

m = Map(layers=(watercolor, ), center=(50.9,-1), zoom=9)

draw_control.on_draw(handle_draw)

m.add_control(draw_control)

button = Button(description="Shape co-ords")

feature_collection = {
    'type': 'FeatureCollection',
    'features': []
}

def on_button_clicked2(b):
    ''' Get the crimes back within the location. '''
    crimes_df = setCrimesAsGeoDataFrame( api.get_crimes_area(drawn_boundary, date='2019-01') )
    m.add_layer(MarkerCluster(
        markers=[Marker(location=geolocation.coords[0][::-1]) for geolocation in crimes_df.geometry])
    )

button.on_click(on_button_clicked2)

widget_control = WidgetControl(widget=button, position='bottomright')
m.add_control(widget_control)
m

Here’s how it looks:

For some reason it also seems to be pulling reports from just outside the boundary area? But it’s close enough…. and we could always post-filter the data returned from the Police API to select just the crimes reported within the drawn boundary.

Related notebooks are in this repo: https://github.com/psychemedia/crime-data-demo. The above is demoed in the Initial Look At the Data notebook, I think. The repo should be fully Binderised (I think).

PS Getting this far is easy, and demonstrates one of the differences between making “functional” code available to students and “polished” apps. If you draw lots of shapes, or the wrong sort of shape on the map, things will go wrong. The code to handle the complex shapes will have to get more defensive and more complex, either telling students that the shapes profile is not acceptable, or finding ways of handling multiple shapes properly. At the moment, if you want to use the same map to clear the shape and get some new data in, any crime markers on the canvas are not cleared. So that’s a usability thing. And fixing it adds complexity to the code and complexity (possibly) to the UI. Or you could just rerun the cell and you will start with a fresh map canvas.

Author: Tony Hirst

I'm a Senior Lecturer at The Open University, with an interest in #opendata policy and practice, as well as general web tinkering...