论文导读：本软件系统将总的功能实现分解为原始数据输入、理论变差函数计算、变差函数拟合、拟合参数调整、克里格估值、软件操作帮助六个模块，如图1所示。5.1本软件以地学处理软件中克里格估值功能的实现为参考，对点克里格估值的相关模型进行编程，最终实现了区域内任意点的克里格估值。5.3本软件可以对具有空间相关性的地质属性进行估值，实现目标区域地质储量估算。

　　关键词：空间信息统计学，点克里格估值，地质储量估算，软件

　　1引言

　　克立格法是由法国G·马特隆教授以南非矿山地质工程师D·G·克立格的名字命名的，是一种对空间分布数据求最优、线性、无偏内插估计量的方法。克里格法是地质统计学（空间信息统计学）的核心，最早应用于矿山开采。

　　目前地质统计学在理论和软件方面发展迅速。理论方面：Journel将克立格法的估值问题，从一般矢量空间扩展到 n个原始数据的全部可测度函数所形成的矢量空间（希尔伯特空间）进行考察。为地统计学由低维向高维动态发展奠定了基础。Matheron为首的参数地质统计学派以正态假设为前提，在协同区域化理论的基础上，提出多元地质统计学的基本思想。Journel发展了无须对数据分布作任何假设的非参数地质统计学，提出了一些非参数地质统计学克立格方法。同时，多学科相互渗透，拓展了地质统计学的理论体系，使地统计学呈现出以下特点：地质统计学与数学的结合更加紧密；非线性科学思想向地质统计学渗透；稳健统计学与地质统计学的结合更加密切。

　　软件方面：国外的空间信息统计学软件主要有：法国巴黎高等矿院地质统计学研究中心研制出两种大型软件系统：ISATIS系统及HERESIM系统。这两个系统对油气田地质及开采工程等方面考虑得尤为周到。美国斯坦福大学地球科学应用系研制的GEOEAS、GEOSTAT TOOLBOX程序包及1922年、1996年又推出了GSLIB两种版本等。论文发表。另外有法国的BRBM公司、加拿大的GEOSTAT公司与加拿大矿产与能源中心、美国的Mewnont公司等也制作了一些地质统计学软件包，并在生产实践中得到了应用，产生了较大的经济效益。国内的空间信息统计学软件主要有：CGES中文地勘系统软件是原武警黄金指挥部黄金地质研究所在加拿大IGC公司GLS和RCS软件系统基础上全面汉化改造并扩充功能而完成的。北京科技大学侯景儒教授领导下编制的地质统计学程序集，该程序集是为地质统计学理论、方法研究及地质勘探使用而研制的。另外还有李裕伟教授领导下研制的 KPX2.0 固体矿产勘查评价自动化系统，王仁铎、胡光道教授领导下编制的一些地质统计学软件包等。

　　2地统计学原理

　　地统计分析的理论基础包括前提假设、区域化变量、变异分析和空间估值四个部分。地统计学分析的核心是假设区域化变量满足二阶平稳假设或本征假设，依据采样点确定目标对象随空间位置的变化规律，进而估算待估点的属性值。

　　实验变差函数是反映目标对象空间变化规律的关键，也是地统计学的基本工具。根据观测数据构造变差函数的估计值（见式1）。

　　实验变差函数的基本计算公式：

　　（1）

　　为了估算区域化变量的未知值，需要将实验变差函数拟合成相应的理论变差函数模型参与Kriging 估值。

　　地统计学将变异函数理论模型分为三大类：

　　第一类是有基台值模型，包括球状模型、指数模型、高斯模型、线性有基台值模型和纯块金效应模型；

　　第二类是无基台值模型，包括幂函数模型、线性无基台值模型、抛物线模型；

　　第三类是孔穴效应模型。

　　球状模型是地统计分析中应用最广泛的理论模型，本文采用球状模型进行变差函数拟合。

　　3软件设计

　　本软件系统将总的功能实现分解为原始数据输入、理论变差函数计算、变差函数拟合、拟合参数调整、克里格估值、软件操作帮助六个模块，如图1所示。

　　在理论变差函数计算模块中，将存入预定义数组中的数据进行筛选，满足距离误差限和角度误差限的点对参与函数值的计算。

　　在变差函数拟合部分，通过软件系统中的拟合模型对各方向的实验变差函数值进行拟合，得到拟合参数值；进而绘制对应的理论变差函数拟合图。

　　变差函数拟合所得的拟合参数是球状模型计算的基础，也是对地质储量进行克里格估值的基础。因此，拟合参数直接关系到下一步克里格估值的准确性。本模块是在直接法拟合实验变差函数的基础上，根据实际情况调整拟合参数，参数可以反复比较修改，其对应的拟合图可以在窗口上显示，直到达到满意的效果。

　　克里格估值部分是运用克里格估值模型，实现对研究区目标待估点的地质属性值进行估算，并将所得结果在相应的窗口上显示出来。

　　图 1 系统总体功能结构图

　　4实现过程

　　软件编程实现的关键是变差函数部分和克里格估值部分。

　　4.1实验变差函数计算模块 （2）

　　其中：——距距离为，且在一定距离误差限内的数据对的数目

　　L——可以组成数据对的i的变化范围

　　——在一定距离误差限内的距离权系数

　　（3）

　　在实际工作中，观测点的分布并不完全呈网格分布，在求某一方向的实验变差函数值时，由于观测点不一定完全位于这一方向的同一直线上，本部分编程采用距离误差限和角度误差限技术对其进行限制。如求X方向的实验变差函数，则从某一点出发，位于-～之间的任意点都可以看成该点X方向上的点。论文发表。采用用动态数组存放观测点坐标及变量值。

　　4.2变差函数拟合模块

　　球状模型变差函数：

　　（4）

　　为了确定拟合参数，C，a的值，本模块用直接法计算各个参数。论文发表。具体实现步骤：

　　﹒用实验变差函数开始两点的连线与纵坐标的交点的纵坐标值作为理论变差函数的块金值

　　﹒以实验变差函数倒数第二点的值作为拱高值+ C

　　﹒对应的横坐标作为变差a

　　用直接法求出拟合参数，C，a，代入预定义的过程（球状变差函数模型）中，程序运行后，可将拟合曲线绘制在窗口上。此部分可以实现与用户的良好交互，用户可以根据工作的实际经验，输入新的拟合参数进行拟合调整，直到满意为止。

　　4.3克里格估值模块

　　（1）模型选取

　　①克里格估值模型：

　　（5）

　　其中：——观测信息值

　　——估值权系数

　　——估计值

　　②克里格方程组：

　　（6）

　　其中：n——用于估值的观测点的数目

　　——观测点坐标（）

　　——待估点坐标（x,y）

　　——两点间的协方差函数

　　——求条件机制的拉格朗日乘数

　　——克里格估值权系数

　　③点估值克里格方差：

　　（7）

　　其中：——变差函数基台值

　　将其转换为矩阵形式，令：

　　则普通克立格方程组为：（8）

　　其中：==（9）

　　代表点与点j间的距离

　　解方程组，可得：其估计方差为：（10）

　　将克里格方程组（式6）转换成矩阵形式（式8），求方程组的解就转换为求矩阵的解。本部分用高斯消去法求矩阵的解，将矩阵的解代入式5中，进而得到待估点的克里格值和估值方差。

　　本程序对二维区域内的点进行克里格估值。在对选取的点进行估值时，先对周围点进行扫描，取距离待估点最近的若干个点作为信息样品点参与估值。程序计算待估点的克里格值和估值方差。

　　5结论

　　5.1本软件以地学处理软件中克里格估值功能的实现为参考，对点克里格估值的相关模型进行编程，最终实现了区域内任意点的克里格估值。

　　5.2实现了将空间信息点的克里格估值分解为若干个主要模块，最终在实现点克里格估值的同时，用户可根据需要在不同阶段实现对软件运行产生的数据文档进行调用。

　　5.3本软件可以对具有空间相关性的地质属性进行估值，实现目标区域地质储量估算。

　　5.4本软件具有良好的交互性，界面简洁，有针对性，可以提高用户数据处理的效率。

　　6 参考文献：

　　[1] 侯景儒，尹镇南，李维明，等。空间信息统计学。北京：地质出版社，1998，45

　　[2] Armstrong M.Geostatistics-Proceedings of the International GeostatisticsCongress.France:Kluwer Academic Publishers,1958,Vol,1-2

　　[3] Matheron G.Thetheory of Variables and its Application.France:Centre of GeostatisticsFontaine,1970

　　[4] Journel AG.Nonparameteric estimation of spatial distribution. Math Geol,1983,445-468

　　[5] 肖斌，赵大鹏，侯景儒。地质统计学新进展。地球科学进展，2000,15（3），295

　　[6] 孙洪泉。地质统计学及其应用。 北京：中国矿业大学出版社，1990，39-41

　　[7] 孙洪泉。地质统计学及其应用。 北京：中国矿业大学出版社，1990，66-69