1、地质统计学Geostatistics第四讲 克里金插值Estimation - KrigingMapping is the estimation of the value of anattribute at an unmeasured locationThere are many methods of point estimationPolygons of influenceTriangulationInverse distance They are allLinearUnbiasedEstimation - KrigingEstimation - KrigingClosest (Algori
2、thms)Estimation - KrigingMoving Average (Algorithms)Estimation - KrigingDirectional trend (Algorithms)Estimation - Kriging克里金法是法国G.Matheron教授以南非矿山地质工程师D.G.Krige的名字命名的一种方法,也称为克里金法(Kriging)。从数学角度抽象来说,普通克里金法是一种对空间分布数据求最优、线性、无偏内插估计量(Best Linear Unbiased Estimation,简写为BLUE)的方法。若用矿业上的术语具体说,它是根据一个待估块段邻域内的若
3、干信息样品的数据,在考虑了这些样品的形状、大小及相互位置关系,以及品位的变差函数模型提供的结构信息之后,为了对该块段品位作出一个线性、无偏、最小估计方差的估计而对每个样品值分别赋予一定的权系数,最后加权平均来估计该块段品位的方法。简言之,克里金法就是一种特定的滑动加权平均法。Introduction Several types of kriging: Simple Kriging: known mean m Ordinary Kriging: unknown mean Kriging with a trend model: known mean at every locationm(u) Kr
4、iging with an external drift: trend model scaled from asecondary variable m(u)=a0+a1y(u) Factorial Kriging: RF model Z(u) is splitted intoindependent components (factors) Nonlinear Kriging: Lognornal, MG, Rank Kriging, IK, DK, Indicator Kriging: SIK, OIK, MIK, Probability Kriging: Use indicators and
5、 ranking of the dataEstimation - Kriging Inverse distance 145LinearUnbiasedWeightsKriging7 1 2 12ni11 12 2 2ini11212 1214013512?5 6LinearUnbiasedWeightsw = C-1 . DC11 21C31C41C61C71C12C22C32C42C52C62C72C13C23C33C43C53C63C73ni1C44 C45 C46ini1C14 C15 C16 C17 w1 D12C34 C35 C36 C37 w3 D3C64 C65 C66 C67
6、w6 D6C74 C75 C76 C77w7 D713012560 653 470XY D1i1 D2Z xi wi w 1Di D3 D4 D5 D6 D7CC51Z xi wi w 1C27w2 C24 C25 C26 DC47w4 D4C54 C55 C56 C57 w5 D5 757802 5 61 ?3 4 72 5 61 ?3 4 7Y YEstimation - Kriging14514013513012514514013513012560 65 70 75 80 60 65 70 75 80XInverse distance etcXKrigingGamma(h)Estimat
7、ion - Kriging Kriging involves settingup and solving a set oflinear equations foreach target node to beestimatedC0 - C27C0 - C15Variance C0 All actual distances areconverted to covarianceusing the variogrammodel as a “lookuptable”.C0 - C34d34 d15 d27 Lag (h)克里金法与传统估值法的区别多边形法主要是根据多边形块段内的一个已知资料来估算值,没有
8、考虑周围其它已知信息,可说是一孔之见三角形法所利用的每一个已知数据在估算中的贡献是相同的,即都是等权的,没有区别不同情况给以不同的权距离反比法前进了一步,考虑了周围的样品,而且也对各数据用样品到待估块段中心的距离的导数为权进行了加权平均,但没有考虑样品之间和样品与待估块段之间的空间构形几何因素的影响,同时也没有考虑到所研究变量的空间分布结构信息(即变差函数)克里金法克里金法最大限度地利用地质上提供的上述各种信息,对各样品数据赋以适当的权系数,就可给出更为切合实际的、更精确的估计。克里金法可以避免系统误差)()(* ii uZuZWeighted Linear Estimators The ba
9、sic idea is to estimate the attribute value (say, Augrade) at a location where we do not know the true valueni1where u refers to a location, Z*(u) is an estimate atlocation u, there are n data values Z(ui), i=1,.,n, andirefer to weights.)()(* ii uZuZWeighted Linear Estimators What factors could be c
10、onsideredin assigning the weights? closeness to the location being estimated redundancy between the data values anisotropic continuity (preferential direction) magnitude of continuity / variabilityni1Weighted Linear Estimators Assign all of the weight to the nearest data(polygonal-type estimate) i1i
11、 wdcWeighted Linear Estimators Assign the weights inversely proportional to the distancefrom the location being estimated (inverse distance schemes)1c d iwn 1iwhere di is the distance between data i and the locationbeing estimated, c is a small constant, and is a power(usually between 1 to 3).Weight
12、ed Linear Estimators How about using the variogram? that is krigingSome Definitions Consider the residual data values:Y(ui)= Z(ui) - m(ui), i=1,nwhere m(u) could be constant, locally varying, orconsidered constant but unknown. Variogram is defined as:2(h) = E Y(u) - Y(u + Covariance is defined as:h)
13、2C(0) = sill(h)C(h)C(h) = E Y(u) Y(u + h)Some Definitions Link between the Variogram and Covariance:2(h) = E Y2(u) + E Y2(u + h) - 2 E Y(u) Y(u + h)= VarY(u) + VarY(u + h) - 2 C(h)= 2 C(0) - C(h)So, C(h) = C(0) -(h) C(0) = sill(h)C(h)Y* (u) i Y(u i )Simple Kriging (1) Consider a linear estimator:ni1
14、 where Y(ui) are the residual data (data values minus themean) and Y*(u) is the estimate (add the mean back in) The error variance is defined asE Y * (u) Y (u) 2 EY (u ) Y (u ) 2 EY (u) Y (u ) C (u , u ) 2 C (u, u )Simple Kriging (1)The error variance is defined asE Y * (u) Y (u) 2EY * (u)2 2 EY * (
15、u) Y (u) EY (u)2 A2-2ab+b2n n ni1 j1 i1i j i j i i C (0)n n ni1 j1 i1i j i j i i C (0) 2 jC(u i , u j ) 2 C(u, u i ) , i 1,., n C(u , u ) C(u, u ) ,Simple Kriging (2) Optimal weightsi,i=1,n may be determined by takingpartial derivatives of the error variance w.r.t. the weights inj1and setting them t
16、o zeronj1j i j i i 1,., n This system of n equations with n unknown weights isthe simple kriging (SK) systemProperties of Simple Kriging Kriging variance can be calculated beforegetting the information Kriging takes into account: Geometry of volume being estimated: Distance of the information: Confi
17、guration of the data: Structural continuity of the variable beingconsidered:C(2,1) C(2,2) C(2,3)C(0,2)Simple Kriging (3) There are three equations todetermine the three weights: 1 C(1,1) 2 C(1,2) 3 C(1,3) C(0,1) 1 C(2,1) 2 C(2,2) 3 C(2,3) C(0,2) 1,2 0,1 1,3 0,2 2,3 0,3 1 C(3,1) 2 C(3,2) 3 C(3,3) C(0
18、,3) In matrix notation:(Recall that C(h) = C(0) -( h)C(1,1) C(1,2) C(1,3) 1C(0,1) 2C(3,1) C(3,2) C(3,3) 3C(0,3)Simple Kriging (4)Changing the Rangerange = 1 range = 5 range = 10DistanceSimple kriging with a zero nugget effect and an isotropicspherical variogram with three different ranges:Simple Kri
19、ging (5)Changing the Nugget Effect100%75%nugget = 25%DistanceSimple kriging with an isotropic spherical variogram with a rangeof 10 distance units and three different nugget effects:Simple Kriging (6)Changing the AnisotropySimple kriging with a spherical variogram with a nugget of 25%, aprincipal range of 10 distance units and different “minor” ranges:Sample Redundancy