备注

互动在线版: Binder badge

GeoPandas简介#

本快速教程介绍了GeoPandas的主要概念和基本功能,以帮助您开始您的项目。

概念#

GeoPandas,顾名思义,扩展了科普数据科学库 pandas 通过添加对地理空间数据的支持。如果您不熟悉 pandas, we recommend taking a quick look at its Getting started documentation 在继续之前。

GeoPandas中的核心数据结构是 geopandas.GeoDataFrame ,是的子类 pandas.DataFrame ,它可以存储几何图形列并执行空间操作。这个 geopandas.GeoSeries ,是的子类 pandas.Series 处理几何图形。因此,您的 GeoDataFrame 是一种组合 pandas.Series ,以及传统数据(数字、布尔、文本等),以及 geopandas.GeoSeries ,使用几何图形(点、多边形等)。您可以拥有任意数量的带有几何图形的列;桌面GIS软件没有典型的限制。

地理数据框架方案

每个人 GeoSeries 可以包含任何几何类型(甚至可以在单个数组中混合它们),并且具有 GeoSeries.crs 属性,该属性存储有关投影的信息(CRS代表坐标参考系)。因此,每个 GeoSeries 在一个 GeoDataFrame 可以在不同的投影中,例如,允许您拥有同一几何体的多个版本(不同的投影)。

只有一个 GeoSeries 在一个 GeoDataFrame 被认为是 活动的 几何图形,这意味着所有几何运算都应用于 GeoDataFrame 对此进行操作 活动的 列。

用户指南

请参阅更多信息 data structures in the User Guide 。

让我们看看其中一些概念是如何在实践中发挥作用的。

读取和写入文件#

首先,我们需要读取一些数据。

正在读取文件#

假设您有一个同时包含数据和几何的文件(例如,GeoPackage、GeoJSON、shapefile),则可以使用 geopandas.read_file() ,它会自动检测文件类型并创建 GeoDataFrame 。本教程使用 "nybb" 数据集,纽约区地图,是GeoPandas安装的一部分。因此,我们使用 geopandas.datasets.get_path() 以检索数据集的路径。

[1]:
import geopandas

path_to_data = geopandas.datasets.get_path("nybb")
gdf = geopandas.read_file(path_to_data)

gdf
[1]:
BoroCode BoroName Shape_Leng Shape_Area geometry
0 5 Staten Island 330470.010332 1.623820e+09 MULTIPOLYGON (((970217.022 145643.332, 970227....
1 4 Queens 896344.047763 3.045213e+09 MULTIPOLYGON (((1029606.077 156073.814, 102957...
2 3 Brooklyn 741080.523166 1.937479e+09 MULTIPOLYGON (((1021176.479 151374.797, 102100...
3 1 Manhattan 359299.096471 6.364715e+08 MULTIPOLYGON (((981219.056 188655.316, 980940....
4 2 Bronx 464392.991824 1.186925e+09 MULTIPOLYGON (((1012821.806 229228.265, 101278...

正在写入文件#

要编写一个 GeoDataFrame 返回文件使用 GeoDataFrame.to_file() 。默认文件格式为shapefile,但您可以使用 driver 关键字。

[2]:
gdf.to_file("my_file.geojson", driver="GeoJSON")

用户指南

请参阅更多信息 reading and writing data in the User Guide 。

简单的访问器和方法#

现在我们有了我们的 GeoDataFrame 并且可以开始使用它的几何体。

由于纽约区数据集中只有一个几何列,因此该列自动成为 活动的 使用的几何和空间方法 GeoDataFrame 将应用于 "geometry" 列。

测量面积#

要测量每个多边形(或在此特定情况下为多个多边形)的面积,请访问 GeoDataFrame.area 属性,该属性返回一个 pandas.Series 。请注意, GeoDataFrame.area 只是 GeoSeries.area 应用于 活动的 几何图形列。

但首先,为了使结果更易于阅读,请将各行政区的名称设置为索引:

[3]:
gdf = gdf.set_index("BoroName")
[4]:
gdf["area"] = gdf.area
gdf["area"]
[4]:
BoroName
Staten Island    1.623822e+09
Queens           3.045214e+09
Brooklyn         1.937478e+09
Manhattan        6.364712e+08
Bronx            1.186926e+09
Name: area, dtype: float64

获取多边形边界和质心#

若要获取每个面(线串)的边界,请访问 GeoDataFrame.boundary :

[5]:
gdf['boundary'] = gdf.boundary
gdf['boundary']
[5]:
BoroName
Staten Island    MULTILINESTRING ((970217.022 145643.332, 97022...
Queens           MULTILINESTRING ((1029606.077 156073.814, 1029...
Brooklyn         MULTILINESTRING ((1021176.479 151374.797, 1021...
Manhattan        MULTILINESTRING ((981219.056 188655.316, 98094...
Bronx            MULTILINESTRING ((1012821.806 229228.265, 1012...
Name: boundary, dtype: geometry

由于我们已将边界保存为新列,因此现在同一列中有两个几何图形列 GeoDataFrame 。

我们还可以创建新的几何图形,其可以是例如原始几何图形的缓冲版本(即, GeoDataFrame.buffer(10) )或其质心:

[6]:
gdf['centroid'] = gdf.centroid
gdf['centroid']
[6]:
BoroName
Staten Island     POINT (941639.450 150931.991)
Queens           POINT (1034578.078 197116.604)
Brooklyn          POINT (998769.115 174169.761)
Manhattan         POINT (993336.965 222451.437)
Bronx            POINT (1021174.790 249937.980)
Name: centroid, dtype: geometry

测量距离#

我们还可以测量每个质心距离第一个质心位置有多远。

[7]:
first_point = gdf['centroid'].iloc[0]
gdf['distance'] = gdf['centroid'].distance(first_point)
gdf['distance']
[7]:
BoroName
Staten Island         0.000000
Queens           103781.535276
Brooklyn          61674.893421
Manhattan         88247.742789
Bronx            126996.283623
Name: distance, dtype: float64

请注意, geopandas.GeoDataFrame 是的子类 pandas.DataFrame ,因此我们拥有可在地理空间数据集上使用的所有PANAS功能-我们甚至可以一起使用属性和几何信息执行数据操作。

例如,要计算上面测量的距离的平均值,请访问‘Distance’列并对其调用Mean()方法:

[8]:
gdf['distance'].mean()
[8]:
76140.09102166798

制作地图#

GeoPandas还可以绘制地图,所以我们可以检查几何图形在太空中的显示方式。若要绘制活动几何图形,请调用 GeoDataFrame.plot() 。要按另一列进行颜色编码,请将该列作为第一个参数传入。在下面的示例中,我们绘制激活的几何图形列和颜色代码 "area" 纵队。我们也想展示一个传奇 (legend=True )。

[9]:
gdf.plot("area", legend=True)
[9]:
<AxesSubplot:>
../_images/getting_started_introduction_16_1.png

您还可以使用以交互方式浏览数据 GeoDataFrame.explore() ,它的行为方式与 plot() 而是返回交互式地图。

[10]:
gdf.explore("area", legend=False)
---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
File /usr/local/lib/python3.10/dist-packages/geopandas-0.10.2+79.g3abc6a7-py3.10.egg/geopandas/explore.py:271, in _explore(df, column, cmap, color, m, tiles, attr, tooltip, popup, highlight, categorical, legend, scheme, k, vmin, vmax, width, height, categories, classification_kwds, control_scale, marker_type, marker_kwds, style_kwds, highlight_kwds, missing_kwds, tooltip_kwds, popup_kwds, legend_kwds, map_kwds, **kwargs)
    270     import matplotlib.pyplot as plt
--> 271     from mapclassify import classify
    272 except (ImportError, ModuleNotFoundError):

ModuleNotFoundError: No module named 'mapclassify'

During handling of the above exception, another exception occurred:

ImportError                               Traceback (most recent call last)
Input In [10], in <cell line: 1>()
----> 1 gdf.explore("area", legend=False)

File /usr/local/lib/python3.10/dist-packages/geopandas-0.10.2+79.g3abc6a7-py3.10.egg/geopandas/geodataframe.py:1903, in GeoDataFrame.explore(self, *args, **kwargs)
   1900 @doc(_explore)
   1901 def explore(self, *args, **kwargs):
   1902     """Interactive map based on folium/leaflet.js"""
-> 1903     return _explore(self, *args, **kwargs)

File /usr/local/lib/python3.10/dist-packages/geopandas-0.10.2+79.g3abc6a7-py3.10.egg/geopandas/explore.py:273, in _explore(df, column, cmap, color, m, tiles, attr, tooltip, popup, highlight, categorical, legend, scheme, k, vmin, vmax, width, height, categories, classification_kwds, control_scale, marker_type, marker_kwds, style_kwds, highlight_kwds, missing_kwds, tooltip_kwds, popup_kwds, legend_kwds, map_kwds, **kwargs)
    271     from mapclassify import classify
    272 except (ImportError, ModuleNotFoundError):
--> 273     raise ImportError(
    274         "The 'folium', 'matplotlib' and 'mapclassify' packages are required for "
    275         "'explore()'. You can install them using "
    276         "'conda install -c conda-forge folium matplotlib mapclassify' "
    277         "or 'pip install folium matplotlib mapclassify'."
    278     )
    280 # xyservices is an optional dependency
    281 try:

ImportError: The 'folium', 'matplotlib' and 'mapclassify' packages are required for 'explore()'. You can install them using 'conda install -c conda-forge folium matplotlib mapclassify' or 'pip install folium matplotlib mapclassify'.

切换激活的几何体 (GeoDataFrame.set_geometry )到质心,我们可以使用点几何绘制相同的数据。

[11]:
gdf = gdf.set_geometry("centroid")
gdf.plot("area", legend=True)
[11]:
<AxesSubplot:>
../_images/getting_started_introduction_20_1.png

我们还可以将这两个层次 GeoSeries 在彼此的顶端。我们只需要使用一个地块作为另一个地块的轴。

[12]:
ax = gdf["geometry"].plot()
gdf["centroid"].plot(ax=ax, color="black")
[12]:
<AxesSubplot:>
../_images/getting_started_introduction_22_1.png

现在,我们将活动几何图形设置回原始状态 GeoSeries 。

[13]:
gdf = gdf.set_geometry("geometry")

用户指南

请参阅更多信息 mapping in the User Guide 。

几何图形创建#

我们可以进一步处理几何体,并在已有的基础上创建新的形状。

凸壳#

如果我们对面的凸包感兴趣,我们可以访问 GeoDataFrame.convex_hull 。

[14]:
gdf["convex_hull"] = gdf.convex_hull
[15]:
ax = gdf["convex_hull"].plot(alpha=.5)  # saving the first plot as an axis and setting alpha (transparency) to 0.5
gdf["boundary"].plot(ax=ax, color="white", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5
[15]:
<AxesSubplot:>
../_images/getting_started_introduction_27_1.png

缓冲层#

在其他情况下,我们可能需要使用缓冲几何 GeoDataFrame.buffer() 。几何方法会自动应用于活动几何,但我们可以将它们直接应用于任何 GeoSeries 也是。让我们缓冲这些行政区和它们的质心,并在彼此的顶部进行绘制。

[16]:
# buffering the active geometry by 10 000 feet (geometry is already in feet)
gdf["buffered"] = gdf.buffer(10000)

# buffering the centroid geometry by 10 000 feet (geometry is already in feet)
gdf["buffered_centroid"] = gdf["centroid"].buffer(10000)
[17]:
ax = gdf["buffered"].plot(alpha=.5)  # saving the first plot as an axis and setting alpha (transparency) to 0.5
gdf["buffered_centroid"].plot(ax=ax, color="red", alpha=.5)  # passing the first plot as an axis to the second
gdf["boundary"].plot(ax=ax, color="white", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5
[17]:
<AxesSubplot:>
../_images/getting_started_introduction_30_1.png

用户指南

请参阅更多信息 geometry creation and manipulation in the User Guide 。

几何关系#

我们还可以询问不同几何图形的空间关系。使用上面的几何图形,我们可以检查哪些缓冲的行政区与布鲁克林的原始几何图形相交,即距离布鲁克林不到10,000英尺。

首先,我们得到布鲁克林的一个多边形。

[18]:
brooklyn = gdf.loc["Brooklyn", "geometry"]
brooklyn
[18]:
../_images/getting_started_introduction_32_0.svg

该多边形是一个 shapely geometry object ,与GeoPandas中使用的任何其他几何体一样。

[19]:
type(brooklyn)
[19]:
shapely.geometry.multipolygon.MultiPolygon

然后我们可以检查其中的哪些几何图形 gdf["buffered"] 与其相交。

[20]:
gdf["buffered"].intersects(brooklyn)
[20]:
BoroName
Staten Island     True
Queens            True
Brooklyn          True
Manhattan         True
Bronx            False
dtype: bool

只有布朗克斯(北部)距离布鲁克林超过10000英尺。所有其他的都更近,并与我们的多边形相交。

或者,我们可以检查哪些缓冲质心完全在原始行政区的多边形内。在这种情况下,两者都 GeoSeries 对齐,并对每行执行检查。

[21]:
gdf["within"] = gdf["buffered_centroid"].within(gdf)
gdf["within"]
[21]:
BoroName
Staten Island     True
Queens            True
Brooklyn         False
Manhattan        False
Bronx            False
Name: within, dtype: bool

我们可以把结果画在地图上,以确认这一发现。

[22]:
gdf = gdf.set_geometry("buffered_centroid")
ax = gdf.plot("within", legend=True, categorical=True, legend_kwds={'loc': "upper left"})  # using categorical plot and setting the position of the legend
gdf["boundary"].plot(ax=ax, color="black", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5
[22]:
<AxesSubplot:>
../_images/getting_started_introduction_40_1.png

推算#

每个人 GeoSeries 其坐标参考系(CRS)可在 GeoSeries.crs 。CRS会告诉GeoPandas这些几何图形的坐标在地球表面的位置。在某些情况下,CRS是地理上的,这意味着坐标在纬度和经度上。在这些情况下,其CRS是WGS84,具有授权代码 EPSG:4326 。让我们看看我们纽约各区的投影 GeoDataFrame 。

[23]:
gdf.crs
[23]:
<Derived Projected CRS: EPSG:2263>
Name: NAD83 / New York Long Island (ftUS)
Axis Info [cartesian]:
- X[east]: Easting (US survey foot)
- Y[north]: Northing (US survey foot)
Area of Use:
- name: United States (USA) - New York - counties of Bronx; Kings; Nassau; New York; Queens; Richmond; Suffolk.
- bounds: (-74.26, 40.47, -71.8, 41.3)
Coordinate Operation:
- name: SPCS83 New York Long Island zone (US Survey feet)
- method: Lambert Conic Conformal (2SP)
Datum: North American Datum 1983
- Ellipsoid: GRS 1980
- Prime Meridian: Greenwich

几何图形已加入 EPSG:2263 以英尺为单位的坐标。我们可以很容易地重新投射一个 GeoSeries 到另一个CRS,比如 EPSG:4326 使用 GeoSeries.to_crs() 。

[24]:
gdf = gdf.set_geometry("geometry")
boroughs_4326 = gdf.to_crs("EPSG:4326")
boroughs_4326.plot()
[24]:
<AxesSubplot:>
../_images/getting_started_introduction_44_1.png
[25]:
boroughs_4326.crs
[25]:
<Geographic 2D CRS: EPSG:4326>
Name: WGS 84
Axis Info [ellipsoidal]:
- Lat[north]: Geodetic latitude (degree)
- Lon[east]: Geodetic longitude (degree)
Area of Use:
- name: World.
- bounds: (-180.0, -90.0, 180.0, 90.0)
Datum: World Geodetic System 1984 ensemble
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich

请注意沿曲线图轴线的坐标差异。以前我们有12万到28万英尺,现在我们有40.5到40.9度。在这种情况下, boroughs_4326 有一个 "geometry" WGS84中的列,但其他所有列(带有质心等)保留在原始的CRS中。

警告

对于依赖于距离或面积的操作,您始终需要使用投影CRS(以米、英尺、公里等为单位)不是地理上的(以度计)。GeoPandas的操作是平面的,而度数反映了球体上的位置。因此,使用度的空间运算可能不会产生正确的结果。例如,结果为 gdf.area.sum() (预计CRS)为8 429 911 572 ft2,但结果为 boroughs_4326.area.sum() (地理CRS)为0.083。

用户指南

请参阅更多信息 projections in the User Guide 。

接下来是什么?#

有了GeoPandas,我们可以做的比目前为止引入的要多得多,从 aggregations ,至 spatial joins ,至 geocoding ,以及 much more 。

去那边的 User Guide 要了解有关GeoPandas的不同功能的更多信息,请参阅 Examples 查看如何使用它们,或者到 API reference 了解详情。