-
图像分割是合成孔径雷达(SAR)图像分析和解译的关键步骤之一。近年来,单极化SAR图像的基于水平集的分割方法是研究热点之一,其优点包括简单的展开式、优良的鲁棒性和准确的边界定位等。
由于SAR图像特有的乘性噪声,基于边界信息的水平集分割方法效果并不理想,因此大多数文献采用区域统计特征对SAR图像进行分割。然而已有的方法需要预先指定区域的统计分布[1-4],当假定概率分布与真实的分布偏离较大时,分割结果准确率会降低。另外,由于不同的区域不能由一种特定的概率分布来描述[5-6],因此对于多区域的分割需求,基于指定统计分布的分割方法效果并不令人满意。
针对该问题,本文提出一种无需先验信息基于参数化模型的水平集多区域分割方法,采用改进Edgeworth展开式对不同的区域进行自适应拟合,将其结果嵌入到对应数量的能量泛函模型中,推导数值求解过程最终实现分割。与基于G0分布的水平集分割方法相比,该方法更适用于多区域的分割,分割更精确。
-
基于参数化模型的水平集分割方法如下:首先根据分割区域个数建立对应数量的统计分布函数,各区域的统计分布函数由改进Edgeworth展开式得到;其次将统计分布嵌入到能量泛函,通过变分原理和梯度下降流原理最小化能量泛函,并采用水平集方法求解,最后实现图像分割。
-
在大多数基于统计信息的水平集SAR图像分割方法中,需要事先根据图像信息假定SAR图像的统计分布。通常认为SAR图像服从Gamma分布、Weibull分布或者G0分布。然而,目前还没有一种分布能完全描述不同的场景和地物,因此对于多区域分割,先假定概率分布不能得到满意的分割结果。本文采用改进的Edgeworth展开式估计不同区域的概率分布,克服了需要先假定概率分布的问题,使分割精度提高[7]。
根据Edgeworth展开式规律可知[8-9],如果一个分布不是完全远离正态分布,那么该分布可以由标准正态分布和Chebyshev-Herimite多项式逼近。由于SAR图像某些区域分布不一定满足该条件,因此本文基于中心极限定理先对分布进行变换,再采用Edgeworh展开式对分布进行拟合,该方法称为改进的Edgeworh展开式。其原理如下:
假设独立随机变量序列服从某种相同的分布,根据中心极限定理可知,若n充分大时,其分布近似服从正态分布:
$$Y = \frac{{\sum\limits_{k = 1}^n {{X_k}} - n\mu }}{{\sqrt n \sigma }}$$ (1) 式中,$\mu $和$\sigma $分别为随机变量的均值和方差。
另外,假设SAR图像的像素值由描述,则可以认为总的像素数n足够大。
因此,根据式(1)用变量Y表示图像像素,则图像的分布可以由下式逼近:
$$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;p(y) = \\ G(y)\left[ {1 + \frac{1}{6}{\rho _3}{C_3}(y) + \frac{1}{{24}}{\rho _4}{C_4}(y) + \frac{1}{{72}}\rho _3^2{C_6}(y)} \right] \end{array}$$ (2) 式中,变量y由式(1)得到,G(y)表示标准正态分布;C3、C4、C6表示三、四、六阶Chebyshev-Herimite多项式;${\rho _3}$、${\rho _4}$表示变量y的三阶和四阶累积量。
C3、C4、C6的表达式为:
$$\left\{ {\begin{array}{*{20}{l}} {{C_3}(y) = {y^3} - 3y}\\ {{C_4}(y) = {y^4} - 6{y^2} + 3}\\ {{C_6}(y) = {y^6} - 15{y^4} + 45{y^2} - 15} \end{array}} \right.$$ (3) ${\rho _3}$、${\rho _4}$的表达式为:
$$\begin{array}{*{20}{c}} {{\rho _3} = \frac{{{O_{X;3}}}}{{\sqrt n {\sigma ^3}}},}&{{\rho _4} = \frac{{{O_{X;4}}}}{{n{\sigma ^4}}} - 3} \end{array}$$ (4) 式中,${O_{X;3}}$、${O_{X;4}}$为变量X的三、四阶累积量:
$$\left\{ \begin{array}{l} {O_{X;3}} = {\mu _{X;3}} - 3{\mu _{X;2}}{\mu _{X;1}} + 2\mu _{X;1}^3\\ {O_{X;4}} = {\mu _{X;4}} - 4{\mu _{X;3}}{\mu _{X;1}} - 3\mu _{X;2}^2 + 12{\mu _{X;2}}\mu _{X;1}^2 - 6\mu _{X;1}^4 \end{array} \right.$$ (5) 式中,${\mu _{X;1}}$、${\mu _{X;2}}$、${\mu _{X;3}}$、${\mu _{X;4}}$为变量X的一、二、三、四阶原点矩。
图 1给出了Edgeworth展开式、G0分布、Weibull分布对不同场景分布曲线的拟合程度。两幅SAR图像为均没有经过滤波的强度图。实验场景为复杂度较低的水域场景和复杂度较高的建筑场景。
本文采用Kullback-Leibler (KL)距离对各自模型估计准确度进行评价[10]。距离数值越小表示拟合效果越好。表 1给出了KL距离值。由表 1可以看出,Weibull分布对均匀场景拟合较好,但是对复杂场景效果较差;相反,G0分布不能很好地描述简单场景的统计分布,更适用于复杂场景。Edgeworth展开式虽然在不同场景下不是最好的拟合方法,但是对不同场景都能得到很好的拟合效果,因此更适用于多区域的分割。
表 1 不同模型下的KL距离
方法 水域 建筑 Edgeworth 0.235 9 0.277 8 G0 0.531 5 0.253 8 Weibull 0.210 1 0.576 5 -
水平集函数的定义为曲线等于0,曲线内部大于0,曲线外部小于0。设水平集个数为L,由于水平集函数通常为符号距离函数,因此L个水平集至多可以得到N=2L个分割区域。为了区分不同区域,定义一个特征函数为:
$${\mathit{\boldsymbol{A}}_i} = \left\{ {\begin{array}{*{20}{c}} 1&{属于第i个区域}\\ 0&其他 \end{array}} \right.i = 1,2, \cdots ,N$$ (6) 为了求解过程的计算稳定性,通常水平集方法用Heaviside函数代替1和0[11],常用的函数定义以及其导数Dirac函数为:
$$\left\{ {\begin{array}{*{20}{l}} {H(z) = \frac{1}{2} + \frac{1}{{\rm{ \mathsf{ π} }}}\arctan (z)}\\ {\delta (z) = \frac{{\partial H}}{{\partial z}} = \frac{1}{{\rm{ \mathsf{ π} }}}\left( {\frac{1}{{1 + {z^2}}}} \right)} \end{array}} \right.$$ (7) 假设水平集函数为${u_k}(k = 1,2, \cdots ,L)$,那么其所有的组合方式可以用向量表示为,$i = 1,2, \cdots ,N$,$c_k^i = 0$或1。从而式(6)可以重写为:
$${\mathit{\boldsymbol{A}}_i} = \prod\limits_{k = 1}^L {{\phi _i}({u_k})} \;\;\;\;\;\;\;i = 1,2, \cdots ,N$$ (8) 式中,,其导数为。
由于水平集函数为符号距离函数,因此在求解过程中每次迭代都需要重新初始化水平集函数,为了避免该问题,引入一个惩罚项从而加快结果的收敛速度[12]。设图像由I(x, y)表示,则根据基于区域信息的传统水平集方法可以推导出嵌入Edgeworth展开式的能量泛函为:
$$\begin{array}{*{20}{l}} {\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{E}}\left( {\bigcup\limits_{k = 1}^L {{u_k}} ;\bigcup\limits_{i = 1}^N {{p_i}} } \right) = - \alpha \sum\limits_{i = 1}^N {\int_\mathit{\Omega } {\ln {p_i}{\mathit{\boldsymbol{A}}_i}{\rm{d}}x{\rm{d}}y} } + }\\ {\beta \sum\limits_{k = 1}^L {\int_\mathit{\Omega } {\delta ({u_k})\left| {\nabla {u_k}} \right|{\rm{d}}x{\rm{d}}y} } + \lambda \sum\limits_{k = 1}^L {\int_\mathit{\Omega } {\frac{1}{2}{{(\left| {\nabla {u_k}} \right| - 1)}^2}{\rm{d}}x{\rm{d}}y} } } \end{array}$$ (9) 式中,pi为由式(2)得到的分割区域概率密度;α、β、λ为权重;Ω为图像域;∇为哈密顿算子;第1项为统计区域信息,各个区域的概率密度函数由式(2)估计得到;第2项为保证曲线的光滑,避免过小的分割区域;第3项为惩罚项,用于调整水平集函数为符号距离函数。
-
根据水平集定义可知,将图像界面嵌入到高一维水平集函数中,在演化过程中,零水平集即为所需的分割曲线。根据变分原则和梯度下降流原理,最小化式(9)即可得到其解为:
$$\begin{array}{*{20}{l}} {\;\;\;\;\;\;\frac{{\partial {u_k}}}{{\partial t}} = \alpha \sum\limits_{i = 1}^N {\ln {p_i}} \left( {\prod\limits_{j = 1,j \ne k}^L {{\phi _i}({u_j})} } \right){\psi _i}({u_k}) + }\\ {\beta \delta ({u_k}){\rm{div}}\left( {\frac{{\nabla {u_k}}}{{\left| {\nabla {u_k}} \right|}}} \right) + \lambda \left[ {\Delta {u_k} - {\rm{div}}\left( {\frac{{\nabla {u_k}}}{{\left| {\nabla {u_k}} \right|}}} \right)} \right]}\\ {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;k = 1,2, \cdots ,L} \end{array}$$ (10) 式中,t为引入的时间变量;Δ为拉普拉斯算子;div为散度。
水平集方法中,通常采用迭代方法对偏微分方程进行求解。因此式(10)的解为:
$$\frac{{{u_{k,s + 1}} - {u_{k,s}}}}{{\Delta t}} = \frac{{\partial {u_k}}}{{\partial t}}$$ (11) 式中,s为迭代次数;Δt为时间步长。
为了保证迭代求解的稳定性,通常需要复杂的差分方案,本文采用文献中半隐式的差分方案,即可获得很好的结果。
在迭代过程中需要每次重新计算各自类别的概率密度函数。由1.1节可知,概率密度函数的确定与变量的原点矩有关。第i类中的j阶原点矩$\hat \mu _{_{I;j}}^i$可以由下式估计:
$$\begin{array}{*{20}{l}} {\;\;\;\;\;\;\;\;\;\hat \mu _{I;j}^i = \frac{{\int_\mathit{\Omega } {I{{(x,y)}^j}{\phi _i}({u_k}){\rm{d}}x{\rm{d}}y} }}{{\int_\mathit{\Omega } {{\phi _i}({u_k}){\rm{d}}x{\rm{d}}y} }}}\\ {i = 1,2, \cdots ,N;j = 1,2,3,4;k = 1,2, \cdots ,L} \end{array}$$ (12) 对于分割正确性的验证,本文采用两个水平集函数,由式(10)可得,两个水平集函数分别为:
$$\begin{array}{*{20}{l}} {\;\;\;\;\;\;\;\;\;\;\frac{{\partial {u_1}}}{{\partial t}} = \alpha \delta ({u_1})[(\ln {p_1} - \ln {p_3})H({u_2}) + }\\ {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;(\ln {p_2} - \ln {p_4})(1 - H({u_2}))] + }\\ {\delta ({u_1})\beta {\rm{div}}\left( {\frac{{\nabla {u_1}}}{{\left| {\nabla {u_1}} \right|}}} \right) + \lambda \left[ {\Delta {u_1} - {\rm{div}}\left( {\frac{{\nabla {u_1}}}{{\left| {\nabla {u_1}} \right|}}} \right)} \right]} \end{array}$$ (13) $$\begin{array}{*{20}{l}} {\;\;\;\;\;\;\;\;\;\;\;\;\frac{{\partial {u_2}}}{{\partial t}} = \alpha \delta ({u_2})[(\ln {p_1} - \ln {p_1})H({u_1}) + }\\ {(\ln {p_3} - \ln {p_4})(1 - H({u_1}))] + \delta ({u_2})\beta {\rm{div}}\left( {\frac{{\nabla {u_1}}}{{\left| {\nabla {u_1}} \right|}}} \right) + }\\ {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\lambda \left[ {\Delta {u_2} - {\rm{div}}\left( {\frac{{\nabla {u_2}}}{{\left| {\nabla {u_2}} \right|}}} \right)} \right]} \end{array}$$ (14) 两个水平集可以至多分割4个区域,本文分别给出了仿真和真实SAR图像的例证。
Multi-Region Segmentation Method for SAR Images Based on a Parametric Model
-
摘要: 提出了一种基于参数化模型的水平集合成孔径雷达(SAR)图像多区域分割方法。该方法采用改进的Edgeworth展开式自适应地对SAR图像统计信息进行拟合。由于无需预先估计SAR图像待分割区域的概率密度函数,因此该方法更适用于多区域分割。该方法根据分割区域数量,将改进的Edgeworth展开式嵌入到对应个数的能量泛函模型中,并给出水平集方法求解过程及数值实现方案,最终实现图像多区域分割。实验结果表明,同其他水平集方法相比,该方法能获得更高的分割精度,更适用于多区域分割。Abstract: In this paper, a multi-region segmentation method for synthetic aperture radar (SAR) images based on a parametric model is proposed. The modified Edgeworth expansion series method is employed to fit the statistical information of the SAR image adaptively. Since the estimation of the probability density function of the SAR image is not required, it is more suited for the multi-region segmentation. Based on the number of the segmentation, the modified Edgeworth expansion series is introduced to the corresponding energy functional model and then the solution based on the level set and its numerical method are developed, thus attaining the segmentation. The experimental results indicate that the higher accuracy of the segmentation is achieved by the proposed method compared with those by the other methods, and the proposed method is more suitable for the multi-region segmentation.
-
表 1 不同模型下的KL距离
方法 水域 建筑 Edgeworth 0.235 9 0.277 8 G0 0.531 5 0.253 8 Weibull 0.210 1 0.576 5 -
[1] FENG J L, CAO Z J, PI Y M. Multiphase SAR image segmentation with G (0)-statistical-model-based active contours[J]. IEEE Transactions on Geoscience and Remote Sensing, 2013, 51(7):4190-4199. doi: 10.1109/TGRS.2012.2227754 [2] GALLAND F, NICOLAS J M, SPORTOUCHE H, et al. Unsupervised synthetic aperture radar image segmentation using fisher distributions[J]. IEEE Transactions on Geoscience and Remote Sensing, 2009, 47(8):2966-2972. doi: 10.1109/TGRS.2009.2014364 [3] AYED I B, HENNANE N, MITICHE A. Unsupervised variational image segmentation/classification using a weibull observation model[J]. IEEE Transactions on Image Processing, 2006, 15(11):3431-3439. doi: 10.1109/TIP.2006.881961 [4] SHUAI Y, SUN H, XU G. SAR image segmentation based on level set with stationary global minimum[J]. Geoscience and Remote Sensing Letters, IEEE, 2008, 5(4):644-648. doi: 10.1109/LGRS.2008.2001768 [5] GAO G. Statistical modeling of SAR images:a survey[J]. Sensors, 2010, 10(1):775-795. doi: 10.3390/s100100775 [6] BIAN Y, MERCER B. SAR probability density function estimation using a generalized form of K-distribution[J]. IEEE Transactions on Aerospace and Electronic Systems, 2015, 51(2):1136-1146. doi: 10.1109/TAES.2014.120825 [7] LUO S, TONG L, CHEN Y, et al. SAR image segmentation based on the advanced level set[J]. IOP Conference Series:Earth and Environmental Science, 2014:012246. http://cn.bing.com/academic/profile?id=2134478556&encoded=0&v=paper_preview&mkt=zh-cn [8] PAPOULIS A, PILLAI S U. Probability, random variables, and stochastic processes[M]. New Delh, India:Tata McGraw-Hill Education, 2002. [9] STUARD A, ORD J. Kendall's advanced theory of statistics[M]. Griffin London:Oxford University Press, 1987. [10] JOHNSON D H, SINANOVIC S. Symmetrizing the kullback-leibler distance[J]. IEEE Transactions on Information Theory, 2003, 9(3):96-99. [11] ZHAO H K, CHAN T, MERRIMAN B, et al. A variational level set approach to multiphase motion[J]. Journal of Computational Physics, 1996, 127(1):179-195. doi: 10.1006/jcph.1996.0167 [12] LI C, XU C, GUI C, et al. Level set evolution without re-initialization:a new variational formulation[C]//Computer Vision and Pattern Recognition, 2005.[S.l.]:IEEE, 2005:430-436. [13] CHAN T F, VESE L A. Active contours without edges[J]. IEEE Transactions on Image Processing, 2001, 10(2):266-277. doi: 10.1109/83.902291