《地质统计学》课程读书报告地质统计学读书报告地质统计学包含经典统计学与空间统计学，按其基本原理可定义为：地质统计学是以区域化变量理论为基础，以变异函数为主要工具，研究那些在空间分布上既有随机性，又有结构性的自然现象的科学。

其为数学地质领域中一门发展迅速且有着广泛应用前景的新兴学科。

国内外的生产实践表明，地质统计学除了在异常评价、找矿勘探、矿体圈定、储量计算、采矿设计、矿山生产及地学科研等方面具有明显的优越性外，它在石油地质、第四纪地质、地层学、生物学、生态学、岩石学、地球化学、构造地质、地震地质、海洋地质、农业、水文地质、工程地质、古气候、古地理、环境、林业、医学等许多方面都有成功应用的实例。

地质统计学在不到50年的研究和实践中得到了很大的发展[1]。

一、理论研究及进展经历了数十年的发展，地质统计学的理论与方法研究有了很大的提高[2-3]。

包括：①从初期二维平面分析到三维立体空间的静态估计，发展到今天在时空域内对研究对象进行四维乃至更高维空间的动态估计和模拟。

Journel[4]将克立格法的估值问题，从一般矢量空间扩展到个原始数据的全部可测度函数所形成的矢量空间(希尔伯特空间)进行考察；②在单变量区域化变量理论的基础上，提出了适合多变量的协同区域化理论[4]；③发展了许多计算变异函数(或协方差函数)的方法；④线性地质统计学与非线性地质统计学共同发展；⑤参数地质统计学与非参数地质统计学相互补充。

Matheron[5]为首的参数地质统计学派以正态假设为前提，在协同区域化理论的基础上，提出多元地质统计学的基本思想。

Journel发展了无须对数据分布作任何假设的非参数地质统计学，提出了一些非参数地质统计学克立格方法；⑥由于时空多元地质统计学的研究得到重视，早期空间域静态建模技术的研究逐渐过渡到研究时空域多元动态条件模拟，各种模拟方法得到了发展；⑦早期的等因子模型的因子是埃尔米特多项式，它要求原始数据服从正态分布。

为了拓宽等因子模型的应用，Matheron提出了离散的等因子模型和连续的等因子模型，Rivoirard利用析取克立格技术建立了正交指标剩余模型，Lajauine和La ntuejoul等也提出了建立等因子模型的一些方法；⑧已有的地质统计学方法相互融合。

如指示克立格法与协同克立格法相结合形成指示协同克立格法;指示克立格法与因子克立格法相结合形成主分量指示克立格法；协同克立格法与其它不同的线性地质统计技术相结合形成各种协同克立格技术等[6]。

这里重点介绍一下多点地质统计学[7]。

多点地质统计学是相对于基于变差函数的两点地质统计学而言的。

在两点统计里，储集层相关性通过空间两点协方差( 变差函数) 进行描述。

在多点统计里，则是利用空间多个点组合模式进行描述。

空间多点组合样式称为数据样板，如果在空间点赋予了值，则为一个特定的空间多点组合模式，称为数据事件。

在建模时，对每一个未知点，估计在其处满足给定条件的数据事件出现的概率，随后抽样获得未知点处值或者数据事件，即完成单次模拟。

一旦所有节点得到访问，即完成一次模拟实现。

二、地质统计方法应用及进展克立格估计方法概述[8-9]：克立格法是法国Mathron以南非矿山工程师D.G.Krige的名字命名的一类方法。

简单地说是一种特殊的加权移动平均法。

是一种对空间分布数据求最优线性无偏内插估计量的方法。

根据估计式形式、区域化变量平稳性、分布及涉及的区域化变量个数等克立格方法有如下分类[10]：1.线性克立格法—非线性克立格法2.参数克立格法—非参数克立格法（指示、概率）3.平稳克立格法（普通）—非平稳克立格法（泛克立格）4.单变量克立格法—多变量克立格法(多元克立格)指示克立格法1.基本思路在地质、物探化探数据处理及相关储量计算中往往存在以下情况：(1)特异值的出现。

所谓特异值是指那些比全部数值的均值或中位数高得多的数值。

它既非分析误差所致；也非采样方法等人为误差引起，而是实际存在于所研究的母体之中，这些特异值只占全部数据的极少部分，但却控制了总金属资源量的很大比例，同时局限于一定的空间位置[11]。

(2)在一个研究区或一个矿床中存在几种不同类型的矿化作用，这也影响了品位和储量的精确估计。

(3)在矿山地质工作中经常要知道诸如：在区A及其附近，大于边界品位值Z的可选采矿单元占多大比例全部可选采矿单元的平均品位是多少等等。

为了解决上述问题，地质学家及统计学家进行了许多有成就的研究，而指示克立格法是一种更好的非参数地质统计学方法，它可以在不必去掉重要而实际存在的高值(特异值)数据的条件下处理各种不同的现象，而且给出在一定风险条件下未知量Z(x)的估计量及空间分布。

2.指示函数及其二阶矩指示函数定义：设在矿床D 上得到了某金属的品位值，约定其边界品位为Z ，则在D 上的每一样品点x ∈D 上定义一个Z 的阶梯函数：⎩⎨⎧>≤=Z x Z x Z x Z x Z x I )(,0)(,1);(点上品位值当点上品位值当在D 上的任一区域A ∈D 内，低于边界品位Z 的品位值Z(x)占区域A 的比例表示如下：]1,0[);(1);(∈=⎰dx z x i A Z A A φ);(Z A φ称为废石函数，而高于边界品位Z 的品位值Z(x)所占区域A 的比例为：},)({);(1);(A x Z x Z P Z A Z A ∈>=-=φψ称为矿石回收率。

指示函数的二阶矩。

在边界品位Z 确定的条件下，随机函数I(x;Z)服从二项分布，其期望值是平稳的，与x 无关。

∑∑<=<=====≤Z Z k Z Z k k k p Z x Z prob Z x Z })({F(Z)} )( prob{=Z)}I(x; E{二项分布的概率函数k n k p p k X prob --==)1(C }{kn非中心协方差是})()({)];(),;([);(Z x Z Z h x Z prob Z x I Z h x I E Z h K I ≤≤+=+=同时方差是(Z)F -F(Z)Z)(0;C =Z)}V ar{I(x ; 2I =中心协方差是)();();(211Z F Z h K Z h C -=变异函数是);()();();0(})];();({[21);(2Z h K Z F Z h C Z C Z x I Z h x I E Z h I I I I -=-=-+=γ I(x;Z) 满足平稳假设指示克立格方法废石函数的估计：估计量设计：设在矿床D 上有N 个有效数据，在D 上的一个域A ∈D 内有n 个有效数据},...,2,1,),({n A x x Z =∈ααα，给定边界品位值为Z ，则指示函数空间为：},......3,2,1),;({n Z x i =αα有);(Z A ϕ的估计值可以表示如下：∑==nZ x i Z Z A 1*);()();(αααλϕ若给一系列边界品位值l Z 时：{}L l Z l ,...,2,1,=，这时，估计量可写为： ∑==n l l l Z x i Z Z A 1*);()();(αααλϕ克立格方程组建立：在无偏和估计方差极小的条件下为求);();(**l Z A Z A ϕϕ或可得权系数),,...,2,1)((n Z =αλα满足的方程组，即指示克立格方程组：n Z Z A x Z x x Z n ni i ,...,2,11)();,();,()(11=⎪⎪⎩⎪⎪⎨⎧==+∑∑==αλγμγλβββαβαβ估计方差的表示：其指示克立格方差μγγλσααα+-=∑=n i i KI Z A A Z A x Z 12);,();,()();,(Z x x i βαγ表示在给定的边界品位Z 条件下，矢量h 的两个端点分别在信息域βαx x ,内的所有对点的平均指示变异函数值；);,(Z A A i γ表示在给定边界品位Z 条件下，矢量h 的两个端点在待估域A 内所有对点的平均指示便宜函数值；);,(Z A x i αγ表示在给定边界品位Z 条件下，矢量h 的一个端点在信息域αx 内，另一端点在待估域A 内所有对点的平均指示半变异函数值；μ为拉格朗日乘子。

协方差函数表示的方程及方差：用协方差函数分别表示如下：n Z Z A x C Z x x C Z n ni i ,...,2,11)();,();,()(11=⎪⎪⎩⎪⎪⎨⎧==-∑∑==αλμλααβαβαβμλσααα+-=∑=);,()();,(12Z A x C Z Z A A C n KI 若有l 个边界品位),...,2,1(L l Z l =就应该解l 个指示克立格方程组。

废石函数的意义：相应位置品位值低于Z 的概率废石函数的作用：用于概率估计及品位值估计待估域A 平均品位的指示克立格法估计，待估域A 平均品位的指示克立格法是应用某种克立格法求得);(Z A ϕ的线性估计值);(*Z A ϕ，最后得到待估域A的平均品位及储量。

设待估域A 约定在位置x 上，分两种情况讨论如何用指示克立格法求待估域的估计值*[Z(x )]。

由两种矿化类型组成的矿床中位置x 的平均品位的估计。

指示函数值定义如下： ⎩⎨⎧=类型矿化属于当类型矿化属于当2..........01..........1)(x x x I 则任一位置x 处的I(x)的平均值)(*x I 为：{}{}个数据类型矿化，给定的个数据n x prob n x I prob x I a x I n1|1)()()(***1*∈====∑=ααα{}个数据类型矿化，给定的n x prob x I 2)](1[**∈=-)(*x I ]10[，∈；)](1[*x I -]10[，∈。

位置x 上1类型矿化的品位是根据x 周围1n 个1类型矿化品位值),...,2,1)((111n x Z =αα求得：)(]1|)([11111*αααx Z b x x Z n ∑==∈类型矿化位置x 上2类型矿化的品位是根据x 周围2n 个2类型矿化品位值),...,2,1)((222n x Z =αα求得：∑==∈22221*)(]2|)([n x Z b x x Z ααα类型矿化最后，位置x 的估计品位*)]([x Z 可上述计算结果综合得到：*****]2|)()][(1[]1|)()[()]([型矿化型矿化∈-+∈=x x Z x I x x Z x I x Z 。

由l 个矿化类型)2(>l 组成的矿床中位置x 的平均品位的估计。

在实际工作中，任一位置x 上只能有一种矿化类型（如l 型矿化）占优势，这时：⎩⎨⎧=否则类型矿化属于当..........0..........1)(l x x I l位置x 处矿化属于l 型的概率为：{}个数据型矿化，给定的l n a a l a l n l x prob x I a x I l l ∈=∈=∑=*1*1]1,0[)()]([同上此处l 型矿化品位为：)(]|)([11111*αααx Z b l x x Z n ∑==∈型矿化最后，位置x 的估计品位*)]([x Z 为*1**]|)([)]([)]([型矿化l x x Z x I x Z ll l ∈=∑=概率克立格指示克立格只用了截断值k z 的指示数据，指示协同克立格还考虑了所有的指示数据。