SU里面有克里金插值程序包么

克里金插值法

克里金插值法又称空间局部插值法，是以变异函数理论和结构分析为基础，在有限区域内对区域化变量进行无偏最优估计的一种方法，是地统计学的主要内容之一，由南非矿产工程师D. Matheron于1951年在寻找金矿时首次提出，法国著名统计学家G. Matheron随后将该方法理论化、系统化，并命名为Kriging，即克里金插值法。

1 克里金插值法原理

克里金插值法的适用范围为区域化变量存在空间相关性，即如果变异函数和结构分析的结果表明区域化变量存在空间相关性，则可以利用克里金插值法进行内插或外推。其实质是利用区域化变量的原始数据和变异函数的结构特点，对未知样点进行线性无偏、最优估计，无偏是指偏差的数学期望为0，最优是指估计值与实际值之差的平方和最小[1]。因此，克里金插值法是根据未知样点有限领域内的若干已知样本点数据，在考虑了样本点的形状、大小和空间方位，与未知样点的相互空间关系，以及变异函数提供的结构信息之后，对未知样点进行的一种线性无偏最优估计。

假设研究区域a上研究变量Z（x），在点xiA（i=1，2，„„，n）处属性值为Z（xi），则待插点x0A处的属性值Z（x0）的克里金插值结果Z*（x0）是已知采样点属性值Z（xi）（i=1，2，„„，n）的加权和，即：

Z(x0)iZ(xi) （1）*

i1n

式中i是待定权重系数。

其中Z(xi)之间存在一定的相关关系，这种相关性除与距离有关外，还与其相对方向变化有关，克里金插值方法将研究的对象称“区域化变量”

针对克里金方法无偏、最小方差条件可得到无偏条件可得待定权系数

n)满足关系式： i (i=1，2，„„，



i1ni 1 （2）

以无偏为前提，kriging方差为最小可得到求解待定权系数i的方程组：

根据项目对数据处理的要求，采用了优化的克里金插值算法，将等值线地化数据插值转换为格网数据，以便实现地化数据的三维显示（王家华等，1999）。其主要实现过程如下：

第一步，计算半变异图，用非线性最小二乘拟合半变异函数系数；

第二步，数据点进行四叉树存储；

第三步，对每一格网点搜索邻近数据点；

第四步，由待预测网格点和邻近数据点计算克里金算法中系数矩阵，及右端常数向量；

第五步，对矩阵进行LU分解，回代求解待预测点的预测值。

克里金插值算法主要包括半变异函数和邻近点搜索的计算，实现方法如下。

（1）半变异函数计算

半变异函数是地质统计学中区域化变量理论的基础。地质统计学主要完成2方面的任务：利用半变异函数生成半变异图来量化研究对象的空间结构；通过插值方法利用半变异图中拟合模型和研究对象周围的实测值来对未知值进行预测。

半变异函数是用来描述区域化变量结构性和随机性并存这一空间特征而提出的。在满足假设的条件下，随机函数z（x）和z（x+h）为某一物理参数测定值的一一对应的2组函数，h为每对数之间的距离。半变异函数γ（h）可用下式来计算：

γ（h）＝ 1/2E｛［z（x）-z（x +h）］2｝

4种基本的半变异函数模式（除了这4种基本模式以外，还有很多模式），包括：

1）线形模式（Linear Model）

浙江省农业地质环境GIS设计与实现

2）球面模式（Spherical Model）

浙江省农业地质环境GIS设计与实现

3）指数模式（Exponential Model）

浙江省农业地质环境GIS设计与实现

4）高斯模式（Gaussian Model）

浙江省农业地质环境GIS设计与实现

半变异函数γ（h）会随距离h增大而增大，并逐渐逼近一定值（C0 +C），称为基台值（Sill）；而逼近基台值所对应的距离，称为影响范围（Range），表示空间中两位置间的距离小于影响范围时，是空间相关性的。在线性和球面模式中，影响范围等于a；在指数和高斯模式中，影响范围则分别等于3a和 。而模式于半变异函数轴的截距（C0）成为块金系数（Nugget Effect），产生的原因主要是样本测定的误差和最小采样间距内的变异。在应用上，为探讨说明空间变异在不同方向上的差距，也可利用非等向性的变异函数模式。半变异图拟合半变异函数模式的拟合方法可采用非线性最小二乘法拟合。

（2）邻近点搜索算法

由于矩阵LU分解求解方程的算法会随着矩阵维数的增加计算量增大，所以针对大量采样数据点时不能采用全部数据进行估计，必须采用插值点的临近点数据进行计算，即采用局部数据进行克里金算法进行计算。搜索邻近点可采用四叉树结构存储总数据，以提高搜索邻近点的速度。

对于选取邻近点的数目要有所限制，因该值的大小选择会影响插值的计算结果。若太大，则内插结果过于平滑；太小，则无法反映地表的变化；距离预测点较远的实测点可能与待估样点已经不存在自相关关系，也不能参与插值计算。采取以插值点为圆心，以R为半径的圆来确定取样的范围和参加计算的实测样点数目（如果存在各向异性，则可考虑划定一椭圆作为研究区域）。为了避免方向上的偏差，将圆平均地分为4个扇区，每个扇区内实测点数目在2～5之间，这样总共参与每个待估点预测的实测点数目平均达到8个。

区域内临近点的选择，存在着两种策略。

1）以邻近点的个数为基准。通常情况下，邻近点的个数以8～12个为宜，并且个数不能少于2个。此时计算出来的图像较为光滑。

2）以邻近点的半径尺度为基准。通常情况下，选择5～10 倍栅格间距的距离为宜。此时必须定义选择邻近点的最小和最大个数，当在一定半径内查找的邻近点个数小于最小个数时，应扩大搜索半径，使之达到最小查找个数；反之在一定半径内查找的邻近点个数大于最大个数时，应缩小搜索半径，使之小于最大查找个数。通常情况下最大最小个数分别可以定为20和4。

克里金算法的优点在于它基于一些可被验证的统计假设。根据这些假设，克里金算法产生的栅格节点估计量是最佳的，所有的估计量都依赖于可获得的观测值，并且平均误差最小。克里金算法提供了方差误差分析的表达式，可以表明每一个栅格节点的估计精度。

协同克里金插值设置三个数据集的步骤。

1、在arcmap上方右键点击GeostatisticalAnalyst。

2、点击地统计向导。

3、在输入数据中导入所需的源数据集。

4、选择泛克里金，后面就根据大家的需要进行设置了，如没想法一路下一步到底就可以完成了。

---