This is a modified version of http://nbviewer.ipython.org/github/geopandas/geopandas/blob/master/examples/choropleths.ipynb with the patch to add a legend to those figures.
%matplotlib inline
import geopandas as gp
# we use PySAL for loading a test shapefile
# replace this cell if you have a local shapefile and want to use GeoPandas readers
import pysal as ps
pth = ps.examples.get_path("columbus.shp")
tracts = gp.GeoDataFrame.from_file(pth)
ax = tracts.plot(column='CRIME', scheme='QUANTILES', k=3, colormap='OrRd')
C:\Anaconda\envs\cartopy\lib\site-packages\pandas\computation\expressions.py:184: UserWarning: evaluating in Python space because the '*' operator is not supported by numexpr for the bool dtype, use '&' instead unsupported[op_str]))
tracts.head()
AREA | COLUMBUS_ | COLUMBUS_I | CP | CRIME | DISCBD | EW | HOVAL | INC | NEIG | ... | NSA | NSB | OPEN | PERIMETER | PLUMB | POLYID | THOUS | X | Y | geometry | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | 0.309441 | 2 | 5 | 0 | 15.725980 | 5.03 | 1 | 80.467003 | 19.531 | 5 | ... | 1 | 1 | 2.850747 | 2.440629 | 0.217155 | 1 | 1000 | 38.799999 | 44.070000 | POLYGON ((8.624129295349121 14.23698043823242,... |
1 | 0.259329 | 3 | 1 | 0 | 18.801754 | 4.27 | 0 | 44.567001 | 21.232 | 1 | ... | 1 | 1 | 5.296720 | 2.236939 | 0.320581 | 2 | 1000 | 35.619999 | 42.380001 | POLYGON ((8.252790451049805 14.23694038391113,... |
2 | 0.192468 | 4 | 6 | 0 | 30.626781 | 3.89 | 1 | 26.350000 | 15.956 | 6 | ... | 1 | 1 | 4.534649 | 2.187547 | 0.374404 | 3 | 1000 | 39.820000 | 41.180000 | POLYGON ((8.653305053710938 14.00809001922607,... |
3 | 0.083841 | 5 | 2 | 0 | 32.387760 | 3.70 | 0 | 33.200001 | 4.477 | 2 | ... | 1 | 1 | 0.394427 | 1.427635 | 1.186944 | 4 | 1000 | 36.500000 | 40.520000 | POLYGON ((8.459499359130859 13.82034969329834,... |
4 | 0.488888 | 6 | 7 | 0 | 50.731510 | 2.83 | 1 | 23.225000 | 11.252 | 7 | ... | 1 | 1 | 0.405664 | 2.997133 | 0.624596 | 5 | 1000 | 40.009998 | 38.000000 | POLYGON ((8.685274124145508 13.63951969146729,... |
5 rows × 21 columns
Slightly adapted functions:
def __pysal_choro(values, scheme, k=5):
""" Wrapper for choropleth schemes from PySAL for use with plot_dataframe
Parameters
----------
values
Series to be plotted
scheme
pysal.esda.mapclassify classificatin scheme ['Equal_interval'|'Quantiles'|'Fisher_Jenks']
k
number of classes (2 <= k <=9)
Returns
-------
values
Series with values replaced with class identifier if PySAL is available, otherwise the original values are used
"""
try:
from pysal.esda.mapclassify import Quantiles, Equal_Interval, Fisher_Jenks
schemes = {}
schemes['equal_interval'] = Equal_Interval
schemes['quantiles'] = Quantiles
schemes['fisher_jenks'] = Fisher_Jenks
s0 = scheme
scheme = scheme.lower()
if scheme not in schemes:
scheme = 'quantiles'
print('Unrecognized scheme: ', s0)
print('Using Quantiles instead')
if k < 2 or k > 9:
print('Invalid k: ', k)
print('2<=k<=9, setting k=5 (default)')
k = 5
binning = schemes[scheme](values, k)
values = binning.yb
except ImportError:
print('PySAL not installed, setting map to default')
return binning
# adding linewidth parameter
def plot_polygon(ax, poly, facecolor='red', edgecolor='black', alpha=0.5, linewidth=1):
""" Plot a single Polygon geometry """
from descartes.patch import PolygonPatch
a = np.asarray(poly.exterior)
# without Descartes, we could make a Patch of exterior
ax.add_patch(PolygonPatch(poly, facecolor=facecolor, alpha=alpha))
ax.plot(a[:, 0], a[:, 1], color=edgecolor, linewidth=linewidth)
for p in poly.interiors:
x, y = zip(*p.coords)
ax.plot(x, y, color=edgecolor, linewidth=linewidth)
def plot_multipolygon(ax, geom, facecolor='red', edgecolor='black', alpha=0.5, linewidth=1):
""" Can safely call with either Polygon or Multipolygon geometry
"""
if geom.type == 'Polygon':
plot_polygon(ax, geom, facecolor=facecolor, edgecolor=edgecolor, alpha=alpha, linewidth=linewidth)
elif geom.type == 'MultiPolygon':
for poly in geom.geoms:
plot_polygon(ax, poly, facecolor=facecolor, edgecolor=edgecolor, alpha=alpha, linewidth=linewidth)
import numpy as np
from geopandas.plotting import (plot_linestring, plot_point, norm_cmap)
def plot_dataframe(s, column=None, colormap=None, alpha=0.5,
categorical=False, legend=False, axes=None, scheme=None,
k=5, linewidth=1):
""" Plot a GeoDataFrame
Generate a plot of a GeoDataFrame with matplotlib. If a
column is specified, the plot coloring will be based on values
in that column. Otherwise, a categorical plot of the
geometries in the `geometry` column will be generated.
Parameters
----------
GeoDataFrame
The GeoDataFrame to be plotted. Currently Polygon,
MultiPolygon, LineString, MultiLineString and Point
geometries can be plotted.
column : str (default None)
The name of the column to be plotted.
categorical : bool (default False)
If False, colormap will reflect numerical values of the
column being plotted. For non-numerical columns (or if
column=None), this will be set to True.
colormap : str (default 'Set1')
The name of a colormap recognized by matplotlib.
alpha : float (default 0.5)
Alpha value for polygon fill regions. Has no effect for
lines or points.
legend : bool (default False)
Plot a legend (Experimental; currently for categorical
plots only)
axes : matplotlib.pyplot.Artist (default None)
axes on which to draw the plot
scheme : pysal.esda.mapclassify.Map_Classifier
Choropleth classification schemes
k : int (default 5)
Number of classes (ignored if scheme is None)
Returns
-------
matplotlib axes instance
"""
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D
from matplotlib.colors import Normalize
from matplotlib import cm
if column is None:
return plot_series(s.geometry, colormap=colormap, alpha=alpha, axes=axes)
else:
if s[column].dtype is np.dtype('O'):
categorical = True
if categorical:
if colormap is None:
colormap = 'Set1'
categories = list(set(s[column].values))
categories.sort()
valuemap = dict([(k, v) for (v, k) in enumerate(categories)])
values = [valuemap[k] for k in s[column]]
else:
values = s[column]
if scheme is not None:
binning = __pysal_choro(values, scheme, k=k)
values = binning.yb
# set categorical to True for creating the legend
categorical = True
binedges = [binning.yb.min()] + binning.bins.tolist()
categories = ['{0:.2f} - {1:.2f}'.format(binedges[i], binedges[i+1]) for i in range(len(binedges)-1)]
cmap = norm_cmap(values, colormap, Normalize, cm)
if axes == None:
fig = plt.gcf()
fig.add_subplot(111, aspect='equal')
ax = plt.gca()
else:
ax = axes
for geom, value in zip(s.geometry, values):
if geom.type == 'Polygon' or geom.type == 'MultiPolygon':
plot_multipolygon(ax, geom, facecolor=cmap.to_rgba(value), alpha=alpha, linewidth=linewidth)
elif geom.type == 'LineString' or geom.type == 'MultiLineString':
plot_multilinestring(ax, geom, color=cmap.to_rgba(value))
# TODO: color point geometries
elif geom.type == 'Point':
plot_point(ax, geom, color=cmap.to_rgba(value))
if legend:
if categorical:
patches = []
for value, cat in enumerate(categories):
patches.append(Line2D([0], [0], linestyle="none",
marker="o", alpha=alpha,
markersize=10, markerfacecolor=cmap.to_rgba(value)))
ax.legend(patches, categories, numpoints=1, loc='best')
else:
# TODO: show a colorbar
raise NotImplementedError
plt.draw()
return ax
ax = plot_dataframe(tracts, column='CRIME', scheme='QUANTILES', k=3, colormap='OrRd', legend=True)
# thicker linewidth for polygons
ax = plot_dataframe(tracts, column='CRIME', scheme='QUANTILES', k=3, colormap='OrRd', legend=True, linewidth=2.0)
tracts.plot(column='CRIME', scheme='equal_interval', k=7, colormap='OrRd')
<matplotlib.axes._subplots.AxesSubplot at 0x184964e0>
plot_dataframe(tracts, column='CRIME', scheme='equal_interval', k=7, colormap='OrRd', legend=True)
<matplotlib.axes._subplots.AxesSubplot at 0x179d9208>