（数据科学学习手札82）基于geopandas的空间数据分析——geoplot篇(上)

本文示例代码和数据已上传至我的Github仓库https://github.com/CNFeffery/DataScienceStudyNotes

1 简介

　　在前面的基于geopandas的空间数据分析系列文章中，我们已经对geopandas的基础知识、基础可视化，以及如何科学绘制分层设色地图展开了深入的学习，而利用geopandas+matplotlib进行地理可视化固然能实现常见的地图可视化，且提供了操纵图像的极高自由度，但对使用者matplotlib的熟悉程度要求较高，制作一幅地图可视化作品往往需要编写较多的代码，而geoplot基于geopandas，提供了众多高度封装的绘图API，很大程度上简化了绘图难度，就像seaborn之于matplotlib。

图1

　　本文是基于geopandas的空间数据分析系列文章的第6篇，通过本文你将学习geoplot中的基础绘图API。

2 geoplot基础

　　推荐使用conda install --channel conda-forge geoplot来安装geoplot，可以避免很多恼人的依赖包问题。首先我们从一个简单的例子来初探一下geoplot的基础使用流程：

2.1 从一个简单的例子出发

　　我们下面所使用到的数据：nyc-boroughs.geojson，记录了纽约的行政区域面文件：

import geopandas as gpd

%matplotlib inline

# 读入纽约行政区域面文件
nyc_boroughs = gpd.read_file('geometry/nyc-boroughs.geojson')
nyc_boroughs.head()

图2

　　以及nyc-collision-factors.geojson，包含了纽约所发生的车祸记录点以及相关信息数据：

# 读入纽约车祸记录点文件
nyc_collision_factors = gpd.read_file('geometry/nyc-collision-factors.geojson')
nyc_collision_factors.head()

图3

　　首先我们使用geoplot中的polyplot来绘制纽约行政区划，这里使用geoplot自带的Albers等面积投影作为投影：

import geoplot as gplt
import geoplot.crs as gcrs
import matplotlib.pyplot as plt

ax = gplt.polyplot(df=nyc_boroughs,
                   projection=gcrs.AlbersEqualArea())

plt.savefig("图4.png", bbox_inches='tight', pad_inches=0, dpi=300)

图4

　　接着我们使用geoplot中的pointplot将点叠加到图4上：

ax = gplt.polyplot(df=nyc_boroughs,
                   projection=gcrs.AlbersEqualArea())

ax = gplt.pointplot(df=nyc_collision_factors,
                    s=2,
                    color='grey',
                    alpha=0.2,
                    linewidth=0, # 设置轮廓粗细为0
                    ax=ax)

plt.savefig("图5.png", bbox_inches='tight', pad_inches=0, dpi=300)

图5

　　为了让车祸密集的区域更突出，我们将点图层换成核密度图层：

ax = gplt.polyplot(df=nyc_boroughs,
                   projection=gcrs.AlbersEqualArea())

# 叠加核密度图层
ax = gplt.kdeplot(df=nyc_collision_factors,
                  cmap='Reds',
                  shade=True,
                  shade_lowest=True,
                  clip=nyc_boroughs,
                  ax=ax)

plt.savefig("图6.png", bbox_inches='tight', pad_inches=0, dpi=300)

图6

　　从这个简单的例子中我们可以大致了解到，geoplot在geopandas处理好的数据基础上，针对不同类型图层封装了各自不同的API，由用户自主传入对应类型的矢量数据进行图层叠加，以得到最终结果，且可以兼容matplotlib，譬如上面我们最终使用plt.savaefig()对图片进行保存，下面我们就来详细学习geoplot的基础知识。

2.2 geoplot绘图API

　　在geoplot中内置了功能丰富的绘图API，只需要传入GeoDataFrame格式的矢量数据即可进行绘图（但切记geoplot中传入的数据必须为WGS84地理坐标系，所有的投影转换在geoplot各绘图函数内部传参实现即可！）

2.2.1 Pointplot

　　geoplot中的pointplot即为散点图，其针对点数据进行可视化，其主要参数如下：

df：传入对应的GeoDataFrame对象

projection：用于指定投影坐标系，传入geoplot.crs中的对象

hue：当需要根据df中的某列或外部的其他序列数据来映射散点的色彩时，可传入对应df中指定列名或外部序列数据，默认为None即不进行设色

cmap：和matplotlib中的cmap使用方式一致，用于控制色彩映射方案

scheme：作用类似geopandas中的scheme参数，用于控制分层设色，详见本系列文章的分层设色篇，但不同的是在geoplot0.4.0版本之后此参数不再搭配分层数量k共同使用，而是更新为传入mapclassify分段结果对象，下文中会做具体演示

scale：用于设定映射散点大小的序列数据，格式同hue，默认为None即每个散点等大小

limits：元组型，当scale不为None时，用于设定散点大小尺寸范围，格式为(min, max)

s：当scale设置为None时，用于控制散点的尺寸大小

color：当hue设置为None时，用于控制散点的填充色彩

marker：用于设定散点的形状

alpha：控制全局色彩透明度

linewidths：控制散点轮廓宽度

edgecolors：控制散点轮廓颜色

legend：bool型，用于控制是否显示图例

legend_var：传入'hue'或scale，当设定为hue时图例显示色彩映射信息，当设定为'scale'时图例显示大小映射信息

legend_values：list型，用于自定义图例显示的各个具体数值

legend_labels：list型，用于自定义图例显示的各个具体数值对应的文字标签，与legend_values搭配使用

legend_kwargs：字典，在legend参数设置为True时来传入更多微调图例属性的参数

extent：元组型，用于传入左下角和右上角经纬度信息来设置地图空间范围，格式为(min_longitude, min_latitude, max_longitude, max_latitude)

figsize：元组型，用于控制画幅大小，格式为(x, y)

ax：matplotlib坐标轴对象，如果需要在同一个坐标轴内叠加多个图层就需要用这个参数传入先前待叠加的ax

　　知晓了上述主要参数之后，下面我们通过实际案例来学习修改各个参数得到的效果，使用到的数据为波士顿区划面数据以及波士顿部分地区Airbnb房源点数据：

图7

• 普通散点分布

  首先我们来简单绘制房源分布散点图，对大小、色彩、透明度等基础属性进行简单调整：

# 简单绘制波士顿行政区划
ax = gplt.polyplot(df=boston_zip_codes,
                   projection=gcrs.AlbersEqualArea(),
                   edgecolor='lightgrey',
                   linewidths=0.5)

gplt.pointplot(df=boston_airbnb_listings,
               ax=ax, # 叠加图层
               s=1,
               linewidths=0.1,
               color='grey',
               alpha=0.4)

plt.savefig("图8.png", bbox_inches='tight', pad_inches=0, dpi=300)

图8

　　通过这样一张简单的图我们是看不出太多信息的，只能大致看出哪些地方房源分布较多。

• 映射房源价格到色彩上

　　将房源价格列作为色彩映射列，使用mapclassify中的分位数法将价格区间等分成五段，并使用其他的视觉参数和自定义图例参数：

import mapclassify as mc

#解决中文显示问题
plt.rcParams['font.sans-serif'] = ['SimHei'] # 指定默认字体
plt.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题

# 简单绘制波士顿行政区划
ax = gplt.polyplot(df=boston_zip_codes,
                   projection=gcrs.AlbersEqualArea(),
                   edgecolor='lightgrey',
                   linewidths=0.5)

scheme = mc.Quantiles(boston_airbnb_listings['price'], k=5)
gplt.pointplot(df=boston_airbnb_listings,
               ax=ax, # 叠加图层
               s=1, # 散点大小
               linewidths=0.1, # 散点轮廓宽度
               hue='price', # 以price作为色彩映射列
               cmap='Reds', # 色彩方案为Reds
               scheme=scheme, # 传入mapclassify对象
               legend=True, # 开启图例
               legend_kwargs={
                   'loc': 'lower right', # 图例位置
                   'title': '价格区间', # 图例标题
                   'title_fontsize': 8, # 图例标题字体大小
                   'fontsize': 6, # 图例非标题外字体大小
                   'shadow': True, # 添加图例阴影
               },
               legend_labels=['80%-100%价格房源',
                              '60%-80%价格房源',
                              '40%-60%价格房源',
                              '20%-40%价格房源',
                              '前20%价格房源'])

plt.savefig("图9.png", bbox_inches='tight', pad_inches=0, dpi=300)

图9

• 映射房源价格到尺寸上

　　看完了如何映射颜色，下面我们来看看如何将值映射到散点大小上，使用scale='price'来将房源价格映射到散点大小上，再配合一些相关参数进行绘图：

import numpy as np

# 简单绘制波士顿行政区划
ax = gplt.polyplot(df=boston_zip_codes,
                   projection=gcrs.AlbersEqualArea(),
                   edgecolor='lightgrey',
                   linewidths=0.5)

ax = gplt.pointplot(df=boston_airbnb_listings,
                    ax=ax, # 叠加图层
                    linewidths=0.2, # 散点轮廓宽度
                    scale='price', # 以price作为色彩映射列
                    color=np.array([0., 0., 0., 0.]), # 设置填充色为透明
                    edgecolor='grey', # 设置轮廓颜色
                    limits=(1, 16), # 设置散点的尺寸范围
                    legend=True, # 开启图例
                    legend_kwargs={
                        'loc': 'lower right', # 图例位置
                        'title': '价格区间', # 图例标题
                        'title_fontsize': 8, # 图例标题字体大小
                        'fontsize': 6, # 图例非标题外字体大小
                        'shadow': True, # 添加图例阴影
                        'markeredgecolor': 'grey', # 图例标记的轮廓色彩
                        'markeredgewidth': 0.2 # 图例标记的轮廓粗细
                    })

plt.savefig("图10.png", bbox_inches='tight', pad_inches=0, dpi=300)

图10

　　现在我们可以一眼看出那些半径较大的圆圈对应着价格较高的房源，值得注意的是在我们映射值到散点大小上时，默认条件下会自动在图例中按照等间距法分出5段，这样得到的图例各个圆圈大小过渡保证了均匀，当然你也可以自由地通过legend_values和legeng_labels这两个参数自定义图例内容。

• 同时映射颜色与尺寸

　　geoplot允许用户同时映射色彩和尺寸，但同一张图中的图例只能显示色彩或尺寸其中之一的信息，使用legend_var参数来选择让哪一种映射信息显示在图例上：

# 简单绘制波士顿行政区划
ax = gplt.polyplot(df=boston_zip_codes,
                   projection=gcrs.AlbersEqualArea(),
                   edgecolor='lightgrey',
                   linewidths=0.5)

scheme = mc.Quantiles(boston_airbnb_listings['price'], k=5)
gplt.pointplot(df=boston_airbnb_listings,
               ax=ax, # 叠加图层
               scale='price', # 以price作为尺寸映射列
               limits=(1, 16), # 设置散点的尺寸范围
               alpha=0.6, # 设置散点透明度
               linewidths=0.1, # 散点轮廓宽度
               hue='price', # 以price作为色彩映射列
               cmap='Reds', # 色彩方案为Reds
               scheme=scheme, # 传入mapclassify对象
               legend=True, # 开启图例
               legend_kwargs={
                   'loc': 'lower right', # 图例位置
                   'title': '价格区间', # 图例标题
                   'title_fontsize': 8, # 图例标题字体大小
                   'fontsize': 6, # 图例非标题外字体大小
                   'shadow': True, # 添加图例阴影
               },
               legend_labels=['80%-100%价格房源',
                              '60%-80%价格房源',
                              '40%-60%价格房源',
                              '20%-40%价格房源',
                              '前20%价格房源'])

plt.savefig("图11.png", bbox_inches='tight', pad_inches=0, dpi=300)

图11

2.2.2 Polyplot

　　geoplot中的polyplot用于绘制基础的面数据，并不像pointplot那样带有值映射功能，其主要参数如下：

df：传入对应的GeoDataFrame对象

projection：用于指定投影坐标系，传入geoplot.crs中的对象

extent：元组型，用于传入左下角和右上角经纬度信息来设置地图空间范围，格式为(min_longitude, min_latitude, max_longitude, max_latitude)

figsize：元组型，用于控制画幅大小，格式为(x, y)

ax：matplotlib坐标轴对象，如果需要在同一个坐标轴内叠加多个图层就需要用这个参数传入先前待叠加的ax

alpha：控制全局色彩透明度

linewidths：控制线宽度

edgecolors：控制线颜色

facecolor：控制填充颜色

linestyle：控制线样式，详情见本系列文章前作基础可视化篇图5

hatch：控制填充阴影纹路，详情见本系列文章前作基础可视化篇图7

　　下面我们就对纽约区划面数据进行举例说明：

gplt.polyplot(df=nyc_boroughs,
              projection=gcrs.AlbersEqualArea(),
              figsize=(10, 10),
              linewidths=0.5,
              linestyle='-.',
              edgecolors='grey',
              facecolor='#d9c09e',
              hatch='--')

plt.savefig("图12.png", bbox_inches='tight', pad_inches=0, dpi=300)

图12

2.2.3 Webmap

　　geoplot中的webmap用来添加在线瓦片地图底图，使得我们可以在在线地图上图层，但目前暂时只支持叠加基于点要素的图层。值得注意的是，因为常见在线地图如谷歌地图、OpenStreetMap、高德地图等的投影均为EPSG:3857也就是我们常说的Web Mercator，所以一旦要使用webmap，则投影锁死为EPSG:3857，其主要参数如下：

df：传入对应的GeoDataFrame对象

extent：元组型，用于传入左下角和右上角经纬度信息来设置地图空间范围，格式为(min_longitude, min_latitude, max_longitude, max_latitude)

figsize：元组型，用于控制画幅大小，格式为(x, y)

ax：matplotlib坐标轴对象，如果需要在同一个坐标轴内叠加多个图层就需要用这个参数传入先前待叠加的ax

zoom：int型，控制在线地图底图的缩放级别，越大越清楚，同时获取瓦片地图资源从而渲染地图所耗费的时间也越多，上限由具体所使用的在线地图所决定，通常情况最大缩放级别为18

provider：str型，用于指定在线地图底图的类型，下面会举例说明

　　下面我们将纽约车祸点数据叠加到在线地图上，这里我们选择provider参数为ST_TERRAIN_LINES，并设置缩放级别为11级：

ax = gplt.webmap(df=nyc_boroughs,
                 provider='ST_TERRAIN_LINES',
                 zoom=10)

ax = gplt.pointplot(df=nyc_collision_factors,
                    color='lightyellow',
                    edgecolor='yellow',
                    alpha=0.4,
                    s=1,
                    ax=ax)

plt.savefig("图13.png", bbox_inches='tight', pad_inches=0, dpi=300)

图13

　　如果想要切换底图样式，可以修改provider参数的输入，目前为止所有可用的地图如下图所示：

图14

2.3 在模仿中学习

　　在本系列文章基础可视化篇的最后我们对数据可视化专家用R绘制的澳大利亚火灾影响地图进行了模仿，从而加深对geopandas数据可视化的融会贯通，而本文作为geoplot篇的上半篇，介绍了geoplot中最基本的几种数据可视化API，使得我们足以完成较为基础的数据可视化作品，而同样为了加深对上文所介绍知识的理解掌握，接下来我们再次对其他优秀的数据可视化作品进行模仿。

　　这次我们要模仿的作品来自Github仓库https://github.com/Z3tt/30DayMapChallenge，是利用R进行地理空间数据可视化的一个集锦仓库，要用geoplot来模仿复现的作品如图15所示，展示了柏林所有电动汽车充电桩的分布情况：

图15

　　我们主要复现的是图15中柏林地图以及内部元素部分，使用到的数据在我的Github仓库对应本文路径下的Berlin文件夹中，其中ladesaeulen_bnetza_und_be_emobil.xlsx记录了EPSG:25833投影坐标系格式下的充电桩经纬度点信息，gis_osm_roads_free_1.shp记录了柏林市OSM路网信息，Bezirke__Berlin.shp记录了柏林行政区划信息。

　　在分析了原图的R代码之后，我们将整幅图拆解分为四个图层，1是柏林最边缘的灰色轮廓，这其实是整个柏林区域面数据向外生成缓冲区之后的效果；2是柏林各行政区区划，3是柏林内部的部分OSM路网，构成了图中依稀可见的类似纹路的要素，4是所有的充电桩点数据，即图中黄色的半透明散点，其中除路网线数据可视化以外的其他图层我们均使用geoplot来实现。

　　数据预处理部分分步骤代码较多，不便在文章中全盘放出，你可以到文章开头的Github仓库中对应路径下查看和下载，下面只贴出绘图部分的代码以方便理解思想：

# 绘制最底层柏林缓冲区
ax = gplt.polyplot(df=gpd.GeoSeries([GeometryCollection(berlin_area.geometry.tolist())], crs='EPSG:4326') 
                         .buffer(0.005, resolution=100),
                   projection=gcrs.WebMercator(),
                   facecolor='grey',
                   edgecolor='None',
                   alpha=0.6)

# 绘制柏林区划
ax = gplt.polyplot(df=berlin_area,
                   facecolor='black',
                   edgecolor='white',
                   alpha=0.65,
                   ax=ax)

# 利用geopandas绘制内部OSM路网
ax = intersect_roads.plot(ax=ax, 
                          linewidth=0.1,
                          edgecolor='black',
                          alpha=0.25)

# 绘制充电桩散点
ax = gplt.pointplot(df=df_emobil,
                    ax=ax,
                    extent=berlin_area.total_bounds,
                    color='#edc00d',
                    alpha=0.3,
                    linewidth=0.2,
                    s=4.5)

# 绘制充电桩中心点
ax = gplt.pointplot(df=df_emobil,
                    ax=ax,
                    extent=berlin_area.buffer(0.01).total_bounds,
                    color='#d29c14',
                    edgecolor='#6f603e',
                    linewidth=0.01,
                    s=0.7)

plt.savefig("图16.png", bbox_inches='tight', pad_inches=0, dpi=600)

　　最终得到的效果如图16所示：

图16

　　

　　以上就是本文的全部内容，我将在下一篇文章中继续与大家一起探讨学习geoplot中更高级的绘图API。如有疑问和意见，欢迎留言或在我的Github仓库中发起issues与我交流。