一种自适应选取参考点的Kriging 插值方法研究

邓显石

（湖南环境生物职业技术学院园林学院，湖南 衡阳421005）

空间插值是利用已知的样本点来估算其它未知点数值的过程，常用的空间插值方法主要有IDW、Kriging、趋势面法、自然邻域法等。其中Kriging 以变异函数理论和结构分析为基础，理论上能够获得更高的插值精度，在地质、水文、土壤等领域获得广泛应用[1]。

Kriging 插值算法精度的提高主要受到变异函数、样点的分布、参考点的数量等因素的影响[2]。作为直接影响Kriging 插值误差的变异函数，关于其研究较多。如李莎与徐爱萍等在2012 年提出采用反复实验人工拟合变异函数[3,4]，2016 年徐武平提出一种变异函数的自动拟合方法[5]。作为影响程度次于变异函数的参考点的分布及数量选取，目前却基本是通过实验对比方式来进行调优。如张仁铎在2005 年认为选择远距离的数据点可能会更容易违反二阶平稳的假设，一般选择邻近的10 个点左右即可[6]；李君在2010 年通过实验发现参考点数在40-60 个点最为合适[7]；olea 认为邻近点的选择一般应介于3-25 点之间[8]。关于参考样本点分布均匀能够提高插值精度已得到多数学者的研究证实，如孙益权通过实验发现样本点分布均匀越能够提高kriging 的插值精度[9]，但目前顾及样本点分布不均匀性的自适应选取参考样点的Kriging 插值方法研究较少。

针对此种现状，本文提出一种自适应选取参考样点的Kriging 方法--OOK。该方法将待插值点的“一阶邻近点”（见后文定义1）作为最终插值的参考点，能够有效降低

图7 选择特征随吸合次数增加的变化曲线样本点采样分布不均匀对最终插值结果造成的影响，同时不需要用户指定参考点的数目。

1 普通Kriging 模型

在Kriging 插值中，待插值点的属性由周边已知样本点变量值的线性组合求得[10]，其估计值为：

式中，Y(xi,xj) 为邻近点xi，xj之间的半方差Y(xi,x0)为邻近点xi与待插值点x0之间的半方差，可从拟合的半变异函数模型（球形、高斯、线性模型）中计算得到，φ 为拉格朗日算子。

2 OOK 插值方法

2.1 定义1：将待插值点与研究区域内已知的样本点建立不规则三角网，将与待插值点形成邻接拓扑关系的样本点称为“一阶邻近点”，如图（1）所示。由图可以发现，采用“一阶邻近点”作为待插值点的参考点，能够在一定程度上降低样本点在一侧发生聚集对参考点选择造成的影响。

图1 一阶邻近点

2.2 技术路线

OOK 插值算法主要包括数据清洗，参考点的自适应选取，半变异函计算，待求点插值4 个过程，如图2 所示。

2.2.1 数据清洗：现实世界复杂多变，无论采用何种方式获取的空间数据，均存在一些不可避免的问题[12]，如重复的空间数据点。由于需要将空间点数据生成不规则三角网，因此首先应清洗空间数据中存在的重复点。

2.2.2 邻近点的自适应选取：①首先利用刘永和于2008 年提出的“一种快速生成平面Delaunay 三角网的横向扩展法”[13]，将n 个已知样本点建立Delaunay 三角网。对于每个待插值点，不将其与全部样点重新建立不规则三角网，而是采取逐步加点法构建新的Delaunay 三角网，如此将提高Delaunay 三角网的构网效率；②将待插值点q 的“一阶邻近点”作为初始参考点集合；③计算“一阶邻近点”与待插值点之间的距离edgei，当edgei＞R(变程[2])，则将此点从初始参考点集合中去除。

图2 OOK 插值流程

2.2.3 半变异函数值的计算：本文中半方差函数模型的拟合采用GS+软件，其提供的半变异函数有球状、高斯等模型，在拟合各类模型时都需设置步长与“最大滞后距”两个关键参数。根据相关研究指出，“最大滞后距”一般取研究区域样点间最大长度的1/2，步长值的设置应使得各组步长内存在不低于30 个样本点对[7]。

2.2.4 计算待求点的属性值：根据步骤2.2.3 求得参考点与待插值点以及参考点之间的半变异值，代入上式（3），计算得到待插值点对应的权重系数，从而代入公式（1）求得待插值点的属性值，并对其进行精度评估。

3 实验对比分析

本文采用“Spatial Interpolation Comparison 97 (SIC97)”数据集，收集了瑞士1986 年5 月8 日467 个降水站点观测值，随机从467 个站点中抽取100 个站点作为训练集，另367 个站点的数据作为检验数据集，467 个站点的分布如图3 所示，数据来源于https://www.researchgate.net/profile/Gregoire_Dubois/publication/281292076_Spatial_Interpolation_Comparison_97_ (SIC97) _dataset）。

图3 气象站点分布图

为检验OOK 方法的可行性，将其与OK、AIDW[14]进行比较。对比方法为交叉验证法，比较OK 在不同的选点条件下与OOK、AIDW 的MAE、RMSE、MAX、MIN 值，实验包含如下两种情形：

①将上述100 个站点的数据作为已知样本数据（单位km2存在2‰个采样点），估计剩余非边界点外站点的降水值，下文简称为实验1。

②从467 个站点除边界点外每次取出一个点作为待插值点（单位km2有12‰个点，实验2 样本数据的密度较实验1 大），从剩余样本点中抽取OK、AIDW 与OOK 插值方法需要的参考点来预测待插值点的降水值，简称为实验2。

实验1：研究区内最大降雨量为58.5mm，最低值为1mm，均值为18mm。首先利用GS+拟合出该区域的半变异函数模型，在给定“步长”与“最大滞后距”后，对样本数据的Z 值进行对数转换，选择RSS 值最低的球状模型（C0=0.071，C0+C=0.554，A0=71800，r2=0.847，RSS=0.0477），可得块金系数为12%，低于25%[15]。将OOK、OK（4-16 个参考样点）、AIDW 应用于该数据集的插值实验，得到如图4。

图4（左- 右，上- 下为a,b,c,d）

从图（4-a）可以发现，OK 插值算法的MAE 误差随着参考点的增加而逐步降低，在参考点数达到12 以上时变化趋于平缓，在14 个参考点时处于全局最低点。OK 在选择4-11 个参考点时的MAE 误差大于OOK，在12 以上略低于OOK 算法，约为0.06mm。

从图（4-b）可以发现，OK 插值算法的RMSE 误差同样随着参考点的增多而逐步减少，在参考点数为12 后呈现平缓的趋势，在点数为13 时全局最低。在OK 插值方法选择5 个参考点后，AIDW 的RMSE 误差大于OK，同时大于OOK。OOK 插值算法的误差比OK 选择4-10 个参考点要低。

从图（4-c）可以发现，OK 插值算法的MAX 误差基本随参考点数的增多而下降，在参考点大于5 以后，OK 与AIDW 插值算法的MAX 误差要优于OOK 算法。经过分析发现，OOK、AIDW与OK 插值算法MAX 值最大的待插值点的真实降水量为研究区域最大-51.8mm，其“一阶邻近点”中存在着近距离的局部最低点-19.2mm，可见其空间异质性显著。由于OK 与OOK 使用的半变异函数模型不能反映此局部变化，因此造成此点插值误差较大。但OK 算法随着搜索范围的扩大，通过引入更多与待插值点具有弱相关性的参考点，对该点的插值误差有轻微的抑制作用，OK 与OOK 的MAX 误差的极差为1.2mm，发生在OK 算法选取14 个参考样点。

从图（4-d）可知，OK 插值方法的最小值呈现不同振幅的震荡，其中MIN 值总体上最优的插值方法为OOK，其次AIDW，OK最后。

从上述实验可以看出：1）OOK 算法在OK 算法选取参考点数为4-11 点时占据绝对的优势，OK 算法在参考点大于12 时相较于其余两种插值算法精度更高；2）OOK 插值算法的精度总体上优于AIDW。

实验2：首先将研究区内样本数据降雨量进行均方根转换，拟合的半变异函数模型为球状模型（C0=0.001，C0+C=1.99，A0=83700，r2=0.959，RSS=0.21），可知块金系数为0.05%，相较于实验1 的12%，空间自相关性更强。将OOK、OK（4-16 个参考样点）、AIDW 应用于插值实验，得到如图5。

从图（5-a）中可以看出，OK 插值算法参考点为4-6 点时，MAE 误差逐步降低，其中在参考点为6 时，在3 种插值算法中最低。但随着参考点的数量逐渐增多，MAE 值呈现类对数方式上升，最后在参考点为13 时呈现平稳的趋势，在16 个点时，MAE 在3 种插值算法中误差最大。OOK 插值方法在OK 参考点数量为6 时，其MAE 误差略大于OK。AIDW 在OK 的参考点为9 后其MAE 值低于OK，同时从图可以发现AIDW 比OOK 的MAE 误差大。

图5（左-- 右，上-- 下为a,b,c,d）

从图中（5-b）可以发现，OK 插值算法的RMSE 误差在插值点数为6-9 个时呈现"J" 形增长趋势，随后趋于稳定，在16 个点时误差在3 种算法中最高，在参考点为6 时，RMSE 值最低。AIDW 在参考点为9 后RMSE 值低于OK。OOK 的误差比AIDW低，仅在参考点为6 时大于OK。

从图（5-c）中可以发现，OK 插值算法的MAX 值5、6 点时全局最低，此后总体趋势呈现逐步增大并最终趋于稳定，AIDW 的MAX 误差在3 种插值算法中最大。OOK 在OK 选择参考点为5、6 时，其MAX 值略大，约为0.3mm。

从图（5-d）可以看出，OK 插值算法的MIN 值随着点数的增多呈现上下震荡，OOK 的MIN 比AIDW 略低，在参考点为8--15 时较OK 算法略低，但总体而言，3 种方法的MIN 相差都较小，极差约为0.03mm。

由可知：a.OOK 插值方法仅在OK 选择6 个参考点时，其插值精度略低于OK，其余在3 种插值方法中均占据绝对优势；b.OOK 插值算法的插值精度比AIDW 要高，这是由于该数据集空间自相关性强烈且该区域的变异函数模型较好的反映了空间数据的分布特征[16]。

综合以上2 个实验可以发现：a. 没有哪种算法在任何情形下都占据绝对的优势；b.在空间自相关性强的数据中，顾及样本分布不均匀性的OOK 插值方法的插值精度（MAE、RMSE、MAX、MIN）较高，尤其在样本数据密度较大时。

4 结论

本文针对空间样本点数据的分布存在不均匀性的情形，提出了一种自适应选取参考样点的插值方法OOK，通过与经典的OK、AIDW 对比发现：OOK 插值算法无需用户指定参考样点的数量，同时能够顾及样点分布的方向不均匀性，能够提高插值精度。

但自适应选取参考点的Kriging 插值方法也存在一定的局限性：①从本文实验1 的MAX 值误差分析发现，方法的抗异质性影响的能力弱；②目前Kriging 插值算法的半变异函数普遍是通过不断试验比较得来，加大普通用户使用Kriging算法的难度。因此未来关于Kriging 插值方法的研究将聚焦于将ML（机器学习）用于空间异质性的判别与半变异函数模型的自动求取，进一步提高kriging 插值算法的精准度与自动化程度。