IDW空间插值法matlab,基于IDW对PM2.5进行空间插值及可视化

前面几篇推文我们分辨介绍了使用Python和R绘制了二维核密度空间插值方法，并使用了Python可视化库plotnine、Basemap以及R的ggplot2完成了相关可视化教程的绘制推文，详细内容如下：

接下来，我们将继续介绍空间插值的其他方法，本期推文，我们将介绍IDW(反距离加权法(Inverse Distance Weighted)) 插值的Python计算方法及插值结果的可视化绘制过程。主要涉及的知识点如下：

IDW简介

自定义Python代码计算空间IDW

分别使用plotnine、Basemap进行IDW插值结果可视化绘制

IDW简介

反距离权重 (IDW) 插值假设：彼此距离较近的事物要比彼此距离较远的事物更相似。当为任何未测量的位置预测值时，反距离权重法会采用预测位置周围的测量值与距离预测位置较远的测量值相比，距离预测位置最近的测量值对预测值的影响更大。反距离权重法假定每个测量点都有一种局部影响，而这种影响会随着距离的增大而减小。由于这种方法为距离预测位置最近的点分配的权重较大，而权重却作为距离的函数而减小，因此称之为反距离权重法。(解释来源于网络)，繁琐的公式也没放，这里我们给出几张示意图即可，原理不解的小伙伴可自行百度。

(基于采样点距离的IDW插值(左)从高程矢量点插值的IDW曲面(右))

自定义Python代码计算空间IDW

我们免去了了繁琐的IDW插值原理部分，这节我们直接根据原理自定义IDW函数，根据已有样例站点位置及对应值，计算IDW结果。在这之前，我们给出所需样例的预览及地图文件的范围(构建插值网格所需)，结果如下：

样例点：

地图文件范围信息：

js_box = js.geometry.total_bounds

js_box

#array([116.36196 ,  30.757975, 121.975185,  35.122924])

小伙伴们对上述计算结果有疑惑的地方可以详细阅读之前的插值文章(文前链接)，或者等我将这系列做完会推出详细的源码及解释文档(目前在整理中)

定义IDW计算函数

这里主要涉及两个计算函数，计算经纬度点转实际距离(km)的haversine方法和计算IDW的函数，定义函数如下：

haversine方法:

import math

import numpy as np

#更换求距离的函数

from math import radians, cos, sin, asin, sqrt

def haversine(lon1, lat1, lon2, lat2):

R =  6372.8

dLon = radians(lon2 - lon1)

dLat = radians(lat2 - lat1)

lat1 = radians(lat1)

lat2 = radians(lat2)

a = sin(dLat/2)**2 + cos(lat1)*cos(lat2)*sin(dLon/2)**2

c = 2*asin(sqrt(a))

d = R * c

return d

IDW

def IDW(x, y, z, xi, yi):

lstxyzi = []

for p in range(len(xi)):

lstdist = []

for s i