k均值聚類法在識別低弱地球化學異常中的應用：以內蒙古阿木牛林場地區(qū)為例

周吉晨，李隨民，張詩敏，焦 建

河北地質大學，河北 石家莊 050031

0 前言

數(shù)據(jù)處理與解釋是勘查地球化學工作的重要環(huán)節(jié)，在化探數(shù)據(jù)的處理中，背景值及異常下限的確定對地球化學異常的圈定與解釋具有重要的意義。土壤是在已風化基巖之上巖石(礦石)風化作用的殘留疏松物[1]，對原巖具有繼承性，土壤數(shù)據(jù)會因受到原巖的影響而出現(xiàn)異常，若不能準確確定背景值，就無法客觀的反映地質信息。許多學者對消除背景值的影響進行了研究，提出了多種解決方法[2-7]，這些方法對于強異常的發(fā)現(xiàn)和提取有較好的效果，但對低弱異常的識別效果不佳。目前針對低弱異常識別的方法有：子區(qū)中位數(shù)襯值濾波法[8]、趨勢面法[9-10]、等。低弱異常的識別對找礦具有重要的意義。

論文利用內蒙古阿木牛林場地區(qū)1∶5萬化探數(shù)據(jù)，通過k均值聚類方法實現(xiàn)對研究區(qū)樣品的分類，每一類可視為一個單獨的個體。然后對每個子類中的數(shù)據(jù)進行標準化和分析，從而達到突出低弱地球化學異常的目的，為進一步找礦提供參考。

1 研究區(qū)地質概況

調查區(qū)位于大興安嶺北段東坡，大興安嶺主脊在西北部通過。地勢西高東低，山脈走向以近南北向、北西向、北東向為主，總體構造線方向為北東向、北東東向，其中北東向構造控制著地層和侵入巖的展布方向，次之為北西向構造，褶皺構造不太發(fā)育。地層自下而上劃分為古生界下奧陶系銅山組、侏羅系滿克頭鄂博組、瑪尼吐組、白音高老組，第四系全新統(tǒng)沼澤堆積物和沖洪積物等。其中侏羅系分布最廣，奧陶系少量分布，第四系沿溝谷、河床分布。古生界下奧陶系銅山組的巖石組合以變質砂巖、變泥質粉砂巖、變細砂巖、變砂質泥巖為主，中間夾有變中酸性火山熔巖和大理巖。滿克頭鄂博組為一套酸性火山熔巖、火山碎屑巖噴發(fā)沉積組合，多數(shù)直接噴發(fā)不整合覆蓋與古生界地層和中侏羅世侵入巖之上?，斈嵬陆M是以中性火山巖—火山碎屑巖噴發(fā)沉積為主的一套地層，巖石組合以安山巖、安山質含角礫熔結凝灰?guī)r、英安巖、多斑角閃安山巖為主，局部以安山巖與火山碎屑巖互層出現(xiàn)。白音高老組地層是一套偏堿性酸性火山巖—火山碎屑巖夾少量中酸性火山巖等特征的巖石組合，其巖石組合顏色整體偏淺，多呈現(xiàn)灰白色，為灰白色流紋巖與灰白色流紋質凝灰?guī)r互層產(chǎn)出為特征。區(qū)內侵入巖十分發(fā)育，主要有中三疊紀二長花崗巖，中侏羅世二長花崗巖，早白堊紀二長花崗巖、花崗斑巖、石英二長斑巖、石英正長巖等，其中中侏羅世二長花崗巖在研究區(qū)發(fā)育大面積出露，主要分布在研究區(qū)中部和東部(圖1)。

圖1 研究區(qū)地質概略

2 基本原理及方法

聚類是數(shù)據(jù)挖掘的一種重要的手段，k均值聚類算法是一種得到最廣泛使用的基于劃分的聚類算法，該算法簡單、快速，對處理大數(shù)據(jù)集，該算法是相對可伸縮和高效的[11]。論文利用該方法對數(shù)據(jù)進行分類，然后進行正態(tài)轉換，對異常樣品進行剔除，將剩余數(shù)據(jù)進行標準化。該方法可形成一個突出低弱地球化學異常的標準化數(shù)據(jù)集，作為進一步計算和編圖的基礎數(shù)據(jù)。

2.1 k均值聚類算法的基本思想

在整個數(shù)據(jù)集中，隨機選擇k個數(shù)據(jù)，每個數(shù)據(jù)作為一個簇的中心，然后計算其它數(shù)據(jù)與這個中心的距離，跟據(jù)距離將這些數(shù)據(jù)劃分到最近的簇。所有數(shù)據(jù)劃分好后，再計算每個簇的平均值，并以該平均值為中心重新聚類，以此往復，直至滿足收斂函數(shù)的要求，即：

2.2 K均值聚類方法描述

輸入：n個數(shù)據(jù)對象，數(shù)量為k

輸出：k個滿足收斂函數(shù)的簇

(1)從n個對象中任意選擇k個對象，每個對象作為初始聚類中心；

(2)計算每個對象與聚類中心的距離，根據(jù)距離把每個數(shù)據(jù)劃分到距離最近的簇；

(3)重新計算每個簇中的平均值，重新選定聚類中心；

(4)計算收斂函數(shù)，滿足收斂條件則結束，否則回到第(2)步。

3 內蒙古阿木牛林場地區(qū)低弱地球化學異常的識別

本次土壤測量采樣密度為8.90點/km2，采集樣品2 518個，所有樣品均采集為基巖風化產(chǎn)物中的殘積層或巖屑，且為多點采集，組合成為一個樣品。采樣時避開了各種污染、廢石堆和河床堆積物，不能取樣時棄點在記錄中進行了注明；采樣后留有明顯標志。土壤測量以-4～+20目為采樣粒級。通過對數(shù)據(jù)的預處理，剔除特異值后，本次分類的有效數(shù)據(jù)為1 675個。

3.1 分類數(shù)的確定

分類主要依據(jù)來源為內蒙古阿木牛林場地區(qū)地質單元的數(shù)目和k均值聚類法的輪廓系數(shù)。

最優(yōu)分類數(shù)k的選取對于k均值聚類方法至關重要，最優(yōu)分類數(shù)可以通過輪廓系數(shù)確定[12]。輪廓系數(shù)的下公式計算如下[13]：

其中，Si代表第i個變量的輪廓系數(shù)；ai代表一個簇中第i個變量到其它所有變量的平均距離；bi代表到其它所有變量的最小距離。通過給定的初始分類數(shù)可以計算出相應的平均輪廓系數(shù)，輪廓系數(shù)的值在-1至1之間，值越大則代表聚類效果越好。同時選取的類數(shù)不應過小。

根據(jù)內蒙古阿木牛林場地區(qū)的地質單元數(shù)目和K均值聚類的輪廓系數(shù)，該地區(qū)可分為5大類。(圖2)

圖2 研究區(qū)分類圖

3.2 分類標準化

首先檢測每個子類內元素是否符合正態(tài)分布，對于近似正態(tài)分布的數(shù)據(jù)，采用3S法剔除異常值，直至符合正態(tài)分布。剔除異常值后，通過檢測，近似正態(tài)分布的元素按照以下公式進行均值和標準差的運算，即：

呈近似對數(shù)正態(tài)分布的元素按照以下公式進行均值和標準差的運算，即：

其中，μ為數(shù)據(jù)的自然對數(shù)值；σ為數(shù)據(jù)的自然對數(shù)標準離差。

根據(jù)各類數(shù)據(jù)的均值和標準差，計算異常下限值和襯度值，然后將所有子類數(shù)據(jù)合并為一個數(shù)據(jù)集。

內蒙古阿木牛林場地區(qū)1：5萬化探數(shù)據(jù)可分為5類，論文選取研究區(qū)潛在成礦元素Pb、Zn元素與傳統(tǒng)方法圈定的異常進行對比。首先，對每類數(shù)據(jù)正態(tài)化后進行3S檢驗，將異常值剔除，然后計算標準差、方差和異常下限值，如表1、表2所示：

表1 k均值聚類法Pb元素異常下限值

根據(jù)各類的異常下限值，利用Excel表分別計算5類中Pb、Zn元素的襯度值，即Pb、Zn元素與該元素異常下限的比值。將得到的襯度值合成為一個數(shù)據(jù)集，并用surfer軟件制作等值線圖。

3.3 分類數(shù)據(jù)與原始數(shù)據(jù)進行對比

將未分類原始數(shù)據(jù)中的Pb、Zn元素進行特異值剔除，使其符合或近似符合正態(tài)分布，然后采用平均值加2倍標準差的方法確定Pb、Zn元素的異常下限值，作為傳統(tǒng)方法圈定元素異常的依據(jù)(表3)。

表3 傳統(tǒng)方法Pb、Zn元素異常下限值

從傳統(tǒng)方法和k均值分類法(圖3、圖4)圈定的異?？梢钥闯?，Pb、Zn元素的分布受到了巖性地層分布的控制，大部分異常區(qū)與中侏羅世花崗巖相對應。若研究區(qū)內巖性變化較大，按統(tǒng)一的異常下限確定元素異?？赡軆H是對巖性的反映，不能真實反映異常情況。因此，傳統(tǒng)方法圈定異常并不是適用地質條件較復雜的地區(qū)。而k均值聚類法不僅使原有的北北西方向的Pb異常得到了加強，而且顯示出Pb、Zn的低弱異常(實線圈定區(qū))。若使用傳統(tǒng)方法，這些異常往往被忽略。

圖3 Pb元素異常圈定對比圖

圖4 Zn元素異常圈定對比圖

4 結束語

論文利用k均值聚類法對研究區(qū)元素進行分類，并對每個子類的數(shù)據(jù)進行標準化，計算其均值、方差和襯度值。以Pb、Zn元素為例，在與傳統(tǒng)方法圈定地球化學異常的對比中發(fā)現(xiàn)，分類處理可以有效的強化原有異常，識別低弱地球化學異常。與傳統(tǒng)方法(全地區(qū)統(tǒng)一異常下限值)圈定異常不同的是，分類方法和傳統(tǒng)方法圈定異常的區(qū)域有所差異，分類方法可以消除巖性對異常的影響，顯示出一些低弱異常的區(qū)域，而傳統(tǒng)方法圈定的異?？赡芤驇r性的影響而出現(xiàn)偏差。