地質統計學原理及其在礦床建模與儲量估算中的應用

會計學1地質統計學原理及其在礦床建模與儲量估算中的應用礦床品位建模及其應用需求礦體表面模型（礦化邊界）礦床品位模型勘探線剖面品位分析品位－噸位曲線分析第1頁/共98頁礦床品位建模及儲量估算流程組合樣品分析樣品確定礦床塊體模型參數選擇插值類型設置插值參數等確定搜索鄰域精度驗證滿意估值礦床品位模型否是第2頁/共98頁回顧：地理學第一定律及應用地理學第一定律：距離越近，兩點的地理現象相似性越大逐點移面內插：以待插點為中心，確定一個鄰域范圍，用該鄰域內的采樣點計算內插點的高程值。反距離加權平均法第3頁/共98頁內容介紹地質統計學簡介區(qū)域化變量變差函數建?？死锔癫逯邓惴ǖV體儲量估算應用第4頁/共98頁歷史背景與產生為解決礦床從普查勘探、礦山設計到礦山開發(fā)整個過程中各種儲量計算和誤差估計問題發(fā)展起來的。地質統計學是數學地質的重要分支，它首先由

D·G·克立格（Krige）工程師在南非的金屬礦產儲量計算中使用，后由法國馬特隆（G·Mathreon）教授領導的小組對此作了深入的研究并系統地總結出地質統計學的理論和方法。第5頁/共98頁地質統計學定義地質統計學(Geostatistics)是以區(qū)域化變量理論作為理論基礎，以變差函數作為主要工具，對既具有隨機性又具有結構性的變量（如品位值）進行研究的科學。其核心即“克里格法”，它是一種無偏的最小誤差的儲量計算方法。區(qū)域化變量變差函數克里格估值第6頁/共98頁與傳統儲量估算方法相比從傳統方法把部分鉆孔品位當作一個塊段的品位，從而使高品位估計偏高，低品位估計偏低，而且沒有考慮礦石品位的空間變異性，在計算塊段平均品位時，每個樣品的貢獻僅僅是若干個幾何因素。地質統計學方法避免了傳統方法的兩個缺陷。其加權因子是以礦床的各個方向變差函數的參數為基礎計算出來的,這種加權方法充分考慮了礦體形態(tài)的空間變化及其品位空間變化特征,并且采用了無偏的、誤差最小的數理統計方法計算樣品的加權因子和塊段的品位。第7頁/共98頁地質統計學的發(fā)展完善的理論基礎基本概念—區(qū)域化變量基本工具—變差函數基本假設—本征假設基本方法—克里格法方法與技巧不斷涌出析取克里格、多元高斯克里格和各種條件模擬技術的應用和發(fā)展地質統計學的軟件包及應用軟件不斷推出

美國斯坦福大學的GSLIB軟件包

挪威ODEN公司的STORM隨機建模軟件

加拿大的Geostat地質統計學軟件

澳大利亞的SurpacVision\Micromine礦山工程軟件第8頁/共98頁內容介紹地質統計學簡介區(qū)域化變量變差函數建?？死锔衿肺还乐档V體儲量估算應用第9頁/共98頁區(qū)域化變量G.馬特隆定義區(qū)域化變量是：一種在空間上具有數值的實函數，它在空間的每一個點取一個確定的數值，即當由一個點移到下一個點時，函數值是變化的.特征：隨機性和結構性

隨機性

結構性第10頁/共98頁區(qū)域化變量從地質及礦業(yè)角度來看，區(qū)域化變量具有如下性質：

（1）空間局限性：即它被限制在一個特定的空間（如一個礦體內）；該空間稱為區(qū)域化的幾何域；區(qū)域化變量是按幾何支撐定義的。

（2）連續(xù)性：不同的區(qū)域化變量具有不同的連續(xù)性，這種連續(xù)性是通過相鄰樣品之間的變差函數來描述的。

（3）異向性：當區(qū)域化變量在各個方向上具有相同的性質時稱各向同性，否則稱各向異性。

（4）相關性：一定范圍內、一定程度上的空間相關性，當超出這一范圍后相關性減弱以至消失。 （5）對于任一區(qū)域化變量而言，特殊的變異性是疊加在一般規(guī)律之上。第11頁/共98頁內容介紹地質統計學簡介區(qū)域化變量變差函數建模克里格插值算法礦體儲量估算應用第12頁/共98頁變差函數建模為表征一個礦床金屬品位等特征量的變化，經典統計學通常采用均值、方差等一類參數，這些統計量只能概括該礦床中金屬品位等特征量的全貌，卻無法反映局部范圍和特定方向上地質特征的變化。地質統計學引入變差函數這一工具，它能夠反映區(qū)域化變量的空間變化特征——相關性和隨機性，特別是透過隨機性反映區(qū)域化變量的結構性，故變差函數又稱結構函數。第13頁/共98頁變差函數定義我們可以把一個礦床看成是空間中的一個域，如圖中

為沿方向被矢量分割的兩個點，其觀測值分別為及，該兩者的差值就是一個有明確物理意義的結構信息，因而可以看成是一個變量。區(qū)域化變量在空間相距的任意兩點和處的值與差的方差之半定義為區(qū)域化變量的變差函數，記為

第14頁/共98頁變差函數定義定義：在任一方向，相距的兩個區(qū)域化變量和的增量的方差的一半。公式：變差函數值與區(qū)域化變量位置無關二階平穩(wěn)假設和本征假設第15頁/共98頁二階平穩(wěn)假設當區(qū)域化變量滿足下列兩個條件時，稱該區(qū)域化變量滿足二階平穩(wěn)：（１）在整個研究區(qū)內，區(qū)域化變量的期望存在且等于常數：

（常數）（２）在整個研究區(qū)內，區(qū)域化變量的空間協方差函數存在且平穩(wěn)：

當時，上式變成：

即它有有限先驗方差。第16頁/共98頁本征假設當區(qū)域化變量的增量滿足下列兩個條件時，稱該區(qū)域化變量滿足本征假設：（１）在整個研究區(qū)內，區(qū)域化變量的增量的期望為０：

（２）對于所有區(qū)域化變量的增量的方差函數存在且平穩(wěn)：

即要求的變差函數存在且平穩(wěn)第17頁/共98頁實驗變差函數計算其中:

=兩個樣本點間的距離

=樣本點屬性值(位置)=樣本點屬性值(位置)=樣本點數變差函數計算公式：第18頁/共98頁變差函數計算實例某地區(qū)規(guī)則采樣數據，數據為屬性值，樣本間距為100米。第19頁/共98頁實驗變差函數計算實例圖中表示的是東西方向，相距為100米的樣本點對。第20頁/共98頁實驗變差函數計算實例通過變差函數計算公式得到東西方向上，滯后距為100米的變差函數值。第21頁/共98頁實驗變差函數計算實例變差函數圖：滯后距100米的變差函數點024681012141618200100200300400500滯后距變差函數第22頁/共98頁實驗變差函數計算實例相距為200米的樣本點對。第23頁/共98頁實驗變差函數計算實例滯后距為200米的變差函數值。第24頁/共98頁變差函數計算實例變差函數圖：滯后距200米的變差函數點024681012141618200100200300400500滯后距變差函數第25頁/共98頁變差函數計算實例變差函數圖：滯后距300米、400米的變差函數點024681012141618200100200300400500滯后距變差函數第26頁/共98頁變差函數計算實例計算南北方向滯后距為100米、200米和300米的變差函數。第27頁/共98頁實驗變差函數計算實例南北方向400m點數過少，不參與計算。滯后距東西方向南北方向1001.465.352003.39.873004.3118.884006.7變差函數值第28頁/共98頁實驗變差函數計算實例變差函數圖：東西方向和南北方向024681012141618200100200300400500滯后距變差函數東西方向南北方向第29頁/共98頁實驗變差函數計算--距離和角度容差對于不規(guī)則采樣點：沿某一特定方向和特定滯后距上并沒有足夠的樣本點采用距離和角度容差解決該問題第30頁/共98頁實驗變差函數計算步長：4m步長容差：2m方位角：60傾角：0方位容差：22.5傾角容差：22.5水平帶寬：5m垂直帶寬：5m第31頁/共98頁實驗變差函數計算(3D)第32頁/共98頁變差函數的計算過程是由系統自行完成的，而合適的參數大小將直接影響計算結果的好壞。關于參數的選取第33頁/共98頁實驗變差函數參數選擇步長大小的選擇：步長間距太小步長間距較合適第34頁/共98頁實驗變差函數參數選擇步長個數的選擇：原則:步長大小*步長個數=研究區(qū)域長度的一半步長總間距第35頁/共98頁理論變差函數實驗變差函數并不能定量的反映數據空間相關性，需要對實驗變差函數進行擬合得到理論變差函數。理論變差函數三參數：塊金值/基臺值/變程(基臺值=先驗方差)SamplesnotspatiallycorrelatedSamplesSpatiallyCorrelated基臺值變程塊金值0.......(h)gh樣本空間相關樣本空間不相關第36頁/共98頁理論變差函數模型SamplesnotspatiallycorrelatedSamplesSpatiallyCorrelated球狀模型線性模型指數模型高斯模型第37頁/共98頁球狀模型球狀模型公式：接近原點處，變差函數呈線性形狀，在變程處達到基臺值。原點處變差函數的切線在變程的2/3處與基臺值相交。實驗變差函數在大多數情況下可以擬合成球狀模型。因此，球狀模型是應用最廣的一種變差函數模型。第38頁/共98頁指數模型指數模型公式：變差函數漸近地逼近基臺值，在實際變程處，變差函數為0.95，模型在原點處為直線。在原點處連續(xù)性最好，是一種較穩(wěn)定的模型。第39頁/共98頁高斯模型高斯模型公式：變差函數漸近地逼近基臺值，在實際變程處，變差函數為0.95，模型在原點處為拋物線。為一種連續(xù)性好但穩(wěn)定性較差的模型。第40頁/共98頁變差函數擬合用球狀模型、指示模型或高斯模型對實驗變差函數進行擬合。得到塊金值、基臺值和變程三個參數。第41頁/共98頁變差函數擬合球狀模型變程為4141m，指數模型變程為5823m，高斯模型變程為2884m觀察圖形：高斯模型擬合最好，其次是球狀模型。根據實際情況確定變差函數類型，結果因人而異。第42頁/共98頁變差函數擬合—過程第43頁/共98頁幾何各向異性基臺值相同變程不同在不同的方向具有相同的變異程度（基臺值相同）但具有不同的連續(xù)程度（變程不同）為幾何各向異性。第44頁/共98頁帶狀各向異性基臺值不同變程可同可不同

在一些不同的方向上具有不同的變異程度（基臺值不同）連續(xù)程度（變程）可以相同也可不同為帶狀各向異性。第45頁/共98頁變差函數結構套合不同方向結構套合幾何各向異性基本思路為通過線性變換將各向異性的坐標向量

轉化為各向同性的新坐標向量設這個線性變換為

，其中對于各向同性模型，，其中對于幾何各向異性變差函數，變化為矩陣形式第46頁/共98頁變差函數結構套合不同方向結構套合帶狀各向異性對于帶狀各向異性，采用分塊處理的方法。具體的變差函數模型公式為，其中對于

做和幾何各向異性相同的處理，對于

做如下處理

，對于做如下處理：

總的來說，對于帶狀各向異性的處理方法是將其看作是幾何各向異性進行坐標變換后，再分別對次軸和垂直軸方向上多出的基臺值進行疊加處理。第47頁/共98頁各向異性橢球各向異性橢球：主軸變程次軸變程垂直軸變程方位角傾角旋轉角度第48頁/共98頁幾何各向異性結構套合第49頁/共98頁變差函數表面圖Allpointsthatfallintheblockarepairedwiththepointat(x,y)tocreatethevariogrammaps.Thesizeoftheblockisthelagsize.第50頁/共98頁變差函數建模變差函數是區(qū)域化變量空間變異性的一種度量反映了空間變異程度隨距離而變化的特征可定量的描述區(qū)域化變量的空間相關性地理學第一定律：距離越近的點相似性越大第51頁/共98頁內容介紹地質統計學簡介區(qū)域化變量變差函數建?？死锔癫逯邓惴ǖV體儲量估算應用第52頁/共98頁如果要估算的值，一般情況下的值應該是的平均值，并且隨著距離的增大而減小。第53頁/共98頁克里格插值算法克里格插值算法建立在變差函數及結構分析理論之上適用條件是變差函數及相關分析的結果表明樣品間存在空間相關性其實質是利用區(qū)域化變量的原始數據和變差函數的結構特點，對未采樣點的區(qū)域化變量的取值進行線性、無偏、最優(yōu)估計。第54頁/共98頁克里格插值過程組合樣品分析樣品確定塊體模型參數選擇克里格類型計算實驗變差函數并擬合確定搜索鄰域交叉驗證滿意估值品位模型否是第55頁/共98頁塊體模型定義將整個研究區(qū)域劃分為多個規(guī)則小塊，分別對小塊屬性進行估值。起始點坐標塊大小塊個數第56頁/共98頁搜索鄰域確定搜索橢圓直接定義點數第57頁/共98頁交叉驗證交叉驗證的原理為將原始的樣品點去除，然后采用原始樣品點周圍的點來進行克里格估值得到原始樣品的估計值，最后做出原始樣品和估值樣品的散點圖，并對估值誤差進行統計。第58頁/共98頁克里格插值算法從礦業(yè)上的術語具體來說，它是根據一個待估塊段鄰域內的若干信息樣品的品位數據，在考慮了這些樣品的形狀、大小及相互位置關系，它們與待估塊段相互之間的空間位置等幾何特征，以及品位的變差函數模型所提供的結構信息之后，為了對該塊段品位作出一種線性、無偏、最小估計方差的估計而對每個樣品值分別賦予一定的權系數，最后進行加權平均來估計該塊段品位的方法。第59頁/共98頁克里格插值算法克里格插值算法:

B.L.U.E

-best,linear,unbiasedestimator

best=

最小估計誤差

linear=

線性估值方式（同距離反比估值）

unbiased=無偏估計，估計誤差之和為0

estimator=

估值方法第60頁/共98頁克里格插值算法定義其中

=待估點位置和其中一個鄰接點位置

=估算未知點用到的鄰接點個數

=和對應的預測平均值

=對應的克里格權重

克里格插值公式第61頁/共98頁克里格插值算法基礎關鍵在于確定鄰接權重最小方差限制條件無偏估計限制條件克里格插值公式第62頁/共98頁克里格插值類型最常用克里格的三種類型簡單克里格普通克里格泛克里格其區(qū)別在于的確定方式不同

非線性克里格

指示克里格（IndicatorKriging）

多元高斯克里格（Multi-Gausskriging)協克里格（Cokriging）塊克里格(Blockkriging)第63頁/共98頁簡單克里格插值應用實例六個樣本點數據，給出樣本點間距、樣本點屬性值和變差函數變差函數模型：球狀模型，塊金值0，基臺值0.78，變程4141m

Pnt1Pnt2Pnt3Pnt4Pnt5Pnt6Pnt1018973130244114001265Pnt2189701281145619702280Pnt3313012810152328003206Pnt4244114561523015231970Pnt514401970280015230447Pnt612652280320619704470第64頁/共98頁簡單克里格插值應用實例由簡單克里格插值公式得

第65頁/共98頁簡單克里格插值應用實例由簡單克里格插值公式得

第66頁/共98頁簡單克里格插值應用實例已知該采樣數據平均值為14.70六個采樣點數據的屬性值分別為

13.84,

12.15,

12.87,

12.68,

14.41,

14.59

由得

第67頁/共98頁簡單克里格插值應用實例第68頁/共98頁普通克里格插值普通克里格插值：

未知普通克里格估值公式為由無偏最優(yōu)估計限制條件，構建拉格朗日函數其中第69頁/共98頁拉格朗日乘數法用“拉格朗日乘數法”求函數f(x,y,z)在條件φ(x,y,z)=0下的極值方法（步驟）是：

1.做拉格朗日函數L=f(x,y,z)+λφ(x,y,z),λ稱拉格朗日乘數

2.求L分別對x,y,z,λ求偏導,得方程組,求出駐點P(x,y,z)

如果這個實際問題的最大或最小值存在,一般說來駐點唯一,于是最值可求.第70頁/共98頁普通克里格插值求偏導分別得到下列公式得到求的方程組估值方差計算公式第71頁/共98頁指示克里格估值在地質、物化探數據處理及礦產儲量計算中影響計算精度的因素有很多，但主要有以下幾個問題：(1)特異值的出現，所謂特異值是指那些比全部數值的平均值或中位數高得多的數值，它既非分析誤差所致，也非采樣方法等人為誤差引起。而是實際存在于所研究的母體之中。這些特異值只占全部數據的極少部分，但卻控制了總金屬資源量的很大比例。(2)在一個研究區(qū)域或一個礦床中存在幾個不同類型的礦化作用，這也影響了品位和儲量的精確估計。為解決上述問題，指示克里格法應運而生，它是在不必去掉重要而實際存在的高值數據的條件下來處理不同的現象，而且給出在一定風險概率條件下未知量的估計值及空間分布。第72頁/共98頁指示克里格估值指示克里格是一種非參數地質統計學方法。它是根據一系列的臨界值，例如邊界品位，先對原始數據如下進行轉換然后對轉換后的數值求變差函數、進行克里格估值。第73頁/共98頁總結簡單克里格（整體平穩(wěn)）普通克里格（局部平穩(wěn)）泛克里格（整體存在趨勢）指示克里格（數據不平穩(wěn)存在極值）協克里格（存在輔助變量）已知常量但未知其中常量但未知第74頁/共98頁幾點注意內容變差函數參數

塊金值：塊金值越小，距離越近的點越重要，這樣會導致權值的變化范圍變大（從負值到大于1的值變化），使數據出現異常。塊金值越大，估值結果越平滑。

變程：變程小于任意兩點之間的距離，變差函數為純塊金模型，隨著變程增大，已知點的位置以及叢聚性變的重要。如果變程很大存在基臺值，則變差函數相當于純塊金模型，如果變程很大但不存在基臺值，則變差函數為線性模型。

比例：比例大小只與克里格估值方差有關，與估值結果無關。

形狀：球狀模型與指數模型在原點位置接近線性關系，指數模型在原點處更加陡峭一些，與變程較小的球狀模型相似。高斯模型在原點處為拋物線形式，這種模型要求原始變量要有很高的連續(xù)性，否則會出現大量負的權值情況，使估值結果極不穩(wěn)定。

各向異性：以各個方向變程不同反映出來，可以理解為對原始數據坐標位置的變化，通過幾何校正轉化為各向同性情況。第75頁/共98頁幾點注意內容屏蔽效應第一個點在任何情況下獲得的權值更大一些。第二個點的權值小于第一個點，在克里格估值過程中，第二個點的權值可能變?yōu)樨撝?，這樣第二個點的信息變得重要。例如：當一個負值權重可能隱含某種趨勢（如果第一個點比第二個點要小，那么估值點可能小與這兩個點。隨著塊金值的增大屏蔽效應減弱。所用樣品點對估值結果作用基本相同。0.220.20.18第76頁/共98頁幾點注意內容搜索策略最小值：作用是限制克里格估值過程中數據個數，如果在搜索橢球范圍內，數據個數過少則不進行估值，否則估值結果可能存在很大誤差。

最大值：估值過程中用到的數據越多，估值結果越平滑。由于數據的局部平穩(wěn)性，數據個數不能太多，否則會使不存在相關性的數據仍然參與到估值過程中。

搜索半徑：搜索半徑越大估值結果越平穩(wěn)，搜索半徑與最大值相互影響。搜索半徑設置與變差函數變程相關。第77頁/共98頁幾點注意內容克里格估值類型

簡單克里格：簡單克里格需要一個固定的平均值。在估值過程中平均值被賦予一個相當大的權值，從而導致估值結果過渡平滑。

普通克里格：與簡單克里格類似，但平滑性減弱。普通克里格估值權值之和為1，可以合理的減少叢聚效應。

泛克里格：與普通克里格類似，但估值結果邊緣會存在一個奇異值，這是由趨勢面構造的不合理造成的。

塊克里格：由于塊克里格估值的平滑效應會使克里格估值結果方差減小。隨著塊的增大，估值結果越來越平滑。第78頁/共98頁內容介紹地質統計學簡介區(qū)域化變量變差函數建模克里格品位估值礦體儲量估算應用第79頁/共98頁礦床品位建模及儲量估算流程組合樣品分析樣品確定礦床塊體模型參數選擇克里格類型計算變差函數并擬合確定搜索鄰域精度驗證滿意估值礦床品位模型否是第80頁/共98頁克里格法資源儲量估算流程－關鍵步驟否否是是數據分析數據變換選擇泛克里格估值選擇其它克里格估值剔出趨勢導入組合樣品數據是否服從正態(tài)分布是否存在趨勢數據分析與變換第81頁/共98頁克里格法資源儲量估算流程－關鍵步驟計算實驗變差函數并繪制變差函數表面圖是否服從各向異性單方向變差函數擬合確定最大變程方向多方向變差函數擬合是否球狀模型擬合變差函數擬合單方向變差函數擬合多方向變差函數擬合指數模型擬合高斯模型擬合幾何各向異性帶狀

各向異性變差函數計算及擬合數據分析與變換結