本文主要介紹GeoPandas結合matplotlib實現地圖的基礎視覺化。GeoPandas是一個Python開源專案,旨在提供豐富而簡單的地理空間資料處理介面。GeoPandas擴充套件了Pandas的資料型別,並使用matplotlib進行繪圖。GeoPandas官方倉庫地址為:GeoPandas。GeoPandas的官方檔案地址為:GeoPandas-doc。關於GeoPandas的使用見[資料分析與視覺化] Python繪製資料地圖1-GeoPandas入門指北。
本文所有程式碼見:Python-Study-Notes
GeoPandas推薦使用Python3.7版本及以上,執行環境最好是linux系統。GeoPandas安裝命令如下:
pip install geopandas
如果上述命令安裝出問題,則推薦使用conda安裝GeoPandas,命令如下:
conda install geopandas
或:
conda install --channel conda-forge geopandas
# jupyter notebook環境去除warning
import warnings
warnings.filterwarnings("ignore")
# 檢視geopandas版本
import geopandas as gpd
gpd.__version__
'0.10.2'
GeoPandas基於matplotlib庫封裝高階介面用於製作地圖,GeoSeries或GeoDataFrame物件都提供了plot函數以進行物件視覺化。與GeoSeries物件相比,GeoDataFrame物件提供的plot函數在引數上更為複雜,也更為常用。
GeoDataFrame物件提供的plot函數的常用輸入引數如下:
def GeoDataFrame.plot(
column: str, np.array, pd.Series (default None), # 用於繪圖的列名或資料列
kind: str, # 繪圖型別
cmap: str, # matplotlib的顏色圖Colormaps
color: str, np.array, pd.Series (default None), # 指定所有繪圖物件的統一顏色
ax: matplotlib.pyplot.Artist (default None), # 指定matplotlib的繪圖軸
cax: matplotlib.pyplot Artist (default None), # 設定圖例的座標軸
categorical: bool (default False), # 是否按照類別設定物件顏色
legend: bool (default False), # 是否顯示圖例,如果column或color引數未賦值,則此引數無效
scheme: str (default None), # 用於控制分層設色
k:int (default 5), # scheme的層次數
vmin:None or float (default None), # 圖例cmap的最小值
vmax:None or float (default None), # 圖例cmap的最大值
markersize:str or float or sequence (default None), # 繪圖點的大小
figsize: tuple of integers (default None), # 用於控制matplotlib.figure.Figure
legend_kwds: dict (default None), # matplotlib圖例引數
missing_kw: dsdict (default None), # 缺失值區域繪製引數
aspect:‘auto’, ‘equal’, None or float (default ‘auto’), # 設定繪圖比例
**style_kwds: dict, # 其他引數,如物件邊界色edgecolor, 物件填充色facecolor, 邊界寬linewidth,透明度alpha
)->ax: matplotlib axes instance
GeoSeries物件提供的plot函數的常用輸入引數如下:
def GeoSeries.plot(
s: Series, # GeoSeries物件
cmap: str (default None),
color: str, np.array, pd.Series, List (default None),
ax: matplotlib.pyplot.Artist (default None),
figsize: pair of floats (default None),
aspect: ‘auto’, ‘equal’, None or float (default ‘auto’),
**style_kwds: dict
)->ax: matplotlib axes instance
想要更好使用以上引數最好熟悉matplotlib,如有不懂可以查閱matplotlib檔案。下面簡要介紹GeoPandas中繪圖函數的使用。
讀取資料集
import geopandas as gpd
# 讀取自帶世界地圖資料集
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))
# 人口,大洲,國名,國家縮寫,ISO國家程式碼,gdp,地理位置資料
world.head()
pop_est | continent | name | iso_a3 | gdp_md_est | geometry | |
---|---|---|---|---|---|---|
0 | 920938 | Oceania | Fiji | FJI | 8374.0 | MULTIPOLYGON (((180.00000 -16.06713, 180.00000... |
1 | 53950935 | Africa | Tanzania | TZA | 150600.0 | POLYGON ((33.90371 -0.95000, 34.07262 -1.05982... |
2 | 603253 | Africa | W. Sahara | ESH | 906.5 | POLYGON ((-8.66559 27.65643, -8.66512 27.58948... |
3 | 35623680 | North America | Canada | CAN | 1674000.0 | MULTIPOLYGON (((-122.84000 49.00000, -122.9742... |
4 | 326625791 | North America | United States of America | USA | 18560000.0 | MULTIPOLYGON (((-122.84000 49.00000, -120.0000... |
world.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b270fe10>
# 讀取自帶世界各大城市資料集
cities = gpd.read_file(gpd.datasets.get_path('naturalearth_cities'))
# 名字,地理位置資料
cities.head()
name | geometry | |
---|---|---|
0 | Vatican City | POINT (12.45339 41.90328) |
1 | San Marino | POINT (12.44177 43.93610) |
2 | Vaduz | POINT (9.51667 47.13372) |
3 | Luxembourg | POINT (6.13000 49.61166) |
4 | Palikir | POINT (158.14997 6.91664) |
分割區統計圖
GeoPandas可以輕鬆建立分割區統計圖(每個形狀的顏色基於關聯變數值的地圖)。只需在plot函數中將引數column設定為要用於指定顏色的列。
# 去除南極地區
world = world[(world.pop_est>0) & (world.name!="Antarctica")]
# 計算人均gpd
world['gdp_per_cap'] = world.gdp_md_est / world.pop_est
# 繪圖
world.plot(column='gdp_per_cap')
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b241f390>
圖例設定
legend=True即可展示圖例。
import matplotlib.pyplot as plt
fig, ax = plt.subplots(1, 1)
# 根據人口繪圖數繪製
world.plot(column='pop_est', ax=ax, legend=True)
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b23585d0>
但是我們會發現繪圖區域和圖例區域不對齊,我們通過mpl_toolkits設定圖例軸,並通過cax引數直接控制plot函數的圖例繪製。
from mpl_toolkits.axes_grid1 import make_axes_locatable
fig, ax = plt.subplots(1, 1)
divider = make_axes_locatable(ax)
# 設定圖例,right表示位置,size=5%表示圖例寬度,pad表示圖例離圖片間距
cax = divider.append_axes("right", size="5%", pad=0.2)
world.plot(column='pop_est', ax=ax, legend=True, cax=cax)
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b21de650>
此外我們可以通過legend_kwds設定圖例資訊和方向,就像matplotlib那樣。
import matplotlib.pyplot as plt
fig, ax = plt.subplots(1, 1)
world.plot(column='pop_est',
ax=ax,
legend=True,
legend_kwds={'label': "Population by Country and Area",
'orientation': "horizontal"})
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b20ec110>
顏色設定
GeoPandas的繪圖函數中通過cmap引數設定顏色圖,具體cmap引數的設定可以見matplotlib顏色圖選擇教學。
# 設定顏色圖為tab20
world.plot(column='gdp_per_cap', cmap='tab20')
<matplotlib.axes._subplots.AxesSubplot at 0x7f17b205a490>
如果只想繪製邊界,可以呼叫boundary屬性進行繪製。
world.boundary.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f15ac1cecd0>
缺失資料繪製
當繪製項某些地圖的資料缺失時,GeoPandas會自動放棄這些區域資料的繪製。如下所示。將Africa地區的gdp_per_cap資料設定為確實,那麼將不在地圖上繪製Africa地區的圖形。
import numpy as np
world.loc[world[world.continent=="Africa"].index, 'gdp_per_cap'] = np.nan
world.plot(column='gdp_per_cap')
<matplotlib.axes._subplots.AxesSubplot at 0x7f15ac0afe50>
如果我們想要在地圖上展示缺失值的地區,可以通過missing_kwds設定相應的引數。如下所示:
# 將缺失值地區顏色設定為lightgrey
world.plot(column='gdp_per_cap', missing_kwds={'color': 'lightgrey'})
<matplotlib.axes._subplots.AxesSubplot at 0x7f15ac094250>
# 更復雜的例子,hatch表示內部影象樣式
ax = world.plot(column='gdp_per_cap', missing_kwds={'color': 'lightgrey','edgecolor':'red',"hatch": "///","label": "Missing values"})
# 不顯示座標軸
ax.set_axis_off();
圖層設定
當我們需要疊加多張地圖或多個資料的結果,則需要用到地圖的圖層。GeoPandas提供了兩種方式進行圖層設定,但是要注意的是在合併地圖之前,始終確保它們共用一個共同的座標參考系crs。如將world地圖和cities地圖合併繪製。
# 以world資料集為準,對齊crs
cities.plot(marker='*', color='green', markersize=5)
cities = cities.to_crs(world.crs)
第一種方式最為簡單,直接通過ax引數控制繪圖軸來實現圖層疊加。
# base為繪圖軸
base = world.plot(color='white', edgecolor='black')
cities.plot(ax=base, marker='o', color='red', markersize=5);
第二種方式最為靈活,即通過matplotlib控制繪圖。如果想要繪製更加美觀的地圖,這種方式更加推薦,但是需要了解matplotlib的使用。
import matplotlib.pyplot as plt
# 設定繪圖軸
fig, ax = plt.subplots()
ax.set_aspect('equal')
world.plot(ax=ax, color='white', edgecolor='black')
cities.plot(ax=ax, marker='o', color='red', markersize=5)
plt.show()
圖層順序控制
如果我們需要控制地圖圖層的展示順序,即哪張圖片顯示在前。一種辦法是通過繪圖順序設定,後繪製的資料的圖層順序越靠上。一種是通過zorder引數設定圖層順序,zorder值越大,表示圖層順序越靠上。如下所示:
# cities先繪製,則圖層順序更靠下,導致一些資料點被world圖層遮蓋。
ax = cities.plot(color='k')
world.plot(ax=ax);
# 調整繪圖順序,完整顯示cities資料點
ax = world.plot()
cities.plot(ax=ax,color='k')
<matplotlib.axes._subplots.AxesSubplot at 0x7f159d535f10>
# 通過zorder引數設定圖層順序,完整顯示cities資料點
ax = cities.plot(color='k',zorder=2)
world.plot(ax=ax,zorder=1)
<matplotlib.axes._subplots.AxesSubplot at 0x7f159d4e8a10>
Pandas Plots
GeoPandas的Plot函數還允許設定kind引數繪製不同型別的圖形,kind預設為geo(地圖),其他可選引數與pandas預設提供的繪圖函數一致。pandas預設提供的繪圖函數見:Chart visualization。
gdf = world.head(10)
gdf.plot(kind='scatter', x="pop_est", y="gdp_md_est")
<matplotlib.axes._subplots.AxesSubplot at 0x7f159d5d69d0>
該範例主要參考:
1 讀取地圖資料
import geopandas as gpd
import matplotlib.pyplot as plt
# 讀取中國地圖資料,資料來自DataV.GeoAtlas,將其投影到EPSG:4573
data = gpd.read_file('https://geo.datav.aliyun.com/areas_v3/bound/100000_full.json').to_crs('EPSG:4573')
# 中國有34個省級行政區,最後一條資料為南海九段線
data.shape
(35, 8)
# 儲存各個省級行政區的面積,單位萬平方公里
data['area'] = data.area/1e6/1e4
# 中國疆域總面積,和官方資料會有差距
data['area'].sum()
989.1500769770855
# 檢視最後五條資料
data.tail()
adcode | name | adchar | childrenNum | level | parent | subFeatureIndex | geometry | area | |
---|---|---|---|---|---|---|---|---|---|
30 | 650000 | 新疆維吾爾自治區 | NaN | 24.0 | province | {'adcode': 100000} | 30.0 | MULTIPOLYGON (((17794320.693 4768675.631, 1779... | 175.796459 |
31 | 710000 | 臺灣省 | NaN | 0.0 | province | {'adcode': 100000} | 31.0 | MULTIPOLYGON (((20103629.026 2566628.808, 2011... | 4.090259 |
32 | 810000 | 香港特別行政區 | NaN | 18.0 | province | {'adcode': 100000} | 32.0 | MULTIPOLYGON (((19432072.340 2517924.724, 1942... | 0.125928 |
33 | 820000 | 澳門特別行政區 | NaN | 8.0 | province | {'adcode': 100000} | 33.0 | MULTIPOLYGON (((19385068.278 2470740.132, 1939... | 0.004702 |
34 | 100000_JD | JD | NaN | NaN | NaN | NaN | MULTIPOLYGON (((20309263.208 2708129.155, 2032... | 0.280824 |
# 拆分資料
nine_dotted_line = data.iloc[-1]
data = data[:-1]
nine_dotted_line
adcode 100000_JD
name
adchar JD
childrenNum NaN
level NaN
parent NaN
subFeatureIndex NaN
geometry MULTIPOLYGON (((20309263.208229765 2708129.154...
area 0.280824
Name: 34, dtype: object
2 地圖繪製
# 建立畫布matplotlib
fig, ax = plt.subplots(figsize=(12, 9))
# 繪製主要區域
ax = data.plot(ax=ax)
# 繪製九段線
ax = gpd.GeoSeries(nine_dotted_line.geometry).plot(ax=ax,edgecolor='red',linewidth=3)
# 儲存結果
fig.savefig('res.png', dpi=300, bbox_inches='tight')
3 繪圖自定義
# 建立畫布matplotlib
fig, ax = plt.subplots(figsize=(12, 9))
# 繪製主要區域
ax = data.plot(ax=ax,facecolor='grey',edgecolor='lightgrey',alpha=0.8,linewidth=1)
# 繪製九段線
ax = gpd.GeoSeries(nine_dotted_line.geometry).plot(ax=ax,edgecolor='grey',linewidth=3)
# 強調首都
ax = data[data.name=="北京市"].representative_point().plot(ax=ax, facecolor='red',marker='*', markersize=150)
# 強調港澳臺
# 設定字型
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "blue",'weight': 'bold'}
for index in data[data.adcode.isin(['710000','810000','820000'])].index:
if data.iloc[index]['name'] == "臺灣省":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "臺灣省"
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
elif data.iloc[index]['name'] == "香港特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "香港"
ax.text(x, y, name, ha="left", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
elif data.iloc[index]['name'] == "澳門特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "澳門"
ax.text(x, y, name, ha="right", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
# 移除座標軸
ax.axis('off')
(15508157.566510478, 20969512.421140186, 92287.20109282521, 6384866.68627745)
4 圖例設定
這裡採用單獨繪製圖例的方式來建立圖例,需要對matplotlib使用有一定的瞭解。
# 建立畫布matplotlib
fig, ax = plt.subplots(figsize=(12, 9))
# 繪製主要區域
ax = data.plot(ax=ax,facecolor='grey',edgecolor='lightgrey',alpha=0.8,linewidth=1)
# 繪製九段線
ax = gpd.GeoSeries(nine_dotted_line.geometry).plot(ax=ax,edgecolor='grey',linewidth=3)
# 強調首都
ax = data[data.name=="北京市"].representative_point().plot(ax=ax, facecolor='red',marker='*', markersize=150)
# 強調港澳臺
# 設定字型
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "blue",'weight': 'bold'}
# 這一段程式碼可能因為不同matplotlib版本出現不同結果
for index in data[data.adcode.isin(['710000','810000','820000'])].index:
if data.iloc[index]['name'] == "臺灣省":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "臺灣省"
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
elif data.iloc[index]['name'] == "香港特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "香港"
ax.text(x, y, name, ha="left", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
elif data.iloc[index]['name'] == "澳門特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "澳門"
ax.text(x, y, name, ha="right", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
# 移除座標軸
ax.axis('off')
# 單獨繪製圖例
plt.rcParams["font.family"] = 'FZSongYi-Z13S'
ax.scatter([], [], c='red', s=80, marker='*', label='首都')
ax.plot([], [], c='grey',linewidth=3, label='南海九段線')
# 設定圖例順序
handles, labels = ax.get_legend_handles_labels()
ax.legend(handles[::-1], labels[::-1], title="圖例",frameon=True, shadow=True, loc="lower left",fontsize=10)
# ax.legend(title="圖例",frameon=True, shadow=True, loc="lower left",fontsize=10)
<matplotlib.legend.Legend at 0x7f159cbe7510>
5 小地圖繪製
在很多中國地圖中,南海諸島區域都是在右下角的小地圖單獨繪製。在GeoPandas中想要實現這一功能,可以通過matplotlib的add_axes函數實現。add_axes主要功能為為新增繪圖子區域,該區域可以位於畫布中任意區域,且可設定任意大小。add_axes輸入引數為(left, bottom, width, height),left, bottom表示相對畫布的比例座標,width和height表示相對畫布的比例長寬。
這種繪圖方式非常不專業,也不推薦,建議只是學習使用思路。具體使用看如下程式碼。
# 建立畫布
fig = plt.figure(figsize=(6,3))
# 建立一個填充整個畫布的子圖
ax = fig.add_axes((0,0,1,1))
# 從畫布寬40%,高40%處繪製子圖
ax_child = fig.add_axes((0.4,0.4,0.2,0.2))
具體新增小地圖首先確定中國大陸區域範圍和南海範圍,然後分開繪製。
from shapely.geometry import Point
# 設中國大陸區域範圍和南海範圍,估計得到
bound = gpd.GeoDataFrame({
'x': [80, 140, 106.5, 123],
'y': [15, 50, 2.8, 24.5]
})
bound.geometry = bound.apply(lambda row: Point([row['x'], row['y']]), axis=1)
# 初始化CRS
bound.crs = 'EPSG:4326'
bound = bound.to_crs('EPSG:4573')
bound
x | y | geometry | |
---|---|---|---|
0 | 80.0 | 15.0 | POINT (15734050.166 1822879.627) |
1 | 140.0 | 50.0 | POINT (20972709.719 6154280.305) |
2 | 106.5 | 2.8 | POINT (18666803.134 309722.667) |
3 | 123.0 | 24.5 | POINT (20344370.888 2833592.614) |
# 建立畫布matplotlib
fig, ax = plt.subplots(figsize=(12, 9))
# 繪製主要區域
ax = data.plot(ax=ax,facecolor='grey',edgecolor='lightgrey',alpha=0.8,linewidth=1)
# 繪製九段線
ax = gpd.GeoSeries(nine_dotted_line.geometry).plot(ax=ax,edgecolor='grey',linewidth=3)
# 強調首都
ax = data[data.name=="北京市"].representative_point().plot(ax=ax, facecolor='red',marker='*', markersize=150)
# 強調港澳臺
# 設定字型
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "blue",'weight': 'bold'}
# 這一段程式碼可能因為不同matplotlib版本出現不同結果
for index in data[data.adcode.isin(['710000','810000','820000'])].index:
if data.iloc[index]['name'] == "臺灣省":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "臺灣省"
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
elif data.iloc[index]['name'] == "香港特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "香港"
ax.text(x, y, name, ha="left", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
elif data.iloc[index]['name'] == "澳門特別行政區":
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = "澳門"
ax.text(x, y, name, ha="right", va="top", fontdict=fontdict)
gpd.GeoSeries(data.iloc[index].geometry.centroid).plot(ax=ax, facecolor='black', markersize=5)
# 設定大陸區域範圍
ax.set_xlim(bound.geometry[0].x, bound.geometry[1].x)
ax.set_ylim(bound.geometry[0].y, bound.geometry[1].y)
# 移除座標軸
ax.axis('off')
# 單獨繪製圖例
plt.rcParams["font.family"] = 'FZSongYi-Z13S'
ax.scatter([], [], c='red', s=80, marker='*', label='首都')
ax.plot([], [], c='grey',linewidth=3, label='南海九段線')
# 設定圖例順序
handles, labels = ax.get_legend_handles_labels()
ax.legend(handles[::-1], labels[::-1], title="圖例",frameon=True, shadow=True, loc="lower left",fontsize=10)
# ax.legend(title="圖例",frameon=True, shadow=True, loc="lower left",fontsize=10)
# 建立南海插圖對應的子圖,調整這些引數以調整地圖位置
ax_child = fig.add_axes([0.75, 0.15, 0.2, 0.2])
ax_child = data.plot(ax=ax_child,facecolor='grey',edgecolor='lightgrey',alpha=0.8,linewidth=1)
# 繪製九段線
ax_child = gpd.GeoSeries(nine_dotted_line.geometry).plot(ax=ax_child,edgecolor='grey',linewidth=3)
# 設定子圖顯示範圍
ax_child.set_xlim(bound.geometry[2].x, bound.geometry[3].x)
ax_child.set_ylim(bound.geometry[2].y, bound.geometry[3].y)
ax_child.text(0.98,0.02,'Produced by luohenyueji',transform = ax.transAxes,
ha='center', va='center',fontsize = 12,color='black')
# 移除子圖座標軸
ax_child.set_xticks([])
ax_child.set_yticks([])
fig.savefig('res.png', dpi=300, bbox_inches='tight')
如下程式碼所示,繪製江蘇省地級市GDP地圖。
# 讀取2019江蘇省各市GDP資料
import geopandas as gpd
import matplotlib.pyplot as plt
import pandas as pd
plt.rcParams["font.family"] = 'FZSongYi-Z13S'
# 資料來自網際網路
gdp = pd.read_csv("2022江蘇省各市GDP.csv")
gdp
排行 | 地級市 | 2022年GDP(億元) | |
---|---|---|---|
0 | 1 | 蘇州市 | 23958.3 |
1 | 2 | 南京市 | 16907.9 |
2 | 3 | 無錫市 | 14850.8 |
3 | 4 | 南通市 | 11379.6 |
4 | 5 | 常州市 | 9550.1 |
5 | 6 | 徐州市 | 8457.8 |
6 | 7 | 鹽城市 | 7079.8 |
7 | 8 | 揚州市 | 6696.4 |
8 | 9 | 泰州市 | 6401.8 |
9 | 10 | 鎮江市 | 5017.0 |
10 | 11 | 淮安市 | 4742.4 |
11 | 12 | 宿遷市 | 4112.0 |
12 | 13 | 連雲港市 | 4005.0 |
# 讀取江蘇地圖資料,資料來自DataV.GeoAtlas,將其投影到EPSG:4573
data = gpd.read_file('https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json').to_crs('EPSG:4573')
# 合併資料
data = data.join(gdp.set_index('地級市')["2022年GDP(億元)"],on='name')
# 修改列名
data.rename(columns={'2022年GDP(億元)':'GDP'},inplace=True)
data.head()
adcode | name | childrenNum | level | parent | subFeatureIndex | geometry | GDP | |
---|---|---|---|---|---|---|---|---|
0 | 320100 | 南京市 | 11 | city | {'adcode': 320000} | 0 | MULTIPOLYGON (((19828216.260 3681802.361, 1982... | 16907.9 |
1 | 320200 | 無錫市 | 7 | city | {'adcode': 320000} | 1 | MULTIPOLYGON (((19892555.472 3541293.638, 1989... | 14850.8 |
2 | 320300 | 徐州市 | 10 | city | {'adcode': 320000} | 2 | MULTIPOLYGON (((19736418.457 3894748.096, 1973... | 8457.8 |
3 | 320400 | 常州市 | 6 | city | {'adcode': 320000} | 3 | MULTIPOLYGON (((19927182.917 3638819.801, 1992... | 9550.1 |
4 | 320500 | 蘇州市 | 9 | city | {'adcode': 320000} | 4 | MULTIPOLYGON (((19929872.531 3547724.116, 1993... | 23958.3 |
fig, ax = plt.subplots(figsize=(9, 9))
# legend_kwds設定matplotlib的legend引數
data.plot(ax=ax,column='GDP', cmap='coolwarm', legend=True,legend_kwds={'label': "GDP(億元)", 'shrink':0.5})
ax.axis('off')
# 設定
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "black",'weight': 'bold'}
# 設定標題
ax.set_title('江蘇省地級市2022年GDP資料視覺化', fontsize=24)
for index in data.index:
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = data.iloc[index]["name"]
if name in ["蘇州市","無錫市"]:
x = x*1.001
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
# 儲存圖片
fig.savefig('res.png', dpi=300, bbox_inches='tight')
可以看到在上述地圖中,由於蘇州市的數值太大,其他資料被壓縮到淺色區域,無法有效展示資料分佈。需要使用地圖分層設色來更好地展示資料。地圖分層設色是一種常見的地圖視覺化方式,它可以將地圖上的資料按照不同的分類進行分層,並對每一層資料進行不同的顏色設定,以便更加直觀地展現地理空間資料的分佈情況和特徵。
在本文通過Python模組mapclassify用於分層設色和資料視覺化。使用mapclassify之前需要輸入以下命令安裝相關模組:
pip install mapclassify
mapclassify官方倉庫見:mapclassify。mapclassify提供了多種分組方法,可以幫助我們更好地理解資料的分佈情況,mapclassify提供的方法包括:
BoxPlot
: 基於箱線圖的分類方法。這種分類方法適用於資料分佈比較規律的情況。EqualInterval
: 等距離分類方法。這種分類方法將資料劃分為等距離的若干區間。適用於資料分佈比較均勻的情況。FisherJenks
: 基於Fisher-Jenks演演算法的分類方法。這種分類方法將資料劃分為若干區間,使得每個區間內部的差異最小,不同區間之間的差異最大。適用於資料分佈比較不規律的情況。HeadTailBreaks
: 基於Head-Tail演演算法的分類方法。這種分類方法將給定的資料集分為兩部分:頭部和尾部。頭部通常包含出現頻率最高的值,而尾部包含出現頻率較低的值。適用於識別資料集中的異常值和離群值。JenksCaspall
: 基於Jenks-Caspall演演算法的分類方法。這種分類方法根據資料中發現的自然分組將資料集劃分為類。適用於需要將資料分類為幾個具有明顯含義的區間的情況。JenksCaspallForced
: 強制基於Jenks-Caspall演演算法的分類方法。與JenksCaspall演演算法類似,但是它對區間的數量和大小有更強的控制力。適用於需要精確控制區間數量和大小的情況。JenksCaspallSampled
: 取樣基於Jenks-Caspall演演算法的分類方法。該方法對資料進行取樣,然後使用Jenks-Caspall演演算法對取樣後的資料進行分類,適用於資料量比較大的情況。MaxP
: 基於最大界限的分類方法。這種分類方法將資料劃分為幾個區間,使得不同區間之間的差異最大。適用於需要將資料分類為幾個具有明顯差異的區間的情況。MaximumBreaks
: 基於最大間隔的分類方法。這種分類方法與MaxP演演算法類似,但是它更加註重區間的可理解性。適用於需要將資料分類為幾個具有明顯含義的區間的情況。NaturalBreaks
: 基於自然間隔的分類方法。這種分類方法將資料劃分為幾個區間,使得每個區間內部的差異最小,不同區間之間的差異最大。適用於資料分佈比較不規律的情況Quantiles
: 基於分位數的分類方法。Percentiles
: 基於百分位數的分類方法。StdMean
: 基於標準差分組的分類方法。UserDefined
: 基於自定義分組的分類方法。關於以上常用方法的具體介紹可以看看基於geopandas的空間資料分析——深入淺出分層設色。mapclassify對資料進行分類簡單使用方法如下:
範例1
import mapclassify
# 匯入範例資料
y = mapclassify.load_example()
print(type(y))
print(y.mean(),y.min(), y.max())
y.head()
<class 'pandas.core.series.Series'>
125.92810344827588 0.13 4111.45
0 329.92
1 0.42
2 5.90
3 14.03
4 2.78
Name: emp/sq km, dtype: float64
mapclassify.EqualInterval(y)
EqualInterval
Interval Count
--------------------------
[ 0.13, 822.39] | 57
( 822.39, 1644.66] | 0
(1644.66, 2466.92] | 0
(2466.92, 3289.19] | 0
(3289.19, 4111.45] | 1
範例2
y = [1,2,3,4,5,6,7,8,9,0]
# 分為四個區間
mapclassify.JenksCaspall(y, k=4)
JenksCaspall
Interval Count
--------------------
[0.00, 2.00] | 3
(2.00, 4.00] | 2
(4.00, 6.00] | 2
(6.00, 9.00] | 3
範例3
y = [1,2,3,4,5,6,7,8,9,0]
# 自定義區間
mapclassify.UserDefined(y, bins=[5, 8, 9])
UserDefined
Interval Count
--------------------
[0.00, 5.00] | 6
(5.00, 8.00] | 3
(8.00, 9.00] | 1
範例4
import mapclassify
import pandas
from numpy import linspace as lsp
demo = [lsp(3,8,num=10), lsp(10, 0, num=10), lsp(-5, 15, num=10)]
demo = pandas.DataFrame(demo).T
demo.head()
0 | 1 | 2 | |
---|---|---|---|
0 | 3.000000 | 10.000000 | -5.000000 |
1 | 3.555556 | 8.888889 | -2.777778 |
2 | 4.111111 | 7.777778 | -0.555556 |
3 | 4.666667 | 6.666667 | 1.666667 |
4 | 5.222222 | 5.555556 | 3.888889 |
# 使用apply函數應用分層,rolling表示是否進行滑動視窗計算以消除隨機波動
demo.apply(mapclassify.Quantiles.make(rolling=True)).head()
0 | 1 | 2 | |
---|---|---|---|
0 | 0 | 4 | 0 |
1 | 0 | 4 | 0 |
2 | 1 | 4 | 0 |
3 | 1 | 3 | 0 |
4 | 2 | 2 | 1 |
方法1
GeoPandas中分層設色可以通過plot函數中的scheme引數和k引數設定資料分層方式和分層類別數。如下所示,通過JenksCaspall將GDP資料分為4級,可以直觀看到GDP資料在第一梯隊的城市有哪些。
fig, ax = plt.subplots(figsize=(9, 9))
# 使用分層設色後,legend_kwds要進行相應修改
data.plot(ax=ax,column='GDP', cmap='coolwarm', legend=True,scheme='JenksCaspall',k=4,legend_kwds={
'loc': 'lower left',
'title': 'GDP資料分級(億元)',
})
ax.axis('off')
# 設定
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "black",'weight': 'bold'}
# 設定標題
ax.set_title('江蘇省地級市2022年GDP資料視覺化', fontsize=24)
for index in data.index:
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = data.iloc[index]["name"]
if name in ["蘇州市","無錫市"]:
x = x*1.001
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
# 儲存圖片
fig.savefig('res.png', dpi=300, bbox_inches='tight')
方法2
在GeoPandas中也可以通過mapclassify直接處理資料,生成新的資料列進行展示。通過該種方式可以看到,蘇州和南京的GDP領先於其他地級市。
# 建立分類器
classifier = mapclassify.HeadTailBreaks(data['GDP'])
classifier
HeadTailBreaks
Interval Count
----------------------------
[ 4005.00, 9473.76] | 8
( 9473.76, 15329.34] | 3
(15329.34, 20433.10] | 1
(20433.10, 23958.30] | 1
# 賦值資料
data['GDP_class'] = data['GDP'].apply(classifier)
data['GDP_class'] = data['GDP_class'].apply(lambda x : int(x))
data.head()
adcode | name | childrenNum | level | parent | subFeatureIndex | geometry | GDP | GDP_class | |
---|---|---|---|---|---|---|---|---|---|
0 | 320100 | 南京市 | 11 | city | {'adcode': 320000} | 0 | MULTIPOLYGON (((19828216.260 3681802.361, 1982... | 16907.9 | 2 |
1 | 320200 | 無錫市 | 7 | city | {'adcode': 320000} | 1 | MULTIPOLYGON (((19892555.472 3541293.638, 1989... | 14850.8 | 1 |
2 | 320300 | 徐州市 | 10 | city | {'adcode': 320000} | 2 | MULTIPOLYGON (((19736418.457 3894748.096, 1973... | 8457.8 | 0 |
3 | 320400 | 常州市 | 6 | city | {'adcode': 320000} | 3 | MULTIPOLYGON (((19927182.917 3638819.801, 1992... | 9550.1 | 1 |
4 | 320500 | 蘇州市 | 9 | city | {'adcode': 320000} | 4 | MULTIPOLYGON (((19929872.531 3547724.116, 1993... | 23958.3 | 3 |
import matplotlib.colors as colors
fig, ax = plt.subplots(figsize=(9, 9))
# 設定分層顏色條
cmap = plt.cm.get_cmap('coolwarm', len(set(data['GDP_class'])))
# vmax和vmin設定是為了讓等級值居中
data.plot(ax=ax,column='GDP_class', cmap=cmap, legend=False,vmin=-0.5,vmax=3.5)
ax.axis('off')
# 設定Colorbar的刻度
cbar = ax.get_figure().colorbar(ax.collections[0],shrink=0.5)
cbar.set_ticks([0,1,2,3])
cbar.set_label('GDP資料分級')
cbar.set_ticklabels(['等級0','等級1','等級2','等級3'])
# 隱藏刻度線
ticks = cbar.ax.get_yaxis().get_major_ticks()
for tick in ticks:
tick.tick1line.set_visible(False)
tick.tick2line.set_visible(False)
# 設定
fontdict = {'family':'FZSongYi-Z13S', 'size':8, 'color': "black",'weight': 'bold'}
# 設定標題
ax.set_title('江蘇省地級市2022年GDP資料視覺化', fontsize=24)
for index in data.index:
x = data.iloc[index].geometry.centroid.x
y = data.iloc[index].geometry.centroid.y
name = data.iloc[index]["name"]
if name in ["蘇州市","無錫市"]:
x = x*1.001
ax.text(x, y, name, ha="center", va="center", fontdict=fontdict)
# 儲存圖片
fig.savefig('res.png', dpi=300, bbox_inches='tight')
方法3
分層設色不僅可以設定各區域的顏色,也可以設定各區域的填充圖案
# 建立分類器
classifier = mapclassify.MaximumBreaks(data['GDP'], k=3)
classifier
MaximumBreaks
Interval Count
----------------------------
[ 4005.00, 13115.20] | 10
(13115.20, 20433.10] | 2
(20433.10, 23958.30] | 1
# 賦值資料
data['GDP_class'] = data['GDP'].apply(classifier)
data['GDP_class'] = data['GDP_class'].apply(lambda x : int(x))
data.head()
adcode | name | childrenNum | level | parent | subFeatureIndex | geometry | GDP | GDP_class | |
---|---|---|---|---|---|---|---|---|---|
0 | 320100 | 南京市 | 11 | city | {'adcode': 320000} | 0 | MULTIPOLYGON (((19828216.260 3681802.361, 1982... | 16907.9 | 1 |
1 | 320200 | 無錫市 | 7 | city | {'adcode': 320000} | 1 | MULTIPOLYGON (((19892555.472 3541293.638, 1989... | 14850.8 | 1 |
2 | 320300 | 徐州市 | 10 | city | {'adcode': 320000} | 2 | MULTIPOLYGON (((19736418.457 3894748.096, 1973... | 8457.8 | 0 |
3 | 320400 | 常州市 | 6 | city | {'adcode': 320000} | 3 | MULTIPOLYGON (((19927182.917 3638819.801, 1992... | 9550.1 | 0 |
4 | 320500 | 蘇州市 | 9 | city | {'adcode': 320000} | 4 | MULTIPOLYGON (((19929872.531 3547724.116, 1993... | 23958.3 | 2 |
import matplotlib.patches as mpatches
# 設定圖案列表
patterns = ["///", "", "*", "\\\\",".", "o", "O",]
cmap = plt.cm.get_cmap('coolwarm', len(set(data['GDP_class'])))
color_list = cmap([0,1,2])
fig, ax = plt.subplots(figsize=(9, 9))
# 自定義圖示
legend_list = []
# 按層次設定legend
for i in set(data['GDP_class']):
tmp = data[data['GDP_class']==i]
tmp.plot(ax=ax,column='GDP_class', legend=False,hatch=patterns[i],edgecolor='black',color=color_list[i], linestyle='-',linewidth=2)
legend_list.append(
mpatches.Patch(facecolor=color_list[i], edgecolor='black',linestyle='-', linewidth=2,hatch=patterns[i], label='等級{}'.format(i))
)
ax.axis('off')
# 設定標題
ax.set_title('江蘇省地級市2022年GDP資料視覺化', fontsize=24)
# 自定義圖示
ax.legend(handles = legend_list, loc='best', fontsize=12, title='GDP資料分級', shadow=True)
# 儲存圖片
fig.savefig('res.png', dpi=300, bbox_inches='tight')