CN111985144B - 一种地学数据多参数协同优化的idw插值方法 - Google Patents

Abstract

本发明公开了一种地学数据多参数协同优化的IDW插值方法，包括以下步骤：A、给定粒子群算法的初始化条件，搜索维度D＝4(α，λ，θ，N)，α为距离衰减系数，λ为距离调节参数，θ为各向异性方向，N为最邻近点数，采用随机初始化种群的方法产生一组包含(α，λ，θ，N)的解集合；B、计算粒子群的适应度值；C、获取每个粒子的个体最优值pbest_i；D、获取群体的全局最优值；E、对粒子的速度和位置进行更新；F、如未满足粒子群算法的终止条件，则返回步骤B，当满足粒子群算法的终止条件，则终止更新；G、输出最佳位置；H、将优化后的参数组合解(α，λ，θ，N)代入IDW插值模型中获取待插值点属性值。本发明能够解决现有技术的不足，获得IDW插值效果在全局意义下的满意解。

Description

一种地学数据多参数协同优化的IDW插值方法

技术领域

本发明属于地学数据处理技术领域，具体是一种地学数据多参数协同优化的IDW插值方法。

背景技术

在地学研究中，通常需要具有一定规则的面状数据。然而，多数可获得的地学数据都是站点采集的，呈现离散和不规则的分布。这就经常需要利用空间插值技术把不规则的离散数据转换成规则的面状空间数据。空间插值是通过有限的离散采样点来建立某种插值函数关系，并把待插值点一定范围内的已知采样点代入该函数表达式来获得该插值点的属性值。常用的空间插值方法包括克里金插值方法、自然邻域法、趋势面法及反距离加权插值(Inverse distance weighting,IDW)等。由于IDW插值原理简单、计算简便且符合地理学第一定律，被认为是GIS领域中一种标准的空间插值方法，在各学科领域都得到了广泛应用。然而，目前的IDW插值算法大多仅考虑单参数的调优，或对各参数独立调优，难实现插值模型的整体优化。此外，传统的IDW插值算法没有顾及各向异性对空间邻近度的影响。

发明内容

本发明要解决的技术问题是提供一种地学数据多参数协同优化的IDW插值方法，能够解决现有技术的不足，获得IDW插值效果在全局意义下的满意解。

本发明的内容包括以下步骤，

A、给定粒子群算法的初始化条件，搜索维度D＝4(α，λ，θ，N)，α为距离衰减系数，λ为距离调节参数，θ为各向异性方向，N为最邻近点数，采用随机初始化种群的方法产生一组包含(α，λ，θ，N)的解集合；

B、计算粒子群的适应度值；

C、获取每个粒子的个体最优值pbest_i；

D、获取群体的全局最优值；

E、对粒子的速度和位置进行更新；

F、如未满足粒子群算法的终止条件，则返回步骤B，当满足粒子群算法的终止条件，则终止更新；

G、输出最佳位置；

H、将优化后的参数组合解(α，λ，θ，N)代入IDW插值模型中获取待插值点属性值。

作为优选，步骤A中，α、λ、θ和N的搜索区间分别为[1，5]，[0.01，100]，[0，2π]和[4，30]。

作为优选，步骤B中，计算粒子的适应度值的方法为，

将每个粒子的x_i代入得到每个粒子的f(x_i)，f(x_i)为适应度函数，x_i表示第i个粒子的位置(α_i，N_i，λ_i，θ_i)，M为研究区域内已知样本点的个数，k代表第k个已知样本点，Z_k表示第k个样本点的值，/>表示估计待插值点的属性时不包括待插值点本身在内，即只根据待插值点周边的样本点进行插值；

把不同的参数组合x_i代入适应度函数绘制成曲面，得到arg min f(x_i)所对应的参数组合，即最邻近点数，距离衰减系数、距离调节及各向异性的方向参数，x_i＝arg min f(x_i)，其中argmin为使其后面式子达到最小值时的变量取值。

作为优选，步骤C中，对每个粒子，用其f(x_i)和值和个体历史最优位置P_best(i)比较，如果f(x_i)优于P_best(i)，则用当前待优化参数的位置x_i更新个体历史最佳位置P_best(i)。

作为优选，步骤D中，对每个粒子，将其当前f(x_i)与全局最佳位置G_best对应的适应度值做比较，如果当前的适应值更佳，则用当前粒子的位置更新全局最佳位置G_best。

作为优选，步骤E中，在D维空间中，D＝4，设定粒子群数为Q，粒子i的位置表示为x_i＝(x_i1，x_i2，...，x_iD)，粒子i的速度为v_i＝(v_i1，v_i2，...，v_iD)，pbest_i＝(p_i1，p_i2，...，p_iD)表示粒子i经历过的最佳位置，G_best＝(g₁，g₂，...，g_D)表示所有粒子的全局最佳位置，粒子i的第d维位置和速度更新方法分别为，

表示粒子i在第t+1次迭代的第d维分量的速度/>表示粒子i在第t+1次迭代的第d维分量的位置；C₁，C₂为学习因子，调节学习步长，C₁和C₂均优选为1.48；r₁，r₂为两个随机变量，取值范围为[0,1]，以增加搜索随机性；W表示惯性权重，调节对解空间的搜索范围；

惯性权重W采用线性递减权值，计算方式如下，

w_max表示最大惯性权重，取值为1.5，w_min为最小惯性权重，取值为0.5，iterations表示当前迭代次数，iterations_max表示最大迭代次数，iterations_max设置为200*D，其中D＝4。

作为优选，步骤F中，粒子群算法的终止条件为，迭代次数达到200或者最佳适应度值的增量比值少于1e^-6。

作为优选，步骤G中，最佳位置x_i为使得适应度函数值最佳时的距离衰减系数、最邻近点数、距离调节及各向异性方向参数。

作为优选，步骤H中，IDW插值模型为，

i代表任意待插值点，Z_j为点j处的实测值，S_ij是未知点i与已知点j之间的距离，N为最邻近点数，α为距离衰减系数。

作为优选，(X_i，Y_i)和(X_j，Y_j)表示点i与j的笛卡尔坐标；

当λ等于1，那么S_ij将等价于欧氏空间距离；当λ＞1，如果两点之间ΔY越大，那么两点Y方向之间的距离将更远，即Y方向的空间邻近性减弱；反之，当λ＜1，那么两点之间X方向的距离相对将增大，即两点的X方向的空间邻近性将减弱；

进一步引入方向参数θ来描述一般情形下的各向异性，建立原始坐标与反映各向异性坐标的转换关系，

E_ij表示待插值点i与样本点j之间顾及各向异性的距离，λ用于调节坐标系统下不同坐标轴的距离，θ∈(0，2π)表示原始坐标与各向异性坐标系的旋转角度，ΔX与ΔY为原始数据中待插值点与已知样本点之间的欧氏距离；

在步骤H中使用IDW插值模型时，使用E_ij替换原有的S_ij。

本发明的有益效果是：本发明首先引入具有物理意义的距离及方向调节参数来反映地理现象的各向异性，把典型的欧氏距离扩展为各向异性的“椭圆”距离，随后基于CV交叉检验方法构造了适合于本发明插值算法(PIDW插值算法)的适应度函数，并采用粒子群算法对IDW插值模型的多个参数进行协同优化，实现IDW插值效果最佳的目的。

附图说明

图1为本发明的流程图。

图2为年均气温插值实验中本发插值算法和现有6中插值算法在MAE指标下的插值精度的对比图。

图3为年均气温插值实验中本发插值算法和现有6中插值算法在RMSE指标下的插值精度的对比图。

图4为年均气温插值实验中本发插值算法和现有6中插值算法在MAX指标下的插值精度的对比图。

图5为DEM插值数据下本发插值算法和现有6中插值算法在MAE指标下的插值精度的对比图。

图6为DEM插值数据下本发插值算法和现有6中插值算法在RMSE指标下的插值精度的对比图。

图7为DEM插值数据下本发插值算法和现有6中插值算法在MAX指标下的插值精度的对比图。

图8为各向同性环境下的空间邻近度示意图。

图9为各向异性环境下的空间邻近度示意图。

图10为各向异性方向与坐标轴重合情况下的空间邻近度示意图。

图11为各向异性方向与坐标轴存在夹角情况下的空间邻近度示意图。

图12为本发明在5组实验中RMSE和MAE指标的对比图。

图13为本发明在5组实验中CV残差的对比图。

图14为本发明各向同性环境下最邻近点的选择。

图15为本发明各向异性环境下最邻近点的选择。

具体实施方式

参照图1，一种地学数据多参数协同优化的IDW插值方法，包括以下步骤，

A、给定粒子群算法的初始化条件，搜索维度D＝4(α，λ，θ，N)，α为距离衰减系数，λ为距离调节参数，θ为各向异性方向，N为最邻近点数，采用随机初始化种群的方法产生一组包含(α，λ，θ，N)的解集合；

B、计算粒子群的适应度值；

C、获取每个粒子的个体最优值pbest_i；

D、获取群体的全局最优值；

E、对粒子的速度和位置进行更新；

F、如未满足粒子群算法的终止条件，则返回步骤B，当满足粒子群算法的终止条件，则终止更新；

G、输出最佳位置；

H、将优化后的参数组合解(α，λ，θ，N)代入IDW插值模型中获取待插值点属性值。

步骤A中，α、λ、θ和N的搜索区间分别为[1，5]，[0.01，100]，[0，2π]和[4，30]。

步骤B中，计算粒子的适应度值的方法为，

将每个粒子的x_i代入得到每个粒子的f(x_i)，f(x_i)为适应度函数，x_i表示第i个粒子的位置(α_i，N_i，λ_i，θ_i)，M为研究区域内已知样本点的个数，k代表第k个已知样本点，Z_k表示第k个样本点的值，/>表示估计待插值点的属性时不包括待插值点本身在内，即只根据待插值点周边的样本点进行插值；

把不同的参数组合x_i代入适应度函数绘制成曲面，得到arg min f(x_i)所对应的参数组合，即最邻近点数，距离衰减系数、距离调节及各向异性的方向参数，x_i＝arg min f(x_i)，其中argmin为使其后面式子达到最小值时的变量取值。

步骤C中，对每个粒子，用其f(x_i)和值和个体历史最优位置P_best(i)比较，如果f(x_i)优于P_best(i)，则用当前待优化参数的位置x_i更新个体历史最佳位置P_best(i)。

步骤D中，对每个粒子，将其当前f(x_i)与全局最佳位置G_best对应的适应度值做比较，如果当前的适应值更佳，则用当前粒子的位置更新全局最佳位置G_best。

步骤E中，在D维空间中，D＝4，设定粒子群数为Q，粒子i的位置表示为x_i＝(x_i1，x_i2，...，x_iD)，粒子i的速度为v_i＝(v_i1，v_i2，...，v_iD)，pbest_i＝(p_i1，p_i2，...，p_iD)表示粒子i经历过的最佳位置，G_best＝(g₁，g₂，...，g_D)表示所有粒子的全局最佳位置，粒子i的第d维位置和速度更新方法分别为，

表示粒子i在第t+1次迭代的第d维分量的速度，/>表示粒子i在第t+1次迭代的第d维分量的位置；C₁，C₂为学习因子，调节学习步长，C₁和C₂均优选为1.48；r₁，r₂为两个随机变量，取值范围为[0,1]，以增加搜索随机性；W表示惯性权重，调节对解空间的搜索范围；

惯性权重W采用线性递减权值，计算方式如下，

w_max表示最大惯性权重，取值为1.5，w_min为最小惯性权重，取值为0.5，iterations表示当前迭代次数，iterations_max表示最大迭代次数，iterations_max设置为200*D，其中D＝4。

步骤F中，粒子群算法的终止条件为，迭代次数达到200或者最佳适应度值的增量比值少于1e^-6。

步骤G中，最佳位置x_i为使得适应度函数值最佳时的距离衰减系数、最邻近点数、距离调节及各向异性方向参数。

步骤H中，IDW插值模型为，

i代表任意待插值点，z_j为点j处的实测值，S_ij是未知点i与已知点j之间的距离，N为最邻近点数，α为距离衰减系数。

(X_i，Y_i)和(X_j，Y_j)表示点i与j的笛卡尔坐标；

当λ等于1，那么S_ij将等价于欧氏空间距离；当λ＞1，如果两点之间ΔY越大，那么两点Y方向之间的距离将更远，即Y方向的空间邻近性减弱；反之，当λ＜1，那么两点之间X方向的距离相对将增大，即两点的X方向的空间邻近性将减弱；如图14和15所示，距离调节参数λ能够改变最邻近点的选择，其中图14为均质空间下所选择的空间邻近点，图15为顾及各向异性时选择的空间邻近点。

进一步引入方向参数θ来描述一般情形下的各向异性，建立原始坐标与反映各向异性坐标的转换关系，

E_ij表示待插值点i与样本点j之间顾及各向异性的距离，λ用于调节坐标系统下不同坐标轴的距离，θ∈(0，2π)表示原始坐标与各向异性坐标系的旋转角度，ΔX与ΔY为原始数据中待插值点与已知样本点之间的欧氏距离；

在步骤H中使用IDW插值模型时，使用E_ij替换原有的S_ij。

实验验证

使用不同分辨率的全国气温站点与黄河堤坝高程数据来验证本发明PIDW算法的可行性，并将PIDW插值算法与经典的反距离平方加权(CIDW,Classical IDW)、四象限IDW(FIDW，Four quadrant IDW)、自适应的反距离平方权重(Adaptive-IDW,AIDW)、k近邻反距离平方加权(k-nearest neighbor Adaptive IDW,kAIDW)、普通克里金(OrdinaryKriging，OK)以及顾及各向异性的克里金(Anisotropic Kriging,AnisOK)进行对比实验。

实验中涉及的各种方法参数设置如下：CIDW、FIDW以及OK参考点数量分别设置为4～20个，KAIDW与AIDW待插值点的参考点数为待插值点的一阶邻近点。CIDW、FIDW、AIDW、KAIDW的距离衰减系数均取值为2，OK与AnisOK插值算法的变异函数采用最优化方式拟合得到块金系数、基台值、变程，其中AnisOK分别在0°、45°、90°以及135°四个方向计算空间结构的各向异性。利用平均绝对误差(Mean Absolute Error，MAE)、均方根误差(Root MeanSquare Error，RMSE)及最大误差(MAX)3个定量指标对插值精度进行评价，其中MAE、RMSE与MAX的定义如下：

式中T_i为测试点的真实属性值，为测试点的属性估计值，n代表测试点的个数。

年均气温插值实验

数据来源于中国气象中心2005年的682个气象站点的观测数据，通过前期的数据预处理，剔除少量观测台站仪器出现故障的数据，对每个台站的日平均气温求和，最终获取各个台站的年平均气温。由于气温与海拔存在约0.6℃/100m的比例关系，本文根据气象站台的海拔数据，将每个气象站台的年平均气温归算到0海拔面。这样处理的理由是剔除不同高程对年均气温的影响，着重讨论经纬度变化对气温差异的影响。

为减少样本选择的影响，本实验重复进行了3次随机抽样，其中每次从总体中抽取30个作为测试集，其余作为训练集。经过PSO优化，PIDW的距离衰减系数α为1，最邻近点数N为7，距离调节参数λ为0.26，各向异性方向θ为0.02°。把计算结果的平均值作为插值精度的最终评价值，结果详见图2-4。

从图2可以发现：CIDW与FIDW插值算法的MAE随着插值点数的增多呈现震荡式变化，并且FIDW的MAE显著低于CIDW。KAIDW、AIDW与FIDW的MAE大致保持一致，但OK和PIDW的MAE显著低于KAIDW、AIDW、FIDW和CIDW。OK的MAE低于PIDW，但最大偏差仅为0.007℃。AnisOK的MAE误差为0.52，低于PIDW的0.54。从图3中可以得出：CIDW、FIDW和OK的RMSE随着点数的增多而下降，在参考点数量为13以上时趋于稳定。OK的RMSE稍好于AIDW与KAIDW方法，但劣于CIDW和FIDW算法。AnisOK的RMSE误差显著低于CIDW、FIDW、OK、AIDW、KAIDW。PIDW的RMSE为0.78℃，比次优的AnisOK方法低0.02℃。图4中的MAX反映了类似的规律。CIDW与FIDW插值算法随着点数的增多稳步下降，其中FIDW的MAX略低于CIDW，而AIDW与KAIDW的MAX表现最差。OK的MAX总体上优于AIDW和KAIDW，但明显差于CIDW和FIDW算法。PIDW的MAX显著低于CIDW、FIDW、AIDW、KAIDW以及OK。AnisOK的MAX具备一定的优势，相较于次优的PIDW方法，其最大值误差还要低0.4℃。

上述结果表明：1)顾及各向异的PIDW和AnisOK都显著优于CIDW、FIDW、OK、AIDW和KAIDW算法；2)不存在哪种插值算法，其MAE、RMSE与MAX的3个指标上都占据一致性的优势；3)顾及样点分布均匀性的FIDW插值算法精度普遍优于CIDW；4)采用一阶邻近点作为参考样点的AIDW与KAIDW插值算法并未体现出较大的优势，可能是由测试数据集中一些位于边界区域的测试点所导致。顾及空间异质性的KAIDW未呈现比AIDW插值算法的优势，可能是由于KAIDW对于局部空间异质性的分类存在较多的误判。5)OK相比于CIDW、FIDW、AIDW、KAIDW总体上并未呈现出明显优势，这可能是由于研究区域不能较好的满足OK插值算法所需的平稳性假设条件，导致变异函数未能较好的反映气温变化的变异规律。6)PIDW插值算法除在MAE指标稍劣于OK插值算法，在其余指标上相较于未顾及各向异性影响插值算法具备明显的优势。7)除RMSE指标外，AnisOK在MAE和MAX指标上低于PIDW插值算法。

DEM插值

DEM插值数据来源于黄河山东济南段大堤的1:500大比例尺地形图，共采集了639个高程点。黄河大堤为“梯形”构造，其中此段大堤南北方向起伏大，而东西向高程则变化平缓，其中平均高程为43.03m，最大高程为49.30m，位于大堤的顶部，最低高程36.31m。

实验中各算法的初始值与气温数据保持一致，PIDW算法得到最终的优化参数结果为[1.28，6，99，0]，其中距离衰减系数为1.28，最邻近点数为6，距离调节系数为99，各向异性方向值为0°，各种算法的结果对比如图5-7所示。

从图5中可以发现：随着插值点数的增多，CIDW、FIDW与OK的MAE逐渐下降，在采样点为12时达到最低，分别为1.43m，1.37m，1.42m。AIDW与KAIDW的MAE低于CIDW算法，但大于其余三种算法。AnisOK的MAE误差明显低于CIDW、FIDW、OK、AIDW以及KAIDW。PIDW表现出明显的优势，其MAE在七种插值算法中的误差最小，比次优AnisOK的最小MAE值再少27.3cm。图6中的RMSE表现出类似的结果。CIDW、FIDW和OK的RMSE随着点数的增多呈现下降，在样点为14以上时基本达到稳定状态，并且OK的RMSE值最小。AIDW与KAIDW的RMSE的结果最差。AnisOK的RMSE低于除PIDW外的5种插值方法。PIDW的RMSE指标在所有算法中同样呈现出明显的优势，比次优的AnisOK低27cm左右。图7表明CIDW、FIDW与OK的MAX随着点数增多而呈现下降趋势，在参考点达到17以上基本保持稳定。AIDW和KAIDW比CIDW与FIDW能够取得更优的效果，但PIDW的MAX显著优于其它插值方法。

总结以上两个实验结果，可得出以下结论：(1)顾及各向异性的AnisOK与PIDW几乎在所有指标中都取得最小的误差，表明在计算过程中顾及到各向异性的空间属性，将有助于提高模型的插值精度；(2)总体上PIDW的插值效果稍好于AnisOK方法，而且PIDW方法对数据无需假设前提，但地统计的AnisOK方法要求数据具有空间的二阶平稳性。

参数分析

为了解释PIDW插值算法有效性的原因，下面对模型的参数优化作进一步实验分析。

1)各向异性距离调节系数对空间邻近性的影响

由于各向异性能够引起空间邻近度的改变，因此设计对照实验来分析各向异性对空间邻近度的影响。图8和图9分别现实未顾及各向异性和顾及到各向异性时对空间邻近度的影响。由图8可知，当不顾及x或y方向的各向异性时，样本点与参考点之间的距离为同心圆模式，而当顾及各向异性的影响时，“同心圆”距离模式被压缩为“椭圆”，位于东西方向的样本点具有更强的空间邻近性，这表明PIDW的距离调节参数符合实验2中堤坝高程南北方向起伏大，而东西向变化平缓的事实。

2)各向异性的方向对空间邻近性的影响

考虑到各向异性的方向并不一定正好沿着坐标轴方向，而是存在一定的角度θ。本文进一步设计对照实验来分析各向异性的方向对空间邻近性的影响。将数据逆时针旋转45°，即变换后的各向异性方向X_New与Y_Old之间的夹角约为135°，如图10-11所示。

3)多参数优化对比分析

影响PIDW插值精度的参数主要有各向异性(距离调节参数λ与各向异性方向参数θ)，距离衰减系数α和最邻近点数N。本文基于实验1的数据，设计如下5个实验来表明PIDW的插值精度由各参数共同决定，而单参数或少量参数的优化难以取得插值的整体满意效果。(a)固定距离调节参数为1，方向值为0，优化最邻近点数与距离衰减系数；(b)固定距离衰减系数为2(ArcMap软件默认值)，优化最邻近点数、距离调节及各向异性方向参数；(c)固定最邻近点数为12(ArcMap软件默认值)，优化距离调节、各向异性方向参数及距离衰减系数；(d)固定方向值为0，同时优化距离衰减系数、最邻近点数及距离调节参数；(e)距离衰减系数、最邻近点数、距离调节与各向异性方向参数协同优化。

实验采用MAE、RMSE和CV残差作为精度的评价指标，结果如图12-13所示。从图中可以看出方案(a)的效果最差，其MAE、RMSE以及训练样本的CV残差值都最大。这是由于气温的分布存在较强的各向异性特征，如果插值过程中将距离调节参数固定为1，导致重要的参数不能正确地反映真实情况。方案(b)和(c)的MAE、RMSE和CV相较于方案(a)都有一定程度的提升，且方案(c)略优于方案(b)，表明本例中最邻近点数对影响插值精度的重要性要低于距离衰减系数。方案(d)的CV残差低于方案(a)、(b)和(c)，表明顾及各向异性影响能够得到较好的插值效果。从图中还可以发现方案(e)的预测结果几乎与方案(d)完全一致，仅有CV的残差和相差0.1。这是由于气温变化沿着纬度方向变化速率较之于经度更快，而方案(e)优化的各向异性方向与方案(d)所固定的方向相同，因此其插值结果一致。由于粒子群算法的随机性，导致它们的CV残差和存在微小的差异。以上实验表明如果仅优化部分参数均难以获得插值效果的整体满意解，采取多参数整体优化方法将有效增强的模型插值精度。由于PIDW算法也可能存在过拟合的现象，导致预测精度的降低，因而在实际应用中需要注意避免PIDW算法的过拟合问题。

Claims (8)

1.一种地学数据多参数协同优化的IDW插值方法，其特征在于包括以下步骤，

A、采集不同分辨率的全国气温站点与黄河堤坝高程数据；全国气温站点数据来源于中国气象中心2005年的682个气象站点的观测数据，通过数据预处理，剔除观测台站仪器出现故障的数据，对每个台站的日平均气温求和，最终获取各个台站的年平均气温；根据气象站台的海拔数据，将每个气象站台的年平均气温归算到0海拔面；黄河堤坝高程数据来源于黄河山东济南段大堤的1:500大比例尺地形图，共采集639个高程点，初始值与气温数据保持一致；给定粒子群算法的初始化条件，搜索维度D＝4(α，λ，θ，N)，α为距离衰减系数，λ为距离调节参数，θ为各向异性方向，N为最邻近点数，采用随机初始化种群的方法产生一组包含(α，λ，θ，N)的解集合；

B、计算粒子群的适应度值；

C、获取每个粒子的个体最优值pbest_i；

D、获取群体的全局最优值；

E、对粒子的速度和位置进行更新；

F、如未满足粒子群算法的终止条件，则返回步骤B，当满足粒子群算法的终止条件，则终止更新；

G、输出最佳位置；

H、将优化后的参数组合解(α，λ，θ，N)代入IDW插值模型中获取待插值点属性值；

IDW插值模型为，

i代表任意待插值点，Z_j为点j处的实测值，S_ij是未知点i与已知点j之间的距离，N为最邻近点数，α为距离衰减系数；

(X_i，Y_i)和(X_j，Y_j)表示点i与j的笛卡尔坐标；

当λ等于1，那么S_ij将等价于欧氏空间距离；当λ＞1，如果两点之间ΔY越大，那么两点Y方向之间的距离将更远，即Y方向的空间邻近性减弱；反之，当λ＜1，那么两点之间X方向的距离将增大，即两点的X方向的空间邻近性将减弱；

进一步引入方向参数θ来描述各向异性，建立原始坐标与反映各向异性坐标的转换关系，

E_ij表示待插值点i与样本点j之间顾及各向异性的距离，λ用于调节坐标系统下不同坐标轴的距离，θ∈(0，2π)表示原始坐标与各向异性坐标系的旋转角度，ΔX与ΔY为原始数据中待插值点与已知样本点之间的欧氏距离；

在步骤H中使用IDW插值模型时，使用E_ij替换原有的S_ij。

2.如权利要求1所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤A中，α、λ、θ和N的搜索区间分别为[1，5]，[0.01，100]，[0，2π]和[1，30]。

3.如权利要求1所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤B中，计算粒子的适应度值的方法为，

将每个粒子的x_i代入得到每个粒子的f(x_i)，f(x_i)为适应度函数，x_i表示第i个粒子的位置(α_i，N_i，λ_i，θ_i)，M为研究区域内已知样本点的个数，k代表第k个已知样本点，Z_k表示第k个样本点的值，/>表示估计待插值点的属性时不包括待插值点本身在内，即只根据待插值点周边的样本点进行插值；

把不同的参数组合x_i代入适应度函数绘制成曲面，得到arg min f(x_i)所对应的参数组合，即最邻近点数，距离衰减系数、距离调节及各向异性的方向参数，x_i＝arg min f(x_i)，其中argmin为使其后面式子达到最小值时的变量取值。

4.如权利要求3所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤C中，对每个粒子，用其f(x_i)和值和个体历史最优位置P_best(i)比较，如果f(x_i)优于P_best(i)，则用当前待优化参数的位置x_i更新个体历史最佳位置P_best(i)。

5.如权利要求4所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤D中，对每个粒子，将其当前f(x_i)与全局最佳位置G_best对应的适应度值做比较，如果当前的适应值更佳，则用当前粒子的位置更新全局最佳位置G_best。

6.如权利要求5所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤E中，在D维空间中，D＝4，设定粒子群数为Q，粒子i的位置表示为x_i＝(x_i1，x_i2，...，x_iD)，粒子i的速度为v_i＝(v_i1，v_i2，...，v_iD)，pbest_i＝(p_i1，p_i2，...，p_iD)表示粒子i经历过的最佳位置，G_best＝(g₁，g₂，...，g_D)表示所有粒子的全局最佳位置，粒子i的第d维位置和速度更新方法分别为，

表示粒子i在第t+1次迭代的第d维分量的速度，/>表示粒子i在第t+1次迭代的第d维分量的位置；C₁，C₂为学习因子，调节学习步长，C₁和C₂均为1.48；r₁，r₂为两个随机变量，取值范围为[0,1]，以增加搜索随机性；W表示惯性权重，调节对解空间的搜索范围；

惯性权重W采用线性递减权值，计算方式如下，

w_max表示最大惯性权重，取值为1.5，w_min为最小惯性权重，取值为0.5，iterations表示当前迭代次数，iterations_max表示最大迭代次数，iterations_max设置为200*D，其中D＝4。

7.如权利要求1所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤F中，粒子群算法的终止条件为，迭代次数达到200或者最佳适应度值的增量比值少于1e^-6。

8.如权利要求1所述的地学数据多参数协同优化的IDW插值方法，其特征在于：步骤G中，最佳位置x_i为使得适应度函数值最佳时的距离衰减系数、最邻近点数、距离调节及各向异性方向参数。