Python 空间显著性打点图绘制

使用python Basemap绘制GPCP降水趋势以及显著性打点：

数据处理：

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import xarray as xr
from scipy.stats import linregress
import matplotlib as mpl
import pandas as pd

def linregress_params(time, values):
    """
    计算时间序列的线性回归参数
    """
    result = linregress(time, values)
    return result.slope, result.intercept, result.rvalue, result.pvalue, result.stderr

#xarray读取GPCP月尺度降水数据
path = 'Pr_GPCP_180_360_200301_202212.nc'
ds = xr.open_dataset(path)
ds = ds.assign_coords(time=pd.date_range("2001-01-01", "2020-12-31", freq="M"))
lons = ds ['lon'][:] 
lats = ds ['lat'][:]
time = ds ['time'][:]

# 按年份分组，计算年平均降水
annual_mean = ds.groupby("time.year").mean("time")
annual_mean = annual_mean.rename(year="time")  

#计算趋势、P值
results = xr.apply_ufunc(
    linregress_params,                  # 函数
    annual_mean["time"],            # 时间变量
    annual_mean['Pr'],                         # 数据变量
    input_core_dims=[["time"], ["time"]], 
    vectorize=True,                     # 自动向量化
    dask="parallelized",                # 支持并行计算
    output_core_dims=[[], [], [], [], []],  # 每个输出单独一个维度
    output_dtypes=[float, float, float, float, float],  # 每个输出的数据类型
)

slope = results[0] #趋势 
pvalue = results[3] #p值
p_point = np.where(pvalue.values<0.1,1,np.nan) #p值小于0.1设置为1

画图

#绘图
plt.figure(figsize=(12,6),dpi=200)
world = Basemap(projection='mill', lat_0=0, lon_0=0, resolution='l', area_thresh=1000.0,
                llcrnrlon=-180, urcrnrlon=180, llcrnrlat=-60, urcrnrlat=85)
lon1,lat1=np.meshgrid(lons,lats)
x, y = world(lon1,lat1)
world.drawparallels(np.arange(-90.,91.,30.),labels=[1,0,0,0],linewidth=0.5,fontsize=10)
world.drawmeridians(np.arange(-180.,181.,60.),labels=[0,0,0,1],linewidth=0.5,fontsize=10)
world.drawcoastlines()

#设置绘图数值范围
norm = mpl.colors.BoundaryNorm(boundaries=np.linspace(-0.05,0.05,25),ncolors=256)
#画图
pc = world.pcolormesh(x,y,slope,cmap='RdBu',norm=norm,shading='nearest')
cb=world.colorbar(pc,label='slope (mm/day/year)')
num=2 #设置打点多少个点打一个点，调整密集度
significance_mask = p_point==1
sampled_mask = np.zeros_like(significance_mask, dtype=bool)
sampled_mask[::num, ::num] = significance_mask[::num, ::num]
x_point, y_point = world(lon1[sampled_mask], lat1[sampled_mask])
world.scatter(x_point, y_point,s=0.1, color='black',alpha=1,edgecolors='black',linewidths=0.6 )
plt.title('GPCP Pr ',fontsize=16)

可设置更密集，（num=1）：