赵 园,潘 斌,邰建豪
(武汉大学 遥感信息工程学院,湖北 武汉 430079)
SAR图像的分割是实现SAR图像目标识别的关键步骤,所以,快速、准确地分割图像,是SAR图像解译的关键技术之一。但是,由于SAR图像特有的斑点噪声的影响,使SAR图像目标的分割问题较为复杂,是SAR图像分割研究中的难点之一[1]。
基于图论的图像分割方法的核心是对图像建模,生成描述图的能量函数并使其最小化,再通过能量最小化,求取图像的最大流进行图像分割。常见的算法有归一化分割(Normalized cuts)、随机游走(Random work)、图割(Graph cuts)、迭代图割(Grab Cut)等。GREID等提出了Graph Cuts算法,该算法指定背景区域,通过直方图建模生成能量函数,完成对图像的分割[2-3],交互简洁,处理速度快,但仅限于对灰度图像的分割;Boykov提出了GrabCut算法,该算法利用高斯混合模型(GMM)代替灰度直方图,对前景和背景区域建立GMM模型,构造能量函数,并使用迭代的Graph Cuts算法完成图像的分割[4]。迭代求解使得分割结果更为准确,但同时也会减缓处理速度。
针对Grab Cut算法,目前已有不少改进。在提升算法效率和准确度上,Tomoyuki Nagahashi 等在 Grab cut 方法上引入多尺度高斯平滑;翟玲提出一种基于超像素与特征改进的 Grab cut 算法[5];有学者提出一种结合分水岭和Grab cut基于前景形态的多尺度图像分割方法,以少量的像素点迭代估计GMM参数,提高算法稳定性与普适性[6]。在减少人工交互操作上,D.Khattab以聚类算法为基础,提出了Automatic Grab cut算法的目标分割[7];张林对整幅图像进行建模,利用全局最大流邻域生长算法进行SAR目标的分割等[8]。
本文结合SAR图像的乘性噪声干扰、目标较小等特点,从提高算法准确度和减少人工交互操作出发,通过预分割的方式,利用FCM算法提供两组准确度较高的GMM参数初始值,迭代求得能量函数的最小值,得到含目标、阴影、噪声点的分割图;对分割图做滤波处理后,再结合二维熵滤去孤立的噪声点和阴影,最终分割出目标区域。
Grab cut是在Graph cuts的基础上改进得来的,该算法利用高斯混合模型(GMM)取代原有的直方图,并通过多次迭代求得能量函数的最小值,精简了交互过程,使图像分割准确度更高[10]。其核心算法:
1)矩形框标定。对于待分割图像Z=(z1,z2,…,zn),如图1所示,人工画出一个包含目标的矩形框。框内表示未知区域,即“可能为目标的像素”,初始化时,该区域一般包含了目标和一部分背景区域,标记为TU;框外表示背景区域,标记为TB;前景区域标记为TF,即“确定为目标的像素”,初始化时TF=Ø。图像像素点集合α=(α1,α2,…,αn),α∈(0,1),0表示TB区域,1表示TU区域。
图1 图像区域示意图
2)建立GMM模型及能量函数。由式1可以看出,Grab cut分割过程即为能量最小化的过程:
α=argminEα.
(1)
该算法的能量函数E表示:
E(α,θ,k,Z)=U(α,θ,k,Z)+V(α,Z),
(2)
(3)
其中,U表示能量函数的数据项,即某个像素属于目标或者背景的概率的负对数;V表示能量函数的边界项;D为表示GMM的函数,θ为选用的GMM模型参数,α为像素的标记,k表示GMM个数。用来表示高斯混合模型的函数D包含多个高斯分量,每一个高斯分量又都是相互独立的,因此,引入一组新的向量K={k1,k2,k3,…,kN},k=1,2,3,…,N,用来标记像素Zn属于第Kn个高斯分量。高斯混合模型的函数定义为
D(αn,kn,θ,Zn)=-logP(Zn|αn,kn,θ,Zn)-
logω(αn,kn).
(4)
化简后得到
D(αn,kn,θ,Zn)=-logω(αn,kn)+
μ(αn,kn)]T∑(αn,kn)-1[Zn-μ(αn,kn)].
(5)
其中,GMM的参数θ有三个:ω(αn,kn)为GMM的权重,μ(αn,kn)为GMM的均值向量,∑(αn,kn)为GMM的协方差矩阵。一旦确定了这三个参数,带入目标的GMM和背景的GMM,确定能量函数的数据项。故GMM的参数可表示为
θ={ω(α,k),μ(α,k),∑(α,k)}.
(6)
同时,对于边界项V:
V(α,Z)=γ∑(m,n)∈cexp(-β‖Zm-
Zn‖2)×δ(Zm,Zn).
(7)
边界项主要体现邻域像素之间不连续的惩罚,像素之间的相似性用欧式距离来表达,两邻域像素差别越大,能量越小。γ,β是常数,Rother[10]等人在实验中取得γ=50,β根据图像对比度来确定。
3)迭代能量函数最小化。首先利用人工标记的TU和TB两个像素集合初始化GMM模型,建立能量函数,并通过最大流-最小割算法进行初始分割;根据初始分割结果更新前景区域TF和背景区域TB集合,更改前景和背景区域对应的GMM参数θ,并通过最大流算法[11-12],使得能量函数E(α,k,θ,Z)越来越小,重复此过程,直到迭代终止。
Grab cut算法主要缺点有两个:一是需要进行人工交互操作,不够方便;二是算法稳定性不够高,GMM初始化参数的好坏会直接影响到分割结果[8]。Grab cut进行图像分割时,利用背景集合和前景集合初始化了两个GMM,相当于将图像分为了两类,因此利用无监督的聚类方法对图像自动预分割,并利用统计结果代替原有的GMM初始化参数,既可以将算法改为自动化,又可以得到更为准确的GMM初始化参数。
GMM的参数估计由EM算法得到,GMM的初始化依赖于EM算法的初始化[13-14],则该问题转化为EM算法的初始值选择。针对SAR图像含噪的特点,考虑用FCM算法代替原有的随机参数的方法来选取EM算法的初始值,即利用FCM算法将SAR图像粗略的划分为两类TU和TB,根据统计结果得到两组参数作为EM算法的初始值。
SAR图像中的乘性噪声模糊了相邻区域的对比度,使得相邻区域的跃变趋于平缓,不同地物边缘处出现过渡带,难以确定边缘的准确位置,因此SAR图像的地物边缘具有模糊性。FCM 算法属于模糊聚类的一种,该算法利用隶属度来表示数据的真实分布,尤其适用于处理各类间有重叠的情况。它的目标函数:
(8)
(9)
其中,U为uki组成的矩阵集合,uki表示某个像素的隶属度,其和为1;V为的集合,vk表示第k个聚类中心;m=2,表示模糊指数,c表示聚类数,‖·‖表示欧氏距离。
由式(6)可知GMM的参数θ是由μ,∑,w组成的。在改进的算法中,采用FCM算法将特征向量归为k类中,并将各类的方差和均值作为∑和μ,w是各个类中所含的特征向量百分比。
由于SAR目标图像中目标较小,目标的阴影较大且含噪,因此在使用本文改进的Grab cut算法进行分割后,图像中不仅包含了目标,还有目标的阴影及斑点噪声,如图2(b)所示。
图2 原图、含阴影的分割图及灰度直方图对比
为了将SAR目标图像中的目标分割出来,把阴影和噪声滤去,根据图2中显示的灰度信息,阈值方法是比较简单实用的方法。二维熵法同时考虑了像素灰度分布信息和其邻域空间相关信息,具有较强的抗干扰能力[15],且二维熵法可以产生多个阈值,有利于后续进行SAR图像多目标分割的研究。故通过建立目标及背景区域的概率密度函数迭代求出各自的二维熵,选择合适的阈值,分割出SAR图像中的目标。
1)改进的Grab cut分割。通过FCM聚类将原图像分为两组数据,初始化两个更准确的GMM参数;建立能量函数,迭代更新至收敛阈值,得到目标、阴影、噪声分割图。
2)采用中值滤波,消除图像中孤立的噪声点,减少分割干扰,考虑到实验所用SAR图像尺寸,窗口大小选择3×3;
3)利用二维熵算法确定分割阈值,滤除阴影及噪声,提取出目标区域;
4)对分割后的图像进行二值化处理,完成SAR目标的分割。
本文实验数据采用Sandia 实验室 MSTAR 数据库的 SAR 图像,MSTAR成像参数为 X波段、单极化(HH),分辨率为 0.3 m×0.3 m,大小为128×128,包含BMP-2、BTR-70、T-72 3类目标。
根据本文算法,分别对数据库中的3类目标进行实验,结果如图3、图4、图5所示。
由图3(b)、图4(b)、图5(b)可以看出,利用本文改进的自动Grab cut算法进行SAR图像分割后,可以比较明显的分割出目标及阴影区域,且目标及阴影的外部轮廓都比较完整,证明该改进算法的可行性;但是还存在两个问题:一是有残余的SAR目标图像的斑点噪声,尤其是当待分割的原图中噪声较为严重时,结果会更明显,这点可以从图4(a)、图4(b)看出;二是目标与阴影相接区域的灰度值变化不明显,轮廓显得较模糊。
图3 BMP-2数据分割算法结果
图4 T-72数据分割算法结果
图5 BTR-70数据分割算法结果
对分割图进行滤波可进一步抑制噪声及杂波,若是直接对原图进行滤波,可能会丢失一部分边缘信息,导致分割结果不准确。由图3(c)、图4(c)、图5(c)可以看出,中值滤波对于脉冲噪声具有较好的抑制作用,图中的孤立噪声明显减少,且图像的轮廓边缘变形较小。
由图3(d)、图4(d)、图5(d)可以看出,经过二维熵算法分割后,SAR目标图像中的噪声点进一步减少,这是因为该算法利用迭代技术,抗干扰能力也有所提高;另外,通过设定合适的阈值,将阴影部分滤去,只保留了目标区域,实现了SAR目标图像的分割。通过观察可以看出,分割后的目标轮廓较为完整,但是目标内部有部分缺失,这是由于目标内部灰度分布不均匀,导致部分区域随着阴影被一起滤除掉。
为了更加客观地评价本文算法,此处选用错误概率(probability of error)和无监督的两种评价指标,分别将改进后的Grab cut算法与改进前的Grab cut算法进行比较;同时将本文所提出的经过后处理的改进算法与传统的ICM_MRF分割方法、双参数CFAR分割方法进行对比。
错误概率[16](probability of error)属于有监督评价指标的一种。其定义:
P(error)=P(O)P(B|O)+P(B)P(O|B).
(10)
其中,P(R)表示由人工分割图计算得到的区域R的先验分割概率,P(R1|R2)表示把区域R2的像素划分到区域R1的概率。错误概率常用于二类分割算法的评价,而本文所用图像都为单目标图像,因此该评价方法具有适用性。
区域内部均匀性度量[17]属于无监督评价的一种,不需要先验知识和标准参考分割图,其定义:
(11)
该方法利用同类内部均匀性的程度来表示图像的分割质量,即PP值越接近1说明图像分割质量越高。
从MSTAR数据库中随机挑选50张图像,评价指标取结果均值,结果如表1所示。
表1 改进前后Grabcut算法的性能指标比较
由表1可以看出,使用改进后Grab cut分割出的图像质量更好,但是准确率却低于改进前的Grab cut 分割方法。主要原因是人工交互可以准确的识别简单目标,而利用函数模型表达分割问题则会有一定的出错率。
由表2可以看出,使用本文方法得到的分割结果在区域内部均匀性度量和错误概率两项指标上均优于传统的ICM-MRF和双参数CFAR分割方法,具有更高的精度和一定的适用性。
表2 不同分割算法的性能指标比较
本文通过FCM聚类进行初始化分割,得到两组质量较好的GMM参数,把传统需要交互的Grab cut算法自动化,得到了初步的分割图。同时为了滤去SAR目标图像中的阴影及噪声,引入二维熵算法,通过选定合适的阈值,最终分割出SAR图像中的目标,改善了SAR图像目标分割方法对于噪声敏感且需要人工交互的问题。根据实验结果,本文方法取得的效果良好,有一定的适用性。
本文使用的实验图像都为单一目标SAR图像,除了特有的噪声外,SAR图像背景也较简单,没有其他地物等的干扰,因此在利用聚类算法初始化时能取得较好的预分割效果,后续进行阈值分割时也比较容易;但是在多目标图像或者背景分布比较复杂的SAR图像上时,算法是否适用以及该怎样改进,是后续研究重点之一。另外, Grab cut算法中进行能量函数最小化的过程以及二维熵算法都使用了迭代技术,增加了计算量,因此,在进行大场景数据分割时,还需要做出下一步的改进提高算法的效率。
[1] 焦李成, 张向荣, 侯彪, 等. 智能SAR图像处理与解译[M]. 北京:科学出版社,2008:21-27.
[2] GREIG D,PORTEOUS B,SEHEULT A. Exact maximum a posteriori estimation for binary images[J]. Journal of the Royal Society, Series B,1989,51(2):271-279.
[3] SZELISKI R, ZABIN R, SCHARSTEIN D, et al. A Comparative Study of Energy Minimization Methods for Markov Random Fields with Smoothness Based Priors [J]. IEEE Transaction on Pattern Analysis and Machin Intelligence(S0162-8828),2008,30(6): 1068-1080.
[4] BOYKOV Y,KOLOMOGOROY V. An experimental comparison of min-cut/maxflow algorithms for energy Minimization in Computer Vision and Pattern Recognition[M]. Springer Verlag, 2001:359-374.
[5] 翟玲, 朱敏, 戴李君. 基于超像素与特征改进的Grab cut前景分割[J].微型电脑应用,2015,31(11):48-50.
[6] 伊聪聪, 吴斌, 张红英. 一种改进的Grab cut 图像分割方法[J].小型微型计算机系统,2014(35):1164-1168.
[7] KHATTAB D, EBIED H M, HUSSIEN A S,et al. Automatic Grab Cut based on unsupervised clustering for image segmentation[J]. Intelligent Systems and Computing, 2015( 323):579-592.
[8] 张林, 朱兆达. 基于全局Maxflow领域生长算法的SAR图像目标分割[J]. 南京航空航天大学学报, 2010(12):15-16.
[9] BLAKE A, ROTHER C, BROWN M, et al.Interactive Image Segmentation using an adaptive GMMRF model[C]. In Proc. European Conf. Computer Vision,2004, 30(21):428-441.
[10] ROTHER C, KOLMOGOROV V, BLAKE A. Grabcut: Interactive foreground extraction using iterated graph cuts[C]. Transactions on Graphics (TOG). ACM, 2004.23(3): 309-314.
[11] 陈国良, 孙广中, 徐云, 等. 并行算法研究方法学[J]. 计算机学报,2008, 31(9):1493-1502.
[12] YUAN J,BAE E,TAI X C. A study on continuous max-flow and min-cut approaches[C].Computer Vision and Pattern Recognition (CVPR), 2010 IEEE Conference on. IEEE, 2010: 2217-2224.
[13] 王鑫.基于高斯混合模型的聚类算法及其在图像分割中的应用[D]. 太原:中北大学,2013.
[14] WU W Y, WANG M J. Elliptical object detection by using geometric properties[J]. Pattern Recognition. 1993. 26(10):1499-1509.
[15] YIMIT A. A New Two-Dimensional Direction Histogram Based Entropic Thresholding[C].International Conference on Image & Graphics, 2011:106-110.
[16] 张晓敏. 一类马氏样本下假设检验问题错误概率的估计[J]. 应用数学, 2008,21(1):179-184.
[17] ZELNIO E G, GARBER F D. Algorithms for Synthetic Aperture Radar Imagery X[J]. Society of Photo-optical Instrumentation Engineers Conference Series, 2012:8394.