前几天有个朋友在微信群里面问,如何画一个这样的图:

其实这种图还是很简单的,我记得我之前在做美赛的时候,就用python的basemap画过。当时也没公众号,代码也丢了,现在basemap也基本上过时了,都是使用cartopy画图。
这种图其实我也不知道叫什么名字,描述一下就是:map add 3d bar
⚠️:代码在文末
把这个描述在谷歌里面搜索,其实可以找到很多类似的图:










上面的若干截图只是为了展示图像样式!!!!
那么本文,将介绍,如何使用python来画,基于中国地图做一些创作。
# 导入包
from shapely.geometry import Polygon, MultiPolygon
from matplotlib import cm
import matplotlib
from getchinamap.getchinamap import DownloadChmap
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import platform
from mpl_toolkits.mplot3d.art3d import Poly3DCollection
import geopandas as gpd
import cartopy.feature
from cartopy.mpl.patch import geos_to_path
import itertools
from mpl_toolkits.mplot3d import Axes3D
# import matplotlib.pyplot as plt
from matplotlib.collections import LineCollection,PolyCollection
# matplotlib 显示中文的问题
if platform.system() == 'Darwin':
plt.rcParams["font.family"] = 'Arial Unicode MS'
elif platform.system() == 'Windows':
plt.rcParams["font.family"] = 'SimHei'
else:
pass
chinamapdata = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/geojson?code=100000_full")
chinamapdata.head(3)

sampledata = pd.DataFrame({'lon':chinamapdata.centroid.x, 'lat':chinamapdata.centroid.y, 'name':chinamapdata['name']})
sampledata['value'] = np.random.randint(0, 1000, sampledata.shape[0])
sampledata.head(4)
lon、lat、value,分别是经纬度和对于的值,一般这个值都是正数。name 其实要不要无所谓,你想展示就加,不想展示就算了。反正我这里有(主要是想分享如何在3d空间 添加文本)
bounds_box = chinamapdata.bounds
minx = bounds_box['minx'].min()
miny = bounds_box['miny'].min()
maxx = bounds_box['maxx'].max()
maxy = bounds_box['maxy'].max()
这个步骤是核心,因为现成的geopands、cartopy、matplotlib组合,是没办法直接画出来这种特定需求的图的。那么就需要对数据做一系列的转换。
主要步骤有:
geoms = chinamapdata.geometry
# target_projection = ccrs.PlateCarree()
# geoms = [target_projection.project_geometry(geom, target_projection)
# for geom in geoms]
paths = list(itertools.chain.from_iterable(geos_to_path(geom) for geom in geoms))
# At this point, we start working around mpl3d's slightly broken interfaces.
# So we produce a LineCollection rather than a PathCollection.
segments = []
for path in paths:
vertices = [vertex for vertex, _ in path.iter_segments()]
vertices = np.asarray(vertices)
segments.append(vertices)
lc = LineCollection(segments, color='black')
%matplotlib widget
# part 1
segments = []
for path in paths:
vertices = [vertex for vertex, _ in path.iter_segments()]
vertices = np.asarray(vertices)
segments.append(vertices)
# with plt.style.context('fivethirtyeight'):
fig = plt.figure()
ax = Axes3D(fig, xlim=[minx, maxx], ylim=[miny, maxy])
# ax.set_zlim(bottom=0)
lc = LineCollection(segments, color='black',linewidths=1)
ax.add_collection3d(lc)
ax.bar3d(x=sampledata['lon'], y=sampledata['lat'], z=np.zeros_like(sampledata['value']),
dx=np.ones_like(sampledata['value']),
dy=np.ones_like(sampledata['value']),
dz=sampledata['value'],alpha=0.8)
for index, iterrow in sampledata.iterrows():
ax.text(iterrow['lon'], iterrow['lat'], iterrow['value']+2,iterrow['name'], color='green',size=9)
ax.text(80, 30, 1000, '公众号: world of statistics', size=20, color='gray')
ax.set_xlabel('经度')
ax.set_ylabel('维度')
# ax.set_xlim([minx, maxx])
# ax.set_ylim([miny, maxy])
ax.set_zlabel('value')
# ax.set_title("公众号: world of statistics")
plt.show()
效果如下:

如果你希望可以对区域填充,只需要稍微更改几行代码即可.
# geoms
# part 2
fig = plt.figure()
ax = Axes3D(fig, xlim=[minx, maxx], ylim=[miny, maxy])
concat = lambda iterable: list(itertools.chain.from_iterable(iterable))
polys = concat(path.to_polygons() for path in paths)
lc = PolyCollection(polys, edgecolor='black',
facecolor='green', closed=False, alpha=0.4)
ax.add_collection3d(lc)
ax.bar3d(x=sampledata['lon'], y=sampledata['lat'], z=np.zeros_like(sampledata['value']),
dx=np.ones_like(sampledata['value']),
dy=np.ones_like(sampledata['value']),
dz=sampledata['value'],alpha=0.8)
ax.text(80, 30, 1000, '公众号: world of statistics', size=20, color='gray')
ax.set_xlabel('经度')
ax.set_ylabel('维度')
# ax.set_xlim([minx, maxx])
# ax.set_ylim([miny, maxy])
ax.set_zlabel('value')
# ax.set_title("公众号: world of statistics")
plt.show()
效果如下:

整体的感觉是:如果地图填充了颜色,其实还好看一点。
本文章的代码全部在GitHub上免费分享。
01开头的ipynb文件。