空间插值是通过有限的离散采样点来建立某种插值函数关系f(x)，并将已知采样点范围内的任意位置代入函数关系式计算未知点属性值，常用的空间插值方法有反距离加权(inverse distance weighted，IDW)、克里金(Kriging)、趋势面、样条函数、自然邻域法。其中，IDW以其插值原理简单、计算简便且符合地理学第一定律而被广泛使用^[1]。

近年来关于IDW算法的研究较多，主要集中在幂指数^[2-3]、搜索点数^[4]、搜索方向^[5]等方面。通过分析现有的IDW插值方法，发现大多IDW插值算法均假设空间过程的平稳性^[6]，而由于地表空间异质性的普适存在^[7]，现有IDW及其改进算法的精度仍然不够令人满意。基于此，李佳霖等^[8]与樊子德等^[6]分别在2015年与2016年提出基于空间分层异质性的IDW插值方法，该算法的核心步骤为根据采样点的位置和属性进行分区。此类方法能够使得插值精度明显提高，但也存在一些问题：(1)分层数目的确定需要先验知识，而且通常会比较困难；(2)位于异质分层边界区域的插值精度会下降，尤其是点数较少的层内；(3)在空间异质分层后的子区域内部并不能消除空间异质性的影响。

针对现有基于空间分层异质性的IDW插值算法存在的问题，本文将机器学习的监督分类思想应用于IDW插值中待插值点的空间局部异质性判别，提出一种新的k近邻反距离加权(k-nearest neighbor adaptive IDW，kAIDW)插值算法。该算法首先根据样点属性数据的数理统计特征对样点进行自动分类，随后利用机器学习中的k近邻法在未得到插值结果的前提下判定待插值点的所属类别，从而自适应地为周边一阶邻近样点再设置不同的权重调和因子。该算法在无需人为干预的情况下能同时顾及空间相关与空间异质性对插值结果的影响。

1 IDW与k‐NN算法 1.1 IDW插值算法

IDW插值算法以其简便、高效等诸多优点被广泛应用于各大GIS软件中^[2]，如ArcGIS、MapInfo、MapGIS、SuperMap等，计算模型如下：

$ z\left(x\right)=\stackrel{m}{\sum\limits_{i=1}}\left[\right(1/{D}_{i}^{n})\cdot {z}_{i}]/\stackrel{m}{\sum\limits_{i=1}}\left(1/{D}_{i}^{n}\right) $

(1)

式中，z(x)为待插值点的属性值；z_i为样本点的属性值；D_i为样本点与待插值点之间的欧氏距离；n为幂指数，一般取1~3^[9]，但多数学者认为n为2能获得更好的实验结果^[10]；m为选择的样点个数。

1.2 k-NN分类算法

k-NN算法简单，易于实现，无需估计参数，作为数据挖掘的经典算法被广泛使用^[11-12]。经典的k-NN算法^[13]选取测试样本周边的k个近邻对象，应用“投票”规则确定测试样本的类别，即待判定点在邻近空间中的k个最相似的样本中多数属于同一类别，则该点将被预测为此类别，如此可实现算法误分类的经验风险最小化^[14]。

本文中k值为待插值点一阶邻近点数，从而使得选取的k个样点既考虑就近原则，同时又避免由于样点数据发生聚集现象而导致的分类结果严重依赖于某个方向样点的问题。并且由于k值的选取不需要训练，因此能大幅提高算法效率。

2 kAIDW算法 2.1 基本原理

经典IDW(classical IDW，CIDW)算法以假设空间过程的平稳性为前提，通常采用距离平方的反比进行定权，但忽略了待插值点与样点之间空间异质性的存在。本文提出的kAIDW算法以地理学第一定律为核心思想^[15]，同时通过机器学习算法k-NN判定参考样点与待插值点之间的空间异质性，最终构造出一个集空间相关与异质性于一体的插值模型，降低异质性参考样点的权重，从而提高插值的精度。

2.2 算法的设计与分析

kAIDW插值算法主要包括分类、定权、插值3个步骤，详细流程如图 1所示。

2.2.1 分类

k-NN算法确定待插值点属性类别的前提是获得样点属性的类别，该过程的主要任务是将样点连续的属性数据离散化，获取样点待插值属性数据离散化后的分类结果。

1) 数据预处理。由于空间点可能存在压盖的情况，影响后续Delaunay三角网的生成，因此将压盖位置的空间点保留其中一个。

2) 样点分类。根据预处理得到的样点，计算得到全部样点属性值的均值Z_mean与中误差Z_std，从而可以获取每个样点属性值分类的阈值H_i和L_i，计算公式如下：

$ {H}_{i}={Z}_{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}}+\left({Z}_{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}}/{z}_{i}\right)\times {Z}_{\mathrm{s}\mathrm{t}\mathrm{d}} $

(2)

$ {L}_{i}={Z}_{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}}-\left({Z}_{\mathrm{m}\mathrm{e}\mathrm{a}\mathrm{n}}/{z}_{i}\right)\times {Z}_{\mathrm{s}\mathrm{t}\mathrm{d}} $

(3)

从式(2)、式(3)中可以看出，H_i与L_i具备自适应特征，对于样点中较大的属性值，则易被式(2)检测，相反属性值愈低，则易被式(3)发现。通过如上方法，可以将全部样点按照属性值分成3类值：

$ {t}_{i}=\left\{\begin{array}{l}\mathrm{高}, {z}_{i}>{H}_{i}\\ \mathrm{中}, {L}_{i}\le {z}_{i}\le {H}_{i}\\ \mathrm{低}, {z}_{i} <{L}_{i}\end{array}\right. $

(4)

3) 一阶邻近样点的选择。Delaunay三角网广泛应用于空间邻近分析，如空间数据聚类、空间邻近关系^[16]，故本文应用Delaunay三角网获取待插值点周边的一阶邻近点^[1]，使得参考样点较为均匀地分布在待插值点的周围。

kAIDW插值算法的具体分类步骤如下：(1)移除原始样点中压盖的空间点；(2)计算全部样点属性值的均值Z_mean与中误差Z_std；(3)选取样点中任意点P_i，计算P_i对应属性值z_i的H_i与L_i；(4)如果z_i > H_i，P_i的类别赋为高，如果z_i < L_i，P_i的类别赋为低，如果 L_i$ \le $z_i$ \le $H_i，P_i的类别赋为中；(5)重复步骤(3)、(4)，直至所有样点均分类完成。

通过上述分类步骤，将§3.2实验2中的降雨数值离散化得到如图 2所示的类别数据。从图 2中可以发现，该地区降雨量主要可以分成3大区域，经过地理探测器^[17]软件探测得到q-统计值为0.65，且通过了显著性水平为0.05的分层异质性假设检验，这说明该地区降雨存在较明显的空间分层异质性^[18]，同时由图 2还可以发现该区存在空间局部异质性。因此，需要算法能够针对层内与层间做出不同的定权处理，降低异质性样点数据对插值结果造成的影响。

2.2.2 定权

常规的IDW插值算法一般采用距离倒数的平方定权，本文中不仅要考虑参考样点与待插值点之间的空间自相关，而且还需顾及一阶邻近样点与待插值点之间的空间异质性。因此定权过程的主要任务是在反距离平方定权的基础上自适应确定权重调和因子，使得最终的权重不仅考虑到样点与待插值点之间的空间相关性，还要顾及空间异质性。算法过程如下：(1)获取待插值点一阶邻近点的分类结果，统计一阶邻近点属性值分类为高、中、低各自的数量k₁、k₂、k₃。(2)定权，根据待插值点一阶邻近点中高、中、低3种类别点的数量，利用k‐NN算法“少数服从多数”的原则确定待插值点的属性类别。

k‐NN算法得到待插值点的预测值可能与真值不一致，如图 3所示，未知点实为“低”类，但预测值为“高”类，这将造成参考样点与待插值点异质性的误判。为了削弱误分类对插值结果造成的偏差，kAIDW算法将3类点各自的权重调和因子a_i设置成与样点数量成正比例的关系，而不是将异质参考样点的权重设置为0。即：

$ {a}_{i}=\left({k}_{i}/W\right)/{D}_{i}^{2} $

(5)

式中，W为一阶邻近点中的总参考样点数；D_i表示待插值点与参考样点之间的欧氏距离。从式(5)中可以发现，权重调和因子$ {\alpha }_{i}$的确定同样是自适应的，无需人工干预，在空间相关性的基础上顾及了空间异质性对插值结果的影响。

2.2.3 插值

1) 模型表达式

由于参考样点中可能存在空间异质性数据，特别是距离待插值点较近的空间异质点会使常规IDW插值算法得到的结果偏离真值较大。基于此，本文在反距离平方定权的基础上，提出了一种顾及空间异质性的IDW插值模型:

$ z\prime \left(x\right)=\frac{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{{k}_{1}}{{d}_{i}^{n}}\cdot {z}_{i}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{{k}_{2}}{{d}_{j}^{n}}\cdot {z}_{j}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{{k}_{3}}{{d}_{k}^{n}}\cdot {z}_{k}}{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{{k}_{1}}{{d}_{i}^{n}}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{{k}_{2}}{{d}_{j}^{n}}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{{k}_{3}}{{d}_{k}^{n}}} $

(6)

式中，$ z\prime \left(x\right) $为kAIDW算法计算得到的待插值点属性值; z_i、z_j、z_k分别为高、中、低3类点各自对应的属性值；d_i、d_j、d_k分别为高、中、低3类点各自对应的欧氏距离。

2) 算法特点分析

由于kAIDW插值模型不仅考虑了待插值点与样点之间的空间相关性，而且还顾及了空间异质性对插值结果的影响，因此能够得到更高的插值精度。

首先将式(1)按照式(6)的形式改写，则常规IDW插值算法得到的待插值点属性值为：

$ z\left(x\right)=\frac{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{1}{{d}_{i}^{n}}\cdot {z}_{i}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{1}{{d}_{j}^{n}}\cdot {z}_{j}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{1}{{d}_{k}^{n}}\cdot {z}_{k}}{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{1}{{d}_{i}^{n}}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{1}{{d}_{j}^{n}}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{1}{{d}_{k}^{n}}} $

(7)

(1) 若$ {k}_{1}>{k}_{2}>{k}_{3}\ne 0, $且$ {d}_{i} <{d}_{j} <{d}_{k} $，根据k-NN算法的“少数服从多数”的原则，那么待插值点的类别即为高类。将式(6)分子与分母同时除以k₁，得到kAIDW插值模型的变形式$ z\prime \left(x\right)$，即：

$ \begin{array}{c}z\prime \left(x\right)=\\ \mathrm{ }\frac{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{1}{{d}_{i}^{n}}\cdot {z}_{i}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{1}{{d}_{j}^{n}}\cdot {z}_{j}\cdot \frac{{k}_{2}}{{k}_{1}}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{1}{{d}_{k}^{n}}\cdot {z}_{k}\cdot \frac{{k}_{3}}{{k}_{1}}}{\stackrel{{k}_{1}}{\sum\limits_{i=1}}\frac{1}{{d}_{i}^{n}}+\stackrel{{k}_{2}}{\sum\limits_{j=1}}\frac{1}{{d}_{j}^{n}}\cdot \frac{{k}_{2}}{{k}_{1}}+\stackrel{{k}_{3}}{\sum\limits_{k=1}}\frac{1}{{d}_{k}^{n}}\cdot \frac{{k}_{3}}{{k}_{1}}}\end{array} $

(8)

对比式(7)与式(8)，可以明显看出式(8)中高类点的权重比式(7)中高类点的权重要大。同时从式(8)还可以看出，高类点的权重比中、低类点在等距离的前提下权重都要大。反之，即使d_i > d_j > d_k，由于kAIDW插值模型考虑了空间数据存在异质性的情况，通过调和因子k₁、k₂、k₃，同样可以使插值结果的精度高于传统IDW插值算法，从而可以发现，顾及空间异质性的kAIDW插值算法优于选择同样参考点的IDW插值算法，如传统的IDW、自适应的反距离权重(adaptive-IDW，AIDW)、四方向反距离权重(four-IDW，FIDW)等。

(2) 若k₁、k₂、k₃中存在两类点的数量为0，如图 2中的同一层间(待插值点位于“中”或“高”区域的内部)，则待插值点周边一阶邻近点中不存在空间异质性，那么式(6)将等价于AIDW插值算法，其插值精度仍高于常规的IDW插值算法^[1]。

(3) 若$ {k}_{1}={k}_{2}={k}_{3}\ne 0 $，式(6)同样等价于AIDW插值算法，通过(2)、(3)可以发现，AIDW插值算法在一定程度上属于kAIDW插值算法的一种特例。

(4) 若k₁、k₂、k₃存在其中一种类别为0，即图 2的不同层间(即待插值点位于“中”与“高”区域之间)，其分析结果与(1)一样。

通过以上分析可以发现：①针对样本空间点数据中存在空间异质性的情形，kAIDW插值算法在同一层间仅考虑空间相关性，在不同层间或局部存在空间异质性区域能够同时顾及空间异质与相关性对插值结果的影响。②kAIDW插值精度并不是在任何插值点上都高于传统的IDW插值算法，kAIDW插值算法的精度能否高于传统的关键在于k-NN算法对待插值点属性类别的分类是否准确，这也成为阻碍该IDW插值算法精度提高的关键因素，也是下一步重点研究的区域。

3 实验对比分析

本文设计两个不同应用领域的实验来验证kAIDW算法的可行性，将kAIDW算法与CIDW、FIDW、AIDW算法进行对比分析，其中CIDW与FIDW参考点数量均设置为4~12个，kAIDW与AIDW待插值点的参考点数量为一阶邻近点，4种IDW插值算法幂指数均取值为2，检查点为除了研究区域边界点以外的样点，实验方案采用留一交叉验证法^[10]，选取平均绝对误差(mean absolute error，MAE)、均方根误差(root mean square error，RMSE)、最大误差(MAX)、最小误差(MIN)4个指标对实验结果进行分析。

3.1 实验1：DEM插值

将上述4种算法应用于高程的插值研究当中，实验数据来源于山东黄河大堤的大比例尺地形图中，研究区域呈带状分布，形成中高两低的阶梯，通过野外采集了308个高程点，平均高程为42.90 m，最大高程为49.01 m，最低处为36.30 m，高程中误差为3 m，相对中误差为7%，实验结果如图 4所示。

从图 4(a)中可以发现，CIDW、FIDW插值算法的MAE基本都随着插值点数的增多而逐渐下降，在选择12个采样点时达到最小，分别为1.776 m、1.718 m；AIDW插值精度比CIDW、FIDW基本都要高，为1.721 m，而顾及空间异质性的kAIDW插值算法的MAE在4种算法中最低，为1.689 m。

从图 4(b)中可以发现，CIDW、FIDW插值算法的RMSE随着插值点数的增多基本呈对数方式下降，其中FIDW下降幅度更快；在参考点数量达到10个以上时，FIDW的RMSE在4种方法中最低，这是由于该区域高程值变化较为稳定，高程相对误差仅为7%，待插值点周边区域中不存在变化较大的高程点数据，而FIDW插值算法在4个方向的样点数量较为均匀，其RMSE随参考点数量增多而稳定下降直至稳定；AIDW与kAIDW插值算法在RMSE方面相差较小，约为2 mm。

从图 4(c)中可以发现，CIDW与FIDW插值算法在参考点数量大于7时，MAX值围绕AIDW与kAIDW呈类正弦波动，在参考点达到11个以上时，FIDW略优于其余3种插值方法，AIDW与kAIDW插值算法得到的MAX值结果基本一致。经分析发现，kAIDW的MAX值发生的待插值点的真实高程，通过图 1所示的分类步骤结果为高类点，但该点周边一阶邻近点中高类点的权重调和因子仅占1/6(即5个中类点，1个高类点)，从而造成该点kAIDW算法的MAX值相比于FIDW算法要略大一些。

从图 4(d)中可以发现，CIDW与FIDW的MIN值随着样点数的增多呈震荡形式，且CIDW较FIDW小；AIDW与kAIDW插值算法的MIN值保持一致，这是由于待插值点一阶邻近点中不存在异质点，kAIDW与AIDW定权一致；kAIDW与CIDW插值算法相比，当样点的数量为7、12，其MIN值超过CIDW约1 mm，这是由于该插值点附近地形平坦，空间自相关性较强，当参考样点的个数达到合适的条件，将使得CIDW法的MIN值较小。

通过对图 4进行分析可知：(1)不存在一种算法的插值精度在任何数据分布上MAE、RMSE、MAX、MIN的4个指标都占绝对的优势。(2)顾及样点分布均匀性的FIDW插值算法精度优于CIDW。(3)AIDW插值算法的精度通常要优于FIDW。(4)待插值点一阶邻近点中不存在空间异质性，kAIDW插值结果与AIDW算法保持一致，反之，在待插值点属性类别分类准确的前提下，kAIDW插值精度将明显高于AIDW。

3.2 实验2：降雨数据的插值

为了进一步验证kAIDW插值算法的适用性，将其应用于气象数据的插值研究中。实验采用SIC97数据中192个站点1986年5月8日降雨观测数据。192个站点的平均降雨量为19.8 mm，中误差为12.1 mm，相对中误差约为61.1%，最大降雨量为51.7 mm，最小降雨量为1.7 mm，实验结果如图 5所示。

从图 5中可以发现，kAIDW插值算法在4~12个点中的MAE、RMSE、MAX值均完全优于CIDW、FIDW、AIDW 3种插值算法，其中kAIDW插值算法与AIDW的MIN值相等，CIDW与FIDW的MIN值围绕kAIDW上下波动，且4种插值方法的MIN值相差较小，极差为0.79 mm。基于该实验结果的分析，可以发现kAIDW得到的实验结果优于其余3种插值算法。

除此之外，通过对比DEM与降水实验的数据，瑞士192个站点的降水值相对中误差为61.1%，而DEM实验高程数据的相对中误差仅为7%，根据离散方差大则其空间异质性强的原理^[19]，实验2数据的空间异质性更为明显，但kAIDW算法在实验2的实验结果优于实验1，从而说明顾及空间异质性的kAIDW算法在一定程度上较另3种IDW插值算法的抗差能力更强。

4 结语

本文提出了一种顾及空间异质性的自适应IDW插值算法，该算法简便、易于程序实现，通过与其他3种IDW插值算法进行对比分析后发现：(1)kAIDW算法具备良好的自适应特征，用户不需指定任何初始化参数，降低用户使用IDW算法的知识储备要求。(2)应用机器学习监督分类思想的kAIDW算法能够有效地对空间异质性数据进行处理，使得其抗大粗差能力更强，能够有效地提高插值精准度。然而kAIDW算法并不是在任何插值点上都是有效的，算法的适用条件在于k-NN算法在插值场景中对待插值点属性类别分类的准确率，因此未来关于kAIDW算法的研究工作主要集中在两个方面：(1)将多种机器学习算法引入到插值算法的研究当中，进一步提高待插点属性分类的准确率。(2)根据可靠性更高的判别结果进一步加大同类样点间的权重，从而降低异质样点对最终插值结果的影响。