原生晕-因子分析地球化学异常识别<br/>——以内蒙古查干花地区1:5万水系沉积物Mo异常提取为例

原生晕-因子分析地球化学异常识别
——以内蒙古查干花地区1:5万水系沉积物Mo异常提取为例

2023-03-18赖联新魏俊浩徐崇文曾红红

地质与勘探 2023年6期

赖联新，魏俊浩，高强，李斌，徐崇文，曾红红

（中国地质大学（武汉）资源学院，湖北武汉 430074）

0 引言

复杂成矿地质背景对地球化学元素间隐藏变量与相互关系的识别和解释提出了挑战。降维被广泛用于揭示地球化学元素之间的相互依赖性及其与所研究问题的相关性（石文杰等，2011；2019；Zuo et al.，2021； Parsa et al.，2022）。目前主流的处理方法有主成分分析（PCA）、因子分析（FA）、独立成分分析、聚类分析、多维缩放以及基于知识驱动的专家根据成矿地质条件对每个元素赋予权重的元素组合（Wold et al.，1987； Hyvärinen and Oja，2000；Filzmoser et al.，2009； Kohonen，2012； Wang et al.，2017；张七道等，2021；冶建虎等，2022；徐龙生等，2022；黄啸坤等，2023）。

水系沉积物地球化学异常的识别是找矿技术方法中最直观的一种。通过对水系沉积物地球化学数据进行处理，可以分析元素之间的相互关系和空间分布，从而了解区域内成矿元素的富集规律，以指导预找矿工作。该方法在找矿应用中取得了显著的成果，表明水系沉积物中的次生晕和矿化原生晕在成因和空间上密切相关。次生晕指的是矿化原生晕在表生风化作用下发生活化、迁移和再富集，在水系沉积物中形成次生异常。因此，水系沉积物中的次生晕反映了地表上矿化原生晕（钻孔）的特征，而矿化原生晕则是水系沉积物地球化学异常形成的来源之一，它们两者具有成因上的直接联系。

但是水系沉积物次生晕异常与矿化原生晕并不是一种完全吻合的关系。由于不同元素由原生晕经过表生风化作用发生活化、迁移和再富集的过程中，不同元素组成的岩石矿物物理化学性质具有很大的差异，在不同的Ph-Eh 环境下的活动性表现并不相同，从而在表生环境中迁移存在差异，其在找矿的实际运用过程中具有一定的局限性（陈静生，1958；蒋敬业，2006；詹小飞等，2021）。

在水系沉积物地球化学异常的识别过程中，元素间的权重组合是一个重要的步骤。然而，传统方法通常基于整个图幅的地球化学数据进行权重分析，将背景值和异常值都视为影响元素间权重的因素。这种方法忽视了不同元素在表生环境下活动性的差异，从而稀释了与成矿密切相关的几个成矿元素之间的相关关系。同时，对于研究区的地球化学异常识别，该方法缺乏元素的针对性，忽略了研究区内成矿原生晕的指导作用。因此，在水系沉积物地球化学异常的识别过程中，有必要加强有利异常的考虑，以减少伪异常所带来的不确定性。

本研究利用内蒙古查干花地区的1∶5 万水系沉积物地球化学测量数据以及查干花矿床钻孔地球化学测量数据资料，通过研究岩石钻孔原生晕数据反向指导刻画水系沉积物地球化学元素之间相关关系，使元素之间的关联性更接近区域内已有成矿事件的元素关联性。首先通过中心对数（clr）变换对地球化学数据进行预处理，打开数据的闭合效应，减少数据间的伪相关性；利用R 型聚类分析确定水系沉积物地球化学元素中与Mo 最相关的元素；并进一步分析原生晕对应元素钻孔地球化学数据，计算与Mo 矿化相关的元素的权重得分，并将其与FA 的Mo 矿化的因子载荷替换，对各元素数据进行权重耦合，将结果作为PH-FA 得分，通过该方法以针对性统计学的方式强化元素之间隐藏关系。最后通过浓度-面积分形（C-A）确认Mo 的地球化学异常下限进行异常提取。通过对比传统因子分析法（FA）的异常提取结果，探讨原生晕-因子分析法（PH-FA）的有效性，指导矿床外围矿化和与矿床具有共同成矿期次矿化的找矿工作。

1 地质背景与化探数据

1.1.地质背景

狼山成矿带位于兴蒙造山带西段（图1），是中蒙边境巨型成矿带的重要组成部分。研究区位于该成矿带北部，近年来发现多个钼多金属矿床（点），成为了大型钼矿产地（刘玉堂等，2005；聂凤军等，2007；彭润民等，2007；刘翼飞，2013；王树庆等，2020）。

图1 内蒙古查干花地区大地构造位置（a）和地质图（b）（据刘翼飞， 2013改）Fig.1 Geotectonic location （a） and geological map （b） of the Chaganhua area in Inner Mongolia （modified from Liu， 2013）

研究区地层主要包括新元古代宝音图群、白垩系固阳组（图1）。岩浆岩沿北东-南西走向发育，记录了长期复杂的构造-岩浆活动历史。构造表现为北东向主断裂与北西向次级断裂。此外，该地区经历了多次岩浆活动，包括新元古代花岗岩、早石炭世黑云二长花岗岩和印支期复式岩体，展示了岩浆岩的多样性，时代跨度从新元古代到侏罗纪（刘翼飞，2013；张若然，2015；李光耀等，2021）。

研究区以Mo 矿化为主，伴生有Bi、W、Cu 等元素，矿化集中分布于东北以及西南区域，东南区域也有少量分布，矿床类型以浅层低温热液脉型和斑岩型矿床（点）为主。典型矿床主要有：查干花斑岩型钼矿床、查干德尔斯斑岩型钼矿床、查干敖包低温热液脉型铜钼矿床。该3个矿床集中分布于研究区的东北区域，从区域地质背景及矿化特征分析，其有着相同的成矿地质背景和找矿标志，其中成矿期岩体为晚二叠世花岗岩，SHRIMP 锆石U-Pb 测定的形成年龄为253 Ma（刘翼飞，2013），北东向和北西向断裂构造同样也与成矿具有重要的控制作用。

查干花斑岩型钼矿床是研究区内规模最大、成矿最为典型的钼矿床，其能很好代表该地区的钼成矿事件。矿区及外围主要有三种类型的花岗质侵入体（图1、图2），均为二叠纪侵入岩，早二叠世花岗岩位于矿区西部，岩性为黑云母二长花岗岩和花岗闪长岩。晚二叠世花岗岩为成矿期岩体，同时也是矿床的主要赋矿岩体。这些岩体均侵入于宝音图群内，宝音图群分布于矿区的东边，覆盖于晚二叠世成矿岩体上方。矿区构造总体较为简单，但对于成矿具有重要的控制作用，主要发育北东向和北西向的两组断裂构造，成矿岩体主要发育于两组断裂形成的“X”型交会部位（图1）。钼矿体基本未出露地表，埋深13.7～160.2 m，多呈近水平脉状和凸镜状产出，矿体厚度1.5～173 m，最厚可达356.73 m，Mo的平均品位0.12%，最高可达4.44%。

图2 查干花钼矿床5号勘探线剖面图（据刘翼飞， 2013修改）Fig.2 Profile of exploration line No.5 of the Chaganhua Mo deposit （modified from Liu， 2013）

1.2 样品数据来源

研究数据来自4个1∶5万图幅，共计6534个水系沉积物地球化学样本，平均采样密度为3.5件/km2，资料由内蒙古自治区地质矿产勘查院提供。分析元素包括Cu、Pb、Zn、W、Mo、Au、Sn、Ag、As、Sb、Bi、Hg、Ni共13种元素，通过原子吸收光谱（AAS）测定以下9种元素：Zn、As、Bi、Cu、Hg、Mo、Sb、Pb、W、Ni。其他元素Ag、Au和Sn通过电感耦合等离子体原子发射光谱法（AES）测定。样品采集在有利于冲积物汇集的河道变宽、河道内测、河道变缓等处。采样介质粒径为0.18～2.00 mm。样品制备、样品分析和质量控制均以《地球化学普查规范（1∶50000）》（DZ/T0011-91）、《岩矿分析试样制备规程》（DZ/T0130.3-94）相关规定执行。通过计算分析值和标准参考样品的推荐值的对数差（ΔlogC），用以衡量偏差的大小：

其中C计算所确定的值；是标准样本的平均值；Cs是标准样品的推荐值（Xie et al.，1997； Zhao et al.，2016）。

1.3 数据统计特征

为了突出水系沉积物测试样品的含量特征及成矿潜力，对数据进行统计处理（表1）。结果显示，W、Mo、As、Sb、Bi的变异系数均大于1，一定程度上表明它们具有局部富集的趋势，高值点可能更趋近于矿体的近端。对于水系沉积物中的元素含量来说，其应服从对数正态分布，即在Q-Q图中表现为一条直线；但由于成矿元素存在不同母体以及不同因子（如地质体背景、成矿作用）的叠加影响，往往造成成矿元素偏离对数正态分布。图3显示Mo、Bi、Hg、Zn、W、Sn元素呈现为曲线，表明变异程度较高，这些元素在本研究区具有一定的成矿潜力（Reimann et al.， 2002； Zuo et al.，2009； Zhao et al.，2016）。

表1 查干花地区水系沉积物地球化学描述性统计Table 1 Descriptive geochemical statistics of stream sediments in the Chaganhua area

图3 查干花地区水系沉积物地球化学Mo、Bi、Hg、Zn、W、Sn对数转换Q-Q图Fig.3 Geochemical Q-Q diagrams of Mo， Bi， Hg， Zn， W， Sn logarithmic conversion of stream sediments in the Chaganhua area

本次研究收集了区内钼矿化规模最大的查干花矿床的21 条钻孔资料，共计11476 件钻孔地球化学测试样品，每块样品分别测试了Mo、W、Bi地球化学元素含量。通过统计学方法对数据进行处理（表2），数据显示Mo、W、Bi 平均含量分别为0.042%、0.008%、3.901%，为水系沉积物对应元素含量的100～1000 倍。其中Mo、W 的变异系数大于1，表明在矿区内仍有局部富集的趋势；然而Bi 相对于水系沉积物Bi 的变异系数来说，变异系数更低，在一定程度上反映Bi 的集中富集趋势比Mo 弱，但相对于整个研究区来说矿床范围仍有较强的Bi 的富集趋势。

表2 查干花矿床钻孔地球化学数据描述性统计（n=11476）Table 2 Descriptive statistics of geochemical data from drilling holes in the Chaganhua deposit（n=11476）

2 方法

2.1 clr成分数据变换

成分数据是一种比例数据，只携带相对信息。地球化学数据是一种典型的成分数据，其组分间观察到的相关性中至少有一部分是由闭合效应带来的元素间的伪相关（Aitchison，1982；周蒂，1998；Yin et al.，2021）。因此在数据处理前需要进行数据变换，打开数据闭合效应（Zhao et al.，2016）。此次研究采用clr变换进行数据的预处理。

clr变换的基本思想是将数据中每个样品的各组成部分之间的相对含量关系转换为对数比值，再通过对每个样品的对数比值取中心化对数变换，得到一个新的数据矩阵。这个新的矩阵中的每个元素表示一个样品的某个组成部分与其他所有组成部分相对含量的比值的对数值，这样可以消除数据的闭合效应，使得数据更加符合多元统计分析的假设。

设X=(x1，x2，…，xD)'，则有：

其中X为某一观测成分数列，x为成分分量，D为成分分量的个数，为D个成分分量的几何均值。

2.2 原生晕-因子分析

矿体中不同元素的富集程度代表了与成矿期成矿元素的紧密关系，同时也是水系沉积物次生晕的直接影响因素。如何将水系沉积物次生晕与矿体中的元素富集程度联系起来进行地球化学异常识别是一个值得研究的课题。

因子分析（FA）是地球化学数据聚类和降维应用最为广泛的方法，是一种重要的研究地球化学数据间相关性和组合特征的方法。该方法可以降低数据的维数并检测地球化学元素间的隐藏结构，并解释这些结构对变量之间的相关性的影响，从而把错综复杂关系的变量归结为少数几个综合因子。（Filzmoser et al.，2009；李申等，2018；Zhao et al.，2022）。

原生晕-因子分析（PH-FA）是基于FA 改进的分析聚类方法。采用了R型聚类分析来确定水系沉积物中与Mo 矿化最相关的地球化学元素。随后，进一步分析这些元素与原生晕相对应的钻孔元素的关联程度（公式3），并通过关联程度计算与矿化相关的元素的权重得分。接着，我们将这些权重得分用于替代因子分析（FA）中Mo 矿化的因子载荷，以实现对各元素数据的权重耦合并以该值作为PH-FA 得分（公式4）。这样的处理方式旨在通过针对性统计学的方法增强元素之间的潜在关联性。

原生晕地球化学元素与水系沉积物对应地球化学元素的关联程度得分计算表示为：设X=(x1，x2，…，xD)' ∈SD

X为某一水系沉积物地球化学观测成分数列，x为成分分量，D是成分分量的总个数，xi表示这组数据第i个点的数据值。表示j元素原生晕地球化学数据的平均值，μ是j元素水系沉积物地球化学数据的平均值。

对σ(X j)进行归一化处理以获得权重系数，PH - FA因子得分函数表示如下：

2.3 浓度-面积分形（C-A）

C-A 分形技术被广泛用于揭示地球化学异常（Cheng et al.， 1994； Zhao et al.， 2022）。该方法可以较好地将背景值和异常值区分，通常将数据分成三个（背景，低异常和高异常）或四个（背景，高背景、低异常和高异常）级别。C-A模型表示为：

其中A( ≥C)表示浓度大于或等于C的面积，β是分形维数。分形维数β可以从A(C)与C的对数-对数图的斜率中获得。C的值对应于拟合对数图中线条的交点，并用作将地球化学图划分为不同水平区域的阈值。

3 结果和讨论

对研究区内1：5万水系沉积物地球化学数据分别进行传统FA 和PH-FA 分析，对比这两种方法的结果，判断PH-FA分析的优劣性。首先对6534个原始数据进行clr变换减少数据间的伪相关，选择Ni作为clr变换的除数，得出的结果数据进行以下处理。

3.1 传统因子分析（FA）

采用SPSS 软件对基于clr 变换的水系沉积物地球化学数据进行因子分析，通过使用Z-Score 标准差法对数据进行标准化处理，消除不同元素之间的量纲，对标准化数据通过主成分分析法采用特征值大于1 的标准提取公因子，使用Kaiser 标准化最大方差法对因子载荷矩阵作正交旋转，最终计算得出5个主因子（表3）。

表3 旋转后成分得分系数矩阵Table 3 Component score coefficient matrix after rotation

其中F1为Cu、Zn、Ni元素组合，代表研究区与Cu-Ni矿化相关的高温成矿元素组合；F2为As、Sb元素组合，代表研究区低温成矿流体元素组合；F3为Mo、Bi、W、Hg元素组合，是研究区最重要的成矿元素组合，指示查干花、查干德尔斯、查干敖包钼矿床的Mo-Bi矿化；F4为Pb、Sn、Ag元素组合，代表该地区与Ag-Pb矿化相关的中温成矿元素组合；F5为Au、Hg元素组合，代表了研究区Au矿化，且Au相对于W其得分系数相差巨大，表明Au元素地球化学行为具有一定的独立性。

3.2 原生晕-因子分析（PH-FA）

为了针对性识别矿区外围与成矿密切相关元素水系沉积物地球化学异常，对水系沉积物数据进行R 型聚类分析，组间采用平方欧式（Euclidean）距离，结果如图4 所示。从谱系图中可以看出（图4），全区13种元素可分为8类，5个双元素组合和3个单元素组合。结合实际地质情况、成矿特征研究，研究区与Mo 矿化元素最密切相关的水系沉积物地球化学元素有Bi，Mo、Bi显示出强烈的相关性，在原生晕特征上表现为相互伴生元素。

图4 水系沉积物地球化学数据R型聚类谱系图Fig.4 R-type clustering spectrum of geochemical data of stream sediments

对于原生晕地球化学特征，对查干花钻孔的Mo和Bi含量按照深度进行做图，探讨这两个元素的相关性。将0.1 作为Mo 元素的异常下限，将CZK703、CZK406、CZK1306 这3 条钻孔数据归一化后进行投图（图5）。结果显示其具有以下三个特征：（1）钻孔原生晕地球化学异常和水系沉积物地球化学异常在Mo、Bi 异常特征上表现出明显的相似性，表明原生晕异常和次生晕水系沉积物异常具有一定的吻合性，两元素在表生环境下有着相似的活动性，具有相似的迁移距离；（2）Mo、Bi含量具有明显的正相关性，在Mo 高值异常区内同样有Bi 的高值异常的异常区比例约占95.28%；（3）Mo、Bi异常呈现不相关的个别区域中主要分布在深度50 m 范围内，这可能和后期浅部的流体活动导致元素的迁移有关。

图5 查干花矿床钻孔归一化地球化学数据Mo、Bi异常折线图（a：CZK703； b：CZK406； c：CZK1306）Fig.5 Mo and Bi anomalies lines of normalized geochemical data from drilling holes in the Chaganhua deposit（a：CZK703； b：CZK406； c：CZK1306）

Mo 的钻孔地球化学数据高值异常与查干花、查干德尔斯矿床形成的成矿作用有关，在钻孔中其代表Mo 元素富集成矿形成的矿脉，其是矿床周边形成水系沉积物地球化学数据Mo 及其伴生元素异常的直接来源。将Mo 的高值异常钻孔地球化学数据提取出来，共计1885 块钻孔化学样数据，对其中Mo、Bi 元素进行PH-FA 的富集程度的相关性分析，结果显示Mo 的PH-FA 权重系数为0.717607，Bi 的PH-FA 权重系数为0.282393；将两元素进行Z-Score 标准化，把两元素数据变成均值为0，标准差为1 的标准化数据，并通过系数进行权重组合。

3.3 异常识别与对比研究

在水系沉积物地球化学异常图的制作过程中，异常提取结果与异常下限的取值密切相关。我们首先利用浓度-面积分形（C-A）方法来区分异常和背景，然后结合传统因子分析（FA）和原生晕-因子分析（PH-FA）进行异常识别，得分情况如图6所示。结果显示，FA 和PH-FA 两种方法均显示出四种数据类别，即背景区、高背景区、低异常区和异常区。这些结果不仅有效地利用了阈值来区分背景和异常，而且显示改进方法所确定的阈值更为合理。在图6a 中，阈值类别相对较小，导致背景区域的双对数值分散且数量较少，从而使背景区域相对较小。而在图6b 中，阈值相对较大，使背景区域的双对数值更为集中且数量相对更多，进一步细分背景和异常，有效地增强了地球化学异常，并增加了背景与异常之间的对比度。

图6 F3 C-A双对数图（a）和PH-FA浓度-面积双对数图（b）Fig.6 F3 C-A log-log graph （a） and PH-FA C-A log-log graph （b）

利用F3（Mo、Bi、W、Hg）因子得分和改进的原生晕-因子分析Mo-Bi 得分两种方法绘制地球化学等值线图（图7）以及提取异常（图8），进一步揭示了Mo 成矿元素系列在局部的成矿特征。两图中高值异常主要分布在研究区的右上区域，集中于查干花矿床、查干德尔斯矿床和查干楚鲁矿床周围，与已知矿床的位置基本套合。根据异常提取结果（图8），F3 因子得分异常图和改进的原生晕-因子分析Mo-Bi 得分异常图的异常与图8 中的4 个Mo 矿床（点）的吻合度均达到75%。这表明两种方法在揭示研究区的Mo 矿床位置上具有较高的准确性和可靠性。

图7 查干花地区因子分析F3因子得分图（a）和原生晕-因子分析Mo-Bi得分图（b）Fig.7 Factor scores of F3 （a） and Mo-Bi scores of primary halo-factor analysis （b） in the Chaganhua area

图8 查干花地区因子分析F3因子得分异常图（a）和原生晕-因子分析Mo-Bi得分异常图（b）Fig.8 Factor analysis of abnormal F3 factor score map （a） and primary halo-factor analysis of abnormal Mo-Bi score map（b） in the Chaganhua area

两种方法所圈定的异常占研究区的比例相似，但两种方法在异常的空间分布的局部细节上有明显的不同。（1）图7b 相对于7a 中高值异常在整体上有北东向、北北西向线状分布的趋势，这与研究区内整体的构造线方向表现一致，也与查干花矿床控矿构造方向表现一致。（2）研究区共有46条断层，F3因子得分图显示的断层与所圈定的异常相交重合率为63.04%，而原生晕-因子分析Mo-Bi 得分图断层与异常相交率达76.09%，显示了PH-FA 相较于FA 有着一定的优势。（3）晚二叠世花岗岩是研究区重要的赋矿岩体，图8中a、b两图高异常区与岩体都有较好的套合性，但是在异常的延伸方面，图8a 表现出往西部延伸，图8b 则相反，具有往东延伸的趋势。由图1 可以看出，图8a 的西部延伸的高值异常对应的岩体是早二叠世黑云母二长花岗岩和花岗闪长岩，两岩体在形成时代上比成矿年龄（晚二叠世）更早，且从查干花矿床的5 号勘探线剖面图（图2）看出，早二叠世黑云母二长花岗岩矿化极少，这些都显示西部延伸的高值异常并不具备Mo 矿床的形成条件。而图8b 的东部延伸高值异常对应宝音图群地层（图1、2），矿床东部的宝音图群地层之下被侵入大量的赋矿岩体，并有非常好的Mo 矿化表现，由此可以推测异常区东部的宝音图群地层之下甚至地层内有着很好的成矿潜力，这与图8b所圈出的异常表现一致。（4）PH-FA 在识别非（弱）异常也有优异的表现，晚元古代花岗岩是研究区内最老的的一套花岗岩，其对于Mo 成矿几乎没有贡献，因此该岩体不存在Mo的高值异常，在图8中，b图能明显看出晚元古代花岗岩非异常的浓度表现，然而在图8a中并不能观察到这种相似的规律。

总体来说，PH-FA与FA相比在识别异常研究区内已知矿化类型相似的异常方面表现出一定的优势，在本区可以有效指示Mo地球化学异常大致分布范围。