《地質統(tǒng)計學》讀書報告

地質統(tǒng)計學》課程讀書報告地質統(tǒng)計學讀書報告地質統(tǒng)計學包含經(jīng)典統(tǒng)計學與空間統(tǒng)計學，按其基本原理可定義為：地質統(tǒng)計學是以區(qū)域化變量理論為基礎，以變異函數(shù)為主要工具，研究那些在空間分布上既有隨機性，又有結構性的自然現(xiàn)象的科學。其為數(shù)學地質領域中一門發(fā)展迅速且有著廣泛應用前景的新興學科。國內(nèi)外的生產(chǎn)實踐表明，地質統(tǒng)計學除了在異常評價、找礦勘探、礦體圈定、儲量計算、采礦設計、礦山生產(chǎn)及地學科研等方面具有明顯的優(yōu)越性外，它在石油地質、第四紀地質、地層學、生物學、生態(tài)學、巖石學、地球化學、構造地質、地震地質、海洋地質、農(nóng)業(yè)、水文地質、工程地質、古氣候、古地理、環(huán)境、林業(yè)、醫(yī)學等許多方面都有成功應用的實例。地質統(tǒng)計學在不到50年的研究和實踐中得到了很大的發(fā)展［1］。一、理論研究及進展經(jīng)歷了數(shù)十年的發(fā)展，地質統(tǒng)計學的理論與方法研究有了很大的提高［2-3］。包括：①從初期二維平面分析到三維立體空間的靜態(tài)估計，發(fā)展到今天在時空域內(nèi)對研究對象進行四維乃至更高維空間的動態(tài)估計和模擬。Journel［4］將克立格法的估值問題，從一般矢量空間擴展到個原始數(shù)據(jù)的全部可測度函數(shù)所形成的矢量空間（希爾伯特空間）進行考察；②在單變量區(qū)域化變量理論的基礎上，提出了適合多變量的協(xié)同區(qū)域化理論⑷；③發(fā)展了許多計算變異函數(shù)（或協(xié)方差函數(shù)）的方法；④線性地質統(tǒng)計學與非線性地質統(tǒng)計學共同發(fā)展；⑤參數(shù)地質統(tǒng)計學與非參數(shù)地質統(tǒng)計學相互補充。Matheron⑸為首的參數(shù)地質統(tǒng)計學派以正態(tài)假設為前提，在協(xié)同區(qū)域化理論的基礎上，提出多元地質統(tǒng)計學的基本思想。Journel發(fā)展了無須對數(shù)據(jù)分布作任何假設的非參數(shù)地質統(tǒng)計學，提出了一些非參數(shù)地質統(tǒng)計學克立格方法；⑥由于時空多元地質統(tǒng)計學的研究得到重視，早期空間域靜態(tài)建模技術的研究逐漸過渡到研究時空域多元動態(tài)條件模擬，各種模擬方法得到了發(fā)展；⑦早期的等因子模型的因子是埃爾米特多項式，它要求原始數(shù)據(jù)服從正態(tài)分布。為了拓寬等因子模型的應用，Matheron提出了離散的等因子模型和連續(xù)的等因子模型，Rivoirard利用析取克立格技術建立了正交指標剩余模型，Lajauine和Lantuejoul等也提出了建立等因子模型的一些方法;⑧已有的地質統(tǒng)計學方法相互融合。如指示克立格法與協(xié)同克立格法相結合形成指示協(xié)同克立格法;指示克立格法與因子克立格法相結合形成主分量指示克立格法；協(xié)同克立格法與其它不同的線性地質統(tǒng)計技術相結合形成各種協(xié)同克立格技術等［6］。這里重點介紹一下多點地質統(tǒng)計學［7］。多點地質統(tǒng)計學是相對于基于變差函數(shù)的兩點地質統(tǒng)計學而言的。在兩點統(tǒng)計里，儲集層相關性通過空間兩點協(xié)方差（變差函數(shù)）進行描述。在多點統(tǒng)計里，則是利用空間多個點組合模式進行描述?？臻g多點組合樣式稱為數(shù)據(jù)樣板，如果在空間點賦予了值，則為一個特定的空間多點組合模式，稱為數(shù)據(jù)事件。在建模時，對每一個未知點，估計在其處滿足給定條件的數(shù)據(jù)事件出現(xiàn)的概率，隨后抽樣獲得未知點處值或者數(shù)據(jù)事件，即完成單次模擬。一旦所有節(jié)點得到訪問，即完成一次模擬實現(xiàn)。二、地質統(tǒng)計方法應用及進展克立格估計方法概述［8-9］：克立格法是法國Mathron以南非礦山工程師D.G.Krige的名字命名的一類方法。簡單地說是一種特殊的加權移動平均法。是一種對空間分布數(shù)據(jù)求最優(yōu)線性無偏內(nèi)插估計量的方法。根據(jù)估計式形式、區(qū)域化變量平穩(wěn)性、分布及涉及的區(qū)域化變量個數(shù)等克立格方法有如下分類[10]：線性克立格法—非線性克立格法參數(shù)克立格法—非參數(shù)克立格法(指示、概率)平穩(wěn)克立格法(普通)—非平穩(wěn)克立格法(泛克立格)單變量克立格法—多變量克立格法(多元克立格)指示克立格法1.基本思路在地質、物探化探數(shù)據(jù)處理及相關儲量計算中往往存在以下情況：特異值的出現(xiàn)。所謂特異值是指那些比全部數(shù)值的均值或中位數(shù)高得多的數(shù)值。它既非分析誤差所致；也非采樣方法等人為誤差引起，而是實際存在于所研究的母體之中，這些特異值只占全部數(shù)據(jù)的極少部分，但卻控制了總金屬資源量的很大比例，同時局限于一定的空間位置[11]。在一個研究區(qū)或一個礦床中存在幾種不同類型的礦化作用，這也影響了品位和儲量的精確估計。在礦山地質工作中經(jīng)常要知道諸如：在區(qū)A及其附近，大于邊界品位值Z的可選采礦單元占多大比例全部可選采礦單元的平均品位是多少等等。為了解決上述問題，地質學家及統(tǒng)計學家進行了許多有成就的研究，而指示克立格法是一種更好的非參數(shù)地質統(tǒng)計學方法，它可以在不必去掉重要而實際存在的高值(特異值)數(shù)據(jù)的條件下處理各種不同的現(xiàn)象，而且給出在一定風險條件下未知量Z(x)的估計量及空間分布。2.指示函數(shù)及其二階矩指示函數(shù)定義:設在礦床D上得到了某金屬的品位值,約定其邊界品位為乙則在D上的每一樣品點xeD上定義一個Z的階梯函數(shù)：一? [1,當x點上品位值Z(x)<Z5 -[0,當x點上品位值Z(x)>Z在D上的任一區(qū)域AeD內(nèi)，低于邊界品位Z的品位值Z(x)占區(qū)域A的比例表示如下：0(A;Z)=Ji(x;z)dxe[0,1]AA0(A;Z)稱為廢石函數(shù)，而高于邊界品位Z的品位值Z(x)所占區(qū)域A的比例為：屮(A;Z)=1-0(A；Z)=P{Z(x)>乙xeA}

稱為礦石回收率。指示函數(shù)的二階矩。在邊界品位Z確定的條件下，隨機函數(shù)I(x;Z)服從二項分布，其期望值是平穩(wěn)的，與x無關。E{I(x;Z)}=prob{Z(x)<Z}=F(Z)=工prob{Z(x)=Z}=工pkk二項分布的概率函數(shù)Zk<=Z Zk<=Z二項分布的概率函數(shù)prob{X二k}二Ckpk(1-p)n-kn非中心協(xié)方差是K(h；Z)二E[I(x+h;Z),I(x;Z)]二prob{Z(x+h)<Z同時Z(x)<Z}、.、.八I—r方差是Var{I(x;Z)}=C(0;Z)二F(Z)-F2(Z)I中心協(xié)方差是C(h;Z)二K(h;Z)-F2(Z)11變異函數(shù)是y(h;Z)=1E{[I(x+h;Z)-1(x;Z)]2}I2=C(0;Z)-C(h;Z)=F(Z)-K(h;Z)II(x;Z)滿足平穩(wěn)假設指示克立格方法廢石函數(shù)的估計：估計量設計：設在礦床D上有N個有效數(shù)據(jù)，在D上的一個域AwD內(nèi)有n個有效數(shù)據(jù){Z(x丿,xJA2二1,2,…,n}，給定邊界品位值為乙則指示函數(shù)空間為：{i(x;Z),a=1,2,3, n}a有申(A;Z)的估計值可以表示如下：申*(A;Z)=Y九(Z)i(x;Z)aaa=1若給一系列邊界品位值Zi時：兀1=1,2,…,L〔這時，估計量可寫為:申*(a；Z)=￡九(Z)i(x;Z)l alala=1克立格方程組建立：在無偏和估計方差極小的條件下為求Q*(A；Z)或9*(a；z)可得權系數(shù)九a(Z)(a=12…,n),滿足的方程組，即指示克立格方程組：工九(Z)丫(x,x;Z)+y=Y(x,A;Z)P iap ia<卩7 a=1,2,...,n工九(Z)=1PP=1估計方差的表示：其指示克立格方差c2=2九(Z)Y(x,A;Z)—Y(A,A;Z)+yKIaia ia=1Yi(xa，xp;Z)表示在給定的邊界品位Z條件下，矢量h的兩個端點分別在信息域"a，xP內(nèi)的所有對點的平均指示變異函數(shù)值；Yi(A,A;Z)表示在給定邊界品位Z條件下，矢量h的兩個端點在待估域A內(nèi)所有對點的平均指示便宜函數(shù)值；Yi(xa,A;Z)表示在給定邊界品位Z條件下，矢量h的一個端點在信息域xa內(nèi)，另一端點在待估域A內(nèi)所有對點的平均指示半變異函數(shù)值；卩為拉格朗日乘子。協(xié)方差函數(shù)表示的方程及方差：用協(xié)方差函數(shù)分別表示如下：廠工九(Z)C(x,x;Z)—卩=C(x,A;Z)P iaP ia<p=1 a=1,2,...,n工九(Z)=1aa=1c2=C(A,A;Z)-2九(Z)C(x,A;Z)+pKI aaa=1若有l(wèi)個邊界品位Zi(l=1,2,…,L)就應該解l個指示克立格方程組。廢石函數(shù)的意義：相應位置品位值低于Z的概率廢石函數(shù)的作用：用于概率估計及品位值估計

待估域A平均品位的指示克立格法估計，待估域A平均品位的指示克立格法是應用某種克立格法求得申(A；Z)的線性估計值申*(A；Z)，最后得到待估域A的平均品位及儲量。設待估域A約定在位置x上，分兩種情況討論如何用指示克立格法求待估域的估計值[Z(x)]*。由兩種礦化類型組成的礦床中位置x的平均品位的估計。指示函數(shù)值定義如下：[1 當x屬于1類型礦化I(x)=<0 當x屬于2類型礦化則任一位置x處的I(x)的平均值1*(x)為:I*(x)=YaI(x)aa=prob*$(x)*=11n個數(shù)據(jù)｝=prob*$xg1類型礦化，給定的n個數(shù)據(jù)｝[1-1*(x)]=prob*4g2類型礦化，給定的n個數(shù)據(jù)｝I*(x)g[0,1]；[1-1*(x)]g[0,1]9位置x上1類型礦化的品位是根據(jù)x周圍n1個1類型礦化品位值a11Z(xa)(a=1,2,...,n1)求得a11[Z(x)lxg1類型礦化=YbZ(x)aa11a=11位置x上2類型礦化的品位是根據(jù)x周圍n2個2類型礦化品位值位置Z(x)(aaZ(x)(aa22=1，2,…,n2)求得:[Z(x)Ixg2類型礦化=YbZ(x)a2 a2a2=1最后，位置x的估計品位[Z(x)]*可上述計算結果綜合得到：[Z(x)]*=I*(x)[Z(x)Ixg1型礦化]*+[1-1*(x)][Z(x)Ixg2型礦化]*由1個礦化類型(l>2)組成的礦床中位置x的平均品位的估計。在實際工作中，任一位置x上只能有一種礦化類型(如l型礦化)占優(yōu)勢，這時：

…、[1 當x屬于l類型礦化I(x)=< ,1 0 否則位置x處礦化屬于1型的概率為：[I(x)]*=tlaI(x)e[0,1]alal=prob*(re1型礦化，給定的n個數(shù)據(jù)｝l同上此處1型礦化品位為：[Z(x)|xe1[Z(x)|xe1型礦化]*a1a=11a1最后，位置x的估計品位[Z(x)]*為[Z(x)]*=￡[I(x)]*[Z(x)Ixe1型礦化]*11=1概率克立格指示克立格只用了截斷值Zk的指示數(shù)據(jù),指示協(xié)同克立格還考慮了所有的指示數(shù)據(jù)。如何能既考慮z-數(shù)據(jù)本身，還考慮它分布于[0,1]之間的標準化次序轉換數(shù)據(jù)，求得Z(x)的分布；概率克立格(PK)估計，實際上是Z(x)的ccdf模型，可以寫為簡單克立格的形式：[i(x;z)]*-F(z)=工九(x;z)-[i(x;z)-F(z)]kPKk akakka=1+工v(x;z)-[p(x)-0.5]akaa=1其中P(xa)=弘①)e[0,1]是數(shù)據(jù)'宀)的cdf變換，它的期望值為0.5，F(xiàn)(z)=Prob｛Z(x)Jz｝是z(xa)的平穩(wěn)cdf，Xa(x;zk)和Va(x;zk)分別是指示數(shù)據(jù)和均勻數(shù)據(jù)的協(xié)同克立格權值。注意這些權值依賴于位置(X)和截斷值(zk)。對應的簡單克立格方程組(PK方程組)要求推導和模擬(2K+1)個(互)協(xié)方差函數(shù)，即K個指示協(xié)方差函數(shù)、K個指示一均勻互協(xié)方差函數(shù)和一個均勻變換數(shù)據(jù)的協(xié)方差函數(shù)。這對于變差函數(shù)模擬的要求仍然較高，在應用中是一個很大的缺陷。軟克里格：馬爾科夫—貝葉斯模型指示克立格法產(chǎn)生后驗條件概率分布(cdf)最主要的優(yōu)勢在于它能夠考慮軟數(shù)據(jù)。只要軟數(shù)據(jù)或模糊數(shù)據(jù)能被編碼成先驗局部概率值，指示克立格就能把那種信息綜合到后驗概率值中。先驗信息可以是以下的某種形式：來自于局部硬數(shù)據(jù)z(xa)的局部硬指示數(shù)據(jù)i(xa;z)：如果z(x)<z,i(x;z)=1,否則i(x;z)=0a a a或者:如果xG類型變量s,i(x;s)=1,否則i(x;s)=0a k ak ak來自于有關局部數(shù)據(jù)z(xa)的硬不等式限制條件的附加信息的局部硬指示數(shù)據(jù)j(xa；z)，如果z(xa)G(aa，ba]，那么：0, 如果z<aaj(x;z)={未定義，如果zg[a,b]a aa1， 如果z>ba3?來自于有關數(shù)值z(xa)的先驗(預后驗)概率值的附加信息的局部軟指示數(shù)據(jù)y(xa;z)：y(x;z)=Prob{Z(x)<zI局部信息}g[0,1],a且y(x;z)豐F(z)：下面定義的局域先驗概率a適用于平穩(wěn)地區(qū)A的所有位置x的區(qū)域先驗信息：F(z)二Prob{Z(x)<z},VxgA對任何位置xGA，有關數(shù)據(jù)z(x)的先驗信息都能表達成以上四種類型中的任何一種。IK過程是根據(jù)鄰域的局部先驗cdf所提供的信息，把局部的先驗cdf通過貝葉斯原理變換成后驗cdf：[Prob{Z(x)<z|(n+n')}]*=九(x)F⑵+工九(x;z)i(x;z)IK0 a aa=1+遲V(x;z)i(x;z)a' aa'=1Xa(x;z)是n個鄰域硬指示數(shù)據(jù)有關的權值，Va'(兀乙)是n'個鄰域軟指示數(shù)據(jù)有關的權值，九0是區(qū)域先驗cdf的權值。為了保證無偏性，九0通常設為：九(x)=1一乙九(x;z)一乙V(x;z)0 a a'a=1 a'=1如果E{Y(x;z)}或E{J(x;z)}不同于F(z)，則需要考慮其他的無偏性條件。上式的ccdf模型能被看成是把不同類型的信息，包括硬的i和j指示數(shù)據(jù)及軟的y-先驗概率信息收集到一起的指示協(xié)同克立格。當不存在或忽略軟信息時(n=0)，表達式變成簡單IK的表達式。參考文獻肖斌，趙鵬大,侯景儒?地質統(tǒng)計學新進展[J].地球科學進展,2000,03:293-29