作者：費(fèi)弗里

博客地址：

https://www.cnblogs.com/feffery/p/11898190.html

說明：本文經(jīng)作者授權(quán)轉(zhuǎn)載，禁止二次轉(zhuǎn)載

全文8500字

本文對(duì)應(yīng)代碼已上傳至我的Github倉(cāng)庫(kù)https://github.com/CNFeffery/DataScienceStudyNotes

1 簡(jiǎn)介

geopandas是建立在GEOS、GDAL、PROJ等開源地理空間計(jì)算相關(guān)框架之上的，類似pandas語法風(fēng)格的空間數(shù)據(jù)分析Python庫(kù)。

其目標(biāo)是盡可能地簡(jiǎn)化Python中的地理空間數(shù)據(jù)處理，減少對(duì)Arcgis、PostGIS等工具的依賴，使得處理地理空間數(shù)據(jù)變得更加高效簡(jiǎn)潔，打造純Python式的空間數(shù)據(jù)處理工作流。

本系列文章就將圍繞geopandas及其使用過程中涉及到的其他包進(jìn)行系統(tǒng)性的介紹說明，每一篇將盡可能全面具體地介紹geopandas對(duì)應(yīng)方面的知識(shí)。

計(jì)劃涵蓋geopandas的數(shù)據(jù)結(jié)構(gòu)、投影坐標(biāo)系管理、文件IO、基礎(chǔ)地圖制作、集合操作、空間連接與聚合。

作為基于geopandas的空間數(shù)據(jù)分析系列文章的第一篇，通過本文你將會(huì)學(xué)習(xí)到geopandas中的數(shù)據(jù)結(jié)構(gòu)。

geopandas的安裝和使用需要若干依賴包，如果不事先妥善安裝好這些依賴包而直接使用pip install geopandas或conda install geopandas，可能會(huì)引發(fā)依賴包相關(guān)錯(cuò)誤導(dǎo)致安裝失敗。

官方文檔中的推薦安裝方式為：

conda install --channel conda-forge geopandas

conda-forge是一個(gè)社區(qū)項(xiàng)目，在conda的基礎(chǔ)上提供了更廣泛更豐富的軟件資源包。通過它我們可以自動(dòng)下載安裝好所有g(shù)eopandas的必要依賴包而無需手動(dòng)繁瑣地去安裝它們。

在完成安裝后，下面我們開始對(duì)geopandas的系統(tǒng)性學(xué)習(xí)之旅。

2 數(shù)據(jù)結(jié)構(gòu)

geopandas作為pandas向地理分析計(jì)算方面的延拓，基礎(chǔ)的數(shù)據(jù)結(jié)構(gòu)延續(xù)了Series和DataFrame的特點(diǎn)，創(chuàng)造出GeoSeries與GeoDataFrame兩種基礎(chǔ)數(shù)據(jù)結(jié)構(gòu)：

2.1 GeoSeries

2.1.1 GeoSeries中的基礎(chǔ)幾何對(duì)象

與Series相似，GeoSeries用來表示一維向量，只不過這里的向量每個(gè)位置上的元素都表示著一個(gè)shapely中的幾何對(duì)象，有如下幾種類型：

• Points

對(duì)應(yīng)shapely.geometry中的Point，用于表示單個(gè)點(diǎn)，下面我們創(chuàng)建一個(gè)由若干Point對(duì)象組成的GeoSeries并像Series一樣定義索引：

from shapely import geometry import geopandas as gpd# 創(chuàng)建存放Point對(duì)象的GeoSeries # 這里shapely.geometry.Point(x, y)用于創(chuàng)建單個(gè)點(diǎn)對(duì)象 gpd.GeoSeries([geometry.Point(0, 0),geometry.Point(0, 1),geometry.Point(1, 1),geometry.Point(1, 0)],index=['a', 'b', 'c', 'd']) 圖1

可以看到創(chuàng)建出的GeoSeries數(shù)據(jù)類型為geometry，即幾何對(duì)象。

• MultiPoint

對(duì)應(yīng)shapely中的MultiPoint，用于表示多個(gè)點(diǎn)的集合，下面我們創(chuàng)建一個(gè)由若干MultiPoint對(duì)象組成的GeoSeries：

# 創(chuàng)建存放MultiPoint對(duì)象的GeoSeries # 這里shapely.geometry.MultiPoint([(x1, y1), (x2, y2), ...])用于創(chuàng)建多點(diǎn)集合 gpd.GeoSeries([geometry.MultiPoint([(0, 1), (1, 0)]),geometry.MultiPoint([(0, 0), (1, 1)])],index=['a', 'b']) 圖2

在jupyter notebook或jupyter lab中可以圖像的形式直接顯示GeoSeries中的單個(gè)元素：

圖3

• LineString

對(duì)應(yīng)shapely中的LineString，用于表示由多個(gè)點(diǎn)按順序連接而成的線。

下面我們創(chuàng)建一個(gè)由若干LineString對(duì)象組成的GeoSeries：

# 創(chuàng)建存放LineString對(duì)象的GeoSeries # 這里shapely.geometry.LineString([(x1, y1), (x2, y2), ...])用于創(chuàng)建多點(diǎn)按順序連接而成的線段 gpd.GeoSeries([geometry.LineString([(0, 0), (1, 1), (1, 0)]),geometry.LineString([(0, 0), (0, 1), (-1, 0)])],index=['a', 'b']) 圖4

同樣地，直接顯示第一個(gè)元素：

圖5

• MultiLineString

對(duì)應(yīng)shapely中的MultiLineString，用于表示多條線段的集合。

下面我們創(chuàng)建一個(gè)由若干MultiLineString對(duì)象組成的GeoSeries：

# 創(chuàng)建存放MultiLineString對(duì)象的GeoSeries # 這里shapely.geometry.MultiLineString([LineString1, LineString2])用于創(chuàng)建多條線段的集合 gpd.GeoSeries([geometry.MultiLineString([[(0, 0), (1, 1), (1, 0)],[(-0.5, 0), (0, 1), (-1, 0)]])],index=['a']) 圖6

同樣地，直接顯示第一個(gè)元素：

圖7

• Polygon(無孔)

geopandas中的Polygon對(duì)應(yīng)shapely中的Polygon，用于表示面，根據(jù)內(nèi)部有無孔洞可繼續(xù)細(xì)分。

下面我們創(chuàng)建一個(gè)由無孔Polygon對(duì)象組成的GeoSeries：

# 創(chuàng)建存放無孔Polygon對(duì)象的GeoSeries # 這里shapely.geometry.Polygon([(x1, y1), (x2, y2),...])用于創(chuàng)建無孔面 gpd.GeoSeries([geometry.Polygon([(0, 0), (0, 1), (1, 1), (1, 0)])],index=['a']) 圖8　　同樣地，直接顯示第一個(gè)元素：圖9

• Polygon(有孔)

區(qū)分于上文中的無孔Polygon，下面我們創(chuàng)建一個(gè)由有孔Polygon對(duì)象組成的GeoSeries：

# 創(chuàng)建存放有孔Polygon對(duì)象的GeoSeries # 這里shapely.geometry.Polygon(polygonExteriors, interiorCoords)用于創(chuàng)建有孔面 # 其中polygonExteriors用于定義整個(gè)有孔Polygon的外圍，是一個(gè)無孔的多邊形 # interiorCoords是用于定義內(nèi)部每個(gè)孔洞（本質(zhì)上是獨(dú)立的多邊形）的序列 gpd.GeoSeries([geometry.Polygon([(0,0),(10,0),(10,10),(0,10)],[((1,3),(5,3),(5,1),(1,1)),((9,9),(9,8),(8,8),(8,9))])])

同樣地，直接顯示第一個(gè)元素：

圖10

• MultiPolygon

對(duì)應(yīng)shapely中的MultiPolygon，用于表示多個(gè)面的集合。

下面我們創(chuàng)建一個(gè)由MultiPolygon對(duì)象組成的GeoSeries：

# 創(chuàng)建存放MultiPolygon對(duì)象的GeoSeries # 這里shapely.geometry.MultiPolygon([Polygon1, Polygon2])用于創(chuàng)建多個(gè)面的集合 gpd.GeoSeries([geometry.MultiPolygon([geometry.Polygon([(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)]),geometry.Polygon([(2, 2), (2, 3), (3, 3), (3, 2), (2, 2)])])],index=['a']) 圖11

顯示第一個(gè)元素：

圖12

• LinearRing

LinearRing對(duì)應(yīng)shapely.geometry中的LinearRing，是一種特殊的幾何對(duì)象。

可以理解為閉合的線或無孔多邊形的邊框，創(chuàng)建時(shí)傳入數(shù)據(jù)的格式與Polygon相同。

下面我們創(chuàng)建一個(gè)由LinearRing對(duì)象組成的GeoSeries：

# 創(chuàng)建存放LinearRing對(duì)象的GeoSeries # 這里shapely.geometry.LinearRing([(x1, y1), (x2, y2),...])用于創(chuàng)建LinearRing gpd.GeoSeries([geometry.LinearRing([(0, 0), (0, 1), (1, 1), (1, 0)])],index=['a']) 圖13

顯示第一個(gè)元素，可以看出LinearRing就是無孔多邊形的邊框線：

圖14

在同一個(gè)GeoSeries可以混合上述類型中的多種幾何對(duì)象，這意味著點(diǎn)線面在概念上相異的幾何對(duì)象可以共存于同一份數(shù)據(jù)中

2.1.2 GeoSeries常用屬性

類似pandas中的Series，GeoSeries在被創(chuàng)建完成之后也擁有很多實(shí)用的地理屬性，下面對(duì)其中較為常用的進(jìn)行列舉：

• area

area屬性返回與GeoSeries中每個(gè)元素一一對(duì)應(yīng)的面積值（這里的面積單位和下文涉及的長(zhǎng)度單位取決于投影坐標(biāo)系。

之后關(guān)于geopandas投影坐標(biāo)系管理的文章將會(huì)詳細(xì)介紹，這里僅做演示）：

# 創(chuàng)建混合點(diǎn)線面的GeoSeries，這里第5個(gè)有孔多邊形內(nèi)部空洞創(chuàng)建時(shí)使用[::-1]顛倒順序 # 是因?yàn)镚eoSeries.plot()方法繪制有孔多邊形的一個(gè)bug，即外部邊框與內(nèi)部孔洞創(chuàng)建時(shí)坐標(biāo) # 方向同為順時(shí)針或順時(shí)針時(shí)內(nèi)部孔洞會(huì)自動(dòng)被填充，如果你對(duì)這個(gè)bug感興趣，可以前往 # https://github.com/geopandas/geopandas/issues/951查看細(xì)節(jié) s = gpd.GeoSeries([geometry.Polygon([(0, 0), (0.5, 0.5), (1, 0), (0.5, -0.5)]),geometry.Polygon([(1, 1), (1.5, 1.5), (2, 1), (1.5, -1.5)]),geometry.Point(3, 3),geometry.LineString([(2, 2), (0, 3)]),geometry.Polygon([(4, 4), (8, 4), (8, 8), (4, 8)],[[(5, 5), (7, 5), (7, 7), (5, 7)][::-1]])])# 在jupyter中開啟matplotlib交互式繪圖模式 %matplotlib widget s.plot() # 對(duì)s進(jìn)行簡(jiǎn)單的可視化 圖15

可以看到，s中包含了多種幾何對(duì)象，下面直接得到s的面積：

圖16 計(jì)算GeoSeries面積

• bounds

bounds屬性返回每個(gè)幾何對(duì)象所在box左下角、右上角的坐標(biāo)信息：

圖17

• length

length屬性返回每個(gè)幾何對(duì)象邊長(zhǎng)：

圖18

• geom_type

geom_type返回每個(gè)幾何對(duì)象類型：

圖19

• exterior與interiors

對(duì)于多邊形對(duì)象，exterior返回LinearRing格式的外邊框線，對(duì)于有孔多邊形，interiors返回所有內(nèi)部孔洞LinearRing格式邊框線集合：

圖20

• is_valid

在shapely中涉及到很多拓?fù)溆?jì)算操作時(shí)，對(duì)幾何對(duì)象的合法性有要求。

譬如定義多邊形時(shí)坐標(biāo)按順序連線時(shí)穿過了之前定義的邊就屬于非法，因?yàn)間eopandas對(duì)矢量對(duì)象的計(jì)算依賴于shapely，于是引進(jìn)了屬性用于判斷每個(gè)幾何對(duì)象是否合法。

下面我們創(chuàng)建兩個(gè)形狀相同的多邊形，其中一個(gè)滿足上述所說的非法情況，另一個(gè)由兩個(gè)多邊形拼接而成：

s_ = gpd.GeoSeries([geometry.Polygon([(4, 0), (6, 1), (4, 1), (6, 0)]),geometry.MultiPolygon([geometry.Polygon([(4, 0), (5, 0.5), (6, 0)]),geometry.Polygon([(5, 0.5), (6, 1), (4, 1)])])])

從形狀上看兩者相同：

圖21

下面我們嘗試用shapely中的interp方法來取得這兩個(gè)幾何對(duì)象的相交部分，出現(xiàn)了拓?fù)溥壿嬪e(cuò)誤：

圖22

查看s_.is_valid，可以看出第一個(gè)自相交的多邊形非法：

圖23

• boundary

boundary返回每個(gè)幾何對(duì)象的低維簡(jiǎn)化表示（點(diǎn)對(duì)象無具體的更低維簡(jiǎn)化，故無返回值）：

圖24

• centroid

centroid返回每個(gè)幾何對(duì)象的重心（幾何中心）：

圖25

• convex_hull

convex_hull返回每個(gè)幾何對(duì)象的凸包，Polygon格式，即恰巧包含對(duì)應(yīng)幾何對(duì)象的凸多邊形：

import numpy as np# 利用獨(dú)立的正態(tài)分布隨機(jī)數(shù)創(chuàng)建兩個(gè)MultiPoint集合 s__ = gpd.GeoSeries([geometry.MultiPoint(np.random.normal(loc=0, scale=2, size=[10, 2]).tolist()),geometry.MultiPoint(np.random.normal(loc=5, scale=2, size=[10, 2]).tolist())])ax = s__.plot(color='red') # 繪制s__ s__.convex_hull.plot(ax=ax, alpha=0.4) # 疊加繪制各自對(duì)應(yīng)凸包，調(diào)低填充透明度以顯示更明顯 圖26

• envelope

envelope屬性返回對(duì)應(yīng)幾何對(duì)象的box范圍，Polygon格式，即包含對(duì)應(yīng)元素中所有點(diǎn)的最小矩形：

import numpy as np# 創(chuàng)建兩團(tuán)獨(dú)立的MultiPoint s__ = gpd.GeoSeries([geometry.MultiPoint(np.random.normal(loc=0, scale=2, size=[10, 2]).tolist()),geometry.MultiPoint(np.random.normal(loc=5, scale=2, size=[10, 2]).tolist())])ax = s__.plot(color='red') # 繪制s__ s__.envelope.plot(ax=ax, alpha=0.4) # 疊加繪制各自對(duì)應(yīng)envelope，調(diào)低填充透明度以顯示更明顯 圖27

2.2 GeoDataFrame

2.2.1 GeoDataFrame基礎(chǔ)

顧名思義，geopandas中的GeoDataFrame是在pandas.DataFrame的基礎(chǔ)上，加入空間分析相關(guān)內(nèi)容進(jìn)行改造而成。

其最大特點(diǎn)在于其在原有數(shù)據(jù)表格基礎(chǔ)上增加了一列GeoSeries使得其具有矢量性，所有對(duì)于GeoDataFrame施加的空間幾何操作也都作用在這列指定的幾何對(duì)象之上。

下面我們舉個(gè)簡(jiǎn)單的例子，基于不同均值和標(biāo)準(zhǔn)差的正態(tài)分布隨機(jī)數(shù)，創(chuàng)建GeoDataFrame來記錄這些信息：

contents = [(loc, 0.5) for loc in range(0, 10, 2)] geo_df = gpd.GeoDataFrame(data=contents,geometry=[geometry.MultiPoint(np.random.normal(loc=loc, scale=scale, size=[10, 2]).tolist())for loc, scale in contents],columns=['均值', '標(biāo)準(zhǔn)差']) geo_df 圖28

其中定義GeoDataFrame時(shí)作為每行所關(guān)聯(lián)幾何對(duì)象的GeoSeries需要通過geometry參數(shù)指定，而除了用上述的方式創(chuàng)建GeoDataFrame，先創(chuàng)建數(shù)據(jù)表，再添加矢量信息列亦可。

這時(shí)幾何對(duì)象列的名稱可以自由設(shè)置，但一定要利用GeoDataFrame.set_geometry()方法將后添加的矢量列指定為矢量主列。因?yàn)槊總€(gè)GeoDataFrame若在定義之處沒有指定矢量列，后將無法進(jìn)行與適量信息掛鉤的所有操作（GeoSeries所有屬性都可同樣作用于GeoDataFrame，因?yàn)樗锌臻g操作實(shí)際上都直接作用于其矢量主列）：

• 添加矢量列但未定義

geo_df = gpd.GeoDataFrame(contents, columns=['均值', '標(biāo)準(zhǔn)差']) geo_df['raw_points'] = [geometry.MultiPoint(np.random.normal(loc=loc, scale=scale, size=[10, 2]).tolist())for loc, scale in contents] # 嘗試查看矢量類型 geo_df.geom_type 圖29

這時(shí)所有直接針對(duì)GeoDataFrame的矢量相關(guān)操作都無法使用。

• 重新為GeoDataFrame指定矢量列

geo_df.set_geometry('raw_points').geom_type

這時(shí)相關(guān)操作可正常使用：

圖30

• 多個(gè)矢量列切換

通過前面的內(nèi)容，我們知道了每個(gè)GeoDataFrame都有一個(gè)矢量主列，相關(guān)操作例如繪圖都基于此列。

實(shí)際上GeoDataFrame允許表中存在多個(gè)矢量列，只要求任意時(shí)刻有且僅有1列為矢量主列即可。

因此我們可以在一個(gè)GeoDataFrame中保存多列矢量，需要用到哪列時(shí)再進(jìn)行切換即可，如下面的例子：

geo_df = gpd.GeoDataFrame(contents, columns=['均值', '標(biāo)準(zhǔn)差']) geo_df['raw_points'] = [geometry.MultiPoint(np.random.normal(loc=loc, scale=scale, size=[10, 2]).tolist())for loc, scale in contents] geo_df.set_geometry('raw_points', inplace=True) # inplace=True表示對(duì)原數(shù)據(jù)進(jìn)行更新# 繪制第一圖層 ax = geo_df.plot(color='red') geo_df['convex_hull'] = geo_df.convex_hull # 切換矢量主列 geo_df.set_geometry('convex_hull', inplace=True)# 繪制第二圖層 geo_df.plot(ax=ax, color='blue', alpha=0.4) 圖31

2.2.2 GeoDataFrame數(shù)據(jù)索引

作為pandas.DataFrame的延伸，GeoDataFrame同樣支持pandas.DataFrame中的.loc以及.iloc對(duì)數(shù)據(jù)在行、列尺度上進(jìn)行索引和篩選。

這里我們以geopandas自帶的世界地圖數(shù)據(jù)為例：

world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres')) world.plot() 圖32 geopandas自帶世界地圖

查看其表格內(nèi)容：

圖33

使用.loc+條件篩選選擇數(shù)據(jù)：

圖34

使用.iloc選擇數(shù)據(jù)：

圖35

而除了這些常規(guī)的數(shù)據(jù)索引方式之外，geopandas為GeoDataFrame添加了.cx索引方式，可以傳入所需的空間范圍，用于索引與傳入范圍相交的對(duì)應(yīng)數(shù)據(jù)：

# 選擇與東經(jīng)80度-110度，北緯0度-30度范圍相交的幾何對(duì)象 part_world = world.cx[80:110, 0:30]# 繪制第一圖層：世界地圖 ax = world.plot(alpha=0.05) # 繪制第二圖層：.cx所選擇的地區(qū) ax = part_world.plot(ax=ax, alpha=0.6) # 繪制第三圖層：.cx條件示意圖 ax = gpd.GeoSeries([geometry.box(minx=80, miny=0, maxx=110, maxy=30).boundary])\.plot(ax=ax, color='red')

示意圖如下：

圖36

放大到所選區(qū)域，可以看出正如前面所說，通過.cx，所有與指定空間范圍有重疊的對(duì)象都被選擇：

圖37

以上就是本文的全部?jī)?nèi)容，如有筆誤望指出，系列文章下一篇將詳細(xì)介紹geopandas中的投影坐標(biāo)系管理，敬請(qǐng)期待。

-END-

往期精彩回顧適合初學(xué)者入門人工智能的路線及資料下載機(jī)器學(xué)習(xí)在線手冊(cè)深度學(xué)習(xí)在線手冊(cè)AI基礎(chǔ)下載（pdf更新到25集）備注：加入本站微信群或者qq群，請(qǐng)回復(fù)“加群”獲取一折本站知識(shí)星球優(yōu)惠券，請(qǐng)回復(fù)“知識(shí)星球”喜歡文章，點(diǎn)個(gè)在看

總結(jié)

以上是生活随笔為你收集整理的Python地信专题 | 基于geopandas的空间数据分析—数据结构篇的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。