基于成分数据分析及S—A分形的地电提取法在川西容须卡地区寻找隐伏锂矿预测研究

石佳磊, 欧阳菲, 罗先熔, 文美兰, 刘攀峰,刘永胜, 梁 鸣, 杨青松, 蔡宇群, 孙 瑶

(1.a桂林理工大学 b地球科学 隐伏矿床研究所,广西 桂林 541004;2.甘肃省有色金属地质勘查局天水矿产勘查院,甘肃 天水 741024)

0 引言

地电提取法是穿透性勘查技术方法之一,具有一定的探测深度,能够有效地寻找隐伏金属矿床。该方法于20世纪80年代由前苏联引进,近50年来在国内得到长足的发展和应用,在不同矿区和不同类型的有色金属矿床均进行已知剖面的找矿试验和未知区找矿预测,取得较多的找矿研究成果[1-10]。但在已有公开发表的成果中,尚未见寻找锂矿的案例。川西容须卡地区是高寒山区景观,野外条件恶劣,导致许多勘查方法的实施受限[11-14]。本文选取川西容须卡地区作为研究区,采用地电提取法寻找隐状锂矿,以期为高寒山区景观下的勘查地球化学方法的应用提供重要参考。

川西分布着大量的伟晶岩型稀有金属矿床,包括甲基卡、可尔因、丹巴白等典型矿床,均与松潘—甘孜造山带有关。其中甲基卡是我国规模最大、品位最高的花岗伟晶岩型锂矿床。在甲基卡矿集区,发育有甲基卡、容须卡和长征3个构造岩浆热穹隆,穹隆内部及周边发育大量的伟晶岩脉。而容须卡的规模仅次于甲基卡,因此该区具有较好的锂矿找矿前景[15-16]。

研究发现地电提取测量结果呈现明显的“闭合效应”,这在一定程度影响异常评价的准确性;因此,需要选用合适的数据处理方法来修正对异常的认识。地电化学数据属于成分数据,也具有“闭合效应”[17]。其数据的空间分布类型归属于艾奇森空间,即成分数据中某一组分的含量发生变化将会影响剩余组分的含量,进而会使原本处于正相关的组分其相关性发生转变。为解决此问题,需要对数据进行对数比变换,从而将数据分布特征从艾奇森空间转换为欧式空间,以便应用统计学方法对数据进行分析[18-20]。经过变换之后的数据更能够体现数据在空间的真实分布模式,进而能准确地识别出元素的背景及异常信息[18-22]。

同时,随着当代计算科学的飞速发展,地球化学数据处理方法已从传统的单一变量发展至多变量分析,从线性分析发展为非线性分析。传统的平均值加标准差求取异常的方法已不足以揭示数据中的非线性特征。受不同成因、多阶段及多期次的地质活动影响,地球化学场往往具有复杂的非线性特征。常规异常圈定方法无法有效分辨由不同地质体导致的地球化学背景差异,往往会丢失或忽视局部的弱异常,从而影响最终的异常划分和靶区圈定[23-24]。为有效识别和分离地质背景和异常,前人提出用多重分形的方法来处理地球化学数据,从而可以更准确地认识元素的地球化学分布特征[25-27],达到有针对性地提取矿致异常信息的目的。

因此,本文首先对研究区的地电化学数据进行成分数据分析,通过与原始数据的特征对比选择正确的数据变换方法,打开数据的“闭合性”,从而发现更明确的地球化学指标。同时,考虑到研究区经历不同期次的地质演化与地质作用的影响,其地球化学分布较复杂,采用S—A异常分解法分离异常与背景,以期能够更准确认识异常特征,识别和筛选出低背景区中有意义的弱异常。

1 地质概况

1.1 区域地质概况

松潘—甘孜造山带位于我国西南部青藏高原东侧,北部与昆仑—秦岭造山带、阿尼玛卿—勉略缝合带相连,西与金沙江缝合带、甘孜—理塘缝合带相连,东与龙门山断裂带相连(图1a)。在印支末期羌塘—昌都陆块与扬子陆块的长期碰撞汇聚及喜马拉雅造山运动时期欧亚板块的碰撞等一系列构造运动叠加的共同作用下,松潘—甘孜造山带整体呈倒三角形,几何形态特殊,构造复杂,蕴含丰富的矿产资源[28-29]。

三叠系为该区分布范围最广的地层,主要岩性为板岩、变质砂岩及千枚岩等。区域构造主要以褶皱穹隆及NW向断裂为主,穹窿构造主要以甲基卡穹窿和容须卡穹窿为典型代表,容须卡穹窿周围发育呈环带状分布的变质带:十字石带→红柱石带→石榴子石带→黑云母带→绢云母→绿泥石带。断裂主要是受印支燕山期NE—SW构造挤压,喜马拉雅造山期NW—SE向剪切作用而呈NW展布(图1b)。岩浆岩则是以印支—燕山早期的中酸性侵入岩为主,无明显分布特征面状分布于全区,与研究区内三叠系地层侵入接触。脉岩沿岩体内外接触带产出,以花岗伟晶岩脉为主,为锂、铍、铯等稀有金属矿体提供物质来源。

图1 容须卡地区大地构造图[26] (a)和区域地质简图[27] (b)Fig.1 Tectonic location map[26] (a) and sketch regional geological map of Rongxuka area[27] (b)1—三叠系新都桥组二段 2—三叠系新都桥组一段 3—三叠系侏倭组 4—三叠系中上统如年组 5—逆断层6—平移断层 7—性质不明断层 8—勘查区 9—穹隆构造

1.2 研究区地质概况

研究区位于容须卡穹隆的SW部。区内出现3条变质带呈NE向分布,分别为石榴石带、红柱石带和堇青石带。地层以上三叠统新都桥组一段为主,岩性主要为深灰色(十字石)红柱石黑云母片岩,其次为石榴石十字石绢云片岩,少许灰、深灰色透镜状变质石英砂岩,局部出露花岗伟晶岩脉和石英岩脉(图2);其次为侏倭组,侏倭组为一套薄-厚层变质石英砂岩、长石石英砂岩、粉砂岩与绢云板岩组成的韵律式互层。

图2 容须卡地区地质简图(据文献[16]修编)Fig.2 Geological sketch map of Rongxuka area [16]1—三叠系新都桥组二段 2—三叠系新都桥组一段 3—三叠系侏倭组 4—红柱石带 5—石榴子石带 6—堇青石带 7—十字石带 8—花岗闪长岩 9—石英脉 10—花岗伟晶岩 11—地电提取测网 12—研究区

2 方法和技术

2.1 地电化学技术及测试分析

2.1.1 地电提取法原理

地电提取法方法通过外加人工电场作用,使地底岩石矿物中离子的动态平衡被破坏,从而使离子向着离子吸收器迁移[3,30]。经过一定时间之后,地底岩石矿物岩石中的离子会重新维持新的动态平衡,同时离子接收器中的离子也会达到饱和状态。通过对离子吸收器进行测试分析指定元素及含量,从而发现深部隐伏矿体有关的元素异常,以此达到找矿目的[31-32]。

2.1.2 测线、测点布设及使用装置

本次地电提取工作设计14条NE向测线, 网度为100 m×20 m,共取样390件。野外采样使用桂林理工大学隐伏矿床预测研究所研制的独立供电偶极地电提取装置,野外工作采用9 V电池供电,提取时间为48 h,采样坑30～40 cm深,电极间距为1 m,提取液为500 ml、15%浓度的稀硝酸;完成野外提取后,取回电极中的泡塑,将其作为样品送往分析测试单位。

2.1.3 样品分析测试

分析测试工作在有色金属桂林矿产地质测试中心进行,选择硝酸-过氧化氢消解-王水提取法作为样品的前处理方法[33],采用X-series电感耦合等离子体质谱仪分析Li、Be、Rb、Cs、Th、U等6种元素的含量。选取合适浓度的标准溶液,使得ICP-MS测定后自动绘制的各元素工作曲线的线性相关系数均大于0.999 8,保证各元素回收率均大于96%,RSD均小于5%。

2.2 成分数据分析

元素之间的关系是研究地球化学数据空间分布特征的基础和前提。通过成分数据分析不仅可以揭示研究区不同元素之间的内在联系,还可以进一步判断元素的空间分布受哪些因素影响。通过讨论判断不同地质过程或地质体对元素空间分布的影响,探讨元素异常与矿体、成矿作用及控矿因素之间的关系,从而达到地质找矿的目的。

地球化学数据作为典型的成分数据,具有明显的“闭合效应”[34-38]。一般而言,地球化学样品中不同元素含量的总和为定值,使得地球化学元素之间具有伪相关现象,数据分析的结果会变得不确定。为消除该效应所引起的负面效果,通常需要相应的方法对参与实验的数据进行转换[35-37]。经历变换后的数据相比而言,通常比未经变换的数据更为具有可解释性[38-39]。

在进行对数比变换的方法中较为常用有加性对数比变换(alr)、中心对数比变换(clr)、等距对数比变换(ilr)3种。其中,前两者无法消除成分数据的闭合效应,只有ilr变换可以。原始成分数据经alr变换和ilr变换之后,减少一个变量;clr变换之后的数据与原始数据的变量数相同,但和为“0”。根据3种变换方法的特点,Filzmoser[36]提出先利用ilr变换打开数据的“闭合效应”,然后利用标准正交基将ilr变换的载荷因子在clr空间中表达出来。这样即克服地球化学数据的“闭合效应”,又保证变换前后变量数目的一致性。通过此变换,能够更加有利于地球化学家们做出合理的解释。以下为各变换方式的理论公式:

以一个原始数据矩阵Xij(m×n,其中m为样品数,n为元素个数)为例,分别介绍2种变换。

(1)clr变换:

(式1)

j=1,2,…,n;g(xi)为第i个样品的几何平均值。

(2)ilr变换:

ilr变换是在clr变换和标准正交基vj的基础上得来的,具体公式如下:

(式2)

j=1,2,…,n-1。

(式3)

j=1,2,…,n-1。

因此clr变换和ilr变换可以通过标准正交基联系起来。

clr(X)=V×ilr(X),ilr(X)=clr(X)×V

(式4)

V=[v1,v2,…,vn-1]

(式5)

2.3 多重分形(S—A)

元素的地球化学(背景—异常)分布模式并非受同一地质过程控制,但勘查地球化学数据的结果却是元素背景和异常复合叠加的结果。在漫长的地质年代中,不同期次、不同种类的地质过程会对元素的空间含量分布造成影响,从某种意义上说明地球化学场并非一种稳定场[40]。在这种影响下,导致地球化学场变得更加复杂,从而加大对元素地球化学空间分布模式的解译难度。能否有效识别并分离出受到不同地质过程影响的背景和异常是找矿勘查工作的关键。随着Cheng等[41-43]用分形方法对勘查地球化学开展研究后发现,绝大多数的地球化学异常都存在分形特征。以元素含量—面积分形模型(C—A)和能谱密度—面积分形模型(S—A)的分形方法开展的地球化学背景—异常分离研究及找矿预测工作已得到广泛应用[44-46]。

元素的异常及背景在频率域中具有广义自相似性。通过将元素地球化学数据从空间域转换至频率域中,并对具有自相似性的背景及异常进行能谱分解,随后再将分解后的背景—异常从频率域重新转换回空间域,实现背景与异常分离的目的。目前S—A法在矿体预测中广泛应用,并取得较好的效果,已发展成为标准的异常分解方法之一。这种自相似性一般表达方程如下:

A(≥S)∝S-β

(式6)

式中:S代表能谱能量密度作为波数向量的函数;A(≥S)表示大于某一阈值(S1)的单位面积;β是各向异性的多维分形指数,即分维数,在双对数图上可以获得,∝表示为成正比。随着S值的增大,相应的A会减小,其变化规律主要取决于β分维数。

为识别研究区元素的背景和异常,本文采用ilr变换方法。首先,对比元素在PCA上的表现,选取合适的元素组合;然后,采用S—A模型提取研究区主成矿元素的异常和背景。其中,S—A 分形模型采用中国地质大学(武汉)GPMR国家重点实验室左仁广团队研发的Arcfractal插件[47]在R语言编程环境中完成数据处理过程。

3 结果分析

3.1 数据统计特征

对地电提取样品中6种金属元素的含量进行分析,基本统计学参数见表1。从元素含量来看,Li、Rb的含量最高,约在10-6的数量级;其次为Cs、Th,约在10-7的数量级;最低的为Be、U,约在10-8的数量级。变异系数是反应元素的相对离散程度的重要指标,较高的变异系数值往往反映该元素在矿区具有较好的成矿潜力或经历较复杂的地质过程[48]。从表1中可以看出,显示除U的变异系数略小以外,其他5种元素都超过0.5,均显示出较高的分散性,故使用这些元素来区分不同地质过程或地质体。箱状图能够直观地展示数据的离散特征,从图3a中可以看出箱体上半部分偏大,说明数据呈现正偏特点;出现大量的离群值和极值,表明数据的离散程度大。

表1 地电提取元素含量参数统计

通过相关转换使数据呈现正态分布,是计算地球化学异常和背景和合理划分异常内、中、外带的前提。通过原始数据、对数变换、ilr变换等3种方法的对比(图3)发现,经过ilr变换的数据(图3c)呈现出最接近理想正态分布的特征。原始数据(图3a)的箱体形态最不符合正态分布,呈现明显的正偏特征,且出现大量离群数据;对数变换(图3b)的箱体接近上下对称,但仍然有少数元素出现偏态,且大量离群数据出现在箱体上方;而经过ilr变换,其箱体形态基本呈上下对称的特点,且离群数据在箱体的上下方都有。

图3 研究区3种变换方法的箱图对比

3.2 成分数据分析

3.2.1 双标图分析

双标图的优势是可以直观地观察到各元素对第一、第二主因子的贡献及不同元素之间的差异程度。为考察ilr变换的优势,分别对原始数据、对数变换、ilr变换等3种处理结果的数据进行对比,经PCA分析后得到PC1/PC2的主成分双标图(图4)。从图中可以看到原始数据、对数变换的双标图(图4a、4b)明显呈“一边倒”的趋势,即数据仍然存在闭合效应[49];而经过ilr变换之后的数据(图4c),原本闭合的变量被打散开,整体呈现放射状,降低了闭合效应(图4c)。

可以看出,原始数据和对数变换的双标图(图4a、图4b)分布相似: Li、Be、Rb、Cs、U、Th均呈现为PC1的正载荷。各元素对PC1、PC2主因子的贡献并没有得到明显的区分,导致大多数元素都聚拢在一起,呈现出典型的闭合效应特征。由于闭合效应导致变量在开展多元统计时产生伪相关,致使产生错误的元素组合解译结果。

而从ilr变换之后的双标图(图4c)来看,元素的矢量指向也更加发散,相较前2种数据双标图,更加平均地分布在双标图中;同时样本的分布比前2种方法更加分散。这表明经过ilr变换之后,已消除数据之间的闭合效应。Li、Rb、Cs、U、Th在PC1的正载荷端;Be元素则单独在第三象限出现,与其余元素之间的相关性较为疏远;其余元素(Li -Rb-Cs-U-Th)之间相关性良好。

根据ilr变换的结果:PC1因子的正载荷端出现Li、Rb、Cs、U、Th等元素,而负载荷端出现Be元素;PC2因子的正载荷端出现Li、U,而负载荷端出现Be、Rb、Cs、Th等4种元素。经统计学计算得到PC1的方差贡献率为70.99%,PC2的方差贡献率为14.85%,PC1足以反映主要的数据特征[50];且Li元素对PC2的正载荷端贡献并不大。故选择PC1作为主要考虑依据,将研究区的元素分为2种组合:①Li -Rb-Cs-U-Th和②Be。由于本区主要寻找隐伏锂矿,故最终选取Li元素所在的PC1正载荷端作为异常划分的主要依据。

图4 研究区3类数据主成分双标图

3.2.2 地球化学平面分布特征

为进一步考察ilr变换的优势,对3种方法得到的PC1主成分的平面分布特征进行对比。图5中可以明显看到原始数据和对数变换的PC1得分图较为相似,而ilr变换与前两者有一定区别。原始数据和对数变换的PC1得分图总体呈NE向展布,两者异常范围和大小也几乎一样(图5a、图5b)。ilr变换的PC1得分图(图5c)与前两者也有一定区别:在研究区南部的异常区较前两者相比被明显分隔开,异常面积也比前两者小;而在研究区北部的异常范围走向更加贴近变质分带界线,异常面积也比前两者要小。

a—原始数据

b—对数变换

c—ilr变换

因此,经过ilr变换的主成分得分图更能反映出研究区地电化学异常整体的NE向分布特征。热接触变质带是川西地区寻找伟晶岩型锂矿的重要标志之一,而红柱石、堇青石变质属于容须卡地区热接触变质带的具体表现[51],矿体往往赋存在这些变质带分界线附近。表明经过ilr变换处理数据后,使地电化学异常的分布与已有的地质资料吻合度较高。

3.3 S—A异常分解

为查明研究区与成矿相关的弱异常信息以便精确限定有利成矿区域,对ilr变换后的PC1组分得分进行计算,用能谱密度—面积分形(S—A)模型来分解异常及背景。异常与背景的滤波值由能谱面积图(图6)确定。

图中各线段的斜率大小代表元素的空间自相似性及分形特征。其中第一段回归方程为y=-0.353x+4.942,代表干扰(噪音);第二段方程为y=-0.758x+6.763,代表Li-Rb-Cs-U-Th元素组合的地球化学异常;第三段方程为y=-1.67x+13.09,代表地球化学背景。分别计算3个方程交点的(S,A),并将其通过傅里叶逆变换得到异常及背景的阈值(图7)。图7a显示研究区的背景呈现东西低中间高的特点。对比图5c和图7b可以发现,通过S—A异常分解之后在西部低背景范围内,异常被明显强化;而在中部高背景的范围内,异常范围缩小。

图6 ilr变换后的能谱面积(S—A)双对数图

3.4 讨论

S—A分解的优势在于强化低背景下的弱异常,缩小高背景下的异常范围。最终划分4个异常(图8),故分别讨论这4个异常在S—A分解前后的变化。Ⅰ号异常主要位于研究区北部的红柱石变质带边缘的高背景区(图7a),经过S—A分解后异常范围变小(图9)。Ⅱ号异常主要位于研究区东部NE走向的堇青石变质带周围的高背景值区(图7a)内,经过S—A分解后异常更精确,范围变小(图8)。Ⅲ号异常位于研究区西部低背景值区内,经过S—A分解后,削弱周围高背景值的影响,而使得该异常得到加强,相较(图5c)该异常区范围变小,使异常更为精确。Ⅳ号异常位于研究区东部,堇青石变质带周围的高背景值区(图7a),经过S—A分解之后,该异常范围变小(图8)。

总体来说,S—A分解后的异常形态更加合理,更加符合地质背景和成矿规律。从图9来看,Ⅰ号异常位于红柱石变质带附近,稍偏西部地区;Ⅱ、Ⅳ号异常位于堇青石变质带附近,异常走向与该变质带走向一致。这一特征符合该研究区的区域成矿规律。区域上形成的几个构造岩浆热穹隆均是岩浆的多期侵入形成围绕隐伏岩体周边的变质分带;而伟晶岩脉则分布于每一次岩浆侵入的内外接触带中,往往也在变质分带的接触带附近形成岩脉。而锂矿也相应赋存在伟晶岩脉中,故形成研究区内位于变质分带接触面附近的带状异常,以NE向为主(少量NW向)。

a—S—A分解后的PC1背景

b—S—A分解后的PC1异常

4 靶区划分

4.1 地质背景分析

由于研究区面积较小,区内仅有上三叠统新都桥组一段存在,其岩性主要为石榴石十字石红柱石二云母片岩。该地层主要是由泥质岩经变质作用后形成,泥质岩密度较大、渗透率较小,具有良好的屏蔽作用,一定程度上限制热量的扩散和成矿有利物质的流失,在成矿作用的“源—运—储”系统中担任重要角色。

其次,受容须卡岩浆热液穹窿影响而形成的热接触变质带是控制矿脉空间展布的重要因素[46],研究区北部,有地表出露与伟晶岩脉关系密切的红柱石带,笔者认为该区域值得重点关注。

4.2 靶区划分

综合分析地电提取异常特征,结合研究区地质背景和成矿规律,共圈定出4个找矿靶区分别是Ⅰ号、Ⅱ号、Ⅲ号、Ⅳ号(图8)。

图8 找矿预测靶区划分图

Ⅰ号靶区

该靶区位于研究区北部,近NE向延伸。原始数据、对数处理及ilr变换等3种数据变换方法的结果都呈现大范围的高异常值;经过S—A分解,异常范围变小,且更贴近红柱石带的分界线。经野外查验发现该靶区主要以红柱石变质带为主。故推测Ⅰ号靶区有良好找矿潜力。

Ⅱ号靶区

该靶区位于研究区北东部。经过S—A分解,该靶区的范围变小,原因可能是受到削弱高背景值的影响。该异常位于红柱石带与堇青石带分界线的西侧红柱石带中,亦具有良好的成矿地质条件。

Ⅲ号靶区

该靶区位于研究区东南部。3种数据变换方法均有异常显示,但经过S—A分解,异常面积有明显变化。其西侧的面积缩小,而东侧的面积扩大;原因是西侧高背景被压制,而东侧低背景区的弱异常被强化。该异常位于红柱石带与堇青石带分界线的西侧红柱石带中,亦具有良好的成矿地质条件。

Ⅳ号靶区

该靶区位于研究区西南部。只有ilr变换才出现异常,另外2种变换方法均为分散的弱异常。该异常位于红柱石带中部,为ilr变换后新发现的异常,值得从地质成矿的角度去探讨其成矿地质前景。

5 结论

本文采用成分数据和多重分形S—A方法的地电提取法在容须卡地区开展地电化学元素组合信息提取与组合异常靶区圈定的工作,得到以下结论:

(1)对比原始数据、对数变换、ilr变换等3种数据的统计分布特征,发现ilr变换后数据更加符合正态分布,能够有效打开原始数据的闭合性,使得元素组合更加符合地质、成矿规律。

(2)利用ilr数据进行主成分分析(PCA),确定研究区与成矿密切相关的元素组合是Li-Rb-Cs-U-Th。

(3)运用能谱—面积(S—A)多重分形模型,对研究区地电提取数据异常、背景的非线性变化特征,使异常形态更加合理和符合地质成矿规律。强化低背景区的弱异常,弱化高背景区的强异常。最终圈定Ⅰ号、Ⅱ号、Ⅲ号、Ⅳ号等4个有利找矿靶区。