地学计算第四章

2023/2/11第四章变异函数和拟合模型资源与环境学院杨勇2023/2/12华中农业大学资源与环境学院目录半方差函数影响半方差函数的因素变异函数的理论拟合模型变异函数（半方差函数）r(h)是地统计学的基石2023/2/13华中农业大学资源与环境学院4.1半方差函数—半方差云图表示两个实测值之间差的一半，称为半方差值，一般意义上，半方差值会随着两点间距离的加大而加大，这是因为距离相近的样品点的性质较为相似。根据点对之间的空间滞后距离h将所有半方差值绘制成的散点图称为半方差云图（如下图，红色圈内为特异值）2023/2/14华中农业大学资源与环境学院4.1半方差函数—经验半方差函数在实践中，一般是将滞后距（h）划分为不同的级别，计算每个滞后级别上的平均半方差值，进而得到经验半方差值和经验半方差图来分析数据的空间关系。表示滞后级别，更为常见的表示方法为：2023/2/15华中农业大学资源与环境学院一般地，在实际计算时，假设步长为lag，当前滞后级别为n(n为正整数)，则h=n*lag，应该这样处理：（1）研究区所有点，找到点对（Pi,Pj），其符合条件：(n-1)*lag<dis(Pi,Pj)<=n*lag，它们之间的距离记为DISi（2）计算[z(pi)-z(pj)]2,记为Si.（3）设找到N(h)个这样的点对，计算平均距离（4）计算，为n滞后级别上的经验半方差值。（5）将各个级别的（havg，r*(havg)），绘制在图上，形成经验半方差图2023/2/16华中农业大学资源与环境学院经验半方差图2023/2/17华中农业大学资源与环境学院课堂练习假设给出某地土壤元素含量采样数据，其格式为（x,y,k），共N个采样点，其中，x和y分别表示采样点的横纵坐标值，k为实验室测得的元素含量值，步长为h，试给出相应的程序步骤，计算该区域内土壤元素的经验半方差图。（需计算的最大滞后据为样点间最大距离的一半）2023/2/18华中农业大学资源与环境学院2023/2/19华中农业大学资源与环境学院变异函数的一些现象1、跃迁现象当h超过某一个范围，例如变程a，变异函数r(h)不再增大，而是在一个极限值r(∞)附近摆动，这种现象称为“跃迁现象”2023/2/110华中农业大学资源与环境学院变异函数的一些现象2、块金效应当h=0时，变异函数r(h)不等于0，而是等于1个常数C0，这种现象称为“块金效应”2023/2/111华中农业大学资源与环境学院变异函数的一些现象3、各向异性2023/2/112华中农业大学资源与环境学院4.2影响半方差函数的因素一、统计方面1、统计分布：特异值的影响，极大或极小值越靠近区域中间，影响越大。采样数据的不均匀分布（h-散点图）。2、屏蔽：特异值可能是错误的值，也可能是真实的值，可区分对待。如果检查发现时错误。的值，必须去除。但如果是真实的值，可以去除，也可以保留（如污染区的污染物含量值）3、偏斜：数据的偏态分布，也会造成方差的不稳定

异常值对变异函数的影响2023/2/113华中农业大学资源与环境学院4.2影响半方差函数的因素二、样品大小和设计1、更多的样品能增强经验半方差值的稳定性，但考虑到工作量和资金，实际取样不可能无限取点，一般要求在变程a以内，各距离上的点对数据不应小于20对，有的学者认为不应小于30对。2、在采样的设计上，一般来讲，规则格网采样更有助于后续的分析。2023/2/114华中农业大学资源与环境学院4.2影响半方差函数的因素三、采样间距：应大致地了解变异的空间尺度以选择一个合适的采样密度。（1）对于完全未知的领域，不知道变异的尺度或模式，应预采样，找到合适的采样间距（2）若先前的采样间距造成经验半方差函数呈纯块金效应，没有明显的空间自相关性，则先前的采样间距过大，应缩小。（3）若依先前的采样间距得到的半方差函数有明显的结构，但研究区内样点分布不均匀，则可进一步优化其分布。对于第（1）（2）情况，可以采取套合采样解决。2023/2/115华中农业大学资源与环境学院4.3变异函数的理论拟合模型理论变异函数用来拟合一些列经验变异函数值，供后续进行插值估计时使用。选用理论变异函数模型是，要根据经验半方差图的形状来选取合适的模型2023/2/116华中农业大学资源与环境学院4.3变异函数的理论拟合模型变异函数在原点处的形状2023/2/117华中农业大学资源与环境学院4.3变异函数的理论拟合模型变异函数的理论模型：有基台值模型无基台值模型2023/2/118华中农业大学资源与环境学院有基台值模型—球状模型C0：块金常数C0+C：基台值C：拱高a：变程应用最广的模型2023/2/119华中农业大学资源与环境学院有基台值模型—指数模型C0：块金常数C0+C：基台值C：拱高3a：变程当C0=0,C=1时，称为标准指数函数模型2023/2/120华中农业大学资源与环境学院有基台值模型—高斯模型C0：块金常数C0+C：基台值C：拱高：变程当C0=0,C=1时，称为标准高斯函数模型2023/2/121华中农业大学资源与环境学院三种常用模型比较0.952023/2/122华中农业大学资源与环境学院有基台值模型—线性有基台值模型C0：块金常数C0+C：基台值C：拱高A：常数，表示直线斜率当C0=0,C=1时，称为标准指数函数模型2023/2/123华中农业大学资源与环境学院有基台值模型—纯块金效应模型2023/2/124华中农业大学资源与环境学院无基台值模型——线性无基台值模型2023/2/125华中农业大学资源与环境学院无基台值模型——幂函数值模型2023/2/126华中农业大学资源与环境学院无基台值模型——对数值模型2023/2/127华中农业大学资源与环境学院套合模型在实际中，有时区域化随机变量Z(x)的变化相当复杂，往往包含各种尺度及各种层次的变化，反映在变异函数r(h)上，就是单一的模型结构不能将其合理表达，而是多层次的结构相互叠加在一起，地统计学上称为套合。所谓套合结构，就是把分别出现在不同距离h上或不同方向上同时起作用的变异性组合起来，对全部有效的结构信息，作定量化的概括，以表示区域化变量的主要特征。2023/2/128华中农业大学资源与环境学院套合模型土壤是一个不均与、具有高度空间异质性的复合体，它与土壤母质、气候、水文、地形和生物等因素有关，分析土壤空间变异的因素，可将其变异分为系统变异（土壤形成因素相互作用造成）和随机变异（可以观测到的，但与土壤形成印务无关且不能直接分析的）两大类。如由h分开的两个点x和x+h的土壤某一性质Z(x)和Z(x+h)。当h趋近于0时，可以认为两点间的差异完全是由取样和测定误差造成，当h逐步增大，如h<1m，差异可能还要加上诸如水分等因素，当h<100m时，在新的变异要考虑地形的作用。2023/2/129华中农业大学资源与环境学院套合模型当h一定时，变异函数r(h)应包含小于h的所有影响因素，因此，绝大多数变异函数都由下面两个变异函数组成：

r(h)=r0(h)+r1(h)，即一个代表纯块金方差，一个代表空间相关的方差。一般情况下，套合模型可以用放映各种不同尺度变化的多个变异函数之和表示，即：ri(h)可以是相同的或不同的理论模型2023/2/130华中农业大学资源与环境学院套合模型如，区域化变量Z(x)的变异性由r0(h),r1(h)和r2(h)组成，其中2023/2/131华中农业大学资源与环境学院套合模型三者组成的套合模型为：2023/2/132华中农业大学资源与环境学院套合模型2023/2/133华中农业大学资源与环境学院套合模型实例某区域镍蕴藏量的变异函数计算值2023/2/134华中农业大学资源与环境学院套合模型实例镍含量的变异函数曲线图可以看出基台值在2.55左右对头两点作切线与纵坐标相交的交点就是块金值，大约在0.4左右0.4/2.55=0.156则大约16%的误差是随机现象造成的，而空间自相关现象造成的误差为84%2023/2/135华中农业大学资源与环境学院套合模型实例C0=0.4C1=1.55a1=14C2=0.6a2=50从图中可看出，理论值与实际值差异较大，尤其是在15到40m之间，因此，需进行反复修改2023/2/136华中农业大学资源与环境学院套合模型实例C0=0.4C1=1.15a1=12C2=1a2=60从图中可看出，理论值与实际值差异拟合较好2023/2/137华中农业大学资源与环境学院变异函数理论模型的最优拟合要使变异函数的理论模型真实地描述变量的变化规律，在建立理论模型过程中，要对模型进行最优拟合，即确定模型中各个参数的值。地统计学中，拟合过程主要包括三个步骤：确定曲线类型、参数最优估计、最优曲线的确定。2023/2/138华中农业大学资源与环境学院最优拟合—确定曲线类型一般来讲，通过(h,r(h))散点图确定曲线的大致类型或初步类型，对这个大致的初步类型进行参数最优估计，确定是否为最优曲线。2023/2/139华中农业大学资源与环境学院最优拟合—参数最优估计变异函数的理论模型主要是曲线模型，将曲线模型经过适当的变换，化为线性模型，然后用最小二乘法原理求未知参数的估计。2023/2/140华中农业大学资源与环境学院最优拟合—参数最优估计上表式地统计学中常用的变异函数理论模型经过适当变换后，化为的线性模型，对于球状模型、指数模型和高斯模型，只讨论0<h<=a时的拟合问题，对变换后的变异函数理论模型，除球状模型为二元线性回归模型外，其余均为一元线性回归模型，根据最小二乘法原理对这两类线性回归模型进行参数估计。2023/2/141华中农业大学资源与环境学院一元线性回归模型的参数估计指数模型，高斯模型，幂函数模型和对数模型，均可以变换成一元线性回归模型根据最小二乘法：通过得到的b0和b1，可以反算理论模型中的C0，C和a等参数值2023/2/142华中农业大学资源与环境学院一元线性回归模型的参数估计但上述一般方法没有注意到变异函数曲线上的头几个点的重要性远大于其他点的重要性（因为后续预测算法中，与预测点距离较近的已知点对预测的结果影响较大），因此，不应把它们与其他实际变异函数曲线上的点平均对待，为了克服这个问题，采用加权回归的方法比较合适，拟合度较高，权重系数主要是采用每一距离上的样本对数N(hi)2023/2/143华中农业大学资源与环境学院加权回归法其中：2023/2/144华中农业大学资源与环境学院二元线性回归模型的参数估计在地统计学的理论模型中，只有球状模型线性化后成为二元线性回归模型，共有三个参数：，三个参数加权最小二乘法的参数最优估计为：2023/2/145华中农业大学资源与环境学院二元线性回归模型的参数估计2023/2/146华中农业大学资源与环境学院计算出的b0,b1,b2后，还要分三种情况讨论（1）b0>0,b1>0,b2<0,此时球状模型中的三个参数C0,C,a为2023/2/147华中农业大学资源与环境学院（2）如b0<0,b1>0,b2<0，此时b0<0，即C0<0，显然不符合球状模型的要求，可设b0=0，这是线性公式为y=b1x1+b2x2，重新根据最小二乘法求出参数b1和b2，在b0=0的条件下，仍可求出C0,C和a2023/2/148华中农业大学资源与环境学院（3）如b0>0,b1>0,b2>=0,此时应分两种情况，一是b2=0，则线性模型变为y=b0+b1x1，为一元线性模型，不是球状模型，可按一元线性回归模型参数估计方法求解参数，另一是b2>0，这时对原始数据进行调整，增加或删除一些不重要的实际变异函数点的数据，反复多次调整，直到b2<0为止，然后代入下述公式中进行计算。2023/2/149华中农业大学资源与环境学院回归模型的检验请回顾第二章的相关内容2023/2/150华中农业大学资源与环境学院实例100个样方土壤表层有机质的含量

2023/2/151华中农业大学资源与环境学院四种变异函数理论模型拟合参数线性有基台值模型的残差最小，决定系数最大，其次是球状模型，基台值基本相同，但线性有基台值模型的变程最小，其次是高斯模型和球状模型。块金值是指数模型的最小，其次是球状模型和线性有基台值模型。对这5个参数，显然最重要的是考虑决定系数R2的大小，其次是考虑残差RSS的大小，然后再考虑变程和块金值的大小，根据这个原则，选择球状模型作为本实例的变异函数理论模型是比较合适的，这个理论模型除了具有较高的拟合精度外，对变程内的模拟可以得到满意的结果。2023/2/152华中农业大学资源与环境学院2023/2/153华中农业大学资源与环境学院基于优化搜索算法的参数拟合对于结构复杂的变异函数理论模型，特别是套合结构模型，参数复杂，难以用一般的通用方法求解出模型中的参数。但一些智能优化算法，如遗传算法、模拟退火算法、蚁群算法能够使用统一的流程求解出接近最优的参数。我们以遗传算法为例，介绍该算法在求解套合结构模型时的流程。2023/2/154华中农业大学资源与环境学院基于遗传算法的变异函数理论模型参数估计1、多尺度套合模型的规范表达2023/2/155华中农业大学资源与环境学院基于遗传算法的变异函数理论模型参数估计从上式可以看出，需求解的参数为2n+1个（因为第一个模型总是纯块金模型）。而在实际计算时，可以令，这样方便从经验半方差图中识别ci取值区间。并有以下约束：2023/2/156华中农业大学资源与环境学院基于遗传算法的变异函数理论模型参数估计编码策略及初始群体产生假设需要顾及m(m<=2n+1)个参数，每个参数的取值范围和估值精度分别是Umin,Umax和Qi，则将m个参数分别以L1,L2,……,Lm为长度进行二进制编码，其中则每条染色体长度为，染色体中每个参数编码对应的解码公式为：

以这种编码方式随机产生T组染色体2023/2/157华中农业大学资源与环境学院基于遗传算法的变异函数理论模型参数估计确定个体适应度评价函数2023/2/158华中农业大学资源与环境学院基于遗传算法的变异函数理论模型参数估计遗传操作遗传算法主要包括3个基本算子，即选择、交叉和变异，为此，需确定交叉概率Pc和变异概率Pm，3个过程执行以后，将产生新一代种群，并记录适应度最高的染色体选择—轮盘赌选择法如下表表示11个个体的适应度函数值、选择概率和累计概率。每一轮产生1个[0,1]随机数，将该随机数作为选择指针来确定被选个体。如随机数为0.81，则第6个个体被选中。2023/2/159华中农业大学资源与环境学院交叉算法单点交叉：考虑如下两个11位变量的父个体父个体1：01110011010父个体2：10101100101交叉点位置为5，则交叉后两个子个体：子个体1：01110100101子个体2：101010110102023/2/160华中农业大学资源与环境学院多点交叉M个交叉位置，可无重复随机地选择，交叉点之间的变量间续地