基于Kriging插值的综合放射性勘探方法在鹿井铀矿找矿中的应用

罗齐彬, 杨亚新, 吴信民, 肖 昆, 付 宸, 吴永鹏, 胡旭东

(1.东华理工大学 核资源与环境国家重点实验室，江西 南昌 330013；2.东华理工大学 地球物理与测控技术学院，江西 南昌 330013)

放射性勘探数据不同于电法、磁法、重力和地震等勘探方法测量的数据，其最大特点就是具有明显的统计涨落。由于地形、采样点环境等因素影响，在进行野外放射性勘探测量时，往往会导致采样点间距疏密不一致，甚至缺失一些测量数据。地面能谱测量、土壤氡测量以及地气测量等勘探方法所测得数据存在采样不均匀、样本规模小等不足(罗齐彬等，2018；罗齐彬，2019；杨亚新等，2018；吴泽民等，2019)。Kriging插值法是一种以变异函数理论和结构分析为基础，对局部区域内数据样本进行无偏最优估计的方法。该插值方法充分利用了数据点之间的空间关联性以及数据空间场性质，适合规模小且不均匀的样品测量或勘探数据，被广泛应用于大气、土壤、地质以及放射性等领域的数据处理(李如仁等，2020；张靖，2014；李华等，2019)。根据鹿井测区内多种放射性勘探数据及其特征，利用Kriging插值法对测区内勘探数据进行插值研究，选取出合适的Kriging插值模型和插值参数，最终得到准确、有效的放射性勘探等值线图，提高综合放射性异常解释科学性和准确度。

1 研究区及测线简介

鹿井铀矿田地处湖南省汝城县与江西省崇义县接壤处，大地构造上位于万洋-诸广山隆起带、常德-安仁断裂、遂川-热水断裂以及仙鹅塘断裂交汇部位(诸广山复式岩体中段)，受遂川断裂及热水断裂组成的地堑式断陷带控制，是南岭EW向铀成矿带的组成部分。其中NE向遂川-热水断裂属走滑剪切带，是主要的控岩、控盆和(铀矿)控矿构造。鹿井地区内铀矿床、矿化点主要分布于丰州盆地边缘，矿床类型为花岗岩型和花岗岩外带型(黄宏业等，2008；周肖华等，2014；孟凡兴等，2016)。鹿井地区控矿构造以NE向断裂构造为主，设计测线共7条，呈NW走向(132°)，线长7 500 m，以1 000 m线距均匀分布，每条测线上测点点距均为100 m，如图1所示。

2 Kriging插值原理

(1)

式中，εi为已知数据点的权重系数，若要使估计值为无偏估计，则需要所有权重系数之和为1。其次，无偏估计需要通过选择合适的权重系数εi使估计偏差期望和估计方差应满足：

(2)

假设区域数据变量Y(xi)满足二阶平稳和本征假设条件，其数学期望为μ，协方差函数C(h)与变异函数γ(h)也仅为数据点之间距离h相关的函数，即有：

E[Y(x0)]=μ

C(h)=E[Y(xi)Y(xi+h)]-μ2

γ(h)=12Var[Y(xi)-Y(xi+h)]

(3)

在二阶平稳条件下，可进一步将协方差函数和变异函数离散化为：

C(h)=1n∑ni=1Y(xi)Y(xi+h)-μ2

γ(h)=12n∑ni=1[Y(xi)-Y(xi+h)]2

(4)

式(2)中的估计方差最优条件则有：

σ2=∑ni=1∑nj=1εiεjC(xi,xj)-2∑ni=1εiC(xi,x0)+C(x0,x0)

(5)

在无偏性估计条件下可构造出有关权重系数的估计方差函数：

f(εi,μ)=σ2+2ξ∑ni=1εi-1

(6)

将式(5)带入式(6)，根据拉格朗日乘数定理，估计方差达到最小值，即为f(εi,μ)对权重系数εi和拉格朗日参数ξ的偏导数均为0。

f(εi,μ)对拉格朗日参数的偏导数为0，即所有权重系数之和为1，f(εi,μ)对权重系数εi的偏导数为：

∂f(εi,ξ)∂εi=∑nj=1εjC(xi,xj)-C(xi,x0)+ξ=0

(7)

由式(7)与条件“所有权重系数εi之和为1”即可解出参数εi和ξ的值，将其带入式(5)中可计算出最小估计方差值：

σ2min=C(x0,x0)-∑ni=1εiC(xi,x0)+ξ

(8)

由此可见，利用Kriging插值法的前提是观测数据的均值和协方差矩阵是已知的。求解出的权重系数值质量取决于所选择的变异函数模型，常见的变异函数模型有指数模型、高斯模型、线性模型、球状模型、幂函数模型等。

3 结果与讨论

(1)测区内放射性勘探数据特征。在本测区内每种放射性勘探数据各获得约530个样本数量，样本量并不大，对这些放射性勘探数据进行统计分析(表1)。结果表明，地面γ能谱铀含量、土壤热释光强度以及210Po活度数据离散程度相对较大、数据分布跨度较长。

表1 鹿井测区放射性勘探数据主要统计参数

(2)Kriging插值变异函数模型的选取。在选择Kriging网格化插值时，可通过拟合出合理的单一或者联合变异函数模型作为放射性勘探区域数据的变异函数模型。笔者分别对鹿井测区内的地面γ放射性元素U/Th值、210Po活度以及土壤热释光勘探数据进行变异函数实验拟合，其中延迟方向为132°(沿测线方向)，方向容差为30°，结果如图2所示。

由图2可以看出，鹿井测区内这几种放射性勘探数据的插值变异函数模型主要为块金效应-线性联合模型，模型函数为：

γ(h)=0h=0

C0+s·hh>0

(9)

式(9)中，C0为块金值，s为斜率，其余设置参数还有各向异性比率(简称异性比，δ)以及各向异性角度(简称异性角，θ)，各个模型函数参数见表2。

(3)Kriging插值法与其他插值法在土壤热释光数据处理中的对比分析。为研究Kriging插值法在放射性勘探数据处理中的插值效果，本研究利用土壤天然热释光强度数据，分别采用Kriging插值法、自然邻点法和最小曲率法等3种插值方法进行处理，并将其与原始数据分布相比较(图3)。

表2 鹿井测区放射性勘探数据插值变异函数参数

由图3可知，在原始土壤热释光强度数据中，由于采样点不均匀性以及部分测点数据的缺失，造成异常分布不连续，这与理论上矿体机械晕、分散晕机制等已有认知是不同的。自然邻点法和最小曲率法对其进行插值后一定程度上形成了“晕状”异常，但仍主要在原始数据单点高值异常位置形成小面积异常。同时在异常分布形态方面，仍然受测线控制，沿测线方向呈“串珠状”异常，这是不利于找矿异常解释的。Kriging插值法则能够克服以上问题，消除数据在空间上分布的不均匀性，同时将异常分布形态信息更好地还原出来，有利于找矿异常解释和成矿范围的圈定。

(4)鹿井测区内地面γ放射性元素异常。鹿井测区内地面γ能谱测量数据中放射性元素U/Th值异常呈圆形集中分布在鹿井矿床附近，且在测区内西南角桥子坑附近伴有中等水平异常；放射性元素Th异常呈片状分布，以丰州盆地边界为限形成内低外高的分布形式(图4，5)。

根据前人研究资料以及已揭露的铀矿床情况来看，鹿井测区铀矿床主要分布在丰州盆地外围梨花开、老虎龙、洞房子、下洞子、高昔以及桥子坑等附近构造比较发育的地区，构造不发育的丰州盆地内亦有部分小型铀矿床揭露(张万良等，2011，2018；周肖华等，2014)。鹿井矿床附近的高U/Th值异常表明该异常区内后期发生过铀的沉淀富集活动，且活动层位较浅。丰州盆地外围高能谱Th异常主要呈东高西低、NE向条状展布的趋势，表明后期铀的活化迁移主要受NE向主断裂控制，从桥子坑一带沿断裂往WS方向进行迁移。

(5)鹿井测区内210Po活度异常。鹿井测区内210Po活度异常平面分布如图6所示。区内3条NE向主断裂从测区内南缘贯穿而过，为氡气迁移提供了良好的通道，从而使得210Po活度高值异常主要集中在桥子坑-高昔南缘一带，并逐渐向界坑-集溪一带展布。另外，由于受到鹿井矿床影响，该附近地区呈现中高值210Po活度异常。

(6)鹿井测区内土壤热释光异常。鹿井测区内土壤热释光异常平面分布如图7所示。土壤热释光高值以上水平的异常主要分布在集溪-界坑一带和金鸡岭-桥子坑-高昔一带，两个异常区呈NW向展布且分界明显。

总的来说，鹿井矿床附近土壤热释光异常属于中低值水平，表明其成矿时间距今相对较短。结合桥子坑-高昔南缘一带出现中低值能谱U/Th值异常、高值210Po活度以及高值土壤热释光异常，表明该区域内可能存在中浅部埋藏的铀矿体。而中低值能谱U/Th值异常、中高值210Po活度以及中高值土壤热释光异常则表明金鸡岭一带和界坑一带很有可能存在中深部埋藏的铀矿体。从丰州盆地内放射性勘探结果异常情况来看，盆地内构造不发育，几种放射性测量结果都属于低值分布，这与矿体、构造控制放射性核素运移的机理是相符合的，也从侧面验证了Kriging插值法在放射性勘探数据处理中应用的适用性和正确性。

4 结论

本研究基于鹿井测区内地面γ能谱放射性元素、210Po活度以及土壤热释光三种勘探数据，针对这些放射性勘探数据特性选用Kriging插值法进行处理，通过研究得到以下几点认识：

(1)Kriging插值法能有效消除采样不均匀性带来的数据误差，使勘探数据成图结果更加接近实际情况，有利于放射性异常解释。

(2)一般来说，大部分放射性勘探数据具有几何各向异性，但其变异函数模型多服从指数函数、线性函数以及块金效应函数中的单一或者联合形式，变异函数模型及参数需要根据测量数据进行实验拟合确定。

(3)受已揭露矿床的影响，鹿井矿床附近具有集中的高值地面γ能谱U/Th值异常；地面γ能谱U/Th值、土壤热释光以及210Po活度在界坑附近均具有明显的中高值异常。鹿井老虎龙附近放射性综合信息异常受揭露的铀矿床影响，丰州盆地内靠近测区边缘一带有小型铀矿床揭露且盆地内构造并不发育。综合考虑，鹿井测区内隐伏铀矿有利区主要在界坑南缘一带、下洞子北部一带以及金鸡岭西侧一带。