GeoPandas - geospatial data in Python made easy

class: center, middle

# GeoPandas

##  Geospatial data in Python made easy

Joris Van den Bossche, EuroScipy, August 30, 2017

https://github.com/jorisvandenbossche/talks/

.affiliations[
  ![:scale 65%](img/logoUPSayPlusCDS_990.png)
  ![:scale 25%](img/inria-logo.png)
]

---
class: center, middle

#  Geospatial data

<!--

# Geospatial data

Geospatial data is everywhere -->

---

background-image: url(img/us_census_density.jpg)

---
class: center

![:scale 100%](img/inbo_bird_tracking.png)

.credits[
Bird tracking. From [slides](https://speakerdeck.com/peterdesmet/fly-my-pretties-fly-exploring-visualizing-and-publishing-bird-tracking-data) by Peter Desmet (INBO)
]

---
class: center

![:scale 90%](img/pickups_dropoffs_narrow.jpg)

.credits[
NYC taxi data (pick-up and drop-off locations), from [blog post](https://medium.com/towards-data-science/if-taxi-trips-were-fireflies-1-3-billion-nyc-taxi-trips-plotted-b34e89f96cfa) by Ravi Shekhar]

---
# Raster vs vector data

![:scale 49%](img/raster_example.png)
![:scale 49%](img/vector_example.png)

--
count: false

.right[
### -> in this talk: focus on vector data
]

???

Two major families of geospatial data

raster: grid based (topology lacking, difficult to link to tabular data)
vector: coordinate based objects, topological

here: vector

vector -> common abstraction model in many software
Open Geospatial consortium standard

---

# OGC Simple Features

Point, LineString, Polygon

.center[
![:scale 60%](img/simple_features.svg)
]

... and MultiPoint, MultiLineString, MultiPolygon, GeometryCollection

--
count: false

Attributes : each vector feature can have a record in attribute table

???

and that is where geopandas comes into play

but before talking about geopandas, first a bit more general about open source geospatial software

WKT / WKB encodings for those simple features

# geospatial software

This presentation: in python

but everything I will present -> builds upon widely used open source libraries

---
class: middle, center

# Open source geospatial software

.center[
![:scale 70%](img/Open_Source_Geospatial_Foundation.svg)
]

???

Open Source Geospatial Foundation

OSGeo was created to support the collaborative development of open source geospatial software, and promote its widespread use.

---

# GDAL / OGR

### Geospatial Data Abstraction Library.

--
count:false

* The swiss army knife for geospatial.
* Read and write Raster (GDAL) and Vector (OGR) datasets
* More than 200 (mainly) geospatial formats and protocols.

.center[
![:scale 100%](img/gdal_formats)
]

.credits[
Slide from "GDAL 2.2 What's new?" by Even Rouault (CC BY-SA)
]

???

GDAL is a translator library for raster and vector geospatial data formats. As a library, it presents a single raster abstract data model and single vector abstract data model to the calling application for all supported formats. It also comes with a variety of useful command line utilities for data translation and processing.

---

# GDAL / OGR

### Widely used (FOSS & proprietary)

.center[
![:scale 100%](img/gdal_users)
]

.credits[
Slide from "GDAL 2.2 What's new?" by Even Rouault (CC BY-SA)
]

---

# GEOS

## Geometry Engine Open Source

* C/C++ port of a subset of Java Topology Suite (JTS)
* Most widely used geospatial C++ geometry library
* Implements geometry objects (simple features), spatial predicate functions and spatial operations

Used under the hood by many applications (QGIS, PostGIS, MapServer, GRASS, GeoDjango, ...)

[geos.osgeo.org](http://geos.osgeo.org)

<!-- ---
# PROJ.4

C library for performing conversions between cartographic projections.

[http://proj4.org/](http://proj4.org/) -->

---

# Python geospatial packages

--
count:false

Interfaces to widely used libraries:

- Python bindings to GDAL/OGR (`from osgeo import gdal, ogr`)
- Python interface to PROJ.4: [`pyproj`](https://jswhit.github.io/pyproj/)
- Pythonic binding to GDAL/OGR:
  - GDAL: [`rasterio`](https://mapbox.github.io/rasterio/)
  - OGR: [`fiona`](http://toblerity.org/fiona/README.html)
- Python package based on GEOS: [`shapely`](https://shapely.readthedocs.io/en/latest/)

---

# Shapely

Python package for the manipulation and analysis of geometric objects

Pythonic interface to GEOS

--
count:false

```python
>>> from shapely.geometry import Point, LineString, Polygon

>>> point = Point(1, 1)
>>> line = LineString([(0, 0), (1, 2), (2, 2)])
>>> poly = line.buffer(1)
```

```python
>>> poly.contains(point)
True
```

???

# Shapely

typical predicates and operations

(images from shapely docs)

---

# GeoPandas

* Started by Kelsey Jordahl in 2013
* Make working with geospatial data in python easier
* Extends the pandas data analysis library to work with geographic objects and spatial operations
* Combines the power of whole ecosystem of (geo) tools (pandas, geos, shapely, gdal, fiona, pyproj, rtree, ...)

Documentation: http://geopandas.readthedocs.io/

???

make working with geospatial data like working with any other kind of data in python
(data stack, numpy, pandas and other tools around those)
analysis for which you otherwise would need desktop GIS applications (QGIS, ArcGIS) or geospatial databases (PostGIS)

makes pandas objects geometry aware

---
class: middle, center

# Demo time!

See [static version](http://nbviewer.jupyter.org/github/jorisvandenbossche/talks/blob/master/2017_EuroScipy_geopandas/geopandas_demo.ipynb)

---

# Summary

* Read and write variety of formats (fiona, GDAL/OGR)
* Familiar manipulation of the attributes (pandas dataframe)
* Element-wise spatial predicates (intersects, within, ...) and operations (intersection, union, difference, ..) (shapely)
* Re-project your data (pyproj)
* Quickly visualize the geometries (matplotlib, descartes)
* More advanced spatial operations: spatial joins and overlays (rtree)

---

# Ecosystem

[geoplot](http://www.residentmar.io/geoplot/index.html) (high-level geospatial visualization), [cartopy](http://scitools.org.uk/cartopy/) (projection aware cartographic library)

[folium](https://github.com/python-visualization/folium) (Leaflet.js maps)

[OSMnx](http://geoffboeing.com/2016/11/osmnx-python-street-networks/) (python for street networks)

[PySAL](http://pysal.readthedocs.io/en/latest/index.html) (Python Spatial Analysis Library)

[rasterio](https://mapbox.github.io/rasterio/) (working with geospatial raster data)

...

---
class: middle, center

# But ... <span style="color: transparent;">GeoPandas is slow</span>

---
count: false
class: middle, center

# But ... GeoPandas is slow

---

# Some timings

Basic `within` and `distance` operation on 100.000 points:

```python
s.within(polygon)
s.distance(polygon)
```

![:scale 49%](img/timings_within.png)
![:scale 49%](img/timings_distance.png)

---

# Comparison with PostGIS

```sql
-- What is the population and racial make-up of the neighborhoods of Manhattan?
SELECT
  neighborhoods.name AS neighborhood_name, Sum(census.popn_total) AS population,
  100.0 * Sum(census.popn_white) / NULLIF(Sum(census.popn_total),0) AS white_pct,
  100.0 * Sum(census.popn_black) / NULLIF(Sum(census.popn_total),0) AS black_pct
FROM nyc_neighborhoods AS neighborhoods
JOIN nyc_census_blocks AS census
ON ST_Intersects(neighborhoods.geom, census.geom)
GROUP BY neighborhoods.name
ORDER BY white_pct DESC;
```

```python
res = geopandas.sjoin(nyc_neighborhoods, nyc_census_blocks, op='intersects')
res = res.groupby('NAME')[['POPN_TOTAL', 'POPN_WHITE', 'POPN_BLACK']].sum()
res['POPN_BLACK'] = res['POPN_BLACK'] / res['POPN_TOTAL'] * 100
res['POPN_WHITE'] = res['POPN_WHITE'] / res['POPN_TOTAL'] * 100
res.sort_values('POPN_WHITE', ascending=False)
```

Disclaimer: dummy benchmark, and I am not a PostGIS expert!

Example from [Boundless tutorial](http://workshops.boundlessgeo.com/postgis-intro/) (CC BY SA)

---

# Comparison with PostGIS

![:scale 49%](img/timings_distance2.png)
![:scale 49%](img/timings_sjoin.png)

Disclaimer: dummy benchmark, and I am not a PostGIS expert!

Example from [Boundless tutorial](http://workshops.boundlessgeo.com/postgis-intro/) (CC BY SA)

---

# Why is GeoPandas slow?

- GeoPandas stores custom Python objects in arrays
- For operations, it iterates through those objects

--
count: false

Assume the following operation: distance from all geometries in a GeoSeries to a point:

.left-column[
```python
geoseries = GeoSeries(...)
point = Point(...)

geoseries.distance(point)
```
]

.right-column[

]

---
count: false

# Why is GeoPandas slow?

- GeoPandas stores custom Python objects in arrays
- For operations, it iterates through those objects

Assume the following operation: distance from all geometries in a GeoSeries to a point:

.left-column[
```python
geoseries = GeoSeries(...)
point = Point(...)

*geoseries.distance(point)
```
]

.right-column[

```python

# List-comprehension under the hood:
[geom.distance(p) for geom in geoseries]
```
]

---

# New developments

- `GeometryArray`: array-like to hold gemetry objects
  - only stores integer pointers to the C GEOS geometry objects, without wrapping in Python geometry object
  - operations iterate through those geometries in cython
  - wrapped into Python shapely object on access
- Custom `GeometryBlock` to integrate this in pandas

Better performance, less memory use, parallelization (nogil), dask?

Many thanks to Matthew Rocklin for his work!

---

# New timings

![:scale 40%](img/timings_within.png)
![:scale 40%](img/timings_distance.png)

![:scale 40%](img/timings_distance2.png)
![:scale 40%](img/timings_sjoin.png)

---
count:false

# New timings

![:scale 40%](img/timings_within_all.png)
![:scale 40%](img/timings_distance_all.png)

![:scale 40%](img/timings_distance2_all.png)
![:scale 40%](img/timings_sjoin_all.png)

---

# Want to try this out?

Feedback and polishing needed!

See https://github.com/geopandas/geopandas/issues/473

Build new version (needs GEOS and shapely already installed):

* clone https://github.com/geopandas/geopandas
* checkout 'geopandas-cython' branch
* `make inplace`

--
count: false

## Friday -> sprints!

---
class: middle

# Thanks for listening!

## Thanks to all contributors!

## Those slides:

- https://github.com/jorisvandenbossche/talks/
- [jorisvandenbossche.github.io/talks/2017_EuroScipy_geopandas](
    http://jorisvandenbossche.github.io/talks/2017_EuroScipy_geopandas)

---
# About me

Joris Van den Bossche

- PhD bio-science engineer, air quality research
- pandas core dev
- Currently working at the Paris-Saclay Center for Data Science (Inria)

https://github.com/jorisvandenbossche

[@jorisvdbossche](https://twitter.com/jorisvdbossche)