基于独立成分分析和形态特征的干涸水系提取

陈军林， 彭润民， 闫岩， 赵伟光

(中国地质大学(北京)地球科学与资源学院，北京 100083)

0 引言

水系提取在地学研究中有很多应用，比如在水系沉积物化探中，需要提取水系网，根据水系网设计采样点；在遥感蚀变矿物异常提取中，河道中的水体和沉积物会对异常提取产生干扰，需要对水系进行提取和掩模；同时水系还是遥感构造分析和岩性解译的重要标志[1]。水系提取根据采用的数据可分为基于数字高程模型(digital elevation model，DEM)的和基于遥感影像的，以及二者结合的方法。基于遥感影像的水系提取大多是依据水系中的水体光谱特征，常用方法有单波段阈值法、比值法、波段差值、谱间分析、模式识别和水体指数法等，其中尤以水体指数法最为常用[2-7]。但是，在冲沟、处于枯水期的季节性河流以及干涸的河道中，因缺乏水体，水系与裸地等其他地物光谱特征相似，故上述基于水体光谱特征的方法并不适用。因此，对于干涸水系，目前主要通过增强水体特征[5]、监督分类以及利用除水体光谱特征以外的其他特征，如形态特征和地形特征等[8-13]进行识别提取。显然，增强水体特征的方法只适合于半干旱的水系；监督分类方法需要人工选择训练样本，存在主观干预大且耗费人力的缺点；而基于水体形态特征及地形特征的方法则需要借助DEM等除遥感影像以外的其他数据。

针对以上问题，本文提出基于独立成分分析(independent component analysis，ICA)和形态特征的干涸水系提取方法。ICA是一种线性分解技术，在目标检测中有很好的应用效果[14-17]。通过优化算法按照统计独立的原则将观测数据分解成若干独立分量，以代表在观测变量上具有不同实际意义的独立信号源。在遥感影像目标检测中，不同地物可以视为不同信号源，它们之间有不同的光谱特征，大多数可以被视为统计独立的实体。因此，通过ICA可以将目标地物和无关地物信号分解到不同分量中，从而分离出目标。利用ICA提取水系的基本思路是通过ICA对包含水系信号的遥感影像进行分离，使得水系独立于其他地物信号，从而识别出水系。ICA是一种非监督方法，不需要选择训练样本，不借助DEM等其他数据，仅利用遥感影像就可进行水系识别，操作简便，能够克服传统方法的诸多不足。

通过ICA提取的水系信息分量图中难免有一些来自于相似地物以及背景的噪声。传统上对于噪声像元通过数学形态学滤波方法[18]进行压制，能够起到一定的作用，但该方法是全局统一尺度的，无法兼顾水系信息的保持和噪声去除，一些较大的噪声图斑无法消除。为了克服这些噪声无法消除的问题，本文提出在数学形态学滤波的基础上再进行基于图斑形态参数的噪声分离。即将经过数学形态学滤波的ICA水系分割二值图像转换为矢量对象，提取每个矢量对象的形态参数，用这些形态参数作为依据进行水系和噪声的分离。

以内蒙古自治区乌拉特后旗获各琦地区为例，应用ICA方法对ASTER影像数据进行线性分解，提取出代表水系信号的独立成分分量，对水系分量进行基于中值滤波的背景抑制，以增强局部信息，减少背景干扰；然后采用数学形态学滤波进行初步去噪；再通过将分割结果转换为矢量对象，提取出对象的形态特征，并用这些特征进行基于规则的噪声消除，最终分割出水系。

1 研究区概况及数据源

1.1 研究区概况

研究区位于内蒙古自治区乌拉特后旗获各琦苏木，阴山山脉狼山段西北侧，处于山脉与高原转换的部位，面积为125 km2。主要地形为山地、低山丘陵以及砂砾戈壁，植被覆盖极少，主要分布在河道及沟谷处。区内水系发育，主要为山洪冲击侵蚀形成，无雨时干涸见底，大雨时则为山洪泄溢通道。该地区属温带大陆性气候区，全年干旱少雨，蒸发量大。绝大多数河床常年暴露，无水体流淌，河道中以第四系河沙、风成沙为主，含有少量砾石及土壤。研究区遥感影像如图1所示。图1(b)为该地区典型的干涸水系地貌，仅少数河道可见地下水形成的持续水体，且规模较小。

(a) 研究区ASTER B3(R),B2(G)，B1(B)合成影像 (b) 典型干涸水系局部放大(必应地图)

1.2 数据源

ASTER影像数据具有光谱和空间分辨率较高、覆盖范围广、数据容易获取的特点，在水系提取中应用广泛[5,19-20]。该数据具体参数见文献[14]。其中热红外波段由于空间分辨率较低，会影响分类精度，故本文仅采用了2005年9月13日的ASTER L1T数据(轨道号130/31)中的B1～B9波段，统一重采样为15 m空间分辨率，用ENVI5.2软件对原始数据进行了辐射定标、大气校正、重采样和裁剪等预处理。研究中截取了水系发育、地物复杂多变的部分区域用来验证本文提出的方法(图1(a))。

2 研究方法

2.1 独立成分分析

ICA是20世纪90年代后期随着盲源信号分离问题发展起来的一种信号处理方法。在不知道源信号和传输通道参数的情况下，从观测信号出发，根据源信号的统计特性找到非高斯数据的线性表示，将其分解成在统计上独立的组分[15,21]，以作为对源信号的一组估计。ICA与主成分分析(principal component analysis，PCA)同属线性变换技术，但PCA只能消除数据之间的二阶相关性，而ICA则考虑到了数据间的高阶统计特性，保证了变换后分量间相互独立[22-24]。

设x=(x1,x2,…,xm)T为观测数据阵，生成该数据阵的独立源为s=(s1,s2,…,sn)T，假定第i个观测信号xi由n个独立成分s线性混合而成，则

xi=ai1s1+ai2s2+…+ainsn,i=1, 2, …,m，

(1)

以向量形式表示为

x=A·s

(2)

式中A=(a1,a2,…,am)是一个n×m的混合矩阵，用来组合叠加信号s。式(2)即ICA的基本模型。

ICA理论认为用来观测的混合数据阵x是由独立源s经过A线性加权获得。A和s都是未知的，x是已知的，ICA的目标就是根据x来求出s的估计y，这个过程也称作为盲源信号分离。求解思路是通过x求得一个分离矩阵W，使得W作用在x上所获得的信号y是独立信号源s的最优逼近，该关系可以表示为

y=W·x=W·A·s

(3)

当估计量与源信号最接近时，有

y=s

(4)

根据式(3)和式(4)可得，W=A-1

(5)

因此，问题求解的思路就是要找到最接近A的逆矩阵A-1的分离矩阵W。ICA模型的估计方法主要有非高斯的最大化、互信息的最小化(互熵最大化)和最大似然函数估计等[25]。本文使用ENVI5.2软件中的ICA模块执行ICA算法。

2.2 水系分量图像增强

通过ICA将数据分解成不同的分量，分别代表对不同独立信号源的近似估计。由于各个分量之间是统计独立的，不同信号源中的信息被集中到不同分量中。对于水系分量来说，水系信息被突出显示，呈现高频特征，而与水系信号在统计上独立性较大的其他信号源的信息在该分量中较少，在图中表现为低频信息，呈现较为连续的低缓背景，即背景值。这些背景值在空间上是有涨落变化的，如果用统一的阈值去分割水系目标和背景，就可能造成错分，导致一些弱的水系信息不能被识别，而一些强的背景信息会被提取为水系。同时，利用ICA分解出的水系分量中可能包含其他与水系信息比较相似的地物信息，会对水系提取造成干扰。综上所述，通过ICA提取的水系分量中可能包含的干扰信号来自于相似信号源和背景信号涨落2方面。对于相似地物造成的干扰较难分离，可以通过后续的基于形态特征的方法进行消除；对于背景造成的干扰，有很多处理方法，其中基于中值滤波的背景抑制方法是一种常用的方法[26]。该方法的具体做法是：选择合适大小的窗口进行中值滤波，用滤波值代表每个像元点的局部背景值，将所有像元点的实际灰度值与滤波值相减，得到的值就能很好地降低背景干扰，增强水系信息。

2.3 基于形态特征的图像去噪

将经过局部信息增强的图像用阈值分割方法分割出水系。这种基于像元灰度值分割出来的结果中必然含有一些噪声，传统上可通过数学形态学滤波将这些噪声与水系信息分离[18,27]。但由于数学形态学滤波是全局统一尺度的滤波，无法兼顾水系信息的保持和噪声去除，仍可能会在消除噪声的同时造成细小水系信息的丢失，且当噪声数据较大时无法消除，因此应用效果比较有限。而将分割出的图斑转换为矢量对象后，图斑即具有形态属性，如面积、长宽比和形状等，这些特征可以帮助区别水系信息和噪声。在转换为矢量的对象中，水系对象往往呈连续的长条形大图斑，而噪声往往呈零散的小图斑，提取出这些矢量对象的形态学属性，从这些属性中挑选出能够区别水系和噪声的特征变量，就可以实现两者的分离。对于提取出的矢量对象，可以用单个形态变量找出合适的阈值分割水系，如果没有能够分离水系和噪声的单个变量，则需要用多个形态特征变量进行多阈值分割或多变量监督分类。

2.4 方法流程

本文提出的水系提取方法的操作流程见图2。其中遥感影像预处理、ICA、图像增强、数学形态学滤波以及监督分类用ENVI5.2软件完成，图斑转矢量及其矢量对象特征提取用eCognition9.0和ArcGIS10.3软件完成。

图2 方法流程

3 结果和讨论

对经过了预处理的研究区ASTER影像数据，以可见光波段和近红外波段(B1～B9)为输入变量，进行ICA运算，共产生9个独立成分分量(图3(a)—(i))，分别记作IC1，IC2，……，IC9。同时，为与常用于特征提取的PCA方法进行对比，将数据利用PCA提取出3个包含水系的主成分分量(图3(j)—(l))，分别记作PC1，PC3和PC4。通过目视解译评价各分量中包含水系信息的情况。

(a) IC1 (b) IC2 (c) IC3

(d) IC4 (e) IC5 (f) IC6

(g) IC7 (h) IC8 (i) IC9

(j) PC1 (k) PC3 (l) PC4

由图3可知，ICA提取的9个独立成分分量中IC7的水系信息最强，所有水系在IC7中都有所显示，水系细节清晰，且其他地物信号的干扰较少。IC7之外的其他分量中干扰信息较多，水系信息较少且都不明显。PCA提取结果显示水系信息分散在多个主成分分量中，而不是主要集中在某一个分量中，且其中非水系信号干扰也较强。

ICA和PCA方法提取水系的结果对比说明，ICA能够将水系信息与其他地物信息分解到不同独立分量中，减少其他地物对水系提取的影响。而PCA不是以分量之间的相互独立为目标的，所以得到的水系信息与其他地物之间独立性差，表现在图上就是在水系分量中存在大量的其他地物信息，与水系信息相互干扰，不利于水系提取。通过与PCA提取结果的对比，可以看出ICA之所以能用于水系提取，关键在于其分解结果是对不同独立信号源的近似估计，独立分量相互之间统计独立，即水系和其他地物之间可视为相互统计独立，因此ICA能够将水系和其他大部分地物分离成相互独立的分量，从而分离出水系。

根据ICA分解结果选择IC7提取水系。IC7中干涸水系集中在图像中较暗的部分，较亮的部分则代表植被。图中水系为高频信息，其余大面积的非水系部分则表现较为连续的低频信息，为背景值。从图中可见，背景值在空间上是有变化的，这种变化会给水系信息的提取带来干扰。除了背景值，还有来自于相似地物的噪声干扰。利用前文提到的背景抑制及去噪方法消除背景及噪声影响，结果见图4。

(d) 滤波图像阈值分割结果 (e) 数学形态学滤波结果 (f) 基于形态特征去噪结果

图4(b)为没有消除背景干扰的情况下通过阈值分割提取的水系，图中大量的背景信息被分割为水系，而一些水系则被分割为背景，在图左上角处最明显。通过基于中值滤波的背景抑制方法来减小背景干扰，具体操作通过ENVI5.2软件执行，经过反复尝试，选择25像元×25像元滤波窗口时背景抑制效果最好(图4(c))。图4(d)为经背景抑制后通过阈值分割出的水系信息，结果显示背景的干扰已大大减少，其中最明显的是图4(d)中左上角的大片背景干扰已被消除，而一些较细小的水系能被正确识别。

去背景干扰之后的分割结果与水系的真实分布符合度很高，但仍然有相当数量的错分像元，噪声像元非常多。通过ENVI5.2软件的数学形态学滤波功能，选择不同尺寸的窗口进行各种数学形态学滤波尝试，发现3像元×3像元窗口下进行先膨胀后腐蚀的操作能最好地消除噪声和保留目标信息。图4(e)为数学形态学滤波结果，从图中可见，虽然通过数学形态学滤波消除了一部分较小的噪声，但仍有很多噪声无法消除。对于这部分噪声，本文通过提取形态特征进行消除，即对经过数学形态学滤波的图像，用eCognition9.0软件进行矢量化，提取出矢量对象，通过矢量对象的形态参数进行噪声分离。以图斑为单位进行矢量化，提取出对象的形态参数，具体包括：面积、长度、宽度、长宽比、形状指数、边界指数、对称度和边界长度。通过观察每个属性的空间变化对于水系分布的响应程度确定出用于分割水系的形态特征集，并经过试验确定出相关变量的阈值，构建水系分类规则，最终分离出水系(图4(f))。从图4(f)可见，噪声大大减少，分割出的水系与实际吻合度很好。

在实际工作中，干涸水系提取更为常用的处理思路是监督分类。通过对比监督分类方法和本文方法的效果来进一步探讨本文方法在水系提取中的表现。通过尝试不同的监督分类方法，本文最终选定效果最好的支持向量机(support vector machine，SVM)算法进行分类和效果比较。SVM算法是一种常用的监督分类方法，在遥感影像分类中应用较多，其原理见文献[28]，本文不再赘述。本文SVM模型采用径向基函数(radial basis function，RBF)作为核函数，通过网格搜索[29]选取最佳模型参数，分别为惩罚系数C=25，RBF核参数γ=0.116。根据目视解译分别选择5 000个水系和非水系的标注样本，进行模型训练和分类。分类结果的总分类精度为93.011%，Kappa系数为0.764 6，输出结果见图5(a)。通过3像元×3像元中值滤波进行分类后处理，结果见图5(b)。将监督分类结果与影像叠加(图5(c))进行对比，发现监督分类结果在影像左上部与实际吻合度较差，这主要是因为左上部一些地层(红色圈内)的光谱特征与水系光谱特征极为相似，不容易区分。图5(d)为本文方法提取结果与影像叠加，对比图5(c)和图5(d)发现监督分类对于较细小的水系识别正确率较差，如图中蓝色圈内。监督分类识别出的部分水系连续性差，而且非水系地物被分为水系的情况较严重(如图5(c)中红色圈内)。综上所述，本文方法整体效果好于监督分类方法，且监督分类方法需要选择训练样本数据，受人主观因素以及模型和参数的影响较大；而本文方法中ICA不需要设置过多参数，也不需要训练样本数据，因此在实际应用中更有优势。

(a) SVM监督分类结果 (b) SVM监督分类结果去噪处理

(c) 处理后的SVM监督分类结果与影像叠加 (d) 本文方法提取结果与影像叠加

4 结论

针对干涸水系与裸地等其他地物光谱特征相似，而无法仅利用光谱信息进行干涸水系提取的问题，本文提出了基于ICA和利用形态特征去噪的干涸水系提取方法。对于ICA分离出的水系信息分量，通过中值滤波降低背景干扰，再通过数学形态学滤波和图斑形态特征消除噪声，最终得到精度较高的水系分布。经过对内蒙古自治区获各琦地区的ASTER数据进行实验，验证了该方法的实际效果。本文研究结果表明：

1)通过ICA可以实现水系和大部分非水系地物的分离，提取的水系信息独立分量中水系信息得到增强和突出，其他地物的干扰被大大减少。

2)通过ICA提取的水系分量中存在少量背景干扰和来自于相似地物的干扰，利用基于中值滤波去背景的方法能够降低背景信息干扰，突出局部水系信息。

3)通过将图斑矢量化提取其形态特征，可以提升噪声分离的效果。

4)本文方法能够很好地用于干涸水系信息的提取，相比于传统的监督分类方法识别精度更高，且无需训练样本数据，操作简便，实用性强。

5)本文提出的水系识别方法识别的水系结果是通过与原图进行目视对比，以及通过与大比例尺谷歌地图进行对比验证的，尚缺少实地验证，这是本文存在的不足之处。