一些常用的空间密度统计方法。主要是cartopy新更新，而matplotlib旧有的hexbin图统计法；matplotlib提供的柱状图统计法；geoplot封装seaborn的kdeplot空间密度估计法。
一、hexbin图
这是一种matplotlib本身携带，但早期cartopy库包不能调用的命令。这个命令会将空间分割为合适大小的六边形，并统计落在该六边形的点的数量，以这个数量作为颜色映射的依据，为每个六边形填上颜色。
在最新的cartopy库包中，已经能够使用这个命令了。传入的值类似散点图命令。
df=pd.read_excel(filename)
lon=df['经度']
lat=df['纬度']
ah=ax.hexbin(lon,lat,gridsize=9,cmap='Blues')
ax.scatter(lon,lat,c='k')
cb=fig.colorbar(ah,ax=ax,pad=0.01)
二、柱状图统计
该命令仅能在PlateCarree投影下使用，通过直方统计出在经向和纬向上，数据的分布样式。可直接查阅matplotlib官网实例Scatter plot with histograms。
df=pd.read_excel(filename)
lon=df['经度']
lat=df['纬度']
data=df['rain']
def scatter_hist(x, y, ax, ax_histx, ax_histy):
ax_histx.tick_params(axis="x", labelbottom=False)
ax_histy.tick_params(axis="y", labelleft=False)
ax.scatter(x, y)
binwidth = 0.05#确定统计间隔，由于地图一度就是100km，所以一般可以取小一点
xymax = max(np.max(np.abs(x)), np.max(np.abs(y)))
lim = (int(xymax/binwidth) + 1) * binwidth
bins = np.arange(-lim, lim + binwidth, binwidth)
ax_histx.hist(x, bins=bins)
ax_histy.hist(y, bins=bins, orientation='horizontal')
left,width=0.1,0.65
bottom,height=0.1,0.65
spacing=0.005
rect_scatter=[left, bottom, width, height]
rect_histx=[left, bottom + height + spacing, width, 0.2]
rect_histy=[left + width + spacing, bottom, 0.2, height]
fig=plt.figure(figsize=(8, 8),dpi=200)
ax=fig.add_axes(rect_scatter)
ax.set_xlim(108.3,109.35)
ax.set_ylim(29.7,30.7)
ax.add_patch(patch)
ax_histx=fig.add_axes(rect_histx, sharex=ax)
ax_histy=fig.add_axes(rect_histy, sharey=ax)
scatter_hist(lon,lat,ax,ax_histx,ax_histy)
三、geoplot封装
这个命令其实是geoplot将seaborn里的KDEPLOT命令封装后的高级命令，使得仅能在sns里使用的二维密度估计能够使用到地图空间里。这时，任意投影都可以使用这个命令。
import geopandas as gpd
import geoplot as gplt
import geoplot.crs as gcrs
proj=gcrs.PlateCarree()
df=pd.read_excel(filename)
lon=df['经度']
lat=df['纬度']
geo=gpd.read_file(shp_path)
df_geo=gpd.GeoDataFrame(data,geometry=gpd.points_from_xy(lon,lat),crs='EPSG:4326')
fig=plt.figure(figsize=(6,6),dpi=200)
ax=plt.axes([0,0,1,0.78],projection=proj)
gplt.kdeplot(df_geo,
cmap='Reds',
projection=proj,
shade=True,
shade_lowest=True,
clip=geo.geometry,
ax=ax,gridsize=100)
gplt.polyplot(geo.geometry,zorder=1,ax=ax)
三种空间密度都可以看出一定的规律，即中部地区密度最高。
