乳腺肿瘤超声图像的反应扩散水平集分割

陈 帅 张 麒

(上海大学通信与信息工程学院，上海 200072)

乳腺肿瘤超声图像的反应扩散水平集分割

陈 帅 张 麒

(上海大学通信与信息工程学院，上海 200072)

针对乳腺肿瘤超声图像分割，提出一种改进的反应扩散(RD)水平集分割算法。先使用Gabor各向异性扩散模型进行滤波，由此构造边界停止函数；再将该函数融入RD水平集演化方程，以控制曲线的演化得到乳腺肿瘤的边界。采用该方法和传统RD方法对77例病人的111幅乳腺超声图像进行分割实验，分割准确率分别为98.5%和98.0%，真阳性率分别为88.2%和82.7%，与金标准之间的均方根误差分别为3.6和4.6像素。结果表明，该改进算法可获得更加准确的乳腺肿瘤分割结果。

边界停止函数 反应扩散 水平集演化 图像分割 乳腺超声图像

0 引言

乳腺癌是女性面临的最常见的恶性肿瘤。基于超声图像的计算机辅助诊断(computer aided diagnosis, CAD)系统是提高乳腺肿瘤诊断客观性的重要手段。肿瘤图像分割是CAD的前提条件和关键步骤，国内外对此提出多种方法，如灰度阈值结合动态规划法[1-2]、一致性直方图法[3]、活动轮廓模型[4-6]。但由于超声图像中斑点噪声严重、对比度低和乳腺肿瘤形状复杂多变等因素[5]，现有方法难以取得理想的分割结果。水平集方法属于几何活动轮廓模型，具有适应拓扑结构变化、稳定性好等优点，广泛用于医学图像分割[7]，成为研究热点。本文结合近年提出的反应扩散(reaction diffusion,RD)水平集[13]和Gabor各向异性扩散模型[16]，提出改进的水平集算法，以提高轮廓演化的效率和细节分割的准确性，实现乳腺肿瘤的精确快速分割。

1 算法设计与实现

1.1 RD水平集原理

原始的水平集方法有两种形式[13]，即基于偏微分方程的水平集方法和基于变分的水平集方法，方程分别如式(1)和式(2)所示。

(1)

(2)

式中：φ为水平集函数(levelsetfunction,LSF)，φ0为初始的水平集函数；F为边界力函数，控制着目标轮廓曲线的演化；为梯度算子；δ(·)为Dirac函数。这两种水平集方法在演化过程中，在过平滑或者过陡峭边缘容易出现严重的数值误差，因此需要重新初始化为一个有符号的距离函数。

反应扩散方程最初用于检测动物皮毛的化学机制[15]，它包括了反应和扩散两个过程，可用于描述纹理分析、自然图像建模与相变建模等领域的动态过程[17]。Zhang等[13]受相变理论启发，提出在传统的水平集演化方程中增加一个扩散项构造成RD方程，以使水平集演化保持稳定且消除重新初始化的过程。该方法不仅适用于变分水平集方法，而且适用于基于偏微分方程的水平集方法。

RD演化方程如下：

(3)

式中：ε为一个小的正常数；Δ为拉普拉斯算子；φ为LSF，φ0为初始的水平集函数；L(φ)为演化方程，在变分水平集方法中，L(φ)=-Fδ(φ)；在基于偏微分方程的水平集方法中，L(φ)= -F|φ|。

式(3)包含了两个动态过程：εΔφ为扩散项，在每个分割区域内逐渐把LSF正则化为一个分段常数；-ε-1L(φ)为反应项，使RD方程趋于稳定，直至L(φ)=0，得到分割区域。在传统的LSM中，由于无此扩散项，因此，不得不通过增加一个额外的过程(如重新初始化)来正则化LSF。

Zhang等[13]在基于相变的VanderWaals-Cahn-Hilliard理论上分析了当ε0+时，基于变分水平集方法的RD方程平衡解，构成与偏微分水平集方法的统一框架，从而，式(3)可看作传统水平集方法的一般形式。采用了两步分裂法(two-stepsplittingmethod,TSSM)来迭代地求解RD方程。

1.2 引入GAD的RD改进算法

在大噪声情况下，传统的RD算法存在分割细节不准确、演化效率低的问题。本文围绕传统RD方程存在的问题，在RD方程基础上增加一个边界停止函数项，构建新的演化方程：

(4)

式中：g为边界停止函数；α为控制边界停止函数权重的系数。

在变分水平集模型下L(φ)=-Fδ(φ)，边界力函数F=μdiv(φ/|φ |)-λ1(I-c1)2+λ2(I-c2)2，其中，μ、λ1、λ2为系数常量；I为目标分割图像；c1、c2分别为曲线内、外部常量值。在进行目标分割演化时。通过曲线内外力值大小的改变，控制曲线向内或者向外演化。当F=0时，表示内外力达到平衡，此时得到最终目标轮廓。增加边界停止函数后，边界力函数F =μdiv(φ/|φ |)-λ1(I-c1)2+λ2(I-c2)2+αg；通过停止函数g降低图像噪声对力的影响，同时改变系数α，增加弱边缘力对演化方程的影响，从而提高分割的准确性。

在弱噪声情况下，边界停止函数可表示为[6]：

(5)

式中：Gσ为方差为σ的高斯函数；I为原始图像；*表示卷积运算。

由式(5)可以看出，当图像梯度较大时，g的取值趋于0，对边界力函数影响相对较弱，曲线演化；反之，当梯度较小时，g的取值趋于1，边界力函数影响相对增强，控制曲线停止演化，此时得到目标的轮廓。

由于超声图像中存在很强的乘性斑点噪声，为克服其对演化曲线的影响，通常需要使用方差很大的高斯滤波器；但滤波器在滤除噪声的同时也使边缘变模糊，容易导致在弱边界处产生“边界泄漏”(即边界不连续、断裂)现象。Gabor各向异性扩散(gabor-basedanisotropicdiffusion,GAD)模型是一种对斑点噪声较为有效的滤波算法，它不仅能够在不同程度斑点噪声下保持较好的降噪鲁棒性，而且具有优越的边缘保持特性[16]。因此，本文在RD方程中引入GAD模型，以提高RD水平集分割的准确性和效率。引入GAD的边界停止函数表示为：

(6)

式中：GAD(I)为经过GAD滤波后的图像。

1.3 算法实现

文献[18]利用TTSM算法实现反应扩散，首先由反应方程生成一个0和1的二值函数;然后将扩散方程应用到二值函数中，得到曲率运动;最后，得到TTSM算法实现RD方程。该方程可以表述为：

① 解反应项φt=-ε-1L(φ0)，φ(x,t=0)=φn，在某时间Tr达到中间解,表示为φ(n+1)/2=φ(x,Tr)。

② 解扩散项φt=εΔφ，φ(x,t=0)=φ(n+1)/2，在某时间Td，最后水平集表示为φn+1=φ(x,Td)。

在以上2个步骤中，通过选择小Tr和Td，φ(n+1)/2和φn+1分别离散地近似为φ(n+1)/2=φn+Δt1[-ε-1L′×(φn)]，φn+1=φn+1/2+ Δt2(εΔφn+1/2)，其中时间步长Δt1和Δt2分别表示时间Tr和Td。将参数ε并入到时间步长Δt1和Δt2中，如Δt1Δt1(-ε-1)和Δt2Δt2ε。因此，类似文献[18]，只要选择两个能够保持数值稳定的时间步长Δt1和Δt2即可。

用n表示迭代次数，φn表示第n次迭代的LSF，基于水平集演化的RD算法实现过程概括如下。

① 初始化：φn=φ0，n=0；

② 计算φ(n+1)/2：

φ(n+1)/2=φn-Δt1L′(φn)

(7)

式中：L′(φn) =L(φn)+αgδ(φn);Δt1为反应项的时间步长，影响着演化的稳定性。

③ 计算φn+1：

φn+1=φn+Δt2Δφn

(8)

式中：φn=φn+1/2；Δt2为扩散项的时间步长，取值范围为0～0.25。

④ 如φn+1满足稳定条件，停止演化；反之，n=n+1，返回步骤②。

2 实验结果与讨论

为验证算法的有效性，采集来自77例患者(年龄17～79岁)的111幅(79幅良性；32幅恶性)乳腺肿瘤B型超声图像进行分割实验。超声仪型号为法国SuperSonic Aixplorer，超声探头频率范围4～15 MHz。

为比较分割算法的性能，将本文提出的改进RD算法与Zhang提出的传统RD 算法及手工分割方法进行实验对比。由于本文算法是在RD算法的基础上提出的，为了客观地对比算法的特性，将两者共有的参数设置为相同的值，并选择相同的初始轮廓(手工描记4～6个点连成多边形)，避免初始轮廓对分割结果的影响。

图1和图2分别为两种算法对良性乳腺肿瘤超声图像和恶性乳腺肿瘤超声图像的分割结果对比。

图1 良性乳腺肿瘤超声图像分割结果

图2 恶性乳腺肿瘤超声图像分割结果

图1是一个形状较为规则的良性乳腺肿瘤的图像。乳腺肿瘤内部灰度不均匀，1点钟方向呈现出弱边缘。由分割结果可看出，在相同参数情况下，本文算法在弱边缘和大噪声下体现出较好的分割特性。图2是一个形状相对复杂、边缘凹凸不平的恶性乳腺肿瘤的图像。在相同参数情况下，本文算法相对传统RD方法能更有效提取9点钟及2点钟方向凸出的边缘，与手工分割结果更加相符。

本文采用真阳性率(true positive rate,TPR)、假阳性率(false positive rate,FPR)、准确率(accuracy,Acc)、特异性(specificity,Speci)、均方根误差(root mean square error,RMSE)共5个定量指标，并以医生手工分割作为金标准，来度量算法的分割质量。计算公式如下：

(9)

TPR代表实际分割肿瘤区域的准确率，FPR是非肿瘤区域被算法分作肿瘤区域的概率，Acc表示综合分隔肿瘤区域的准确率，Speci表示非肿瘤区域被算法正确判别的概率。TPR、Acc、Speci度量了分割算法得到区域与金标准区域之间的值，取值在0～1之间，值越大表示分割越准确。RMSE则反映了分割轮廓坐标点与金标准轮廓坐标点之间的平均距离，值越小表示分割越准确。图3给出GT、TP、TN、FP、FN等参数的物理意义。

对111幅乳腺肿瘤图像分别使用两种算法进行分割，得到的量化指标平均值和标准差如表1所示。从表1可知，RD和本文算法的Acc和TPR分别为98.0%1.0%和98.5%0.7%，82.7%7.7%和88.2%5.7%，RMSE分别为4.6341.373pixels和3.5970.813pixels。本文算法在Acc、TPR和RMSE三个指标上均优于传统RD算法。值得注意的是，本文算法的FPR和Speci指标略逊于RD算法，但这两个指标仅在分割灰度均匀的目标时具有较高参考价值，在分割不均匀的乳腺肿瘤时并不能很准确地表征分割性能。实际乳腺肿瘤病灶内部灰度值的分布往往是很不均匀的，当初始轮廓定义在病灶内部、演化时被病灶内伪边缘所吸引而未到达真实轮廓，Speci和FPR也能达到比较好的值；比如在分割结果完全位于金标准内部时，Speci=1，FPR=0。因此，本文将Acc、TPR和RMSE视为主要指标，FPR和Speci视为次要指标。综合考虑所有指标，可知本文算法对乳腺肿瘤超声图像的分割精度优于传统RD算法。

图3 GT、TP、TN、FP、FN等参数的物理意义

表1 RD算法和本文算法对乳腺超声图像的分割定量指标

Tab.1 The quantitative indicators of breast ultrasound image segmentation(RD vs. improved RD)

指标RD算法均值标准差本文算法均值标准差TPR0．8270．0220．8820．033FPR0．0150．0100．0280．007Acc0．9800．0100．9850．032Speci0．9971．3730．9920．813RMSE/pixel4．6343．597

3 结束语

本文提出了一种改进的RD水平集乳腺肿瘤超声图像分割方法。与传统的RD水平集方法不同的是，本文的方法利用具有斑点噪声强抑制力的GAD算法进行滤波处理，得到边界停止函数，并将该函数作为RD能量函数中的一个增加项，使得改进的RD算法具有更好的演化效率以及捕捉弱边缘的能力。此外，本文的改进RD模型是通过水平集方法实现的，继承了水平集良好的拓扑自适应性优势，可以较好地分割出形状复杂多变的乳腺肿瘤。实验结果表明，本文算法在乳腺肿瘤超声图像分割中取得良好效果。尤其在应用于乳腺超声图像中目标区域灰度不均匀、存在凹凸不平边界及弱边界情况下，与传统RD方法相比，本文算法显示出优越性能。下一步将考虑如何自动提取初始轮廓以进一步提高分割的客观性与自动化水平。

[1] 汪源源，沈嘉琳，王涌，等.基于灰度阈值分割和动态规划的超声图像乳腺肿瘤边缘提取[J].航天医学与医学工程，2005,18(4):281-286.

[2] 汪源源，沈嘉琳，王涌，等.基于形态特征判别超声图像中乳腺肿瘤的良恶性[J].光学精密工程，2006,14(2):333-340.

[3] 黄剑华，于瀛星，张英涛，等.基于一致性直方图的超声乳腺图片分割方法[J].哈尔滨工业大学学报，2008,40(7):1103-1106.

[4] Wang W M,Zhu L,Qin J,et al.Multi-scale geodesic active contours for ultrasound image segmentation using speckle reducing anisotropic diffusion [J].Optics and Lasers in Engineering,2014,54:105-116.

[5] Liu B,Chen H D,Huang J H,el al.Probability density difference based active contour for ultrasound image segmentation [J].Pattern Recognition,2010,43(6):2028-2042.

[6] Chan T F,Vese L A.Active contours without edges [J].IEEE Transactions on Image Processing,2001,10(2):266-277.

[7] 张虎重，康志伟.融合模糊聚类的变分水平集图像分割模型[J].电子测量与仪器学报，2011,25(4):325-330.

[8] 高粱,刘晓云，陈武凡，等.基于相位和GGVF的水平集乳腺超声图像分割[J].仪器仪表学报，2012,33(4)：870-877.

[9] Li L,Kruusmaa M.Ultrasound image segmentation with localized level set based curve evolution and Rayleigh distribution [A],Tuduce R.IEEE International Conference on Systems,Signal and Image Processing[C],Bucharest:IEEE,2013:139-142.

[10]Li C M,Xu C Y,Gui C F,et al.Level set evolution without re-initialization:a new variational formulation [A],Schmid C.IEEE Computer Social on Computer Vision and Pattern Recognition[C],San Diego:IEEE,2005:430-436.

[11]Li C M,Xu C Y,Konwar K M,et al.Fast distance preserving level set evolution for medical image segmentation [A].Processing of the 9th International Conference on Control,Automation,Robotics and Vision[C],Singapore:IEEE,2006:1-7.

[12]Li C M,Xu C Y,Gui C F,et al.Distance regularized level set evolution and its application to image segmentation [J].IEEE Transactions on Image Processing,2010,19(12):3243-3254.

[13]Zhang K H,Zhang L,Song H H,et al.Re-initialization free level set evolution via reaction diffusion [J].IEEE Transactions on Image Processing,2013,22(1):258-271.

[14]虞红伟，厉力华，徐伟栋，等.一种基于水平集的多尺度乳腺肿块分割方法[J].仪器仪表学报，2010,31(6):1418-1423.

[15]Turing A M.The chemical basis of morphogenesis,philosophical transactions of the royal society of London [J].Series B,Biological Sciences,1952,237(641):37-72.

[16]Zhang Q,Han H,Ji C H,et al.Gabor-based anisotropic diffusion for speckle noise reduction in medical ultrasonography [J].Journal of the Optical Social of America A,2014,31(6):1273-1283.

[17]Turk G.Generating textures on arbitrary surfaces using reaction-diffusion[J].Computer Graphics,1991,25(4):289-298.

[18]Merriman B,Bence,Osher S.Motion of multiple junctions:a level set approach [J].Journal of Computational Physics,1994,112(2):334-363.

增强安全意识，提升防爆技能——自仪院防爆技术培训火热招生

随着国民经济的迅速发展，石油、化工、钢铁、煤炭、制药等企业的生产规模不断扩大，防止爆炸性事故的发生已得到社会普遍关注，昆山工厂粉尘爆炸等多起事故以生命和鲜血一次次地敲响了安全生产的警钟。长期实践表明，防爆生产装备相关的管理、设计、安装、使用、维护人员全面掌握防爆电气理论和操作技能对于保障企业安全生产的必要性已毋容置疑。

为增强企业安全意识，提升防爆作业人员技术水平，上海工业自动化仪表研究院继续教育培训中心将于2015年11月10日在上海举办防爆技术培训，现已火热招生，诚邀您的参与。

1.特种作业人员电工(防爆电气)安全技术培训班

(颁发国家安监总局监制《特种作业(电工：防爆电气)操作证》，全国通行)

2.人保部“专业技术人才知识更新工程”防爆技术研修班

(颁发国家人保部“专业技术人才知识更新工程”继续教育证书，全国通行)

招生详情可见www.sipai.com，或询培训中心姜老师、陈老师(021-64368180-324，training@sipai.cn/tccedu@sipai.cn)。

Segmentation of Breast Tumor Ultrasound Images Using Reaction Diffusion Level Set

For segmentation of breast tumor ultrasonic images, the improved reaction diffusion (RD) level set segmentation algorithm is proposed. Firstly, the Gabor anisotropic diffusion model is used for filtering, to construct edge stopping function; then the function is integrated into RD level set evolution equation, for controlling the evolution of curves to obtain the boundary of breast tumor. Segmentations of 111 breast ultrasonic images from 77 patients are conducted with the proposed method and the traditional method, the segmentation accuracies are 98.5% and 98.0% respectively; the true positive rate are 88.2% and 82.7%; and the mean square errors are 3.6 and 4.6 pixels. The results show that the proposed improved algorithm provides more accurate segmentation for breast tumor.

Edge-stopping function Reaction diffusion(RD) Level set evolution Image segmentation Breast ultrasound image

国家自然科学基金青年资助项目(编号：61401267);

上海市自然科学基金资助项目(编号：12ZR1444100);

上海市教委人才计划“晨光计划”资助项目(编号:11CG45)；

陈帅(1987- )，男，现为上海大学生物医学工程专业在读硕士研究生；主要从事医学图像处理方面的研究。

TP391

A

10.16086/j.cnki.issn1000-0380.201509008

上海市教委科研创新项目(编号：12YZ026)。

修改稿收到日期：2015-03-24。