论文部分内容阅读
摘要：本文应用反距离权重法（IDW）、径向基函数（RBF）、局部多项式法（LPI）、普通克里金法（OK）等4种常用的空间插值方法对湖南某研究区土壤重金属Cd的采样数据做插值研究。结果表明：4种常用的插值方法中，局部多项式法和普通克里金法精度优于反距离权重法和径向基函数。插值效果与重金属测量采样点的数量和密度有关，重金属测量采样点的密度越大、数量越多，空间插值结果的精度就越高，研究结果为区域土壤重金属元素污染评价中空间插值方法的优选提供依据和实例。
关键词：土壤重金属；空间分布预计；空间插值；交叉验证；插值精度
中图分类号： X5 文献标识码： A DOI编号： 10.14025/j.cnki.jlny.2015.21.030
土壤重金属污染问题日趋严重，已威胁到农产品质量安全。对污染的农田进行分级评价和区域划分是有效治理土壤重金属污染和保障农产品质量安全的关键措施。
在目前的产地分级划分工作中，仍需要依赖区域采样，监测采样点的重金属含量数据，通过有限采样点的数据进行空间插值分析，估算重金属含量的连续分布情况，从而为后期的污染评价工作提供必要的数据支持。
常用空间插值方法分为以下两类[1]：确定性插值和地统计插值。其中确定性插值是基于样点之间的相似程度来创建一个连续函数。比较有代表性的确定性插值方法包括：反距离权重法（IDW）、样条插值法（Spline）和趋势面分析（TSA）。地统计方法是利用样点的统计特性，量化本点之间的空间自相关性，从而构建样点的空间结构模型，代表性的方法有各种克里金方法（Kriging）。
在空间插值分析过程中，插值方法的选择和参数的优选，极大地影响了插值结果的准确程度。通过交叉验证，比较平均误差和均方根误差，在反距离权重法的研究中，张锦明在研究反距离权重法参数对DEM插值误差影响时，发现权重指数在所有的参数中对插值效果影响最大，他在研究同时发现权重指数大于等于3的时候，对插值精度的影响不明显，2或者3是较好的选择[2]。王宇航在对反距离权重法的研究中也得出类似的结论[3]。
样条插值法（Spline）中，比较常用且有代表性的径向基函数（RBF）。RBF方法是综合5种函数方法，根据情况选择最优的组合产生的插值结果[4]。
与全局趋势面分析（TSA）相比，局部多项式方法对每一个局部区域（即窗口）做模拟，在每一个区域的模拟方法与全局的趋势面类似，精度上要优于全局趋势面分析（TSA）[5]。
克里金方法（Kriging）是一族方法的总称，它包括很多方法，常见的有：简单克里金、普通克里金、泛克里金、漂移克里金、回归克里金、分层克里金、协同克里金。Li综合比较多种克里金插值方法，提出综合精度和易操作性等各项条件，OK表现最好[6]。
本文选取湖南某地部分区域作为研究区域，对土壤重金属Cd数据进行了系统采样。分别采用反距离权重法（IDW）、径向基函数（RBF）、局部多项式法（LPI）、普通克里金法（OK）分别进行实验，根据交叉验证的结果比较各种方法之间的精度差异，优选适合产地分级的插值方法，同时将插值结果和该地区的实际污染情况做简单比较。
1 材料与方法
1.1 研究区域
选择某典型重金属污染区域作为研究区域，区域面积2132.8平方公里，位于中国的中部丘陵地带，气候主要土壤类型，成土母质为常年种植水稻。
1.2 样品采集和分析方法
在研究区域内，采用系统随机采样的方式采集样品942个。2014年，采用仪器取样方式采集土壤。我们选择该区域作为典型区域，采集土壤重金属Cd的含量。
所用仪器采用X射线荧光光谱技术测定土壤样品中Cd含量。
1.3 数据的处理
本文的Cd元素含量数据集的分析采用不同的软件包，数据的处理采用R语言，制图采用ArcGIS10.2完成。
1.4 空间插值方法
本研究中选用反距离权重法（IDW）、径向基函数（RBF）、局部多项式法（LPI）、普通克里金法（OK）这4种常用的空间插值方法对研究区域数据进行插值分析，对插值结果采用标准方法进行精度的比较。
空间插值是应用相同区域中已测样点的测量值，对未抽样点的属性值的一种预测。在预测点处的属性值都可以用下面的公式表示：■ 其中， Z*（X）表示预测点位置的预测属性值， Z（Xi）表示第i个采样点的权重值， 表示第i个采样点位置的属性值。下面分别就4种常用的插值方法展开介绍。
1.4.1 反距离权重法（IDW） 反距离权重法的基础是假设预测值是实测数据的线性组合，同时权指和距离的幂值成反比。公式表示如下：
其中，p是一个任意正实数，通常，p=2；是实测离散点到插值预测点的距离。距离越近，权值越大；距离越远，全值越小。P是用来刻画权值如何随着距离的变化而变变化的。
1.4.2 普通克里金法（OK） 克里金法是基于假设的插值参数可以被视为一个区域化变量。和反距离权重法一样，克里金法的预测值是实测数据的线性组合。根据随机领域的随机特性，可以选用不同类型的克里金方法。其中较常见的包括简单克里金法、普通克立金法、泛克里金法等。这些方法相同之处在于他们都属于无偏估计。
普通克里金的权指由半方差的计算间接得到，半方差是该方法中衡量各点之间空间相关程度的测度值。普通克里金的半方差估计公式，又称作变异函数，如下所示：
使用该方法的关键在于选取合适的变异函数模型。利用做出的半方差图找出与之拟合的最好的理论变异函数模型，可用于拟合的模型包括高斯模型、线性模型、球状模型、指数模型、圆形模型。
1.4.3 径向基函数（RBF） 径向基函数法是精确插值方法中的一种。从概念上讲，径向基函数法是寻求一条穿过实测数据点且总曲率最小的曲线函数。径向基函数的预测值由以下两部分，如下面的公式所示： 其中， ■表示径向基函数， dj表示样点到预测点的距离， fi（x）表示局势函数， 表示权重系数。根据已知点可以计算出ai和dj。
具体表示有以下5种情况：
规则样条函数（CRS）
薄板样条函数（TPS）
高次曲面函数（CRS）
反高次曲面函数（IMQ）
张力样条函数（ST）
其中，d表示样点到预测点的距离。c是一个平滑因子，■表示改进型Bessel函数，表示欧拉常数。
1.4.4 局部多项式 （LPI） 局部多项式插值的实质是搜寻一个多项式公式的过程，这个公式对应的图形经过所有的实测值点。公式表示如下：F（X，Y）=a+bY+cY+dXY+eX2+fY2
全局多项式采用一个多项式来模拟整个区域的表面，而局部多项式则是全局多项式方法和移动窗口方法的结合。与全局多项式方法不同的地方在于，局部多项式方法对每一个局部区域（即窗口）做模拟。确定窗口的大小原则是要保证每个窗口有足够的实测数据点。
1.5 模型评定方法与指标
交叉验证法是评价和比较不同插值方法精度的重要方法。因为样本数量是有限的，所以在本文中也同样选择了交叉验证方法来做插值方法的比较。具体做法是：剔除一个点，用剩余的测量值做插值，得到剔除点位置的预测值，两相比较。
这里介绍评价插值结果精度的两个常用指标： 平均误差（ME）和均方根误差（RMSE）。其中，平均误差反映估值的误差范围，均方根误差反映误差均值的变化范围。这两个值可作为插值方法优劣评价的标准。
其中，Ve，i 、Va，i 分别为样点的实际值和插值预测值，n表示样点的数量。
2 结果与分析
2.1 最优参数确定
插值方法中参数选择的不同会对插值结果产生影响，因此插值方法的比较应以参数的优化为前提。IDW和RBF的参数很少，易于挑选最优参数，而Kriging法的参数较为复杂，主要体现在变异函数的拟合上，通常需反复验证。在局部多项式插值方法应用过程中， 要根据实际需要选择核函数，一般在指数模型、多项式模型、高斯模型和高次模型中选择；本文IDW的权重系数根据最小误差的原则，取值为2；RBF法采用完全规则样条（completely regularized spline）使其实现平滑插值。
2.2 变异函数拟合
采用kolmogorov-Smirnov检验，对研究区的数据集分布特征是否符合正态分布展开讨论。可知两种变换对数据集的表现的差异较大。进一步定量分析可以得到结果。
其中只有Log变换后的结果中p值大于0.05，即满足正态分布。所以，这里的数据采用Log变换来做处理。
OK方法的插值精度受到变异函数影响最大，下面是我们使用研究区域的数据，分别选用4种常用插值模型做出的实验结果，4种模型下的误差，以半球模型为最小，而且4种插值的误差波动相对很接近。因此，对当前重金属Cd数据集来讲，在OK方法的变异函数中，选择半球模型可以取得最好的精度。
2.3 LPI方法参数的选择
在LPI方法的使用中，会涉及核函数的选择，而这个选择会明显影响插值结果的精度，所以，我们对几种常用的函数（包括指数函数、多项式函数、高斯函数和高次曲面函数）分别做了实验，得到结果。
指数模型的平均误差（ME）最小，而其均方根误差与其他的模型不相上下。于是，可知对当前情况而言，指数模型是最合适的选择。
除了这个核函数的选择以外，在局部多项式插值方法应用过程中，要根据实际需要调整搜索椭圆的参数（，R1和 R2 ），其中，表示搜索角度，R1和R2分别表示搜索椭圆长、短半径。
从空间自相关性的概念可知，空间上越靠近，属性就越相似。根据经验或专业背景找出这么一个阈值，作为邻近区域的半径。
2.4 插值结果的验证和分析
分别用4种空间插值方法对土壤重金属空间分布数据进行插值。
对实测数据进行分组，同比例均匀的取样点总数的10%，作为测试组数据。采用交叉验证法对4种插值方法进行精度比较，4种插值方法的MAE、RMSE结果。
从误差分析的结果来看，4种插值方法的精度存在差异。从平均误差ME来看，四种插值方法的误差大小顺序是：IDW>LPI>RBF>OK，即IDW误差最大，LPI和RBF居中，OK最小；从平均误差RMSE来看，四种插值方法的误差变化范围都在0.85到1.0之间，比较接近。因此，四种插值方法的精度优劣程度主要看平均误差，其顺序是OK>RBF>LPI>IDW，即OK最优，RBF和LPI居中，IDW最差。
就OK方法而言，土壤重金属Cd的含量表现出块金效应，表明人为活动对土壤元素的空间变异影响较大。同时就目前的采样数据了解到的信息可知，土壤重金属Cd的污染源分布东北方向较为集中、含量高，西南方向污染源较为分散，含量低，中间地带含量居中。当前的插值结果和实际检测结果反应的情况符合得较为一致。
3 结论
一是OK插值对于刻画区域土壤元素的空间分布具有一定的优势，但该方法操作比较复杂，特别是其半变异函数的模型及参数优选需要反复验证修改。相对而言，RBF的参数设置有简单、易操作的优点，同时方便确定最优参数。IDW精度不高，但算法简单。LPI进度相对较高，但操作相对复杂；二是对于湖南的该研究区，土壤重金属Cd含量的空间分布而言，OK的插值效果最优；三是基于各自的最优插值方法，湖南的该研究区土壤重金属Cd的空间分布可以较好地反映实际的空间相关性，同时插值分布图可以指示其污染分布情况。
参考文献
[1]朱求安，张万昌，余钧辉.基于GIS的空间插值方法研究[J].江西师范大学学报（自然科学版），2004，（02）.
[2]张锦明，郭丽萍，张小丹.反距离加权插值算法中插值参数对DEM插值误差的影响[J].测绘科学技术学报，2012，（01）.
[3]王宇航，缪亚敏，杨昕.采样点数目对反距离加权插值结果的敏感性分析[J].地理信息世界， 2012，（04）.
[4]魏义坤，杨威，刘静.关于径向基函数插值方法及其应用[J]. 沈阳大学学报， 2008，（01）.
[5]汪俊，高金耀，吴招才，张涛.局部多项式插值方法在多源海底沉积厚度数据融合中的应用[J]. 海洋科学，2009，（04）.
作者简介：王超，三峡大学，在读硕士研究生，研究方向：3S技术。
关键词：土壤重金属；空间分布预计；空间插值；交叉验证；插值精度
中图分类号： X5 文献标识码： A DOI编号： 10.14025/j.cnki.jlny.2015.21.030
土壤重金属污染问题日趋严重，已威胁到农产品质量安全。对污染的农田进行分级评价和区域划分是有效治理土壤重金属污染和保障农产品质量安全的关键措施。
在目前的产地分级划分工作中，仍需要依赖区域采样，监测采样点的重金属含量数据，通过有限采样点的数据进行空间插值分析，估算重金属含量的连续分布情况，从而为后期的污染评价工作提供必要的数据支持。
常用空间插值方法分为以下两类[1]：确定性插值和地统计插值。其中确定性插值是基于样点之间的相似程度来创建一个连续函数。比较有代表性的确定性插值方法包括：反距离权重法（IDW）、样条插值法（Spline）和趋势面分析（TSA）。地统计方法是利用样点的统计特性，量化本点之间的空间自相关性，从而构建样点的空间结构模型，代表性的方法有各种克里金方法（Kriging）。
在空间插值分析过程中，插值方法的选择和参数的优选，极大地影响了插值结果的准确程度。通过交叉验证，比较平均误差和均方根误差，在反距离权重法的研究中，张锦明在研究反距离权重法参数对DEM插值误差影响时，发现权重指数在所有的参数中对插值效果影响最大，他在研究同时发现权重指数大于等于3的时候，对插值精度的影响不明显，2或者3是较好的选择[2]。王宇航在对反距离权重法的研究中也得出类似的结论[3]。
样条插值法（Spline）中，比较常用且有代表性的径向基函数（RBF）。RBF方法是综合5种函数方法，根据情况选择最优的组合产生的插值结果[4]。
与全局趋势面分析（TSA）相比，局部多项式方法对每一个局部区域（即窗口）做模拟，在每一个区域的模拟方法与全局的趋势面类似，精度上要优于全局趋势面分析（TSA）[5]。
克里金方法（Kriging）是一族方法的总称，它包括很多方法，常见的有：简单克里金、普通克里金、泛克里金、漂移克里金、回归克里金、分层克里金、协同克里金。Li综合比较多种克里金插值方法，提出综合精度和易操作性等各项条件，OK表现最好[6]。
本文选取湖南某地部分区域作为研究区域，对土壤重金属Cd数据进行了系统采样。分别采用反距离权重法（IDW）、径向基函数（RBF）、局部多项式法（LPI）、普通克里金法（OK）分别进行实验，根据交叉验证的结果比较各种方法之间的精度差异，优选适合产地分级的插值方法，同时将插值结果和该地区的实际污染情况做简单比较。
1 材料与方法
1.1 研究区域
选择某典型重金属污染区域作为研究区域，区域面积2132.8平方公里，位于中国的中部丘陵地带，气候主要土壤类型，成土母质为常年种植水稻。
1.2 样品采集和分析方法
在研究区域内，采用系统随机采样的方式采集样品942个。2014年，采用仪器取样方式采集土壤。我们选择该区域作为典型区域，采集土壤重金属Cd的含量。
所用仪器采用X射线荧光光谱技术测定土壤样品中Cd含量。
1.3 数据的处理
本文的Cd元素含量数据集的分析采用不同的软件包，数据的处理采用R语言，制图采用ArcGIS10.2完成。
1.4 空间插值方法
本研究中选用反距离权重法（IDW）、径向基函数（RBF）、局部多项式法（LPI）、普通克里金法（OK）这4种常用的空间插值方法对研究区域数据进行插值分析，对插值结果采用标准方法进行精度的比较。
空间插值是应用相同区域中已测样点的测量值，对未抽样点的属性值的一种预测。在预测点处的属性值都可以用下面的公式表示：■ 其中， Z*（X）表示预测点位置的预测属性值， Z（Xi）表示第i个采样点的权重值， 表示第i个采样点位置的属性值。下面分别就4种常用的插值方法展开介绍。
1.4.1 反距离权重法（IDW） 反距离权重法的基础是假设预测值是实测数据的线性组合，同时权指和距离的幂值成反比。公式表示如下：
其中，p是一个任意正实数，通常，p=2；是实测离散点到插值预测点的距离。距离越近，权值越大；距离越远，全值越小。P是用来刻画权值如何随着距离的变化而变变化的。
1.4.2 普通克里金法（OK） 克里金法是基于假设的插值参数可以被视为一个区域化变量。和反距离权重法一样，克里金法的预测值是实测数据的线性组合。根据随机领域的随机特性，可以选用不同类型的克里金方法。其中较常见的包括简单克里金法、普通克立金法、泛克里金法等。这些方法相同之处在于他们都属于无偏估计。
普通克里金的权指由半方差的计算间接得到，半方差是该方法中衡量各点之间空间相关程度的测度值。普通克里金的半方差估计公式，又称作变异函数，如下所示：
使用该方法的关键在于选取合适的变异函数模型。利用做出的半方差图找出与之拟合的最好的理论变异函数模型，可用于拟合的模型包括高斯模型、线性模型、球状模型、指数模型、圆形模型。
1.4.3 径向基函数（RBF） 径向基函数法是精确插值方法中的一种。从概念上讲，径向基函数法是寻求一条穿过实测数据点且总曲率最小的曲线函数。径向基函数的预测值由以下两部分，如下面的公式所示： 其中， ■表示径向基函数， dj表示样点到预测点的距离， fi（x）表示局势函数， 表示权重系数。根据已知点可以计算出ai和dj。
具体表示有以下5种情况：
规则样条函数（CRS）
薄板样条函数（TPS）
高次曲面函数（CRS）
反高次曲面函数（IMQ）
张力样条函数（ST）
其中，d表示样点到预测点的距离。c是一个平滑因子，■表示改进型Bessel函数，表示欧拉常数。
1.4.4 局部多项式 （LPI） 局部多项式插值的实质是搜寻一个多项式公式的过程，这个公式对应的图形经过所有的实测值点。公式表示如下：F（X，Y）=a+bY+cY+dXY+eX2+fY2
全局多项式采用一个多项式来模拟整个区域的表面，而局部多项式则是全局多项式方法和移动窗口方法的结合。与全局多项式方法不同的地方在于，局部多项式方法对每一个局部区域（即窗口）做模拟。确定窗口的大小原则是要保证每个窗口有足够的实测数据点。
1.5 模型评定方法与指标
交叉验证法是评价和比较不同插值方法精度的重要方法。因为样本数量是有限的，所以在本文中也同样选择了交叉验证方法来做插值方法的比较。具体做法是：剔除一个点，用剩余的测量值做插值，得到剔除点位置的预测值，两相比较。
这里介绍评价插值结果精度的两个常用指标： 平均误差（ME）和均方根误差（RMSE）。其中，平均误差反映估值的误差范围，均方根误差反映误差均值的变化范围。这两个值可作为插值方法优劣评价的标准。
其中，Ve，i 、Va，i 分别为样点的实际值和插值预测值，n表示样点的数量。
2 结果与分析
2.1 最优参数确定
插值方法中参数选择的不同会对插值结果产生影响，因此插值方法的比较应以参数的优化为前提。IDW和RBF的参数很少，易于挑选最优参数，而Kriging法的参数较为复杂，主要体现在变异函数的拟合上，通常需反复验证。在局部多项式插值方法应用过程中， 要根据实际需要选择核函数，一般在指数模型、多项式模型、高斯模型和高次模型中选择；本文IDW的权重系数根据最小误差的原则，取值为2；RBF法采用完全规则样条（completely regularized spline）使其实现平滑插值。
2.2 变异函数拟合
采用kolmogorov-Smirnov检验，对研究区的数据集分布特征是否符合正态分布展开讨论。可知两种变换对数据集的表现的差异较大。进一步定量分析可以得到结果。
其中只有Log变换后的结果中p值大于0.05，即满足正态分布。所以，这里的数据采用Log变换来做处理。
OK方法的插值精度受到变异函数影响最大，下面是我们使用研究区域的数据，分别选用4种常用插值模型做出的实验结果，4种模型下的误差，以半球模型为最小，而且4种插值的误差波动相对很接近。因此，对当前重金属Cd数据集来讲，在OK方法的变异函数中，选择半球模型可以取得最好的精度。
2.3 LPI方法参数的选择
在LPI方法的使用中，会涉及核函数的选择，而这个选择会明显影响插值结果的精度，所以，我们对几种常用的函数（包括指数函数、多项式函数、高斯函数和高次曲面函数）分别做了实验，得到结果。
指数模型的平均误差（ME）最小，而其均方根误差与其他的模型不相上下。于是，可知对当前情况而言，指数模型是最合适的选择。
除了这个核函数的选择以外，在局部多项式插值方法应用过程中，要根据实际需要调整搜索椭圆的参数（，R1和 R2 ），其中，表示搜索角度，R1和R2分别表示搜索椭圆长、短半径。
从空间自相关性的概念可知，空间上越靠近，属性就越相似。根据经验或专业背景找出这么一个阈值，作为邻近区域的半径。
2.4 插值结果的验证和分析
分别用4种空间插值方法对土壤重金属空间分布数据进行插值。
对实测数据进行分组，同比例均匀的取样点总数的10%，作为测试组数据。采用交叉验证法对4种插值方法进行精度比较，4种插值方法的MAE、RMSE结果。
从误差分析的结果来看，4种插值方法的精度存在差异。从平均误差ME来看，四种插值方法的误差大小顺序是：IDW>LPI>RBF>OK，即IDW误差最大，LPI和RBF居中，OK最小；从平均误差RMSE来看，四种插值方法的误差变化范围都在0.85到1.0之间，比较接近。因此，四种插值方法的精度优劣程度主要看平均误差，其顺序是OK>RBF>LPI>IDW，即OK最优，RBF和LPI居中，IDW最差。
就OK方法而言，土壤重金属Cd的含量表现出块金效应，表明人为活动对土壤元素的空间变异影响较大。同时就目前的采样数据了解到的信息可知，土壤重金属Cd的污染源分布东北方向较为集中、含量高，西南方向污染源较为分散，含量低，中间地带含量居中。当前的插值结果和实际检测结果反应的情况符合得较为一致。
3 结论
一是OK插值对于刻画区域土壤元素的空间分布具有一定的优势，但该方法操作比较复杂，特别是其半变异函数的模型及参数优选需要反复验证修改。相对而言，RBF的参数设置有简单、易操作的优点，同时方便确定最优参数。IDW精度不高，但算法简单。LPI进度相对较高，但操作相对复杂；二是对于湖南的该研究区，土壤重金属Cd含量的空间分布而言，OK的插值效果最优；三是基于各自的最优插值方法，湖南的该研究区土壤重金属Cd的空间分布可以较好地反映实际的空间相关性，同时插值分布图可以指示其污染分布情况。
参考文献
[1]朱求安，张万昌，余钧辉.基于GIS的空间插值方法研究[J].江西师范大学学报（自然科学版），2004，（02）.
[2]张锦明，郭丽萍，张小丹.反距离加权插值算法中插值参数对DEM插值误差的影响[J].测绘科学技术学报，2012，（01）.
[3]王宇航，缪亚敏，杨昕.采样点数目对反距离加权插值结果的敏感性分析[J].地理信息世界， 2012，（04）.
[4]魏义坤，杨威，刘静.关于径向基函数插值方法及其应用[J]. 沈阳大学学报， 2008，（01）.
[5]汪俊，高金耀，吴招才，张涛.局部多项式插值方法在多源海底沉积厚度数据融合中的应用[J]. 海洋科学，2009，（04）.
作者简介：王超，三峡大学，在读硕士研究生，研究方向：3S技术。