This post demonstrates how to
- Display a map from OpenStreetMaps using Folium
- Add custom shape files to define regions of interest
- Color the regions of interest based on data in a pandas dataframe
- (New) Brief introduction to GeoPandas.
The code of this post is available at https://github.com/rsandstroem/IPythonNotebooks/blob/master/GeoMapsFoliumDemo/GeoMapsFoliumDemo.ipynb .
August 1, 2017: This post is almost a year old, but I decided to make some technical improvements to it to improve the viewing experience online. (It is the Swiss National Day after all!) During this renovation, I stumbled upon GeoPandas, which is a great package worth mentioning in this context.
Using the method described here, you do not need any local map data. The maps are automatically pulled from OpenStreetMap.
However, you might want use concepts such as regions, and you need to define them through files containing those geolocations. There are multiple sources of data which is freely available, in this example we are using shape files from this location: http://data.geo.admin.ch.s3.amazonaws.com/ch.swisstopo.swissboundaries3d-kanton-flaeche.fill/data.zip
We have provided the shape files we will use for this demonstration, normally you would need to retrieve and unpack the map data.
There are many file formats commonly used for maps. In this case what we have are a set of shapefiles. Since we will be working with geojson files in this demonstration you need to convert them to geojson.
A sample geojson file has been supplied with this tutorial, but for reference you can create a geojson file from the shape files from a console like this:
ogr2ogr -f GeoJSON -t_srs EPSG:4326 -simplify 1000 switzerland.geojson swissBOUNDARIES3D_1_2_TLM_KANTONSGEBIET.shp
You now have a map in your local directory, formatted as a geojson file. I use the simplify
option to reduce the file size, at the expense of accuracy.
For this demonstration we need "folium". Install it using
pip install folium
or if you are using the Anaconda you can do
sudo conda install --channel https://conda.binstar.org/IOOS folium
where you can replace the channel with your choice from the options listed by
binstar search -t conda folium
import folium
print folium.__version__
First we show the map without any use of the geojson we created from the shapefiles.
kanton_map = folium.Map(location=[46.8, 8.33],
tiles='Mapbox Bright', zoom_start=7)
kanton_map
You should now see an map of Switzerland with surrounding countries and cities. Next, display an OpenStreetMap and overlay the cantons of Switzerland on top from the geojson file created earlier.
kanton_map.choropleth(geo_path='switzerland.geojson')
kanton_map
The purpose of using a color filled map like the one above is usually to have the color represent something in the data, e.g., sales revenue per district. Since this data is separate from the geometrical description we need to retrieve it from an external source and join it with the geometrical shapes from the geojson file.
Read in some data for the canton's of Switzerland from a csv file. We will display this data on the map we just created.
import pandas as pd
kanton_overview = r'Switzerland_overview.csv'
kanton_data = pd.read_csv(kanton_overview)
kanton_data
The KANTONSNUM in the json file is used to ensure correct assignment of the data in the demographics data file. Try using the name of the cantons instead. What happens, and why? (Hint: open the geojson file and compare with the csv file.)
kanton_map2 = folium.Map(location=[46.8, 8.33],
zoom_start=7.5)
kanton_map2.choropleth(geo_path='switzerland.geojson', data=kanton_data,
columns=['CantonNumber', 'Density'],
key_on='feature.properties.KANTONSNUM',
threshold_scale=[100, 200, 300, 500, 1000, 5000],
fill_color='BuPu')
kanton_map2
While making some technical improvements to this post, I stumbled upon GeoPandas. While I cannot say I have used it much yet, I think it is a great way to read in geojson files and work with their contents.
import geopandas as gpd
df = gpd.read_file('switzerland.geojson')
print df.columns
df.head()
Very easy to make a pandas dataframe from this json file! Every row represents one area, and the last column contains the shape of the area as a three dimensional polygon, i.e., the longitude, latitude and altitude of the regional borders. The dataframe also contains data columns, such as number of inhabitants (EINWOHNERZ) and surface area (KANTONSFLA).
A nice feature of using GeoPandas in a Jupyter Notebook is the ease at which we can draw the content of the dataframe:
%matplotlib inline
import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(12,8))
ax.set_aspect('equal')
df = df.to_crs({'init': 'epsg:4326'}) # change projection
df['density'] = df.EINWOHNERZ / df.KANTONSFLA # calculate population density
df.plot(ax=ax, column='density', cmap='OrRd')
Since the geometry of the cantons are represented as Shapely polygons, Jupyter can render them directly. For example, the canton of Valais has index 2 in the dataframe df
:
print type(df.geometry[2])
df.geometry[2]
Let's combined what we learned and add map pins for mountains on the border where the altitude exceeds 3500 m. First we extract the three dimensional points from the dataframe where the altitude exceeds 3500 m.
import numpy as np
peaks = list()
for kanton in df.index:
arr = np.array(df.geometry[kanton].exterior.coords)
peaks.append(arr[np.where(arr[:,2]>3500)].tolist())
peaks = [item for sublist in peaks for item in sublist] # flatten to a single list of 3d points
print peaks[:3]
Then we create a feature group consisting of map markers for the peaks we found in the previous step. The altitude is added as a pop-up that will be displayed if the user clicks on a map marker. Finally we add this new feature group to the populations density map which we created earlier.
feature_group = folium.FeatureGroup("Locations")
for peak in peaks:
feature_group.add_child(folium.Marker(location=[peak[1], peak[0]], popup=str(peak[2])))
kanton_map2.add_child(feature_group)
kanton_map2
The regions with high border mountains are located along the southern border of Switzerland, the Alpes. According to the population density we overlayed earlier, the cantons with high mountains are also the least populated. Hardly surprising, it is hard to commute to work if you live on top of a glacier!
Comments !