结合巴氏系数和灰度共生矩阵的遥感影像分割

杨军，王恒亮

(1.兰州交通大学 电子与信息工程学院，兰州 730070；2.兰州交通大学 测绘与地理信息学院，兰州 730070；3.甘肃省地理国情监测工程实验室，兰州 730070)

0 引言

目前，伴随着遥感技术的飞速进步以及传感器技术的不断发展，遥感影像作为重要数据源已被广泛应用于城市规划、精细农业、林业测量、环境监测、军事目标识别和灾害评估等诸多领域，而图像分割作为遥感影像信息提取和面向对象分类的重要方法，已成为遥感影像处理领域的研究热点和难点之一[1-5]。目前，常见的图像分割算法有：基于阈值的图像分割方法、基于区域的图像分割方法和基于边缘的图像分割算法[6]。

分水岭变换是一种经典的图像分割方法，属于基于区域的分割方法，具有运算简单，边缘定位精度高的优点，且可以得到连通的、封闭的分割轮廓[7]。近年来，国内外研究者针对传统分水岭变换进行了广泛研究，提出了许多改进方法。高丽等[8]提出了一种基于标记的分水岭改进分割算法，该算法使用了一种新的标记提取方法，从原始梯度图像提取出与物体相关的局部极小值，构成标记图像，然后在标记的梯度图像上进行分水岭分割。该算法可以在一定程度上减少普通图像的过分割问题，然而，对于分辨率较高的遥感影像，过分割现象依然比较明显。邵龙等[9]提出了一种结合异质性及分水岭变换的遥感影像分割方法，利用非线性同组滤波(peer group fittering，PGF)去除影像噪声，根据地物的异质性特征，对分水岭分割结果进行合并。该算法分割效果良好，能够用于面向对象的信息提取，但是分割边界比较尖锐，平滑度较差。蔡彩等[10]针对城市高分辨率遥感影像，提出了一种基于分水岭变换和分层区域合并的分割方法。该方法首先使用多通道分水岭方法对影像进行分割，然后分析不同地物光谱变化性的特点，并利用分层区域合并来改进分割结果。Ijitona等人[11]在对SAR海冰影像浮冰尺寸检索的研究中，提出了一种基于分水岭变换和区域强度合并的分割方法。首先采用分水岭变换处理图像，然后基于区域邻接图执行区域强度合并。以上2种基于区域合并的方法分别针对城市影像和SAR影像时分割效果较好，但是对以其他地物为主的遥感影像进行分割时具有一定的局限性。Sahin等[12]提出了一种基于分水岭的高分影像多尺度自动分割算法，首先采用边缘保持平滑滤波器(edge preserving smoothing filter，EPSF)对影像进行平滑处理，然后使用浸入分水岭分割影像，最后进行基于多尺度分割的区域合并优化分割结果，可以有效抑制影像中的噪声，然而却损失了影像的部分纹理信息，对于纹理丰富的区域，分割结果不太理想。针对遥感影像分水岭分割中存在的问题，本文提出了一种基于巴氏系数(Bhattacharyya coefficient，BC)和灰度共生矩阵(gray level co-occurrence matrix，GLCM)的区域合并算法，即BC-GLCM算法，合并分水岭分割后过分割图像中存在的相似区域，改进分水岭分割算法。同时采用遥感影像数据进行实验，以验证本文算法的有效性。

1 分割方法及原理

分水岭变换是一种基于拓扑理论的数学形态学的分割方法，其基本思想是把图像看作是测地学上的拓扑地貌。分水岭变换算法可分为二类，一类是自下而上的模拟泛洪算法，另一类是自上而下的模拟降水算法[13]。二类分水岭算法都是将图像看成立体地形图，分别从地形图的底部和顶部开始算法执行，并按照各自原理执行分水岭变换。根据相关研究，二类算法可以得到大致相同的分割效果。

本文基于巴氏系数和灰度共生矩阵的区域合并方法来改进分水岭算法的分割结果。具体过程为：首先，利用数学形态学的方法提取原始遥感影像的梯度图像；其次，在梯度图像的频率域中，提取图像的低频成分，把提取到的低频成分作为局部极小值标记到梯度图像中，从而得到标记的梯度图像；然后，使用分水岭变换处理标记的梯度图像，并得到初始分割结果；最后，针对初始分割图像中由于过分割出现的同质区域，采用BC-GLCM区域合并方法进行合并。详细流程图如图1所示。

图1 本文方法分割流程图

1.1 提取梯度图像

梯度图像是遥感影像相邻像素的灰度变化情况的反映，相比于原始图像，梯度图像可以更好地表示图像的变化趋势[14]，因此使用梯度图像作为分水岭变换的分割图像。本研究采用数学形态学的方法来提取原始遥感图像的梯度图像。形态学梯度图像的计算过程如下：

假定f(x，y)是待处理遥感影像的灰度图像，b是结构元素，Df表示f的定义域，Db表示b的定义域。

定义1用b对f(x，y)进行膨胀(dilate)，记作f⊕b：

f⊕b=dilate[f(x,y),b]=
max{f(x-x0,y-y0)+b(x0,y0)|
(x-x0,y-y0)∈Df;(x0,y0)∈Db}

(1)

定义2用b对f(x，y)进行腐蚀(erode)，记作f⊙b：

f⊙b=erode[f(x,y),b]=
min{f(x+x0,y+y0)-b(x0,y0)|
(x+x0,y+y0)∈Df;(x0,y0)∈Db)}

(2)

定义3图像f(x，y)形态学梯度图像g(x，y)定义为：

g(x,y)=(f⊕b)-(f⊙b)

(3)

1.2 标记图像

由于遥感影像噪声及地物内部细密纹理的影响，使梯度图像存在大量局部极小值区域，其中一部分局部极小值区域与感兴趣目标无关，称为伪局部极小值区域，会导致感兴趣目标被分割为许多无意义的细小区域[15]。如果直接对梯度图像进行分水岭变换，所得的结果会产生许多零碎区域，造成极其严重的过分割结果，增大了后续区域合并的复杂度。因此，为了减小伪局部极小值区域对分割效果的影响，需要标记梯度图像内部同质性较好的局部极小值区域，把伪局部极小值区域给剔除掉。具体方法是：首先对空间域的梯度图像进行傅里叶变换得到频率域的梯度图像，其次使用巴特沃斯(Butterworth)低通滤波器提取梯度图像的低频成份，然后进行傅里叶逆变换，获取梯度图像的同质区域，最后采用h-极小值变换提取与地物相关的极小值区域，得到二值标记图像。标记图像流程图如图2所示。

图2 标记图像流程图

1.3 遥感影像的初始分割

由于分水岭算法对遥感影像分割具有运算简单、边缘定位精度高的优点，故采用分水岭算法对影像进行初始分割。在1.2节得到只有局部极小值区域的二值标记图像后，根据标记图像对梯度图像进行形态学极小值标定，使局部极小值区域仅出现在被标记的位置，这样就获得了标记的梯度图像。在标记的梯度图像上进行分水岭变换，可以减少分割结果的过分割现象。本研究初始分割采用文献[16]的基于地形学距离的分水岭算法，此算法属于自上而下的模拟降水的分水岭变换。

1.4 BC-GLCM区域合并算法

在使用分水岭变换对遥感影像进行初始分割后，由于影像地物类别复杂，噪声和地物纹理细节等因素，使分割结果往往还存在着过分割，即把同种地物分割成许多零碎区域。为了减少过分割，得到有意义的分割结果，本文提出了一种基于BC和GLCM的区域合并新方法，即BC-GLCM区域合并算法，对分水岭变换后的初始分割结果进行区域合并。

1)巴氏系数。巴氏系数是对2个统计样本的重叠量的近似计算，可以测量2个离散概率分布之间的相似度[17]。因此，本文引入巴氏系数计算初始分割图像中2个区域的灰度概率分布的相似性。假定区域A和区域B为初始分割图像中的2个相邻区域，令pA(x)和pB(x)分别为A和B的像素灰度概率分布，则A和B的巴氏系数BC(pA，pB)可定义为：

(4)

式中：L为图像的灰度级；x=0，1，…，L-1。

BC(pA，pB)的取值为0到1，其值越接近于1，说明2个相邻区域的区域相似度越大，像素灰度概率分布越相似，则2个区域属于同一个地物的可能性就越大。

巴氏系数作为衡量图像区域相似度的一个指标，可以反映区域灰度的相似性，但它仅仅体现了区域像素灰度的概率分布，对于遥感影像来说，地物的纹理信息比较丰富，如果仅仅使用巴氏系数进行区域相似度的判断，会造成区域相似度在纹理方面的缺失，出现一些错误的合并。而灰度共生矩阵作为一种描述纹理的有效方法，可以很好地反映区域的纹理特征。

2)灰度共生矩阵。灰度共生矩阵是分析图像纹理的一种统计方法，在一定程度上反映了图像中各个灰度级在空间上的分布特性，是纹理分析领域中经常采用的特征之一[18]。基于灰度共生矩阵的纹理提取技术作为目前应用范围较广的纹理提取方法，在遥感图像的纹理分析中发挥着十分重要的作用[19]。

灰度共生矩阵Pd(i，j)(i，j=0，1，2，…，L-1)被定义为灰度值为i的像素离开某个固定位置关系d到灰度值为j的像素的概率。其中L是图像的灰度级，i，j分别是像素的灰度，d表示2个像素间的空间位置关系。当2个像素间的位置关系d确定后，就生成一定关系d下的灰度共生矩阵[20]。在灰度共生矩阵中，某个元素代表的是一种灰度组合在一定空间关系下出现的概率。

根据灰度共生矩阵，可以定义出14个纹理特征值，相关研究分析发现只有角二阶矩(angular second moment，ASM)、熵(entropy，ENT)、对比度(contrast，CON)和相关性(correlation，COR)4个纹理特征与人类视觉感知特性有明确的对应关系[21]，因此本研究采用这4个特征量来描述纹理信息。角二阶矩反映图像灰度分布均匀程度和纹理粗细程度。熵反映纹理的复杂程度或非均匀程度。对比度反映了图像的清晰度和纹理沟纹的深浅程度。相关性反映的是纹理的方向性。

本研究利用角二阶矩、熵、对比度和相关性组成一个四维特征向量(ASM、ENT、CON、COR)综合描述图像中区域的纹理属性。欧式距离(Euclidean distance，ED)作为一种相似度计算方法，可以用来描述2个相邻区域纹理特征向量的相似度。设(ASM1、ENT1、CON1、COR1)，(ASM2、ENT2、CON2、COR2)为2个相邻区域的纹理特征归一化向量，其欧式距离的计算公式为：

(5)

式中：M=ASM1-ASM2；T=ENT1-ENT2；N=CON1-CON2；R=COR1-COR2。

2个相邻区域纹理特征向量的欧式距离越接近0，则说明2个区域的纹理属性越相似，那么这2个区域属于同一个地物的可能性就越大。

3)相似区域的合并。基于分水岭变换的图像分割方法往往会得到过分割结果，导致分割区域比较零碎，所以对初始分割结果进行合并是一种常用的改进方法。过分割图像中的相邻区域能否进行合并，取决于2个相邻区域是否属于同一个地物。同一个地物的不同部分具有一定的相似性，因此可以通过计算2个相邻区域的相似度来判断是否属于同一地物，进而确定是否合并。由于遥感影像中地物的灰度信息和纹理信息比较明显，所以本文采用巴氏系数来描述灰度概率分布信息的相似性，采用灰度共生矩阵来描述纹理信息的相似性，结合2种相似性计算结果对过分割的图像区域进行合并。这种区域合并方法充分考虑了相邻区域的灰度信息和纹理信息，能够得到较为准确的区域合并结果。

由1)和2)可知，灰度概率分布相似和纹理信息相似的合并阈值都在0～1之间，但是由于不同的遥感影像所包含的地物类型不一样，且地物特征差别较大，因此对不同的遥感影像，需要经过相应的实验来确定合适的灰度概率分布相似性合并阈值g和纹理特征相似性合并阈值t，不同影像的具体合并阈值g和t等相关参数的设置将在下文的实验结果与分析中进行讨论。

BC-GLCM区域合并方法具体的过程如图3所示。假定区域A和区域B为初始分割图像的相邻区域，首先计算2个区域的巴氏系数(BC)，判断灰度概率分布信息的相似情况，如果BC(pA，pB)t，则纹理特征信息不相似，A和B不合并。

图3 区域合并过程

2 实验结果与分析

为了验证本文分割算法，对比分析了本算法和与基于地形学距离的传统分水岭变换[16]和基于标记的分水岭图像分割算法[8]。

2.1 实验一

选取GeoEye-1高分辨率卫星影像进行实验，如图4(a)所示，为美国迈诺特市洪水时期所拍摄遥感影像的某一部分，拍摄时间为2011年6月，分辨率为0.5 m。可以看出，该影像地物纹理特征明显，轮廓清晰可见，研究区的主要地物类型包括水体、植被、房屋。

图4(b)为基于地形学距离的传统分水岭变换的分割结果，红色曲线为分割线。传统分水岭分割算法对于普通图像有一定的分割效果，但是对高分辨率遥感影像分割效果较差，过分割现象特别严重，基本无法分割出有意义的地物目标对象。图4(c)为基于标记的分水岭图像分割算法的分割结果，相比于传统分水岭分割算法，该算法在一定程度上减少了过分割现象，能够分割出不同的地物，边界较为准确，但是在同一地物内部还是存在着一定的过分割现象，如图4(c)中黄色椭圆标记所示。图4(d)为本文基于BC-GLCM区域合并的分水岭分割算法的结果，经过实验分析，设置区域合并过程中的灰度概率分布相似性合并阈值g为0.95，纹理特征相似性合并阈值t为0.30，在此合并阈值下，分割结果较好。可以看出，与传统的分水岭分割算法和基于标记的分水岭分割算法相比，本算法有效地抑制了过分割现象。如图4(c)中黄色标记的过分割区域在图4(d)中得到了合并，水体、植被、房屋等不同地物都可以清晰地分割出来，并且边界较为准确清晰，同一地物的过分割现象也极少出现，能够得到较为满意的分割结果。

图4 不同算法的分割结果对比

表1所示为3种方法分割结果所得到的分割区域个数，从定量的角度验证本文算法在抑制分水岭算法过分割问题的有效性。本文算法相比于传统分水岭算法，分割结果中区域个数减少了95.8%；相比于基于标记的分水岭算法，分割结果中区域个数减少了51.9%，表明本文算法可以极大降低遥感影像分水岭分割算法的过分割现象，不同地物的分割边界轮廓清晰，分割结果较为准确。

表1 不同算法的分割结果中分割区域的个数

2.2 实验二

选取WorldView-2高分辨率遥感卫星所拍摄的影像进行分割实验，如图5(a)所示，为我国第二艘航母出坞下水时的局部遥感影像，拍摄时间为2017年4月26日，分辨率为0.5 m。该影像包括的主要地物有水域、陆地、航母等。

如图5(b)所示，使用传统分水岭算法对WorldView-2遥感影像进行分割时，过分割现象特别严重，产生许多零碎区域，基本不能分割出航母、水域和陆地。而使用基于标记的分水岭算法相比传统分水岭算法的分割结果能够在一定程度上减少图像的过分割现象，但是还是会不可避免地出现一些过分割区域，如图5(c)中黄色椭圆标记区域所示。实验中，算法设置的灰度概率分布相似性合并阈值g为0.85，纹理特征相似性合并阈值t为0.65，极大地减少分水岭算法的过分割现象，如图5(d)所示，大量的同质区域被合并，航母、水域和陆地都能够清晰地分割出来，能得到较为准确的分割结果。

图5 不同算法的分割结果对比

表2所示为3种方法分割结果所得到的分割区域个数。本文算法相比于传统分水岭方法，分割结果中区域个数减少了99.6%；相比于基于标记的分水岭算法，分割结果中区域个数减少了85.5%，表明本文算法可以极大降低遥感影像分水岭分割算法的过分割现象，所分割出的地物边界轮廓清晰，分割结果较为准确。

表2 不同算法的分割结果中分割区域的个数

2.3 实验三

选取WorldView-2高分辨率遥感卫星所拍摄的影像进行实验，如图6(a)，为江西省新余某地区影像的某部分，拍摄时间为2013年11月。影像由多光谱和全色影像融合而成，分辨率可达0.5 m，主要地物有林地、居民区、厂房、农田、道路等。

如图6(b)，为传统分水岭算法的分割结果，过分割现象极其严重，基本不能分割出有意义的地物对象，有大量零碎同质区域。图6(c)为基于标记的分水岭算法分割结果，相比于传统分水岭算法，能够在一定程度上减少过分割，分割出显著地物，但还存在着同质区域被分割的现象，如图黄色椭圆标记所示。图6(d)为本文算法在灰度概率分布相似性合并阈值g为0.87，纹理特征相似性合并阈值t为0.30时获得的分割结果，可以看出，无论相比于传统分水岭算法，还是相比基于标记的分水岭算法分割结果，本算法都可以极大地减少过分割现象，大量的同质区域被合并，地物对象能够有效地分割出来。

图6 不同算法的分割结果对比

表3所示为3种方法分割结果所得到的分割区域个数。本文算法和传统分水岭算法相比，分割结果中区域个数减少了98.3%；和基于标记的分水岭算法相比，分割结果中区域个数减少了73.1%，表明本文算法可以极大降低遥感影像分水岭分割算法的过分割现象，能够得到比较好的分割结果。

表3 不同算法的分割结果中分割区域的个数

2.4 实验四

选取资源三号卫星所拍摄的影像进行实验，如图7(a)所示，为三亚某地区局部融合影像，成像时间为2012年3月13日，分辨率为2.1 m。该影像比较复杂，地物主要有海域、林地、居民区、草地、裸地、海滩等。如图7(b)所示，为传统分水岭算法分割结果，可以看出过分割极其严重，不能分割出有意义的地物对象。图7(c)所示为使用基于标记的分水岭算法所得的分割结果，相比于传统分水岭算法，过分割现象得到了一定的缓解，但还是在地物内部出现一些过分割区域，如图中黄色椭圆标记所示。图7(d)所示为本文算法的分割结果，合并阈值g取0.88，t取0.28时，分割结果较好，相比于传统分水岭算法和基于标记的分水岭算法，图7(c)中的黄色椭圆标记区域都得到了合并。因此，本文算法能够极大地解决过分割问题，海域、林地、居民区、草地、裸地、海滩等地物都能清晰地分割出来，分割结果较为满意。

图7 不同算法的分割结果对比

表4所示为3种方法分割所得的区域个数。本文算法相比传统分水岭算法，分割结果区域个数减少了97.2%；相比基于标记的分水岭算法，分割区域个数减少了50.9%，表明本文算法可以解决分水岭变换用于遥感影像分割出现的过分割问题，得到较为准确的分割结果。

表4 不同算法的分割结果中分割区域的个数

3 结束语

本文针对分水岭算法应用于遥感影像分割容易出现过分割的问题，提出了一种基于BC-GLCM区域合并的分水岭分割算法。该算法使用分水岭分割结果作为初始分割图像，可以得到准确、连通、封闭的分割边界。此外，BC-GLCM区域合并方法结合遥感影像的灰度和纹理信息，可以较为准确地合并过分割中的同质区域。实验结果表明，本文算法相比于传统分水岭算法和基于标记的分水岭算法，获得了较好的分割效果。然而，对于地物类别比较复杂和地物特征比较相似的遥感影像，利用本文算法进行分割时，可能会产生一些错分割的问题，这也是未来进一步研究的课题。