（数据科学学习手札77）基于geopandas的空间数据分析——文件IO

提示：扫一扫查出行
【扫一扫了解最新限行尾号】
复制提示

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

1 简介

在上一篇文章中我们对 geopandas 中的 坐标参考系 有了较为深入的学习，而在日常空间数据分析工作中矢量文件的读入和写出，是至关重要的环节。

作为 基于geopandas的空间数据分析 系列文章的第三篇，通过本文你将会学习到 geopandas 中的 文件IO 。

2 文件IO

2.1 矢量文件的读入

geopandas 将 fiona 作为操纵矢量数据读写功能的后端，使用 geopandas.read_file() 读取对应类型文件，而在后端实际上是使用 fiona.open 来读入数据，即两者参数是保持一致的，读入的数据自动转换为 GeoDataFrame ，下面是 geopandas.read_file() 主要参数：

filename ：str类型，传入文件对应的路径或url

layer ：str类型，当要读入的数据格式为地理数据库 .gdb 或 QGIS 中的 .gpkg 时，传入对应图层的名称

下面结合上述参数，来介绍一下使用 geopandas.read_file() 在不同情况下读取常见格式矢量数据的方法，使用到的示例数据为中国地图， CRS 为 EPSG:4326 ，本文使用到的所有数据都可以在文章开头提及的 Github 仓库对应本文路径下找到：

图1

2.1.1 shapefile

作为非常常见的一种矢量文件格式， geopandas 对 shapefile 提供了很好的读取和写出支持，下面分为不同情况来介绍：

完整的shapefile

如图2，这是一个完整的 shapefile ：

图2

使用 geopandas 来读取这种形式的 shapefile 很简单：

import geopandas as gpd

data = gpd.read_file('geometry/china_provinces/china_provinces.shp')

print(data.crs) # 查看数据对应的crs

data.head() # 查看前5行

图3

缺少投影的shapefile

当 shapefile 中缺失 .prj 文件时，使用 geopandas 读入后形成的 GeoDataFrame 会缺失 crs 属性：

图4

如果已经知道数据对应的 CRS ，可以在读入数据后补充上 crs 信息以进行其他操作：

import pyproj

data.crs = pyproj.CRS.from_user_input('EPSG:4326')

data.crs

图5

直接读取文件夹

当文件夹下只有单个shapefile时，可以直接读取该文件夹：

图6

读取zip压缩包中的文件

geopandas 通过传入特定语法格式的文件路径信息，以支持直接读取 .zip 格式压缩包中的 shapefile 文件，主要分为两种情况。

当文件在压缩包内的根目录时，使用下面的语法规则来读取数据：

zip://路径/xxx.zip

譬如我们要读取图7所示的压缩包内文件：

图7

按照对应的语法规则，读取该类型数据方式如下：

图8

而当文件在压缩包内的文件夹中时，如图9：

图9

使用下面的语法规则来读取数据：

zip://路径/xxx.zip!压缩包内指定文件路径

将上述语法运用到上述文件：

图10

2.1.2 gdb与gpkg

对于 Arcgis 中的地理数据库 gdb ，以及 QGIS 中的 GeoPackage ，要读取其包含的矢量数据，就要涉及到图层的概念，对应 geopandas.read_file() 的 layer 参数，只需要将 gdb 或 gpkg 文件路径作为 filename 参数，再将对应的图层名称作为 layer 参数传入：

gdb

data = gpd.read_file('geometry/china_provinces.gdb',

                     layer='china_provinces')

print(data.crs) # 查看数据对应的crs

data.head() # 查看前5行

图11

gpkg

类似读入 gdb 文件：

data = gpd.read_file('geometry/china_provinces.gpkg',

                     layer='china_provinces',

                     encoding='utf-8')

print(data.crs) # 查看数据对应的crs

data.head() # 查看前5行

图12

2.1.3 GeoJSON

作为web地图中最常使用的矢量数据格式， GeoJSON 几乎被所有在线地图框架作为数据源格式，在 geopandas 中读取 GeoJSON 非常简单，只需要传入文件路径名称即可，下面我们来读入图13所示的文件：

图13
图14

2.1.4 过滤

geopandas 在0.1.0版本中新增了 bbox过滤 ，在0.7.0版本中新增了 蒙版过滤 和 行过滤 功能，可以辅助我们根据自己的需要读入原始数据中的子集，下面一一进行介绍：

bbox过滤

bbox过滤 允许我们在 read_file() 中传入一个边界框作为参数 bbox ，格式为 (左下角x, 左下角y, 右上角x, 右上角y) ，这样在读入的过程中只会保留几何对象与 bbox 有相交的数据记录，下面我们仍然以上文中使用过的中国地图数据为例，我们在读入的过程中，传入边界框：

from shapely import geometry

data = gpd.read_file('geometry/china_provinces.json',

                    bbox=(100, 20, 110, 30))

%matplotlib widget

ax = data.plot()

# 绘制bbox框示意

ax = gpd.GeoSeries([geometry.box(minx=100,

                                 miny=20,

                                 maxx=110,

                                 maxy=30).boundary]).plot(ax=ax, color='red')

图15

可以看到只有跟红色框有相交的几何对象被读入。

蒙版过滤

蒙版过滤 和 bbox过滤 功能相似，都是筛选与指定区域相交的数据记录，不同的是 蒙版过滤 通过 mask 参数可以传入任意形状的多边形，不再像 bbox过滤 那样只接受矩形：

data = gpd.read_file('geometry/china_provinces.json',

                    mask=geometry.Polygon([(100, 20), (110, 30), (120, 20)]))

ax = data.plot()

# 绘制bbox框示意

ax = gpd.GeoSeries([geometry.Polygon([(100, 20),

                                      (110, 30),

                                      (120, 20)]).boundary]).plot(ax=ax, color='red')

图16

可以看到只有跟红色多边形相交的几何对象被读入。

行过滤

行过滤 的功能就比较简单，通过参数 rows 控制读入原数据的前若干行，可以用于在读取大型数据时先快速查看前几行以了解整个数据的格式：

图17

2.2 矢量文件的写出

在 geopandas 中使用 to_file() 来将 GeoDataFrame 或 GeoSeries 写出为矢量文件，主要支持 shapefile 、 GeoJSON 以及 GeoPackage ，不像 geopandas.read_file() 可以根据传入的文件名称信息自动推断类型，我们在写出矢量数据时就需要使用 driver 参数来声明文件类型：

ESRI Shapefile

我们将上文最后一次读入的 GeoDataFrame 写出为 ESRI Shapefile ，设置 driver 参数为 ESRI Shapefile ，如果你对文件编码有要求，这里可以使用 encoding 参数来指定，譬如这里我们指定为 utf-8 ：

'''在工程根目录下创建output文件夹'''

import os

try:

    os.mkdir('output')

except FileExistsError:

    pass

data.to_file('output/output.shp',

             driver='ESRI Shapefile',

             encoding='utf-8')

可以看到在output文件夹下，成功导出了完整的 shapefile ：

图18

而如果导出的文件名不加后缀扩展名，则会生成包含在新目录下的 shapefile ：

data.to_file('output/output_shapefile',

             driver='ESRI Shapefile',

             encoding='utf-8')

图19

也可以向指定的文件夹下追加图层：

data.to_file('output/output_shapefile_multi_layer',

             driver='ESRI Shapefile',

             layer='layer1',

             encoding='utf-8')

data.to_file('output/output_shapefile_multi_layer',

             driver='ESRI Shapefile',

             layer='layer2',

             encoding='utf-8')

data.to_file('output/output_shapefile_multi_layer',

             driver='ESRI Shapefile',

             layer='layer3',

             encoding='utf-8')

图20

GeoPackage

对于 gdb 文件，由于 ESRI 的限制，暂时无法在开源的 geopandas 中导出，但我们可以用 QGIS 中的 GeoPackage 作为替代方案（开源世界万岁O(∩_∩)O~~），只需要将 driver 参数设置为 GPKG 即可，这里需要注意一个bug：在使用 geopandas 导出 GeoPackage 文件时，可能会出现图21所示错误：

图21