分维变点法在遥感蚀变异常提取中的应用

韩海辉， 王艺霖， 杨敏， 任广利， 杨军录， 李健强， 高婷

(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