python实现概率密度匹配法

Python实现系列：
Python实现贝叶斯优化算法
python实现t-SNE降维
Python实现12种降维算法
Python实现11 种特征选择策略
Python实现10种聚类算法
Python实现8种相似度度量
python实现反距离权重插值(IDW)
Python实现12种概率分布
python实现快速傅里叶变换
python实现隐马尔可夫模型

1 概率密度匹配法作用

在任何测量中都能存在系统误差和随机误差：

• 随机误差呈正态分布，在气候平均处理中可基本被消除。
• 两者的主要差别在于系统误差平均值不为零。

系统误差又有独立误差和非独立误差之分：

• 独立误差，是指误差不随区域、时间或是观测值大小而变化，可通过对整个气象场减去相同的常数来消除，订正方法相对简单。
• 但非独立误差，会随观测量而变化，订正较为困难。

目前概率密度匹配法 (probability density function matching method，简称PDF方法)在订正非独立系统误差方面具有优越性，因此本文准备尝试。

2 PDF方法订正思路

• 对于每个订正格点，根据资料的时空分辨率及误差特征选取适当的时空窗口，收集匹配的地面观测和格点估测资料，样本量足够大时，分别得到二者稳定的累积概率密度分布。
• 相同的累积概率密度值对应不同的地面观测和估测降水量时，根据二者的偏差来订正降水量。例如，当降水累积概率密度为30%时，卫星反演降水量为4.5 mm，而地面观测降水量为5 mm，则此时卫星降水量的误差为-0.5 mm。
• 可采用十折交叉验证的方法检查订正效果。

3 PDF方法订正方案

① 对每个格点选取时空窗口，收集匹配的地面观测和卫星降水格点资料，如选取当前日期起向前共取 $30d$ 作为时间窗口，以目标格点为中心的 $10°\times 10°$ 空间范围作为空间窗口。

② 将上述时空窗口中有观测站点，且地面观测和卫星降水资料均为非缺测的格点值纳入统计，选取的样本数至少为 300 个，如达不到要求，搜索的空间范围适当扩大。

③ 计算每个格点卫星降水量 ($R_s$) 的累积概率密度 ($f_s$)，以及该概率密度所对应的地面降水量 ($R_g$)，得到订正值 (Δr)，则订正后的卫星降水量 ($R_c$) 为$R_c=R_s-\Delta r$。

4 python 实现

1)构建数据集

import scipy.stats as st
import numpy as np
import matplotlib.pyplot  as plt
%matplotlib inline

grid = np.random.randint(0,1000, size=10000)
grid = grid.reshape(100,100)
ix = np.random.randint(0,100,size=300)
iy = np.random.randint(0,100, size=300)
station = np.random.randint(100,800,size=300)
grid_station = grid[iy,ix]
1
2
3
4
5
6
7
8
9
10
11

2)计算概率密度

# 累计概率密度
cdf_station = st.norm.cdf(station,loc=np.mean(station),scale=np.std(station)) # 返回每个数据点累计概率
cdf_grid = st.norm.cdf(grid,loc=np.mean(grid),scale=np.std(grid))
1
2
3

3)概率密度匹配

cdf_grid_station = cdf_grid[iy,ix]
delta = np.zeros(cdf_grid_station.shape)
for i, prob in enumerate(cdf_grid_station):
    for j , prob_true in enumerate(cdf_station):
        if abs(prob-prob_true)<=0.01:
            delta[i]=grid_station[i]-station[j]
            break
1
2
3
4
5
6
7

4)查看效果

蓝色：订正前误差，橙色：订正后误差。从图上来看，有一定效果。

5)格点偏差订正

直接对偏差进行插值即可。

关注公众号我不爱机器学习，后台回复pdf，即可获取完整代码链接。

参考：http://html.rhhz.net/yyqxxb/html/20130504.htm