分維變點法在遙感蝕變異常提取中的應用

韓海輝， 王藝霖， 楊敏， 任廣利， 楊軍錄， 李健強， 高婷

(1.中國地質調查局西安地質調查中心,西安 710054; 2.長安大學地質工程與測繪學院，西安 710054)

分維變點法在遙感蝕變異常提取中的應用

韓海輝1,2， 王藝霖1， 楊敏1， 任廣利1， 楊軍錄1， 李健強1， 高婷1

(1.中國地質調查局西安地質調查中心,西安 710054; 2.長安大學地質工程與測繪學院，西安 710054)

從遙感蝕變異常主分量圖像中提取異常信息，目前主要基于數據的正態分布特征，并未考慮地質異常的非線性特征。針對此問題，提出了分維變點提取算法(fractal dimension-change point method,FDCPM)。首先利用分形模型計算蝕變異常的自相似性參數值，再采用變點模型計算蝕變異常的突變性參數值，然后確定蝕變異常臨界閾值，達到有效區分地質背景和蝕變異常的目的。以甘肅北山方山口地區為研究區，對識別的ASTER蝕變異常進行測試和驗證，并對其提取精度進行初步評價和比較。結果表明： 對于實驗中的褐鐵礦、絹云母和綠泥石3種蝕變礦物來說，分維變點法的總體提取精度略高于門限化方法。利用分維變點法提取3種蝕變礦物的正確率均超過83%，且遙感蝕變異常的分布與化探及重砂異常有較高的吻合度，已知金屬礦(化)點也幾乎都落在遙感蝕變異常區內或其邊緣，表明分維變點法行之有效，可作為今后劃分地質背景和蝕變異常的方法之一。

分維變點法(FDCPM)； 遙感蝕變異常； 自相似性； 突變性； ASTER

0 引言

在遙感地質找礦應用中，遙感圖像識別和提取蝕變礦物信息已成為一個研究熱點[1]，這是因為近礦蝕變作用形成的蝕變礦物能為地質找礦提供重要依據[2]。近年來，大量研究集中于用波段比值、主成分分析、光譜匹配等不同方法識別遙感蝕變礦物異常信息(以下簡稱蝕變異常)[1, 3-4]，但對如何從識別的遙感蝕變分量圖像中有效提取蝕變異常關注較少。目前，最廣泛使用的蝕變異常提取方法大多基于遙感數據的正態分布特征，以均值加n倍標準差作為劃分蝕變異常的閾值，如張玉君等[5]提出的主分量門限化技術。但除正態分布特征外，地質異常還具有自相似性(self-similarity)、突變性(mutability)和不連續性(discontinuity)等非線性特征[6]，為了通過定量特征來有效提取蝕變異常，本文提出分維變點法(fractal dimension-change point method,FDCPM)。該算法首先通過分形模型分析蝕變異常的自相似性特征，再通過變點模型分析蝕變異常的突變性特征，然后確定蝕變異常閾值，區分地質背景和蝕變異常； 并對提取結果以同樣的方法進一步劃分等級。本文還對分維變點法在蝕變異常提取中的有效性進行了分析和探討。

1 分維變點法(FDCPM)

1.1 自相似性和分形模型

分形的本質特性是自相似性[7]，而成礦過程中各種尺度地質異常的大小和個數通常具有自相似性[8]，因而地質異常也具有分形特征。通用的分形模型[9]如下。設有數據集合{Xi，i=1,2，…,N}，則

(1)

式中：r為特征尺度；K為常數；D為分維數；N(r)為尺度大于等于r的Xi的數目。

將式(1)兩邊取對數可得

(2)

由式(2)可以看出，lgr與lgN(r)呈線性關系。計算時，可以遙感蝕變分量圖像的像元亮度值為r值，則N(r)為亮度值大于等于r的Xi的數目。

1.2 突變性和變點模型

地質異常還具有突變性特征[6]，變點是指輸出序列在某未知時刻發生的突然變化，因此變點分析可定量反映地質異常的突變特征。通過變點統計模型可有效檢驗樣本序列是否存在變點，從而判斷變點的位置與個數[10]。采用的均值變點模型[11]如下。

(3)

(4)

(5)

3)計算期望值(E)，即

(6)

由上可見，樣本序列如果存在變點，那么變點會使S和Si的差值增大； 而S是固定值，因此最小的Si值即對應存在的變點。計算時，可先對lgN(r)與lgr的商取對數，以構建樣本序列X0； 再計算X0的統計量Si。

2 研究區概況與數據預處理

2.1 研究區概況

實驗區選擇在甘肅北山方山口地區。區內出露地層主要為中泥盆統三個井組(D2s)，巖性為淺變質中細粒長石石英砂巖、長石砂巖夾粉砂巖和泥質粉砂巖。斷裂構造以近EW向壓扭性斷裂為主，沿走向局部向NEE向偏離，形成弧形構造線。區內巖漿活動強烈，以華力西期為主，次為印支期，從基性―酸性均有產出。其中，輝長巖-輝綠巖呈巖株狀產出，石英閃長巖呈巖基出露于研究區中南部，二長花崗巖呈NE向不規則侵入于石英閃長巖及三個井組地層中； 花崗斑巖呈巖脈狀或不規則巖株狀產出，穿插于石英閃長巖中。區內主要分布熱液型礦床，尤以金礦化為主。區內存在的蝕變類型有硅化、黃鐵礦化(褐鐵礦化)、絹云母化、鐵碳酸鹽巖化和青磐巖化等[12]。

2.2 遙感數據預處理

ASTER傳感器搭載的Terra衛星發射于1999年12月，其數據具有3個可見光和近紅外波段(VNIR，0.52～0.86 μm)、6個短波紅外波段(SWIR，1.6～2.43 μm)和5個熱紅外波段(TIR，8.125～11.65 μm)，較寬的光譜覆蓋范圍可供提取多種遙感蝕變礦物信息。

ASTER數據獲取時間為2000年10月16日，數據級別為L1B，已經過輻射校正和幾何糾正，無云及陰影遮蓋，影像質量較好。本次實驗中對該圖像數據進行了串擾校正[13](串擾是成像光譜儀非理想成像特性中的一種，分為電串擾和光串擾，廣泛存在于成像光譜儀中)、大氣校正和重采樣處理(圖像分辨率重采樣至30 m)，并選擇對河流、植被及第四系等干擾地物敏感的波段，通過掩模處理降低其對異常信息提取的影響。

3 蝕變異常信息提取

3.1 自適應一致估計算法

對于遙感圖像的混合像元分解問題，研究人員提出了大量的光譜解混算法。其中，自適應一致估計法(adaptive coherence estimator,ACE)是一種針對遙感圖像非結構化目標探測的算法，它將像元背景視為符合多元高斯分布的統計模型，將背景和噪聲統一視為背景進行建模，通過貝葉斯統計計算概率得分來識別目標像元[14]，其優點是可以基于較少的已知樣本就能有效檢測目標[15]。

ACE是廣義似然比探測法(generalized likelihood ratio test,GLRT)的一種自適應版本[14]。假設像元背景是一個混合模型，且背景b服從多元正態分布b～N(μ0,Γ)，將均值μ0從觀測值x中移去后，則有

H0：x=b，H1=S·a+σ2·b,

(7)

式中：H0和H1為2種假設地物；H0和H1的分布分別為x～N(0,Γ)和x～N(S·a,σ2·Γ)，即背景在2種假設下具有相同的協方差結構，但有著不同的方差，該方差直接關系到目標物體占像元區域的比例。在上述情況下，ACE算法可表達為

(8)

式中：Γ為圖像協方差矩陣；x為光譜向量；S為目標向量。

3.2 FDCPM蝕變異常提取

本次實驗提取了褐鐵礦、絹云母和綠泥石3種蝕變礦物異常信息。實驗中均采用ACE目標檢測算法，以USGS礦物標準波譜庫為基礎，通過ENVI軟件提取3種蝕變異常圖像，并將提取結果的亮度值拉伸到0～255范圍； 然后應用分形模型計算異常圖像的雙對數圖，以分析其自相似性； 并在此基礎上采用均值變點法分析異常的突變性特征，即通過計算獲取最小的Si值。

對于褐鐵礦異常，在i=59時，Si值最小，對應的像元亮度值r=60，故r=60可作為地質背景和褐鐵礦Ⅲ級異常的臨界點； 繼續利用上述方法對r≥60的數據進行計算，分別得到r=111時的Ⅱ,Ⅲ級異常的臨界點和r=119時的Ⅰ,Ⅱ級異常的臨界點。圖1(a)展示了地質背景和蝕變異常的自相似性分布特征以及各自的自相似性擬合直線，圖1(b)展示了Ⅰ,Ⅱ,Ⅲ級異常的臨界閾值及各自的自相似性擬合直線,不同顏色線段匯交點為各級異常臨界閾值。

(a) 異常分量雙對數圖(b) 蝕變異常分級臨界值雙對數圖

圖1褐鐵礦蝕變分形特征分析圖

Fig.1Diagramoffractalcharacteristicsanalysisforlimonitealteration

從圖1可以看出，褐鐵礦蝕變異常段中3段直線方程的擬合率均較好，表明研究區內褐鐵礦蝕變異常區域具有多重分形特征，不同程度的蝕變異常對應不同的分維值D(圖1(b))。根據確定的臨界閾值對蝕變異常圖像進行密度分割，即可制作出褐鐵礦蝕變異常分級圖。

對于絹云母異常，在i=73時，Si值最小，對應的像元亮度值r=74，故r=74可作為地質背景和絹云母Ⅲ級異常的臨界點； 繼續利用上述方法對r≥74的數據進行計算，分別得到r=125時的Ⅱ和Ⅲ級異常的臨界點以及r=132時的Ⅰ和Ⅱ級異常的臨界點。圖2為絹云母蝕變異常灰度值分形特征分析圖，可以看出異常段中3段直線方程的擬合率均超過0.91，表明研究區內絹云母蝕變異常也具有多重分形特征，不同程度的蝕變異常對應不同的分維值D(圖2(b))，不同顏色線段匯交點為各級異常臨界閾值。根據確定的臨界閾值可制作出絹云母蝕變異常分級圖。

(a) 異常分量雙對數圖(b) 蝕變異常分級臨界值雙對數圖

圖2絹云母蝕變分形特征分析圖

Fig.2Diagramoffractalcharacteristicsanalysisforsericitealteration

對于綠泥石異常，在i=169時，Si值最小，對應的像元亮度值r=170，故r=170可作為地質背景和綠泥石Ⅲ級異常的臨界點； 繼續利用上述方法對r≥170的數據進行計算，分別得到r=221時的Ⅱ,Ⅲ級異常的臨界點和r=229時的Ⅰ,Ⅱ級異常的臨界點。圖3為綠泥石蝕變異常灰度值分形特征分析圖，可以看出異常段中3段直線方程的擬合率均超過0.98，表明研究區綠泥石蝕變異常同樣具有多重分形特征，不同程度的蝕變異常對應不同的分維值D(圖3(b))，不同顏色線段匯交點為各級異常臨界閾值。根據確定臨界閾值可制作出綠泥石蝕變異常分級圖。

(a) 異常分量雙對數圖(b) 蝕變異常分級臨界值雙對數圖

圖3綠泥石蝕變分形特征分析圖

Fig.3Diagramoffractalcharacteristicsanalysisforchloritealteration

圖4 為根據上述臨界閾值制作的實驗區褐鐵礦、絹云母和綠泥石遙感蝕變異常分布圖。

圖4實驗區遙感蝕變異常分布圖

Fig.4Distributionofremotesensingalterationanomaliesinexperimentalarea

從圖 4可以看出，蝕變礦物信息呈近EW向線性展布，與區域構造線方位基本一致，區內已知熱液型金屬礦(化)點幾乎都落在遙感蝕變異常區內或其邊緣； 同時，遙感蝕變異常的分布與已有的化探異常和重砂異常有較好的吻合性。

4 精度評價與分析

在使用蝕變異常提取結果前, 必須進行客觀可靠的精度驗證和分析。目前，點樣本檢驗是遙感分類評估中常用的檢驗方式，它以單個像元作為基本的抽樣單元，隨機選取一定數量的樣本，利用參考數據進行驗證，從而對提取結果進行估計。多項研究也表明，采用點樣本檢驗方法評價遙感分類結果的精度是比較合適的[16]。

本文采用點樣本檢驗法對分維變點法提取結果的精度進行評價。實驗中對128個野外驗證點的異常分布結果(部分驗證點存在同時發育幾種蝕變礦物的情況)以混淆誤差矩陣的形式進行了統計(表1)。雖然有128個野外樣本點，但數量還不足以評價識別率，因此本文僅對提取精度進行討論。

表1 用分維變點模型提取蝕變礦物的混淆矩陣

從表1中的統計結果可以看出，3種蝕變異常提取結果的總體精度較為理想，褐鐵礦和絹云母均超過了80%，綠泥石則超過了90%。不同蝕變異常之間存在誤提，如提取的褐鐵礦在野外觀測為綠泥石或絹云母，對這些抽樣點采集的巖石樣品進行光譜測量后發現其光譜曲線上往往存在多個蝕變礦物的吸收譜帶特征，只是不同礦物特征峰的吸收強弱不同。因此，利用低光譜分辨率的ASTER提取蝕變礦物時，不可避免地會出現這樣的誤提。

為進一步比較3種蝕變異常提取結果，本文實驗還以門限化方法[5]對3種蝕變異常進行了提取，其中N值分別取1.5,2和2.5，計算的臨界閾值見表2； 然后對128個相同驗證點的蝕變異常提取結果進行精度評價(表3)。

表2 遙感蝕變礦物異常信息統計對比表

表3 用門限化法提取蝕變礦物的混淆矩陣

對比表1和表3中的統計結果可以看出，對本文實驗中的3種蝕變異常來說，分維變點法的總體精度略高于門限化方法，表明本文提出的分維變點法取得了較好的效果。

5 結論

1)從地質異常存在非線性特征的因素入手，本文提出了對遙感蝕變異常提取的分維變點法(FDCPM)，包括自相似性分析和突變性分析2個過程。本文實驗獲取的3種蝕變異常的正確率均超過83%，且遙感蝕變異常的分布與化探及重砂異常吻合度較高，證實了FDCPM的可行性。

2)通過野外驗證以及與門限化方法比較，本文提出的FDCPM具有較好的蝕變礦物提取效果，可作為一種劃分地質背景和遙感蝕變異常的有效方法。

3)地質異常還具有不連續性、多樣性等其他非線性特征，如何定量識別這些特征并將其作為參數加入到地質異常識別模型中，還需要進一步研究。

[1] van der Meer F D,van der Werff H M A,van Ruitenbeek F J A,et al.Multi-and hyperspectral geologic remote sensing:A review[J].International Journal of Applied Earth Observation and Geoinformation,2012,14(1):112-128.

[2] 童慶禧,田國良.中國典型地物波譜及其特征分析[M].北京:科學出版社,1990. Tong Q X,Tian G L.Spectral and Analysis of Typical Earth Objects of China[M].Beijing:Science Press,1990.

[3] Sadeghi B,Khalajmasoumi M,Afzal P,et al.Using ETM+and ASTER sensors to identify iron occurrences in the Esfordi 1∶100 000 mapping sheet of central Iran[J].Journal of African Earth Sciences,2013,85:103-114.

[4] Pournamdari M,Hashim M,Pour A B.Spectral transformation of ASTER and Landsat TM bands for lithological mapping of Soghan ophiolite complex,South Iran[J].Advances in Space Research,2014,54(4):694-709.

[5] 張玉君,曾朝銘,陳薇.ETM+(TM)蝕變遙感異常提取方法研究與應用——方法選擇和技術流程[J].國土資源遙感,2003,15(2):44-49,doi:10.6046/gtzyyg.2003.02.11. Zhang Y J,Zeng Z M,Chen W.The methods for extraction of alteration anomalies from the ETM+(TM) data and their application:Method selection and technological flow chart[J].Remote Sensing for Land and Resources,2003,15(2):44-49,doi:10.6046/gtzyyg.2003.02.11.

[6] 陳聆,郭科,柳柄利,等.地球化學礦致異常非線性分析方法研究[J].地球物理學進展,2012,27(4):1701-1707. Chen L,Guo K,Liu B L,et al.The study of non-linear analysis method of geochemical ore-forming anomaly[J].Progress in Geophysics,2012,27(4):1701-1707.

[7] 成秋明,張生元,左仁廣,等.多重分形濾波方法和地球化學信息提取技術研究與進展[J].地學前緣,2009,16(2):185-198. Cheng Q M,Zhang S Y,Zuo R G,et al.The application of fractal dimension-change point method to the extraction of remote sensing alteration anomaly[J].Earth Science Frontiers,2009,16(2):185-198.

[8] 成秋明.成礦過程奇異性與礦床多重分形分布[J].礦物巖石地球化學通報,2008,27(3):298-305. Cheng Q M.Singularity of mineralization and multifractal distribution of mineral deposits[J].Bulletin of Mineralogy,Petrology and Geochemistry,2008,27(3):298-305.

[9] 梁鈺琦,王功文,朱彥彥,等.分形方法在遙感蝕變信息提取中的應用研究[J].遙感技術與應用,2011,26(4):508-511. Liang Y Q,Wang G W,Zhu Y Y,et al.Alteration from ETM+data rating based on fractal technologies[J].Remote Sensing Technology and Application,2011,26(4):508-511.

[10]項靜恬,史久恩.非線性系統中數據處理的統計方法[M].北京:科學出版社,1997. Xiang J T,Shi J E.The Statistical Method of Data Processing in Nonlinear System[M].Beijing:Science Press,1997.

[11]韓海輝,高婷,易歡,等.基于變點分析法提取地勢起伏度——以青藏高原為例[J].地理科學,2012,32(1):101-104. Han H H,Gao T,Yi H,et al.Extraction of relief amplitude based on change point method:A case study on the Tibetan Plateau[J].Scientia Geographica Sinica,2012,32(1):101-104.

[12]任廣利,楊軍錄,楊敏,等.高光譜遙感異常提取在甘肅北山金灘子——明金溝地區成礦預測中的應用[J].大地構造與成礦學,2013,37(4):765-776. Ren G L,Yang J L,Yang M,et al.Application of hyperspectral remote sensing anomaly information on metallogenic prediction in the Jintanzi：Mingjingou area of Beishan,Gansu[J].Geotectonica et Metallogenia,2013,37(4):765-776.

[13]時丕龍,付碧宏,二宮芳樹.基于ASTER VNIR-SWIR多光譜遙感數據識別與提取干旱地區巖性信息——以西南天山柯坪隆起東部為例[J].地質科學,2010,45(1):333-347. Shi P L,Fu B H,Ninomiya Y.Detecting lithologic features from ASTER VNIR-SWIR multispectral data in the arid region:A case study in the eastern Kalpin uplift,southwest Tian Shan[J].Chinese Journal of Geology,2010,45(1):333-347.

[14]張良培.高光譜目標探測的進展與前沿問題[J].武漢大學學報:信息科學版,2014,39(12):1387-1394,1400. Zhang L P.Advance and future challenges in hyperspectral target detection[J].Geomatics and Information Science of Wuhan University,2014,39(12):1387-1394,1400.

[15]Liu W J,Xie W C,Tong H B,et al.Adaptive coherence estimator based on the Krylov subspace technique for airborne radar[J].Journal of Systems Engineering and Electronics,2015,26(4):705-712.

[16]劉旭攏,何春陽,潘耀忠,等.遙感圖像分類精度的點、群樣本檢驗與評估[J].遙感學報,2006,10(3):366-372. Liu X L,He C Y,Pan Y Z,et al.Accuracy assessment of thematic classification based on point and cluster sample[J].Journal of Remote Sensing,2006,10(3):366-372.

(責任編輯：邢宇)

lgN(r)=lgK-Dlgr。

Applicationoffractaldimension-changepointmethodtotheextractionofremotesensingalterationanomaly

HAN Haihui1，2， WANG Yilin1， YANG Min1， REN Guangli1， YANG Junlu1， LI Jianqiang1， GAO Ting1

(1.Xi’anCenterofChinaGeologicalSurvey，Xi’an710054，China; 2.SchoolofGeologicalandSurveying,Chang’anUniversity,Xi’an710054,China)

At present, the extracting method for remote sensing alteration anomalies from principal component image relies mainly on the data’s normal distribution, without considering the nonlinear characteristics of geological anomaly. To tackle this problem, the authors have proposed the fractal dimension-change point method(FDCPM)in this paper. By calculating the self-similarity and mutability of alteration anomalies with fractal dimension-change point model, the critical threshold of an alteration anomaly was acquired quantitatively. The realization theory and access mechanism of the method were elaborated by an experiment with ASTER data in Fangshankou，Beishan，and the results of the proposed method and traditional method (de-interfered anomalous principal component thresholding technique，DIAPCTT) were compared with each other. The results show that the FDCPM has a relatively high extracting precision than the DIAPCTT for three alteration minerals in the experiment. In this experiment, the accuracy of three alteration minerals could reach over 83%. Moreover, the distribution of remote sensing alteration anomalies agrees well with a large amount of evidence from the geochemical anomaly and the heavy sand anomaly. What’s more, the known polymetallic ore spots and mineralized spots fall in the zone of remote sensing alteration anomaly or at its edge. All the results mentioned above show that the FDCPM is one of the effective distinguishing methods for the geological background and the remote sensing alteration anomaly in the future.

fractal dimension-change point method(FDCPM)； remote sensing alteration anomaly； self-similarity； mutability； ASTER

10.6046/gtzyyg.2017.03.20

韓海輝,王藝霖,楊敏,等.分維變點法在遙感蝕變異常提取中的應用[J].國土資源遙感,2017,29(3)：137-142.(Han H H,Wang Y L,Yang M,et al.Application of fractal dimension-change point method to the extraction of remote sensing alteration anomaly[J].Remote Sensing for Land and Resources,2017,29(3)：137-142.)

2016-02-29；

2016-07-11

國家自然科學基金項目“綠泥石礦物近紅外光譜吸收譜帶的位移機理與控制機制研究”(編號： 41502312)、中國地質調查局地質調查項目“天山—北山成礦帶那拉提-營毛沱地區地質礦產調查”(編號： DD20160009)和“十二五”國家科技支撐計劃項目“岔路口―神仙灣鉛鋅成礦帶成礦地質背景及靶區優選”(編號： 41502312)共同資助。

韓海輝(1983-),男,高級工程師,主要從事遙感地質研究。Email： hanhh06@hotmail.com。

TP 79

： A

： 1001-070X(2017)03-0137-06