D0I:10.13374/j.issn1001-053x.1991.02.001 北京科技大学学报 第13卷第2期 Vol.13 No,2 1991年3月 Journal of University of Science and Technology Beijing March 1991 协同区域化与协同克立格法 侯景儒”张树泉·潘汉军·黄竞先· 摘要:描述了协同区域化及协同区城化变量的概念,较详细地讨论了交叉炬及其性质 (例如互协方差及互变异函数)。重点讨论了协同克立格法(包括协同克立格方程组及协 同克立格方差等),给出了根据估计邻域内金和银的数据,并应用协同克立格法实例估计 一个待估块段中金的平均含量。 关键词:克立格法,协同克立格法,协同区城化,互协方差,互变异函数 Coregionalization and Cokriging Hou Jingru Zhang Shuguan'Pan Hanjun'Huang Jingxian' ABSTRACT:The concept of 'coregionalization and coregionalized variable are described.The cross-moments and its properties (such as cross-covariance and cross-variogram)are discussed in more detail.In the view of above analysised, the cokriging (include cokriging systems,cokriging variance and so on)is studied mainly,and the mean grade of gold in a block is estimated.This method makes use of estimated level of Au and Ag in neighborhood. KEY WORDS:kriging,cokriging,coregionalization,cross-covariance,cross- variogram 1 变量相关与协同区域化 在地质、采矿及其他自然现象的研究中,某种区域化现象可以用若干相关的变量进行研 1990一10一09收稿 ,地质系(Department of Geology) ··中国有色金属总公司(Chinese Nonferrous Meta】Company) 95
第 1 , 卷第 么期 1 9 , 1年 8 月 北 京 科 技 大 学 学 报 J o u r n a l o f U n i v e r s i t y o f S e i e n e e a n d T e e h n o l o g y B e i j i n g V o l 。 1 3 N o . 2 M a r e h 1 9 9 1 协同区域化与协 同克立格法 侯景儒 . 张树泉 . 潘汉军 . 黄竞先 .t 摘 要 : 描述 了 协同区城化 及 协同 区城化变是的概念 , 较详细地讨论 了交叉拒及其性质 (例 如互协 方差及互 变异函数 ) 。 重点讨论了协同克立格 法 ( 包括 协同克立格 方程组 及 协 同克立格方差等) , 给 出了根据估计邻域内金和银 的数据 , 并 应用 协同 克立格法实 例估计 一 个持估块 段中金的平均含量 。 关镇饲 : 克立格法 , 协同克立 格法 , 协 同区 城化 , 互协方差 , 互变异 函 数 C o r e g i o n a l i z a t i o n a n d C o k r i g i n g H o “ J i ” g r u . Z h a n 夕 S h u g u a n . P a ” H 。 ” j u n . H “ a 月夕 Ji n 夕x f a 称 二 A B s T RA c T : T h e e o n e e p t o f ` ’ e o r 夕g i o n a l i z a t i o n a n d e o r e g i o n a li z e d v a r i a b l e a r e d e s e r i b e d · T 五e e r o s s 一 m o m e n t s a n d i t s p r o p e r t i e s ( s u e h a s e r o s s 一 e o v a r i a n e e a n d e r o s s 一 v a r i o g r a 爪 ) a r e d i s e u s : e d i n m o r e d e t a i l 。 I n t h e v i e w o f a b o v e a n d a n a l y s i s e d , S `且. :二 t五e e o k r i g i n g s t u d i e d m a i n l y , ( i n e l u d e e o k r i g i n g s y s t e m s , e o k r i g i n g a n d t h e m e a n g r a d e o f g o ld i n a V a r l a n C e 5 0 0 吧) m e t h o d nt a k e s K E Y W O R D S : u s e o f e nt a t e d l e v e l o f A u a n d A g b l o e k 1 5 e s t i m a t e d 。 T h i n n e i g h b o r h o o d 。 k r i g i n g , e 。众 v a r i o g r a m 、 r l g l n g , e o r e g i o n a li乞 。 t i 。 n , e r o s s 一 e o v a r i a n e e 一 e r o s s - 1 变量相关与协同区域化 在地质 、 采矿及其他 自然现象的研究中 , 某种区 域化现象可以用若 干相关 的变量进行研 i 马9 0一 1 0一 0 9收稿 一 地质 系 ( D e P a r t 也 e n t o f G e o l o g y ) 。 . 中国有色金属总公司 ( C h i n e , e N o n f e r r o u , M , t a l C o m P a n y ) 9 5 DOI: 10. 13374 /j . issn1001 -053x. 1991. 02. 001
究。例如,在一个区域的地球化学观测数据中,金、银、砷的含量呈正相关,即在金含量高 的样品中银及砷的含量也程度不同的增高。当把金、银、砷同时进行研究时,一旦一些样品 中由于某种原因而缺少金的含量,就可以用银、砷元素所提供的信息进行金的空间变异性及 统计特征的分析。 在上述的区域地球化学研究中,金、银、砷在样品中的含量是在同一空间域中定义的区 域化变量,它们之间既有空间相关性又具有统计相关性,称金、银、砷是同时区域化的。所 谓协同区域化,是指那些在统计意义上及空间位置上均具有某种程度相关性,并且定义于同 一空间域中的区域化变量,即他们是同时区域化的。 协同区域化变量可用一组K个相关的区城化变量{Z1(x),Z2(x),…Z:(x)}来表 示,或者说可以用一组K个相关的区城化变量组成的向量{Z:(x),Z2(x),·,Z(x)} 来表示。观测前,协同区域化变量是K维区城化变量的向量,即它是一个随机场;观测 后,协同区域化变量就是一个空间点函数,即可以把{Z1(x),Z2(x),…Z:(x)}看成是 上述K维向量的-个现实1?。 2关于协同区域化变量的若干假设 当协同区域化变量满足以下两个条件时称该协同区城化变量满足二阶平稳假设: (1)每一个区城化变量Z:(x)(k=1,2,,K,为变量数)的期望存在且平稳: E{Z:(x)}=m:(常数)青:1,2,…,K,Tx (1) (2)对于每一个区域化变量Z:(x)和Z(x)(k,'=1,2,…,K),其互协方差函 数存在且平稳: EZ(x).Z(x+h)}-m:m=C:(h) (k,=1,2,…,K)Vx (2) 即只依赖于两个样品点x及x+之间的距离h而与样品点的具体位置x无关。式(2)中的 m,m:'分别为Z:(x)及Z:(x)的平均值。显然,当(2)式中的h=0时,该式即为变量 Z(x)与Z:(x)的直接协方差: E{Z:(x)(Zg(x)}-mm:1=C±(0) (3) (k,k'=1,2,…,K) 必须指出的是当卡0时,式(2)中的k,k'的顺序不能随意颠倒。 当协同区域化变量满足以下两个条件时称该协同区域化变量满足内蕴假设: (1)每一个区域化变量Z:(x)(k=1,2,,K)的增量[Z:(x)-Z:(×+h)]的期 望为0: E{Z:(x)-Zk(x+h)}=0 (4) k=1,2,,…,KVx (2)对于每一个区域化变量Z:(x)和Z:(x)(k,k'=1,2,,K),其互变异函数 96
究 。 例如 , 在一个区域的地球化学观侧数据中 , 金 、 银 、 砷的含量呈正 相关 , 即在金含量高 的样品 中银及砷的 含量也程度不同的增高 。 当把金 、 银 、 砷同时进行研究时 , 一旦一些样品 中由于某种原因而缺少金的含 量 , 就可以用银 、 砷元素所提供的信息进行金的空间变异性及 统 计特征的 分析 。 在上述的区域地球化学研究中 , 金 、 银 、 砷在样品 中的 含量是在 同一空间域 中定义 的区 域化变量 , 它 们之 间既有空间相关性又具有统 计相关性 , 称金 、 银 、 砷是 同时区域化的 。 所 谓协同 区域化 , 是指那些在统计意义上及空间位置上均 具有某种栓度相关性 , 并且定义于 同 一空 间域中的 区域化变量 , 即他们是同时区域化的 。 协同 区域化变量 可用一组 K 个相关的区域化变量 丈Z , (x) , Z : ( x) , ” 一 Z * x( )} 来表 示 , 或者说可以 用 一组 K 个相关的区域化变量组成的 向量 毛Z ; ( 二 ) , Z : ( x) , · 一 , Z ; (二 )} 来表示 。 观测前 , 协 同区域化变量是 K 维区域化变量处户量 , 目咋是 一 个萝 机 场 ; 观 侧 后 , 协 同区域化变量就是一个空间点函数 , 即可以把 { Z , (幼 , Z : (幻 , · · … Z * x( )} 看成是 上 述 K 维向量的 一个现实 匕 ’ 〕 。 2 关于 协同区域化变最的 若千假设 当协同 区域化变量满足以下两个条件时称该协同 区域化变量满足二阶平稳假设 : (1 ) 每一个区域化变量Z * ( x) ( 庵= 1 , 2 , “ 一 , K , 为变量数) 的期望存在且平稳 : E { Z * ( 二 ) 卜 m 。 ( 常数 ) 几 二 i , 2 , … , 尤 , v 二 ( i ) (2 ) 对于每一个区域化变量 Z * (川和 Z r (幻 k( , k, = 1 , 2 , … , K ) , 其互协方差函 数存在且平稳 : E { Z 。 ( x ) . 2 * r ( 二 + h ) } 一 m * m : r = C : ` r ( h ) ( 庵 , 吞 , = 1 , 2 , … , K ) V 二 ( 2 ) 即只依赖于两个样 品点 x 及 x + h之间的距离 h 而 与样品点的具体位置 二 无关 。 式 ( 2) 中的 m * , m * , 分 别 为 Z * ( 二 ) 及2 . , ( 二 ) 的平均值 。 显然 , 当 ( 2) 式中的 h = O 时 , 该式即为变 量 Z * x( ) 与 Z * , ( 二 ) 的直接协方差 : E { Z * ( x ) . ( Z * , ( 二 ) } 一 m * 。 : , 二 C 。 : r ( 0 ) ( s ) (吞 , 吞 I = 1 , 2 , … , K ) 必须指 出的是当 h今 。 时 , 式 ( 2) 中的 k , 剐 的顺序不 能随意颠倒 。 当协同 区域化变量满 足以下两个 条件时称 该协同区域化变量满足 内蕴假设 : ( i ) 每一个区域化变量 Z * ( 二 ) ( k = i , 2 , … , K ) 的 增量 〔Z , ( 二 ) 一 Z , ( 二 + h ) 〕 的期 望为 0 : E { Z * ( 二 ) 一 Z * ( 二 + h ) } = o ( 4 ) ` k = 1 , 2 , , 一 , K V 二 (2 ) 对于每一个 区域化变 量 Z 、 ( x) 和 Z * , ( x) (k , ’k = 1 , 2 , … , K ) , 其互变异函数 9 6
y..(h)=0.5E{[Z(x+h)-Z.(x)][Z.(x+h)-Z.(x)]} k,克=1,2,,KVx (5) 存在而且平稳,即它们只依赖于两点x,x+h之间的距离h而与具体位置×无关。 3互协方差与互变异函数的性质 互协方差与互变异函数是研究协同区域化变量的空间变异性及统计特征的最重要的参 数,如同普通克立格法一样,互协方差及互变异函数组成了目的在于对协同区域化变量进行 最优无偏线性估计的协同克立格方程组〔2)。因此,了解它们的性质是十分重要的。互协方 差与互变异函数的主要性质如下: (1)在普通克立格法中的变异函数Y()总是大于等于0,即Y,(h)≥0,但互变异函 数Y'(h)可以有负值。当y'()为负时表示:变量Z:'(x)的增加对应于另一变量Z,(x) 减小,或者Z(x)减小对应于Z()增加,即Z,(x)与Z.(x)在空间分布上呈负相 关。 (2)可以证明:互变异函数对于k'和k对称,即:y'.()=Ysa'(h) (6) 对于h,(-h)也对称,即:Y'(-h)=p'a(h) (7) (3)可以证明:对于互协方差函数而言:C:':(-h)=C:k'(h) (8) 但是,C',(-h)≠C4'()即互协方差函数对于h和(-)无对称性。 (4)在普通克立格法中,当区域化变量服从二阶平稳假设时,其变异函数y()与协方 差函数C(h)及样本方差C(0)之间的关系是: p(h)=C(0)-C(h) (9) 但当协同区域化变量服从二阶平稳假设时,其互协方差函数Cs'()与互变异函数Y:'(h) 之间的关系则如下式所示: (=C:(0)-分(C.(h)+C(h) (k,k'=1,2,…,K)Vx (10) (10)式证明如下: 2y4'.(h)=ECZ(x+h)-Z.'(x)们〔Z.(x+h)-Z.(x)门 =E{〔Zs(x+h)-m门-〔Za(x)-m4门} {CZ.(x+h)-m)-〔Z.(x)-m.]} =ECZ.(x+h)-m门〔Z.(x+h)-m,] -E〔Z.'(x+h)-ma门〔Za(x)-m〕 E[Z(x)-mZ (x+h)-m:3 +ECZ:(x)-m(x)-m 97
, 、 , 、 ( h ) = 0 . S E {〔 Z , r ( 二 + h ) 一 Z , , ( 二 ) ] 一 〔2 . ( 二 + h ) 一 Z , (二 ) 〕} 克 , 秃 , = 1 , 2 , … , K V 二 存在而且平稳 , 即 它们 只依赖于两 点 二 , 二 + h 之间的距离 h 而与具体位置 二 ( 5 ) 无关 。 3 互协方差 与互变异函数的性质 互协方差与互变异函 数是研究协同 区域化变量的空 间变异性及统计特征的最 重 要 的 参 数 , 如 同普通克立格法 一样 , 互协方差及互变异 函数组 成了 目的在于 对协同 区域化变量进行 最优无偏线性估计的协 同克立格方程组 ` 2 ’ 。 因此 , 了解它们 的性质是十分重要 的 。 互 协 方 差与互变异函数的主要性质如下 : ( 1) 在普通克立格法 中的变异 函数 下, ( h) 总是大于等于 。 , 即 , 。 (h) ) 0 , 但互变 异 函 数 ? “ 。 h( ) 可 以有负值 。 当 ? 。 ` , ( h) 为负时表示 : 变量 z “ (幼 的增加 对应于 另 一变量 z 。 (幼 减小 , 或者 Z , , (幻减小对 应于 Z 。 ( x) 增加 , 即 Z 。 , ( 二 ) 与 Z 。 ( x) 在空间 分 布 上 呈 负 相 关 。 ( 2 ) 可以证 明 : 互变异函数对于 壳` 和 k 对称 , 即 : , ` , 。 (h ) 二 , , 、 ` ( h ) ( 6 ) 对于 h , ( 一 几) 也对称 , 即 : ? “ 。 ( 一 h ) = , 。 , , ( h ) ( 7 ) (3 ) 可以证 明 : 对于互协 方差函数而言 : C * , * ( 一 h) = C * * , (h) ( 8 ) 但是 , c , , 。 ( 一 h) 笋 C , , ` h( ) 即互协 方差函数对于 h 和 ( 一 h) 无 对称性 。 (4 ) 在普通克立格法 中 , 当 区域化变量服从二 阶平稳假设时 , 其变异函数 , ( h ) 与协 方 差函数 C (h) 及样本方差 C ( 0) 之间的关 系是 : ? ( h ) = C ( o ) 一 C ( h ) ( , ) 但 当协 同 区域化变量服从二 阶平稳假设时 , 其互协方差函数 c , 。 , h( ) 与互变异 函数 , , , , h( ) 之 间的关系则如下式 所示 : : 。 / 。 ( h ) = e 一 。 ( 。) 一 冬〔e 。 , 。 (、 ) + e * * , ( 。) 〕 乙 ( 寿 , 寿 , = 1 , 2 , … , 尤 ) V 二 ( 1 0 ) 式证明 如下 : 2 ? * ` * ( h ) 二 E 〔Z 。 , (二 + h ) 一 Z 。 z ( x ) 〕〔 Z ` ( x + h ) 一 z * ( x ) 〕 = E { CZ ` , ( x + h ) 一 m ` , 〕 一 〔Z 、 r ( 二 ) 一 。 * ,〕卜 {〔Z 。 ( x + h ) 一 m ` 〕 一 〔Z ` ( 二 ) 一 m ; 〕} = E 〔Z 。 , ( 劣 + h ) 一 杭 , z〕〔Z ` ( x + h ) 一 m , 〕 一 E 〔Z 、 , (盆 + h ) 一 沉 。 , 〕〔Z ` (二 ) 一 m 、 〕 一 E 〔Z , , (劣 ) 一 m 、 , 〕〔 Z ` ( x + h ) 一 m , 〕 + E 〔Z , , ( 二 ) 一 川 * , 〕〔Z 、 (劣 ) 一 拼 、 〕 ( 1 0 )
=ECZ:(x+h)Z:(x+h))-mxim:} -{E〔Z:(x+h)Z:(x)〕-mgm:} -ECZ:(+h)Z:(x)]-m:m +{E〔Z.(x)Z,(x)门-mm,} =C's(0)-C'(h)-C'.(h)+C'(0) =2C'k(0)-[Ck'k(h)+C(h) ,(=C:(0)-2C.h+C(h] 证毕。 (5)式(2)及式(5)中的k'=时,互协方差函数C'()变为协方差函数: C(h)=EZ(x+h).Z(x)}-m=C(h)Vx (11) 互变异函数yk'k(h)变为变异函数: y44(h)=0.5ECZ,(x+h)-Z.(x)〕2=y(h)V¥ (12) 4协同克立格法 设协同区城化变量由K个空间互相关的随机函数{Z.(x),k=1,2,…K}的集合来 表征,它们服从二阶平稳假设,即其数学期望E{Z.(x)}及互协方差C.'()和互变异函数 Y4(h)存在且平稳r8?。设是诸区域化变量(k=1,2,…,K)中某一特定的需要研究 的主变量。 应用估计邻域内定义于支撑{va}上的有效数据{Za,a=1,2,“n,}来估计中心点在 x处的待估域V(x)上变量。的平均值Z~:。的估计量Zo。显然, 2n品2ds (13) 乙=2)1x (14) 而Z0。的估计量Z行,是K个协同区域化变量的全部有效数值的线性组合, (15) 1a=1 为了要使Z了为Z,的最优无偏线性估计量,必须在无偏条件下使其估计方差为最小时来求 式(15)中的诸权系数入a., 无偏性条件定义为Z。与Z5。之差的期望为0,即: E{Z7。-Z。}=0 (16) 98
= {E 〔Z * , ( x + h ) Z 。 ( x + h ) 〕 一 。 、 , m * } 一 { E 〔 Z * z ( 二 + h ) Z ; ( x ) 〕 一 。 * , m * } 一 { E 〔 Z ; ( 二 + h ) Z * , ( 二 ) 〕 一 m : fn ; r } + { E 〔Z , , ( x ) 2 . ( x ) 〕 一 。 ` , m ` } = C , , ` ( 0 ) 一 C . , ` ( h ) 一 C , , , (h ) + C . , ` ( 0 ) = ZC * ` 、 ( o ) 一 〔C 、 ` * ( h ) + C 七 * ` (h ) 〕 。 _ , 。 、 1 , 。 , 。 、 . 。 _ , 二 、 、 一 夕` , ` L“ 少 = “ “ 火U ) 一 下万一 七七 ` , 连气“ 少 十 “ ` 尹 叹“ 少 J ` 证毕 。 (5 ) 式 (2 ) 及式 ( 5) 中的k 尹 = k时 , 互协方差函数 C , 。 ` h( ) 变为协方差函 数 : C , 、 ( h ) = E { Z 。 ( 二 + h ) . 2 ` ( x ) 卜 m 艺= C 。 (盖) V 二 互 变异函数 , 、 ` k h( )变为变异函数 : ? 。 、 (h ) = o 。 S E 〔 Z , ( 二 + h ) 一 Z , ( 二 ) 〕 2 = , 。 ( h ) V 二 ( 1 1 ) ( 1 2 ) 4 协 同 克 立 格 法 设协同 区域化变量由 K 个空间互相关的随机 函数 { Z 。 (幻 , k = 1 , 2 , … K } 的集 合 来 表征 , 它们服 从二阶平稳假 设 , 即其数学期望 E { Z 。 (x )} 及互协方差 C “ , ( h) 和互变异函数 护。 , , h( ) 存在且平稳 〔 “ ’ 。 设 吞 。 是诸区域化变量 庵( 庵二 1 , 2 , … , K ) 中某一特定的需要研究 的主变量 。 应 用估计邻域内定义于支 撑 { 打 。 。 }上 的有效数据{ Z 。 , , a 。 = i , 2 , 一 , 、 } 来 估计中心点 在 “ 。处 的 待估域 v ( , 。 )上变量 k0 的平 均值 Z ; 。 。 的估计量 Z 等 , 。 。 显 然 , Z 二 , 。 = 翩 , 。 Z , 。 ( x ) d x z 。 。 二 壳 .,j 二 。 , Z 。 ( X , d X 而 z , , 。 的估计量 Z导 . 。是 K 个 协 同 区域化变量 的全部有薄数值 的线性组合 , K 月 ` Z导 , 。 二 习 习久 。 。 Z 。 , 几 . 1 仪 为一 1 ( 1 3 ) ( 1 4 ) ( 1 5 ) 为了要使z 专 。 。为 z 令 . 的最优无偏线性估计 量 , 必须 在无偏条件下使其估计方差为最小时来求 式 ( 15 ) 中的诸权系数久。 , 。 无偏性条件定义 为 z 甲 。 。 与 Z备 。 。 之差的 期望为 。 , 即 : E { Z甲 * 。 一 2 7 , 。 } = ` o ( 1 6 )
式(16)可改写为: E(Z。-io}=E{2o}-》a.E2a}-∑∑1aEZ} a4■1 =m〔i-芝…门a.m豆, (17) aho前1 a是■1 要使式(17)等于0必须有如下K个条件成立: ah。=1 (18) 乃 1as=0 atel 程质0 式(18)是无偏性的充分条件组。当”。=0时条件∑a。=1不能成立,也就是说,在 0k0●1 应用协同克立格法估计时,主变量,至少应有一个观测值Z:。 在上述条件下的估计方差是: 品.=E27.-Z} =0)-222e,) ha2p C(mn +宫良点点 (19) 式(19)中C(?。,V)为待估域V。中变量乙,(x)的直接协方差函数的平均值;C (7o,va4)为待估域7。与信息域a4之间区域化变量Z:。(x)与Za(x)的互协方差函数的 平均值;C(vaa,vs)表示信息域VaA与vgs'之间区域化变量Za4(x)与ZBs'(x)的互 协方差函数的平均值。 在式(18)表示的K个无偏条件约束下,要求式(19)所表示的估计方差0品:。达到极小时 来求式(15)中的诸权系数1·,可根据求条件极值的拉格朗日乘数法。 令F=-2(克a-1)-2,(分a) (k牛k。) 将F对诸2a,o,4(k卡k。)求偏导并令其为0, 99
式 ( 16 )可改 写为 : 二、 z二。 。 一 忘 ` 。 } = 二、 z 二 , 。 } 一 誉 ; 。 , 。 E 、z 。 * 。 一 习 习 又 a o E { Z · , } a 几。 一 1 汤专 轰。 汀 玉一 1 月 告。 · m 告。 〔 , 一 习 , a , 。 〕 a 几 。 一 1 要使式 ( 1了) 等于 。 必须有如下 尤 个条件成立 : 习 久 a ` 。 = ; 习 , a 。 = 0 n 乙 一 习 二 , 。 习 久 a , ` 专 白o a 几 一 1 目 ` o ( 1 7 ) ( 1 8 ) 式 ( 18 ) 是无 偏性的 充分条件组 。 当 。 , 。 = 。 时条件 习 而 。 。 = 1 不能成立 , 也 就是说 , 在 J 告 O 一 1 应用协同克立格法 估计时 , 主 变量 寿 。 至少应有一个观测值 Z , 。 。 在上 述条件下 的估计方差是 : 。导 。 。 二 二 { 2 7 。 。 一 z 导 。 。 } K 月 为 = c 。 。 , 。 ( v 。 。 , v * 。 ) 一 2 习 习 几a 。 c , 。 , ( v 。 。 , 。 a * ) 生 J 口 几一 1 尤 f 十 习 习 九 ` ” ` I 习 习 几a , 几, , , C 。 * , ( 。 a 。 , ” , 。 , ) a 。一 1 刀 , I 一 1 ( 1 9 ) ` 一 1 k l . z 式 ( r 。 ) 中乞 , 。 。 。 (、 。 。 , v , 。 ) 为待估域 v ; 。 中变量 Z , 。 ( x ) 的直接协方差函数的平均值; c , 。 。 ( v , 。 , z, a 。 )为待估域 V . 。 与信息域 。 。 * 之’lia 区域化变量 Z , 。 ( 二 ) 与 Z a 。 ( x ) 的 互协 方差函数的 平 均值 ; C , , , ( Z, 。 。 , z, , , , ) 表示信息域 v 。 * 与 。 , , , 之间区域化变量 Z 。 。 (二 ) 与 Z 刀 ` , (二 ) 的互 协方差函数的平均值 。 在式 ( 18 ) 表示的 K 个无偏条件约束下 , 要求式 (1 9 ) 所表示的估计方差 『 令 , 。 达到极小时 来求式 ( 15 ) 中的诸权 系数 只 。 , , 可根据求条件极值的拉格 朗 日乘数法 。 令 二 一 a 导 、 。 一 2 ; 。 。 ( 习 “ a , 。 一 i ) 一 2 ; , ( n ` 习 , a * ) ` k今 k 。 , 口 玉。一 1 将 F 对诸几a , “ . 。 , 拌 。 ( k粉k 。 ) 求偏导并令其为
=2之2i1C4(4,a)-2,-2C(W,7a=0 0入ak k/1 Va&=1,2,…,nik=1,2,…,K -2(空 a。-1)=0 眼.-2(克a)=0 Vk卡k。 经整理得到协同克立格方程组为(20)式所示: ∑Bx1C(wi,a)-4=Ca(。va) Va=1,2,…k5k=1,2,…,K ho at。=1 (20) d南0●1 a=0 (20)式中的 (空:+K)个未知数(2 n个权系数元a及K个拉格朗日乘子4:) 可以由(公”,+K)个线性方程组求解,即用(20)式求出的,才能保证Z.是 Z口。的最优无偏线性估计量。 相应的协同克立格方差为, .a C (m va (21) 在普通克立格法中,变异函数y()与协方差C()之间存在着如(9)式所示的关系,因此可 以很方便的用y(h)代替C(h)表示普通克立格方程组,但在协同克立格法中,y'(h)与 C(h)之间的关系如(10)式所示,从而不能简单地用y4'.(h)来表示协同克立格方程组 (式(20)),只有当 C(h)=C(h) (22) 时,才能使 ys'a(h)=Cs'(0)-C'(h) (23) 成立,这时,协同克立格方程组既可用互协方差函数也可用互变异函数表示,即用(-y (h))代替方程组(20)式中的(C'()导出用互变异函数表示的协同克立格方程组。 要使协同克立格方程组(或(20))具有唯一解的条件是: (1)互协方差矩阵〔Cs(vw',vα)〕严格正定,为此,必须采用正定的点协同区域 100
“ 云 走 , 一 i n ` I 习 人刀 。 , C * , , ( 。 , , , , 。 a , ) 一 2 ; * 一 Z C 。 。 , ( v 。 。 , v a , ) = 口 。 , . 1 V a , = 1 , 2 , … , n , ; 吞= 1 , 2 , … , K n 为。 一 2 ( 习 “ a , 。 一 i 口石 。 一 1 ) 一 刀 舀 一 2 ( 名 “ a 。 ) = 。 V k今掩 。 汀 一 一。 即一、解叭 经整理得到协 同克立格方程组为 (2 0) 式所示 : 兀 打 ` I 习 习 久, 。 ` c , , ` ( v , 。 ` , t, a 。 ) 一 “ * 二 C ` 。 ` ( v , 。 U a . k , 一 l 口 , , 一 i 刀& o 习 久a 肠 “ 1 口 为。 一 1 V a , = 1 , 2 , … k * ; k = 1 , 2 , … , K ( 2 0 ) 久 一 名心 少、厂| | (20 , 式 中的 l(P 一 ` ) 个 ” 数 (恩 二 个权系数 “ 。 及 ` 个拉格 朗 日乘子 动 可以 由 (习 , 。 十 K ) 个线性方程 组 求 解 , 即 用 ( 20 ) 式 求 出 的 肠 , 才 能 保 证 Z认 。 是 曲. 1 Z 甲 , 。 的最优 无偏线性估计量 。 相应的协 同克立格方差为 , a 2 7 吞。 尤 月 舀 = c , 。 * 。 ( v 。 。 , v , 。 ) + ; , 。 一 习 习 “ 。 , c , 。 , ( v ` 。 , v a ; ) 孟七 1 口盖 . 1 ( 2 1 ) 在普通克立格法 中 , 变异函数 夕 (h) 与协方差 C ( h) 之 间存在着如 ( 9) 式所示的关系 , 因 此 可 以很方便的用 夕h( ) 代替 C ( h) 表示普通克立格方程组 , 但在协同 克立格法 中 , 夕, , . (h) 与 C , , , h( ) 之间的关系如 (l 0) 式所示 , 从而 不能简单地 用 夕; , * ( h) 来表示协同 克 立 格 方 程 组 (式 ( 2 0 ) ) , 只 有当 C * , 。 (h ) = C 。 。 , ( h ) ( 2 2 ) 时 , 才能使 夕 , , 。 (h ) = C 、 , , ( o ) 一 C 。 , 。 ( h ) ( 2 3 ) 成 立 , 这时 , 协同克立格方程 组既可 用互协方差函数也可用 互变异函 数表示 , 即用 ( 一 夕。 ` . ( h) ) 代替方程组 ( 2 0) 式 中的 (c , , , ( h ) 导 出用互变异函数表示的协同克立格 方程组 。 要使 协 同克立格方程组 (或 (2 0 ) 具有 唯一解 的条件是 : (1 ) 互协方差矩 阵 〔 C , , , (巧 。 , , 。 a * ) 〕严 格正定 , 为此 , 必须采用 正定的点 协 同 区 域 1 0 0
化模型[C.'s(h)〕。 (2)没有一个数值相对于另一数值是完全多余的,即: Uah'丰VBk' Va,牛β.' (24) (3)待估的主变量Z的观测值的个数不为0,即: ns。卡0 (25) 5协同克立格法计算演示 为了计算某金矿床的金矿储量,对所采样品化验了金和银的品位值,其中A“是所要估 计的主变量,现在待估域又的估计邻域内共有6个信息样品va(a=1,2,…6),其中有3 个样品有主变量Au的品位值ZAu(vα1)(a1=1,2,3),而邻域内6个信息样品均有次要变量 Ag的品位值ZAg(Ua2)(a2=1,2,…,6),这样,待估域7的Au的平均品位值Z7n的估 计量ZA.为: 1 n2 ha1 Z (va+0 Ze(vaa) Gl1 02e1 6 =∑a1ZAu(ea)+∑oa2乙e(wa) (26) a1●1 4201 为了求出上式中的a1及®a2必须解下列协同克立格方程组: 2a,CAaa1,g+2oa,1,ag)-g1=CaAa,a) 6 a2■1 (a1=1,2,3) 6 ∑2a:C,An(eaa,vai)+∑a CAe,Ae(ea2a)-e=CAeW,a) alel a2o1 (a2=1,2,…,6)(27) ∑a1=1 Q11 0a?÷0 a3e1 (27)式中的9个权系数(3个2a1及6个oa2)及两个拉格朗日乘子(41及μ2)可由11个 线性方程组求解,而且用(27)式求出的a1(a1=1,2,3)及@a2(a2=1,2,…,6)才能保证 ZAu是Z了Au的最优无偏线性估计量。 式(27)的矩阵形式是 101
化模 型〔c 。 , , ( h ) 〕 。 (2 ) 没有一个数值相 对于另 一数值是 完全多余 的 , 即 口 a 。 , 等 Z, 刀 , , V a 。今日 。 , (3 ) 待估的主 变量 Z , 。 的观测值的个数不为 。 , 即 : n ` 。今 0 ( 2 4 ) ( 2 5 ) 5 协同克立格法计算演示 为了计算某金矿床的金矿储量 , 对所采样品化验 了金和银的 品位值 , 其 中 A 。是所要佑 计的 主变量 , 现 在待估域 V 的估计邻域 内共 有 6 个信息样品 “ a a( 二 1 , 2 , “ · 6 ) , 其中有 3 个样 品 有主变量 A u 的 品位值 Z A u ( 。 a l ) a( : = 1 , 2 , 3) , 而邻域 内 6 个信息样品 均有次要 变量 A g 的 品位值 Z A : v( 。 : ) ( 。 : 二 1 , 2 , 一 , 6 ) , 这样 , 待估域 v 的 A u 的 平均品 位 值 Z 7 A 。 的佑 计量 Z入 。 为 : 称 1 月 2 ` 甲 A u 一 ` 由 ` . a l O A u 、 . “ 1 , 丁 ` 山 叨 a Z ` A g ` U LZ Z , 仓 2 一 1 3 6 = 习 “ 。 : Z ^ 。 ( 。 a ; ) + 习 。 a , Z A : ( t, a , ) ( 2 6 ) a l 一 1 口 忍一 1 为了 求出上式 中的 久。 : 及 山 a : 必须解下 列协同克立格方程组 : 3 _ 6 习 “ a , C 、 u , A 。 ( 。 。 1 , t, a , , ) + 习 a l 一 l 口 2一 1 3 6 习 ; 。 , 吼 : , A u ( 。 。 , , v a ; ) + 习 。 a : CA u , A : ( 。 a ; , 。 a Z ) 一 拌 ; = C A u , ^ 。 ( v , t, a , ) ( a i , = 1 , 2 , 3 ) 0 。 鉴C ^ : , A : ( “ 。 : , ” a 石) 一 井2 = C ^ 。 , 、 ` (V , ” a , ) a l 一 1 口 2 , 0 1 ( a : = 1 , 2 , … , 6 ) ( 2 7 ) 3 习 几a ; = 1 口 1 一 1 6 名 。 a : = 。 口 习一 工 ( 2 7 ) 式中的 9 个权系数 ( 3 个 久。 , 及 6 个 。 。 : )及两个拉格朗 日乘子 ( 拌: 及 脚 ) 可 由 n 个 线性方程组求解 , 而且用 ( 2 7) 式求 出 的 几。 : a( , = 1 , 2 , 3) 及 。 a Z a( 2 二 1 , 2 , … , 6) 才能 保 证 Z 导 * 。 是 z 二人 。 的最优无偏线性估计量 。 式 ( 2 7 ) 的矩 阵形式是 1 0飞
C,.(w1,v)C,(u1,2)C,.(o1,C,(w1,1)Cb(1,u)C.,(1,u3) C,.(w2,1)C.,.(w2,v2)C,.(w2,s)C.,b(2,v)C,b(w2,2)C,b(2,U) C,w,p)C,(wa,)C,.(e,)C,(,u)C,b(e,2)C,(,3) C,.(w1,v1)C,(o1,u2)C,w1,g)C,(o1,1)C,bo1,v2)C,(u1,w) Cb,.(w2,v1)C,.(w2,u2)C,(e2,v3)C,(w2,v)C,b(w2,v2)C,b(U2,v) C,(wg,v1)C,.(w3,v2)C,.(wg,3)C,b(u3,v)C,b(w3,v2)C,b(w3,s) C,(u4,1)C,.(w4v2)C,.(w4,"g)C,b(w4,v1)C,b(u4,V2)C,b(m4,0g) C,(w6,1)C,.(w6,2)C,(s,U3)C,b(w5,v1)C,b(ws,v2)Cb,b(5,3) C,.(ve,v:)C,.(ve,va)Co,.(vo,v3)C,(v,v1)C,(vo,v2)C,(v0,v3) 1 1 1 0 0 0 0 0 0 1 1 1 C.,b(w,v)C,b(w1,)C.,(w1,vg)10 C,(,1) C,b(w2,w)C,b(w2,v)C.,b(w2,)10 12 C.,.(,2) C.,(v3,0)C.,(v3,v;)C.,(v3,v0)1 0 3 C,.(,s) C,b(w,04)C,b(o1,v)Cb,6(w1,ve)01 ò1 C.,b(7,01) C,b(u2,4)C,b(w2,vs)℃6,b(w2,vg)0 1 02 C.,(,2) C,(v3,v C(v3,v)(v3,v)0 1 03 ℃,b(,3) C,(v1,v)C,(v,v8)Co2b(v,vo)0 1 04 C.,b(,0) Cb,b(w6,v)C,b(u6,"5)C。,b(w6,v)0 1 ①05 ℃.,b(,s) C,(v0,v)C,(v0,v)Co(vo,v0)0 1 C,b(7,g) 0 0 0 0 0 一41 1 1 1 0 0 -42 0 C注):式中的a表示Au,b表示Ag,其他符号同前。 其协同克立格方差是: =CAa,y+-2a1C,a)-∑oaC,a) 721 (28) 6结论 在满足二阶平稳假设或内蕴假设的条件下,对于定义于同一空间域中,并且在统计意义 102
昌 , : ( t, , , U , ) C . , : ( v l , 刀 : ) 厂厄阵闪底除巴陌阮|贬 C 二 , 。 ( 口 : , v : ) C : , : ( 夕 2 , z, : ) C 矶一aC 。 , : ( t, 3 , t, 1 ) C b , . ( t, 1 , t, i ) C b , . ( t, : , 秒 , ) ( v 3 , t, : ) ( v ; , t, : ) 、 , . ( t, 3 , ” ; ) C b , : ( t, 4 , t, 1 ) C 、 , . ( t, 。 , 补 1 ) C 。 , . ( 口。 , t, : ) 1 0 C b , : ( 口 2 , 口 : ) C b , 。 ( 口 3 , 口 : ) C 、 , : ( v 4 , 口: ) C b , 。 ( u s , 口 2 ) C 。 , 。 ( t, 。 , v Z ) 1 0 C 二 , 。 ( t, i , 刀 3 ) C : , : ( 沙 2 , 岁 3 ) C : , : ( 口 : , t, 3 ) C 、 , 。 ( u , , t, 3 ) C 、 , ( t, : , t, 3 ) C b , : ( 口 3 , 口 3 ) C b , 。 ( ` 4 , 口 s ) C b , : (刀 。 , t, s ) C 、 , 。 (刀。 , t, 3 ) 1 0 C : , 、 ( t, 1 , 。 i ) C 。 , b ( p : , p : ) C 。 , b ( 口 3 , 口 i ) C b , 、 ( t, ; , 口 i ) C b , b ( 。 : , v ; ) C 、 , b ( 刀 3 , 口 1 ) C b , b ( 盯 ` , u i ) C 、 , b ( 朴 5 , t, i ) C 、 , 、 ( 口 。 , v i ) 0 1 C 。 , b ( 刀 ; , t, 2 ) C 。 , 。 ( u : , t, : ) C . , 、 ( 口 3 , t, : ) C 、 , b ( 刀 ; , t, : ) C 、 , 、 ( 口 : , 口: ) C 、 , b ( t, 3 , t, : ) C 、 , b ( v 4 , p Z ) C 、 , b ( t, 。 , t, 2 ) C 、 , b ( t, 。 , t, 2 ) 0 1 C . , 、 ( 刀 ; , 口 3 ) C . , 、 ( 岁 : , 2, 3 ) C : , b ( 口 : , t, 3 ) C b , 、 ( 口 : , 朴 : ) C b , 、 ( 刀 : , u s ) C 、 , 、 ( . 3 , 刀 s ) C 。 , 、 ( 刀 ` , p : ) C 。 , b ( p 。 , 口 3 ) C 。 , 、 ( 口 。 , u : ) 0 1 \ elf|esl1| | … | | | es … 1 1 尸 丫 C 01101 。 , 、 ( 口 ; , 口 ; ) C . , b ( t, : , t, ` ) C 。 , 、 ( 刀 3 , 口 ; ) C 、 , 、 ( 口 i , z, ` ) C 、 , b ( u : , 。 ` ) C 、 , b ( u : , u ` ) C b , b ( v ` , u ` ) C 、 , b ( 刀 。 , t, ` ) C b , b ( 刀 。 , 口 ` ) 0 1 。 , 、 ( t, i , 口 。 ) 。 , b ( t, : , z, 。 ) 。 , 。 ( v 3 , t, 。 ) 1 一一c C 。 , 、 ( 沙 i , t, 。 ) C 、 , b ( t, 2 , t, 。 ) C 、 , 、 ( u 3 , v s ) C 卜 , 、 ( t, ` , v 。 ) C b , 、 ( v 。 , u 。 ) C b , b ( v 。 , v s ) 0 1 C 二 , 、 ( t, i , u 6 ) C . , b 扣 : , 。 。 ) C 。 , 、 ( t, 。 , 盯 。 ) C 。 , 、 亡t, : , 。 。 ) C b , b ( t, : , t, 。 ) C b , b ( 刀 3 , u 。 ) C b , b ( v ` , u 。 ) C b , b ( t, 。 , 刀 。 ) C 、 , b ( t, 。 , u 。 ) 0 . , . ( V , t, ; ) . , . (V , 补 : ) . , . (V , t, 3 ) . , 。 (V , t, i ) . , 、 ( V , 口 : ) 。 , b ( V , v 3 ) . , 、 ( V , t, ` ) . , b ( V , v 。 ) . , 、 ( V , 拜。 ) 1 0 〔注〕 : 式中的 a 表示 A u , b 表示 A g , 共他符号同前 。 其协同克立格方差是 : 口 2 7 A u 3 = 瓦 u A 。 ( v , v ) + 。 , 一 云 6 久a : 万 A 。 * 。 ( v , v 。 1 ) 一 习 a , 一 1 a 2 . 1 。 a : C 人 。 人 g ( V , u a : ) (2 8 ) 6 结 论 在满足二阶平稳假 设或内蕴假设的条件下 , 对于定义于 同一空间域中 , 并且在统 计意义 1 0 2
及空间位置上均具一定相关性的区域化变量,应用本文论述的协同克立格法来估计待估域中 主变量的平均含量是提高估计精度的良好方法。在计算诸权系数a(“。=1,2,·,r, k=1,2,…,K)时,要使协同克立格方程组有唯一解的关键问题之一是互协方差矩阵严 格正定,为此,在计算并拟合互变异函数Y'。()时应该进行检验,其检验方法可参阅有 关文献r4,5)。 参老文献 1儒尔奈耳AG等著,侯景儒,黄竞先译。矿业地质统计学。北京:冶金工业出版 社,1982 2 Myers D E.Estimation of Linear Combinations and Cokriging.Mathe- Matical Geology,1983,15 (5):633-637 3 Myers D E.Co-kriging:New developments.NATO Advanced Study Institutes Series.Series C:Mathematical and Physical Sciences,1984, 122:295-305 4侯景儒,黄竞先。地质与勘探,1984,(10):40一43 5侯景儒,黄竞先。地质与勘探,1984,(11):36一39 首钢水厂铁矿矿岩爆破性分区研究 岩石爆破性分区是露天矿制定穿孔、爆破参数、技术经济指标的依据。在充分分析影响 矿岩爆破性诸因素的基础上,在国内首次提出了爆破性分区的概念、分区谁则和测试方法。 该研究成果综合考虑了岩体结构特征、岩石强度和岩石破碎特性做为爆破性分区的主要指 标。在国内首次将岩石动载特性和灰色系统理论运用于爆破性分区上,方法新颖,技术先 进,使用方便,具有较高的理论意义和实用价值,是对传统的岩石分级的一个突破。 该项研究达到的技术指标:爆破效率提高20%以上,根底大块率降低12.3%。在水厂铁 矿推广后年经济效益190万元。该研究是针对水厂铁矿矿岩特征进行的,但是分区原则、指 标的测定和数学模型对其他矿山也是适用的。 该项目由北京科技大学、首钢矿山公司水厂铁矿、长沙矿治研究院共同研制。该项研究 成果处于国内领先地位,其中岩石动载特性测定及应用上达到了国际水平。 103
及空间位置上均具一定相关性 的区 域化变量 , 应用本文论述的协同克立格法来估计待估域 中 主变量的平均含量是提高估计精度的 良好 方法 。 在计算诸权系 数 几。 , ( a 。 二 l , 2 , … , 呱 , 尧= i , 2 , … , K ) 时 , 要使协 同克立格 方程组有唯一解 的关键问题之一是互协方差矩 阵 严 格正 定 , 为此 , 在计算并拟 合互变异函数 ? , ` 。 ( h) 时 应该进行检验 , 其检验方法可参 阅 有 关文献 〔 ` , 。 〕 。 参 老 文 献 1 儒尔奈耳 A G 等著 , 侯景 濡 , 黄竞先译 。 矿业地质统计学 . 北京 : 冶 金工 业 出 版 社 , 1 9 8 2 2 M y e r s D E 。 E s t i m a t i o n o f L i n e a r C o m b i n a t i o n s a n d C o k r i g i n g . M a t h e - M a t i e a l G e o l o g y , 1 9 8 3 , 1 5 ( 5 ) : 6 3 3一 6 3 7 3 M rJ e r s D E 。 C o 一 k r i g i n g : N e w d e v e l o p m e n t s . N A T O A d v a n e e d S t u d y I n s t i t u t e s S e r i e s . S e r i e s C : M a t h e m a t i e a l a n d P h y s i e a l S e i e n e e s , 1 9 8 4 , 1 2 2 : 2 9 5一 3 0 5 4 侯景儒 , 黄竞先 . 地质与勘探 , 1 9 8 4 , ( 1 0 ) : 4 0一 4 3 5 侯 景濡 , 黄竞先 . 地质与勘探 , 1 9 8 4 , ( 1 1 ) : 3 6一3 9 首钢水厂铁矿矿岩爆破性分区研究 岩石爆破性分 区是露天矿制 定穿孔 、 爆破参数 、 技术经济指标的 依据 。 在充分分析影 响 矿岩爆破性诸因素的基础上 , 在 国 内首次提出了 爆破性分区 的概念 、 分区准则 和测试方法 。 该研究成 果综合考虑了岩体结构特征 、 岩石强度和岩石破碎特性做为爆破性分 区的 主 要 指 标 。 在 国内首次 将岩石 动载特性和灰色系统理论运用于爆破性分区 上 , 方 法 新 颖 , 技 术先 进 , 使 用方便 , 具有较高的理论意义和实用 价值 , 是对传统 的岩石分级的一个突破 。 该项研究达 到的技术指标 : 爆破效率提高20 % 以 上 , 根底大块率降低 12 。 3 % 。 在水厂铁 矿推广后年经济效益 1 90 万元 。 该研究是针对水厂铁矿矿岩特征进行的 , 但是分区 原则 、 指 标的测定和数学模型 对其他 矿山也是适用 的 。 该项 目由北京科技大学 、 首钢 矿山公司水厂铁矿 、 长沙 矿冶研究院共同研制 。 该项研究 成果处于国 内领先地位 , 其中岩石动载特性测定及应用 上达到了国际水平 。 1 0 3