克里金插值

Posted UQI-LIUWJ

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了克里金插值相关的知识,希望对你有一定的参考价值。

1 定义

        相比反距离插值反距离插值 IDW_UQI-LIUWJ的博客-CSDN博客,克里金插值公式更加抽象

        其中 是点 (xo,yo)处的估计值

        这里的 λi是权重系数。它同样是用空间上所有已知点的数据加权求和来估计未知点的值。但权重系数并非距离的倒数,而是能够满足点 (xo,yo)处的估计值与真实值 zo的方差最小的一套最优系数,即
        

        同时满足无偏估计的条件

        

2 普通克里金插值

        不同的克里金插值方法的主要差异就是假设条件不同。

        普通克里金插值的假设条件为,空间属性z在各个点有相同的分布。换句话说对于空间任意一点(x,y),z(x,y)都有同样的期望c与方差σ2。即:

        

        换一种说法是,对于任意一点(x,y),该点的值z(x,y),都由区域平均值c和该点的随机偏差R(x,y)组成,即
        

         其中R(x,y)表示点(x,y)处的偏差,其均值为0,方差为常数
 

 2.1 无偏约束条件

        将代入有:

        

         因为对任意的z都有E[z]=c,所以有:

        

          于是有:

              

 2.2 优化目标函数

        我们的目标是

        为了方便起见,我们令

        于是有:

 

         再次简化描述,令:        

        这里:

 

         于是有:

 

2.2.1 优化目标函数的最优解

我们 定义半方差函数,将以及代入J,有:

 

我们的目标是寻找使J最小的一组 λi,且J是 λi的函数,因此直接将J对 λi求偏导数令其为0即可。

但是要注意的是,我们要保证求解出来的最优 λi 满足

这是一个带约束的最优化问题,我们使用拉格朗日乘法进行求解,构造如下的目标函数

求解使这个代价函数最小的参数集 ϕ,λ1,λ2,⋯,λn,则能满足其在 约束下最小化J。即

 

 

由于,而,所以

于是有: 

 

 在以上计算中我们得到了对于求解权重系数 λj的方程组。写成线性方程组的形式就是:
 

 

 写成矩阵形式为:

 对矩阵求逆即可求解。

2.3 半方差函数 

上文中对半方差函数的定义为

它的等价形式为:

 

2.3.1 证明半方差函数的等价

而 

于是有:

所以:

2.3.1 如何计算rij?

 表达了属性的相似度
 

 空间的相似度就用距离来表达,定义i与j之间的几何距离

 3 小节

总的来说,进行克里金插值分为这几个步骤:

1)对于观测数据,两两计算距离与半方差

2)寻找一个拟合曲线拟合距离与半方差的关系,从而能根据任意距离计算出相应的半方差

3)计算出所有已知点之间的半方差 rij

4)对于未知点 zo,计算它到所有已知点 zi的半方差 rio

5)求解方程组,得到最优系数 λi

6)使用最优系数对已知点的属性值进行加权求和,得到未知点 zo的估计值

 

以上是关于克里金插值的主要内容,如果未能解决你的问题,请参考以下文章

克里金插值MATLAB程序

R: Kriging interpolation and cross validation 克里金插值及交叉验证浅析

如何在 Python 中使用克里金法对 2D 空间数据进行插值?

转Kriging插值法

arcgis怎么在指定的研究区范围内做克里金插值,并且把各个取值范围的面积统计出来。

openlayers4 入门开发系列之前端动态渲染克里金插值 kriging 篇(附源码下载)