1、第五章 克里金法提纲1 1克里金法概述克里金法概述2 2线性克里金法性克里金法1.1.简单克里金克里金2.2.普通克里金普通克里金3.3.泛克里金法泛克里金法3 3非非线性克里金法性克里金法1.1.对数正数正态克里金法克里金法2.2.指示克里金法指示克里金法3.3.析取克里金法析取克里金法4 4协同克里金法同克里金法一、克里金法概述一、克里金法概述1 1、克里金法概念及种、克里金法概念及种类概念：概念：又称又称为空空间局部估局部估计或或空空间局部插局部插值法法,克里金法是建立在克里金法是建立在变异函数理异函数理论及及结构分析基构分析基础上，在有限区域内上，在有限区域内对区域化区域化变量的取量的
2、取值进行行线性无偏最性无偏最优估估计的一种方法。的一种方法。主要主要类型：型：简单克里金法克里金法 普通克里金法普通克里金法 Ordinary Kriging 泛克里金法泛克里金法 Universal Kriging 对数正数正态克里金法克里金法 Logistic Normal Kriging 指示克里金法指示克里金法 Indicator Kriging 概率克里金概率克里金 Probability Kriging 析取克里金法析取克里金法 Disjuctive Kriging 协同克里金法同克里金法 Co-Kriging2、克里金估计量设设x为研究区域内任一点为研究区域内任一点待估点的估计待
3、估点的估计值值克里金估计量克里金估计量权重系数权重系数待估点待估点影响范围内影响范围内的的有效样本值有效样本值（1）无偏估计）无偏估计（2）最优估计）最优估计显然，估然，估计的好坏的好坏取决于取决于权重系数重系数i3、克里金法估值过程（1）数据检查（2）模型拟合（3）模型诊断（4）模型比较当区域化变量Z(x)的EZ(x)=m已知，则称为简单克里金法若Z(x)的EZ(x)未知，则称为普通克里金法二、线性克里金法1、简单克里金法设区域化变量设区域化变量Z(x)满足二阶平稳假设，其数学期望为常数满足二阶平稳假设，其数学期望为常数m，协方差函数，协方差函数C(h)和和变异函数变异函数(h h)存在且平
4、稳。存在且平稳。现现要估计中心点在要估计中心点在x0 的待估块段的待估块段V 的均的均值值ZV(x)，ZV(x)表达式表达式为为 由于由于 EZ(x)=m已知已知令令 Y(x)=Z(x)-m则 EY(x)=EZ(x)-m=EZ(x)-m=0待估待估块段新待估段新待估值1、简单克里金法设在待估块段V附近有n个样点xi(i=1,2,n),其观测值为Z(xi)(i=1,2,n)，则观测值新变量为：Y(xi)=Z(xi)-mY(V)的估计值Yv*是Y(xi)(i=1,2,n)的线性组合，则目标：找出一组权重系数目标：找出一组权重系数 ，使得，使得Yv*成为成为Y(V)的线性、无偏、最优估计量的线性、无
5、偏、最优估计量则估计Z(V)的问题转化为估计Y(V)的问题1、简单克里金法在满足以下两个条件时，Yv*是Y(V)的线性、无偏、最优估计量。（1）无偏性）无偏性 由于由于 所以所以 则则 Yv*不需要任何条件即是不需要任何条件即是Y(V)的无偏估计量。（2）最优性）最优性 在满足无偏条件下，可推导估计方差公式为：在满足无偏条件下，可推导估计方差公式为：1、简单克里金法为使估计方差最小，需对上式求为使估计方差最小，需对上式求i的偏导数并令其为0整理得整理得简单简单克里金克里金方程组：方程组：用矩阵表示为：用矩阵表示为：将将简单简单克里金克里金方程组方程组表达式带入估计方差表达式表达式带入估计方差表
6、达式得得简单简单克里金克里金估计估计方差表达式：方差表达式：1、简单克里金法从简单从简单克里金克里金方程组的方程组的n个方程中便可求得个方程中便可求得n个权重系数个权重系数i，则，则YV(x)的简单的简单克里金克里金估计量为：估计量为：简单简单克里金克里金法的估计精度在很大程度上依赖于法的估计精度在很大程度上依赖于m值的准确度，但是通常情值的准确度，但是通常情况下很难正确估计况下很难正确估计m值，从而导致简单值，从而导致简单克里金克里金估计精度降低。估计精度降低。简单克里金法克里金法计算示例：算示例：设某一区域气温数据某一区域气温数据满足二足二阶平平稳假假设，协方差函数和方差函数和变异函数存在
7、所有采异函数存在，所有采样数据的均数据的均值为16.08度，并将均度，并将均值作作为此区域化此区域化变量的数学期望量的数学期望值，将所有，将所有采采样数据剔除数学期望数据剔除数学期望值后后拟合的合的变异函数模型异函数模型为球状模型，如下所示。球状模型，如下所示。现用简单克里金方法根据五个已知点的气温数据来估算0点处的气温值1、简单克里金法2、普通克里金法设区域化区域化变量量Z(x)满足二足二阶平平稳假假设，其数学期望，其数学期望为m，为未知常数，未知常数，协方差方差函数函数C(h)和和变异函数异函数(h)存在且平存在且平稳。现要估要估计中心点在中心点在x0的待估的待估块段段V的均的均值，即，
8、即设待估待估块段段V附近有附近有n个个样点点xi(i=1,2,n),其其观测值为Z(xi)(i=1,2,n)，待，待估估块段段V的真的真值是估是估计邻域内域内n个信息个信息值的的线性性组合，即合，即现要求出要求出权重系数重系数i(i=1,2,n)，使，使Z*V(x)为ZV(x)的无偏估的无偏估计量，且估量，且估计方方差最小。差最小。2、普通克里金法（1）无偏性条件）无偏性条件 由于由于若要满足无偏性条件，需若要满足无偏性条件，需 ，则无偏性条件为：，则无偏性条件为：即在权系数之和为即在权系数之和为1的条件下估计量是无偏的。的条件下估计量是无偏的。（2）最优性条件）最优性条件 即估计方差最小条件
9、在满足无偏性条件下，有如下估计方差公式即估计方差最小条件，在满足无偏性条件下，有如下估计方差公式 要求出在满足无偏性条件要求出在满足无偏性条件 下使得估计方差最小的权系数下使得估计方差最小的权系数i(i=1,2,n)，这这是个求条件极值问题。是个求条件极值问题。2、普通克里金法根据拉格朗日乘数法原理，建立拉格朗日函数F。求出函数F对n个权系数i的偏导数，并令其为0，和无偏性条件联立建立方程组。整理得普通克里金方程组2、普通克里金法将解出的将解出的i(i=1,2,n)带入估计量带入估计量公式得到普通公式得到普通克里金克里金估计量：估计量：从普通从普通克里金克里金方程组可得：方程组可得：将此式带
10、入估计方差公式得将此式带入估计方差公式得普通克里金估计方差，记为：普通克里金方程组和普通克里金估计方差也可用变异函数(h)表示。在在Z(x)满足二阶平稳条件时，可采满足二阶平稳条件时，可采用协方差或变异函数表达的普通用协方差或变异函数表达的普通克克里金里金方程组及方程组及克里金克里金估计方差计算估计方差计算式进行求解计算；但在本证假设条式进行求解计算；但在本证假设条件下，则只可采用变异函数的表达件下，则只可采用变异函数的表达式进行求解计算。式进行求解计算。2、普通克里金法为了书写简便和便于计算，普通为了书写简便和便于计算，普通克克里金里金方程组和普通方程组和普通克里金克里金估计方差估计方差均可
11、用矩阵形式表示。均可用矩阵形式表示。协方差函数表达的普通克里金方程组展开得引入矩阵或普通克里金方程组用矩阵形式表达为：或权重系数或普通克里金估计方差用矩阵表达为：或2、普通克里金法普通克里金计算示例：普通克里金计算示例：设某一区域气温数据满足二阶平稳假设，协方差函数和变异函数存在，拟合的变异函数模型为球状模型，如下所示。数据如下，点的空间分布如图所示。现用普通克里金方法根据已知五个点的气温数据估算0点处的气温值。3、泛克里金法普通克里金法要求区域化变量Z(x)是二阶平稳或本征的，至少是准二阶平稳或准本征的。在此条件下，至少在估计邻域内有EZ(x)=m（常数）。然而实际中，许多区域化变量Z(x)
12、在估计邻域内是非平稳的，即EZ(x)=m(x)，m(x)称为漂移，这时就不能用普通克里金方法进行估计了，而是要采用泛克里金法进行估计。所谓泛克里金法，就是在漂移的形式EZ(x)=m(x)，和非平稳随机函数Z(x)的协方差函数C(h)或变异函数(h)为已知的条件下，一种考虑到有漂移的无偏线性估计量的地统计学方法，这种方法属于线性非平稳地统计学范畴。（1）漂移和涨落漂移漂移：非非平稳区域化变量平稳区域化变量Z(x)的数学的数学期望，在任一点期望，在任一点x上的漂移就是该点上的漂移就是该点上区域化变量上区域化变量Z(x)的数学期望的数学期望。漂移经常用邻域模型来研究。可表达为：在给定的以点x为中心的
13、邻域内的任一点，其漂移m(x)可用如下函数表示。式中，fl(x)为一已知函数；al为未知系数m(x)通常采用多项式形式，在二维条件下，漂移可看成坐标x,y的函数。涨落：涨落：对于有漂移的区域化变量对于有漂移的区域化变量Z(x)，假设可分解为漂移和涨落两，假设可分解为漂移和涨落两部分，部分，式中，m(x)=EZ(x)为点x处的漂移，R(x)称为涨落。（2 2）非平稳区域化变量的协方差函数和变异函数非平稳区域化变量的协方差函数和变异函数1）基本假设假设Z(x)的增量Z(x)Z(y)具有非平稳的数学期望m(x)m(y)和非平稳的方差函数，即假设下式存在：2）协方差函数和变异函数当Z(x)=m(x)+
14、R(x)时，Z(x)的协方差函数C(x,y)为：Z(x)的变异函数(x,y)为：（3）Z(x)的泛克里金法估计设Z(x)为一非平稳区域化变量，其数学期望为m(x)，协方差函数为C(x,y)且已知，则设Z(x)的漂移m(x)可表示为如下k+1个单项式fl(x)(l=0,1,2,k)的线性组合。已知n个样品点xi(i=1,2,n)，其观测值为Z(xi)(i=1,2,n)，现要用这些样品点估计邻域内任一点x的值Z(x)，Z(x)的泛克里金估计量为：为使Z*(x)为Z(x)的无偏最优估计量，需在以下两个条件下求解权重系数i(i=1,2,n)。（3）Z(x)的泛克里金法估计1）无偏性条件）无偏性条件若要
15、满足无偏性条件，需则即对任一组系数a0，a1,ak等式均成立，需成立。这k+1个子式称为无偏性条件。（3）Z(x)的泛克里金法估计2）最优性条件）最优性条件在满足无偏性条件下，用Z*(x)估计Z(x)的泛克里金估计方差为：将无偏性条件带入得要求出在满足无偏性的条件下使得估计方差最小的权系数i(i=1,2,n)，需根据拉格朗日乘数法原理，建立拉格朗日函数F。（3）Z(x)的泛克里金法估计求出函数F对n个权系数i的偏导数，并令其为0，和无偏性条件联立建立如下方程组。整理得估计Z(x)的泛克里金方程组：泛克里金方程组可用矩阵表示为：其中（3）Z(x)的泛克里金法估计从泛克里金方程组可得以下两等式：将
16、等式带入估计方差公式可得泛克里金方差，记为：用变异函数(h)表示如下：（4）泛克里金法计算示例设某一区域气温是非平稳的区域化变量，在南北方向（空间坐标的设某一区域气温是非平稳的区域化变量，在南北方向（空间坐标的y方向）上方向）上存在线性漂移，即存在线性漂移，即。若已知其涨落满足二阶平稳假设，并且拟合的协方差函数。若已知其涨落满足二阶平稳假设，并且拟合的协方差函数模型为球状模型，如下所示。模型为球状模型，如下所示。现用表现用表5 1所示数据，利用泛所示数据，利用泛克里金克里金法根据已知五个点的气温数据来估算法根据已知五个点的气温数据来估算0点处点处的气温值。的气温值。三、非线性克里金法1、对数正
17、态对数正态克里金克里金法法如果区域化变量经对数变换后是正态分布或近正态分布，如果区域化变量经对数变换后是正态分布或近正态分布，则则对区域化变量进对区域化变量进行精确估计的地统计学方法称为对数正态克立格法。行精确估计的地统计学方法称为对数正态克立格法。设区域化变量设区域化变量Z(x)服从对数正态分布，在待估点周围有服从对数正态分布，在待估点周围有n个样点个样点xi(i=1,2,n),其观测值为其观测值为Z(xi)(i=1,2,n)，区域化变量经对数变换后新变量为：，区域化变量经对数变换后新变量为：Y(x)=lnZ(x)，Y(x)为正态分布。假定为正态分布。假定Y(x)满足二阶平稳假设，数学期望为
18、满足二阶平稳假设，数学期望为m，协方，协方差函数差函数C(h)和变异函数和变异函数(h)存在且平稳。存在且平稳。基于对数变换后的采样点数据基于对数变换后的采样点数据Y(xi)(i=1,2,n)，计算实验变异函数并进行变，计算实验变异函数并进行变异函数模型的拟合和选择，然后利用简单克立格或普通克立格估计待估点异函数模型的拟合和选择，然后利用简单克立格或普通克立格估计待估点x处处的值的值Y*(x)。由于估计值由于估计值Y(x)是对数变换后的数值，因此是对数变换后的数值，因此对估计所得对估计所得Y*(x)需进行反变换。需进行反变换。2、指示克里金法实际研究中常常会需要获取研究区内研究对象大于某一给定
19、阈值的概率分布，实际研究中常常会需要获取研究区内研究对象大于某一给定阈值的概率分布，即要获知研究区内任一点即要获知研究区内任一点x处随机变量处随机变量Z(x)的概率分布。的概率分布。还会碰到采样数据中存在特异值的问题。还会碰到采样数据中存在特异值的问题。（特异值是指那些比全部数值的均值特异值是指那些比全部数值的均值或中位数高的多的数值，其既非分析误差所致，也非采样方法等人为误差引起，或中位数高的多的数值，其既非分析误差所致，也非采样方法等人为误差引起，而是实际存在于所研究的总体之中而是实际存在于所研究的总体之中）。指示克立格法就是为解决上述问题而发展起来的一种非参数地统计学方法。指示克立格法就
20、是为解决上述问题而发展起来的一种非参数地统计学方法。指示克立格法不必去掉重要而实际存在的高值数据的条件下处理各种不同现象，指示克立格法不必去掉重要而实际存在的高值数据的条件下处理各种不同现象，并能够给出某点并能够给出某点x处随机变量处随机变量Z(x)的概率分布。的概率分布。2、指示克里金法设一区域化变量Z(x)，对于任意给定的阈值z，引入指示函数I(x,z)，表达式如下：指示克立格法步骤如下：（1）确定一阈值，根据指示函数将原数据转换为0或1；（2）利用转换的数据计算指示变异函数，并进行拟合；（3）建立指示克立格方程组，计算待估点值。若把指示函数看做一普通区域化变量，也可直接由简单或普通克立格
21、方法来计算待估点的值。若选择多个阈值则需重复以上步骤。3、析取克里金法析取克立格析取克立格法法：假设假设已知任意区域化变量（已知任意区域化变量（Z,Z）及（）及（Z0,Z）二维概）二维概率分布条件下，对待估点的值或待估点值超过给定阈值的概率进行估计率分布条件下，对待估点的值或待估点值超过给定阈值的概率进行估计的一种非线性地统计方法的一种非线性地统计方法。估值估值步骤步骤：设设区域化变量区域化变量Z(x)在待估点在待估点x0周围有周围有n个样点个样点xi(i=1,2,n),其其观测观测 值值为为Z(xi)(i=1,2,n)，将原始数据转换为标准正态数据将原始数据转换为标准正态数据对每个新变量对每
22、个新变量Y(xi)(i=1,2,n)计算埃尔米特多项式的值计算埃尔米特多项式的值。计算埃尔米特多项式系数，用埃尔米特多项式来拟合正态变形函数计算埃尔米特多项式系数，用埃尔米特多项式来拟合正态变形函数。计算待估点析取克立格值计算待估点析取克立格值四、协同克里金法1、协同区域化变量理论协同区域化变量理论协同克立格法协同克立格法：是多元地统计学研是多元地统计学研究的基本方法，建立在协同区域化究的基本方法，建立在协同区域化变量理论基础之上，利用多个区域变量理论基础之上，利用多个区域化变量之间的互相关性，通过建立化变量之间的互相关性，通过建立交叉协方差函数和交叉变异函数模交叉协方差函数和交叉变异函数模型
23、型，用易于观测和控制的变量对不用易于观测和控制的变量对不易观测的变量进行局部估计。易观测的变量进行局部估计。协同区域化：协同区域化：在统计意义及空间在统计意义及空间位置上均具有某种程度相关性，并位置上均具有某种程度相关性，并且定义于同一空间域中的区域化变且定义于同一空间域中的区域化变量。量。协同区域化变量可用一组K个相关的区域化变量表示。观测前它是K维区域化变量的向量，即一个随机场，观测后，协同区域化变量是一个空间点函数，可以把看成是上述K维向量的一个实现。1、协同区域化变量理论满足二阶平稳假设的协同区域化变量应满足二阶平稳假设的协同区域化变量应满足：满足：（1）每一个协同区域化变量的数学期望
24、存在且平稳：（2）交叉协方差函数存在，且平稳：满足内蕴假设的协同区域化变量应满足内蕴假设的协同区域化变量应满足：满足：（1）每一个协同区域化变量增量的数学期望为0：（2）对于协同区域化变量，交叉变异函数存在且平稳。即2、交叉协方差函数和交叉变异函数（1）交叉协方差函数和交叉变异函数性质交叉协方差函数和交叉变异函数性质1)当当k=k 时时，交叉协方差函数转化为协方差函数，交叉变异函数转化为变异函，交叉协方差函数转化为协方差函数，交叉变异函数转化为变异函数。数。即即2)交叉交叉变异函数变异函数性质性质交叉变异函数关于交叉变异函数关于k和和k 对称，对称，即即交叉变异函数关于交叉变异函数关于h和和-
25、h对称，即对称，即在普通克立格法中变异函数总是大于等于在普通克立格法中变异函数总是大于等于0，但交叉变异函数可以有负值。，但交叉变异函数可以有负值。3)交叉交叉协方差协方差函数函数性质性质交叉协方差函数交叉协方差函数关于关于h和和-h不对称不对称，即即 ，但，但 当当h0 时时，k和和k 顺序顺序不能随意颠倒，不能随意颠倒，即即当当h=0 时时，交叉，交叉协方差转化为直接协方差协方差转化为直接协方差。2、交叉协方差函数和交叉变异函数4）交叉协方差函数和交叉变异函数具有以下关系：5）同一点两个变量点对点协同区域化变量的相关系数为：2、交叉协方差函数和交叉变异函数（2）交叉协方差函数和交叉变异函数
26、计算公式交叉协方差函数和交叉变异函数计算公式设在点x和x+h处，分别测得两个区域化变量的观测值Zk(x)、Zk(x)、Zk(x+h)、Zk(x+h)，则交叉协方差函数计算公式为:交叉变异函数计算公式为：2、交叉协方差函数和交叉变异函数（3）交叉交叉协方差函数和交叉变异函数计算协方差函数和交叉变异函数计算示例示例采用表51、图51所示的气温和海拔高度数据，以h=0,h1,4为例，交叉协方差和交叉变异计算过程如下：3、协同克立金法估值（1）协同克立协同克立金金估计量估计量（2）协同克立协同克立金法方程组金法方程组（2）协同克立协同克立金法方程组金法方程组1)无偏性条件2)最优性条件（2）协同克立协
27、同克立金法方程组金法方程组对对F求偏导数并令其为零，得协同克立格线性方程组：求偏导数并令其为零，得协同克立格线性方程组：（2）协同克立协同克立金法方程组金法方程组根据协同克立格方程组，协同克立格方差为：根据协同克立格方程组，协同克立格方差为：若有多个变量，则求解若有多个变量，则求解 的协同克立格方程组为：的协同克立格方程组为：协同克立格方差为：协同克立格方差为：（2）协同克立协同克立金法方程组金法方程组要使要使协同克立格方程同克立格方程组具有唯一解的条件是：具有唯一解的条件是：（3）协同克里金法使用条件）协同克里金法使用条件（4）协同克里金）协同克里金法计算示例法计算示例设某一区域有两个某一区
28、域有两个协同区域化同区域化变量：气温量：气温u、海拔高程、海拔高程v，均，均满足二足二阶平平稳假假设和内和内蕴假假设，其中气温是所要估，其中气温是所要估计的主的主变量。量。现在估在估计邻域内共有域内共有5个信息个信息样品，如品，如图5 1所示。假定所示。假定1、2、3号点上有气温号点上有气温值，而，而5个点上均有次要个点上均有次要变量量海拔高程海拔高程值，数据如表，数据如表5 6所示。所示。现拟利用利用协同克立格法估同克立格法估计0号点的气温号点的气温值。（4）协同克里金法计算示例）协同克里金法计算示例估值过程中，气温协方差估值过程中，气温协方差C(ui,uj)根据简单克立根据简单克立金金法计
29、算示例的球状模型法计算示例的球状模型计算，海拔高程协方差计算，海拔高程协方差C(vi,vj)根据式根据式(5.74)计算，气温海拔交叉协方差计算，气温海拔交叉协方差C(ui,vj)根据式根据式(5.75)计算，海拔气温交叉协方差计算，海拔气温交叉协方差C(vi,uj)根据式根据式(5.76)计算。计算。简单简单克立克立金金法计算示例的球状模型法计算示例的球状模型：复习思考题1.请简述克立格法的概念及其种类。请简述克立格法的概念及其种类。2.请简述克立格估值的过程。请简述克立格估值的过程。3.请比较反距离加权法和普通克立格法的异同。请比较反距离加权法和普通克立格法的异同。4.将简单克立格方程组用
30、变异函数进行表达，并写出相应的矩阵表达式，然后利用将简单克立格方程组用变异函数进行表达，并写出相应的矩阵表达式，然后利用变异函数表达的矩阵来计算简单克立格法的计算示例。变异函数表达的矩阵来计算简单克立格法的计算示例。5.用协方差函数表达的矩阵来计算普通克立格法的计算示例。用协方差函数表达的矩阵来计算普通克立格法的计算示例。6.若一区域化变量的拟合模型为纯块金效应模型，请推导说明用普通克立格法进行若一区域化变量的拟合模型为纯块金效应模型，请推导说明用普通克立格法进行估值时，参与插值的采样点的权重特点。估值时，参与插值的采样点的权重特点。7.若普通克立格法的计算示例中，若普通克立格法的计算示例中，
31、1号点的坐标为号点的坐标为(1099904.96,-337837.61)，即与，即与0号点的位置重合，请计算这时号点的位置重合，请计算这时0号点的估计值，并说明克立格插值法的特点。号点的估计值，并说明克立格插值法的特点。复习思考题8、利用普通克立格法和如下球状模型对下表利用普通克立格法和如下球状模型对下表0号点进行估计。号点进行估计。9、根据用变异函数表达的泛克立格方程组写出对应的矩阵表达式，并计算泛克根据用变异函数表达的泛克立格方程组写出对应的矩阵表达式，并计算泛克立格法示例。立格法示例。10、在在ArcGIS软件中实现协同克立格法插值，并对结果进行验证。软件中实现协同克立格法插值，并对结果进行验证。
