一种确定资源储量估算的品位域的方法与流程

本发明属于地质采矿技术领域，尤其涉及一种确定资源储量估算的品位域的方法。

背景技术：

金属矿床的形成受岩性、构造、流体等多个地质要素制约，成矿地质作用往往具有多期多阶段特征，造就了不同地质特征和不同品位分布特征的矿化带。根据控矿地质要素和矿石特征的不同，可以将矿化带划分出不同的矿化域(geological domains or mineralization zones)。例如，不同工业价值共伴生组份的矿化体，不同空间分布的矿化体，层状矿或脉状矿等不同产状特征的矿化体，或者氧化矿与硫化矿等不同矿石类型的矿化体等。资源储量估算就是在这些特定的矿化域内进行。经验表明，矿化元素的分布常呈现出非正态分布，引起这种非正态分布的因素有(侯景儒等，1994)：(1)数据中存在有特异值(Outliers，包括特高值和特低值)；(2)一个矿床存在着多期次矿化或若干个矿化类型，即所谓的混合总体。特异值在全部观测值中只是极少部分，但它对全部观测值的统计结果影响极大。理论上讲，特定矿化域内所有的工程采样数据都应该参与资源储量估算。但在实际操作中，为了最大程度地减少特异值对估值的影响，需要确定资源储量估算的品位域(grade domains)，也即剔除了特异值的样品品位区域。目前，国内外尚没有统一方法来确定该品位域，不过一个好的方法，无疑应当使品位域内的品位数据同时具有地质上和统计上的连续性。在实际操作中，国内勘查规范明确推荐了边界品位、工业品位标准以及特高品位的识别 方式。国外同行在进行资源储量估算时，也习惯于将研究对象划分出废石、低品位域和高品位域(Emery&Ortiz，2005)；总的来讲，特高品位识别和处理方法比较多，相对也较为成熟(丁旭东，1996；GOOVAERTS，1997；Parish，1997；Sinclair and Blackwell，2002；Rossi and Deutsch，2014)，而特低值一般不进行处理，地质解译和矿化体圈定多取决于地质师的地质认识和经验判断。

1、勘查规范推荐指标

国内现行勘查规范根据众多矿山实际生产数据，提出了边界品位、工业品位和矿床平均品位的推荐指标。采用的是3分法：边界品位(划分废石和低品位矿)、工业品位(划分低品位矿和工业矿)、特高品位(划分工业矿和特异值，高于特高品位的数据需要处理)。因此，现行规范推荐的资源储量估算品位域是边界品位‐特高品位之间的区域。具体估算时，边界‐工业品位之间划分为低品位矿，工业‐特高品位之间划分为工业矿。边界品位和工业品位指标确定是根据国内矿山历史生产数据统计给出的推荐值，参考了特定矿种的历史价格因素。《固体矿产地质勘查规范总则(GB/T13908‐2002)》规定，特高品位值一般用矿体平均品位的6‐8倍来衡量，当矿体品位变化系数大时采用上限值，变化系数小时采用下限值。

2、特高品位识别技术(侯景儒等，1994；Marinho et al.，2009；Nowak et al.，2013；Rossi and Deutsch，2014)

(1)根据偏离均值的程度来确定特高品位

①按照样品均方差的倍数来确定是否是特高品位。即特高品位值大于或等于m+3σ，其中m为均值，σ为均方差。

②按照样品品位变化系数(v＝σ/m)来识别特高品位。品位分布很均匀的沉 积矿床，品位变化系数＜20％，特高品位取矿床平均品位的2‐3倍；品位分布均匀的沉积和变质矿床，品位变化系数20‐40％，取4‐5倍；品位分布不均匀的大部分有色金属矿床，品位变化系数40‐100％，取8‐10倍；品位分布很不均匀的有色、稀有、贵金属矿床，品位变化系数100‐150％，取12‐15倍；品位分布极不均匀的稀有、贵金属、放射性元素矿床，品位变化系数＞150％，取值＞15倍。

③估值邻域法。由D.G.克立格和D.M霍金斯把地质统计学的基本思想用于识别和处理特异值的方法，它是把被识别的观测值(称为可疑样品)置于一个空间连续矿化域的背景上进行研究。

(2)根据偏离正常的分布特征来确定特高品位

①影响系数法。该法是通过确定特异值的影响程度进行特异值的识别和处理。

②十分位数分析(Decile Analysis)。Parrish(1997)考察了金属含量的十分位分布以及它们在总金属含量中的比例。如果第10个十分位的金属含量超过了总金属含量的40％，且比前一个十分位含量高出2倍多，那么就存在特高品位。进一步的细分第10个十分位，计算出它的百分位数，如果第100个百分位的金属含量超过了总金属含量的10％，那么移走所有的高品位数据，直至某个百分位金属含量不超过10％为止。这样得到了特高品位的下限值。

(3)图解法

①累计概率分布(CDF)。在分布密度函数曲线上，将拐点对应值作为特高品位的下限值。

②直方图或概率分布图。做出频率直方图和概率分布图(推荐是对数概率分布)，分析直方图中横轴上孤立分布的开始值，或者是概率分布曲线上曲线 梯度变化或者断点处。

③Cutting Curve Plots。Cutoff vs小于cutoff的品位均值，横坐标为cutoff品位，纵坐标为低于cutoff的品位数据均值，均值应趋向于接近水平，那么均值偏离水平线的拐点就是特高品位的下限值。

现有技术存在主要问题如下：

勘查规范推荐标准：边界品位和工业品位的设定过分强调资源储量估算经济性，没有充分考虑品位本身的地质分布和统计分布规律，用较高的圈矿边界品位来剔除低品位样品，无疑会引起矿化域内资源储量的高估。此外，由于矿化元素一般都不符合正态分布，利用不具备平稳性的均值的简单倍数来预测特高品位，可靠性较差。

m+3σ和品位变化系数方法：对于不存在统计均质性的矿化体，均值的不平稳性为特高值预测会带来风险。

直方图或概率分布图：方法的特点是简单快捷，问题是概率曲线的拐点或断点以及频率直方图的孤立分布也可能指示了另一个总体(矿化事件或矿化类型)。统计结果与分组(bin)大小有关，具有一定的人为性。

估值邻域法：在使用该法时，其中的n、σ²以及替代值的确定都是人为确定，没有一致的标准；在面对大量数据时，可操作性不强。估值邻域法只考虑与疑似特高值的相关背景值。

影响系数法：采用此法，过程相对比较繁琐，而且参数K值的确定存在模糊性。

十分位数分析(Decile Analysis)：考虑了金属量且快速，但仅适用于某些特定矿床。

Cutting Curve Plots：这一方法要求至少要有500‐1000个样品数据。

技术实现要素：

本发明针对现有技术存在的问题，提供一种利用分形来确定资源储量估算品位域的方法。

本发明采取如下技术方案：

一种确定资源储量估算的品位域的方法，其特征在于，所述方法包括：

(1)数据整理：将特定矿化域内的所有探矿工程品位数据进行等长组合，消除不同品位数据的权重差异；

(2)构建分形模型：根据Number‐Size模型，若样品品位数据满足分形分布则有，N(≥r)＝Cr^‐D，其中r为品位值，N(r)为品位大于r的样品个数，C为常数，D为分维值；对于给定不同的品位r，统计大于或等于r的品位数据个数N(r)，在双对数图logr～logN(r)上做散点图；用最小二乘法拟合，拟合直线斜率的绝对值就是分形维数D；把具有分形分布特征的品位区间称作分形标度区间；如果品位分布只在一个特定标度区间内具有分形分布特征，叫做单标度分形；如果在多个标度区间内进行拟合，叫做多标度分形，对应多个分维值；所有分形标度区间加起来的变化范围称为分形范围；

(3)计算品位域阈值：分形范围的上下限就是矿化品位域的阈值，用以区分矿化域和特异值；分形分布范围内多标度分形拟合直线的交点为理论上的多标度分形的品位阈值，用于区分低品位域、中品位域和高品位域。

根据上述的方法，其特征在于，所述步骤(1)中，等长样长一般为基本取样样长，组合后的品位数据按照从小到大的顺序排列。

本发明的有益效果：

1、分形方法将地质特征与品位统计特征有机结合。既可以有效区分了资源储量估算的品位域和特异值(特高值和特低值)，又可以对矿化域进行细分出低品位域、中品位域和高品位域，与成矿作用的多期多阶段或多种矿化类型混合的特征相对应，给出了品位分布的地质解释。

2、分形方法可以同时给出品位域的上限和下限值(特高值和特低值)。

3、分形方法是从全局考虑，而不是从局部或者一个邻域去考虑，适用于各种形式特高值的识别。

4、方法的可操作性和重复性强，人为因素相对较少，易于推广。

附图说明

图1是本发明的分形方法确定资源储量估算品位域的理论模型；

图2为K‐氧化矿的lnr～lnN(r)散点图；

图3为K‐氧化矿的品位域分形模拟结果图；

图4为M‐氧化矿的品位域分形模拟结果图；

图5为K‐硫化矿的品位域分形模拟结果图。

具体实施方式

本发明引入的分形方法拟针对传统方法存在的诸多问题：如方法只能针对特殊研究对象；需要大量数据支撑；存在明显人为因素影响；缺少对地质因素的考虑，无法区分多个总体；均为确定特高值的方法，无法确定品位域的下限等。

分形方法是近30多年来被广泛应用于地质地球化学数据处理的一种新方法。大量实践证明，矿化元素在特定的标度区间内具有分形分布特征，指示了特定矿化作用或矿化类型，并且利用分形方法可以有效区分这种多期次或不同类型的矿化叠加作用(Blenkinsop,1991；Cheng,1996；Afzal etal.,2013；Heidari et al.,2013)。基于分形理论，可以获得矿化元素具备分形分布特征的品位范围，其上、下限值作为特异值(特高值和特低值)，进而获得矿化域的品位范围，也即资源储量估算的品位域。如果存在多次成矿作用或不同类型矿化的叠加，那么这一品位域还可以进一步通过多标度分形方法细分，来获得不同成矿作用 或矿化类型的品位区间，借以划分低品位域、中品位域和高品位域。因此，利用品位分布具有分形分布特征的品位范围与地质矿化品位域之间有着良好的对应关系，可以有效确定资源储量估算的品位域。

分形分布模型的建立以及品位域计算的具体过程如下：

1、数据整理。将特定矿化域内的所有探矿工程品位数据进行等长组合，以消除不同品位数据的权重差异。等长样长一般为基本取样样长。组合后的品位数据按照从小到大的顺序排列。

2、构建分形模型。根据Number‐Size(N‐S)模型(Turcotte,1992)，若样品品位数据满足分形分布则有，N(≥r)＝Cr^‐D，其中r为品位值，N(r)为品位大于r的样品个数，C为常数，D为分维值。对于给定不同的品位r，统计大于或等于r的品位数据个数N(r)，在双对数图logr～logN(r)(或lnr～lnN(r))上做散点图。用最小二乘法拟合(拟合优度≥0.950)，拟合直线斜率的绝对值就是分形维数D(Heidari et al.，2013)。并不是所有数据点都符合分形分布特征(Blenkinsop,1991)，把具有分形分布特征的品位区间称作分形标度区间。如果品位分布只在一个特定标度区间内具有分形分布特征，叫做单标度分形；如果在多个标度区间内进行拟合，叫做多标度分形，对应多个分维值。所有分形标度区间加起来的变化范围称为分形范围(Fractal range)。

3、计算品位域阈值。分形范围的上下限就是矿化品位域的阈值，用以区分矿化域和特异值。分形分布范围内多标度分形拟合直线的交点为理论上的多标度分形的品位阈值，用于区分低品位域、中品位域和高品位域。

图1中，横坐标为lnr，纵坐标为lnN(r)，其中r为给定的品位值，N(r)为大于或等于r的品位数据个数。在lnr‐lnN(r)散点图上，散点呈自然弯曲分布，横跨的范 围是测量范围(measured range)，也即原始数据在双对数图上的分布范围。运用最小二乘法拟合，反复调整lnr范围，直至拟合出一条线段或者多段线，A和E为线段或多段线的端点，B、C、D为相邻两条拟合直线的交点。通过反算得到A、B、C、D、E各点对应的品位值r_A、r_B、r_C、r_D和r_E，那么r_A～r_E就是分形范围，r_A为分形范围的下限值，r_E为分形范围的上限值。相应地，r_A～r_E就是资源储量估算的品位域，r_A为品位域的下限值，r_E为品位域的上限值。多段线的各段直线的斜率绝对值就是分维值D_a、D_b、D_c、D_d，其对应的分形标度区间为r_A～r_B、r_B～r_C、r_C～r_D、r_D～r_E，分别为从低品位到高品位的品位域，反映了多期次矿化叠加或者不同类型矿化的混合。

实施例1

本次申请以刚果(金)铜钴矿带某层控型铜矿为例，矿床主要经历了同生沉积、热液改造和表生富集3期成矿作用。该铜矿床分为两个相对独立的矿段K矿段和M矿段，K矿段上部为氧化矿，下部为硫化矿；M矿段均为氧化矿。根据宏观地质特征，将矿床分为3个矿化域，即K‐氧化矿、M‐氧化矿和K‐硫化矿，分别用分形方法确定其资源储量估算的品位域。

1、先计算K‐氧化矿的品位域。

(1)数据整理。将K‐氧化矿原始数据在Micromine软件中进行1m等长组合(也可以在其他矿业三维软件中进行)，组合后的数据按照从小到大排列，共有数据583个，范围从0.01～21.37％。

(2)构建分形模型。令r＝0.01,0.05,0.1,0.2,……,21.1,21.2,21.3,21.3(r的取值取决于对数据精度的要求)；分别统计N(≥r)，并计算lnr和lnN(r)。

表1

做出lnr～lnN(r)的散点图(参见图2)。图中，横坐标为lnr，纵坐标为lnN(r)，其中r为给定的品位值，N(r)为大于或等于r的品位数据个数。在lnr‐lnN(r)散点图上，散点呈自然弯曲分布。

用最小二乘法进行拟合，需要进行多次调整。具体的方法是：从最小的lnr开始，选择至少4个点进行拟合，如果拟合优度R²<0.950，向右平移点，增加下 限值的同时增加上限值，直至R²≥0.950。如果拟合优度R²≥0.950，则保持下限值不变，增加上限值，直至R²<0.950，则该上限值的前一个数值，即为实际的分形标度区间的上限值。计算得到第一个分形标度区间的实际值为0.3～1.2％，其中0.3％为品位域的下限值，分维值D＝0.214。

表2

同理，可以计算出第二个分形标度区间的实际值为1.3～13.8％，其中13.8％为品位域的上限值，分维值D＝1.960。

因此，分形范围(Fractal range)为0.3～13.8％。参见图3，运用最小二乘法拟合，反复调整lnr范围，直至拟合出一条线段或者多段线，A和C为多段线的端点，B为相邻两条拟合直线的交点。通过反算得到A、B、C各点对应的品位值r_A、r_B和r_C，那么r_A～r_C就是分形范围，r_A为分形范围的下限值，r_C为分形范围的上限值。相应地，r_A～r_C就是资源储量估算的品位域，r_A为品位域的下限值，r_C为品位域的上限值。

(3)计算品位域阈值。

通过两条拟合直线求出交点B的lnr值为0.686，反算出r_B＝2.0，即为低品位域和高品位域的拐点理论值。

表3

综上，可以确定K‐氧化矿资源储量估算品位域为0.3～13.8％，进一步细分为0.3～2.0％的低品位域和2.0％～13.8％的高品位域，特异值范围为(等长组合数据中的斜体部分)：0.01～0.29％(特低值)和14.11～21.37％(特高值)。

同理计算M‐氧化矿和K‐硫化矿的资源储量估算品位域。

2、M‐氧化矿的品位域计算。

M‐氧化矿1m等长组合后的数据如下，品位范围为0.02～33.41％。

品位分布的分形模型及品位域计算结果如下。参见图4和下表4。对lnr～lnN(r)散点图运用最小二乘法拟合，反复调整lnr范围，直至拟合出一条线段或者多段 线，A和D为多段线的端点，B、C为相邻两条拟合直线的交点。通过反算得到A、B、C、D各点对应的品位值r_A、r_B、r_C和r_D，那么r_A～r_D就是分形范围，r_A为分形范围的下限值，r_D为分形范围的上限值。相应地，r_A～r_D就是资源储量估算的品位域，r_A为品位域的下限值，r_D为品位域的上限值。

表4

因此，可以确定M‐氧化矿资源储量估算品位域为0.3～17.3％，进一步细分为0.3～3.35％的低品位域，3.35％～9.5％的中品位域，以及9.5％～17.3％的高品位域，特异值范围为(等长组合数据中的斜体部分)：0.02～0.29％(特低值)和17.80～33.41％(特高值)。

3、K‐硫化矿的品位域计算。

K‐硫化矿1m等长组合后的数据如下，品位范围为0.13～31.09％。

品位分布的分形模型及品位域计算结果如下。参见图5和表5，对lnr～lnN(r)散点图运用最小二乘法拟合，反复调整lnr范围，直至拟合出一条线段或者多段线，A和D为多段线的端点，B、C为相邻两条拟合直线的交点。通过反算得到A、B、C、D各点对应的品位值r_A、r_B、r_C和r_D，那么r_A～r_D就是分形范围，r_A为分形范围的下限值，r_D为分形范围的上限值。相应地，r_A～r_D就是资源储量估算的品位域，r_A为品位域的下限值，r_D为品位域的上限值。

表5

因此，可以确定K‐硫化矿资源储量估算品位域为0.1～8.1％，进一步细分为0.1～0.4％的低品位域，0.4％～1.4％的中品位域，以及1.4％～8.1％的高品位域，特异值范围为(等长组合数据中的斜体部分)：18.52～31.09％(特高值)。

4、与给定工业指标和传统方法估算结果的比较

注：传统方法仅给出特高值，分形方法给出品位域；NA‐表示未给出特高值。