王光磊，苑昊，魏帅，时亚松，刘秀玲

(1.河北大学 电子信息工程学院，河北省数字医疗工程重点实验室，河北 保定　071002)



基于区域生长和聚类的主动脉CTA图像序列分割算法

王光磊，苑昊，魏帅，时亚松，刘秀玲

(1.河北大学 电子信息工程学院，河北省数字医疗工程重点实验室，河北 保定071002)

摘要：针对人体主动脉CTA序列图像的特点，提出了一种基于区域生长和聚类的序列分割新算法.在确定好合适的分割阈值范围后，结合基于区域生长的算法对主动脉的目标区域进行轮廓提取，然后对得到的目标轮廓在基于isodata的算法上进行聚类处理，由于主动脉在空间上的连续性，可以将聚类后得到的目标区域的聚类中心作为下一幅CT图像新的种子点再进行区域生长，从而实现仅在选择单幅图像种子点的情况下一次性完成整套CT序列图像的分割.实验结果表明，本文算法计算量小，分割精度高，可以完整准确地将主动脉从CT序列图像中分割出来.

关键词：主动脉； CTA序列图像；空间连续性；区域生长；聚类

第一作者：王光磊(1983-)，男，天津人，河北大学讲师，博士，主要从事医学图像处理等方面研究.

E-mail:windlay@hotmail.com

E-mail:liuxiuling121@hotmail.com

近年来许多主动脉疾病的发病率明显增高，这里包括由动脉粥样硬化引起的腹主动脉瘤和动脉闭塞性疾病等，尤其是腹主动脉瘤，它的发病率已经超过了颅内动脉瘤，成为各种动脉瘤发病率之首，而且还在呈现明显的上升趋势[1].这些疾病严重威胁着人们的身体健康，而如何能快速有效地进行诊断和定量分析就变得十分关键，所以对主动脉进行三维重构是非常有意义的.当然这也是主动脉虚拟手术系统实现的一个重要保障，而前提则是能在二维CT序列图中将主动脉轮廓完整地提取出来，考虑到主动脉的CT序列图像数量过大，以层厚0.6 mm，层间距0.6 mm的螺旋CT为例，要将主动脉完全覆盖到需要大概600多幅CT图像，若采取单一的CT图像提取的方法是十分耗时且无意义的，所以序列分割算法的提出是很必要的.

当前国内外研究者们实现序列分割的方法主要也还是利用连续图像间的相似性.这其中一种实现思想是计算上一张目标区域的质心作为下一张分割的种子点来提取轮廓[2].而另外一种应用较多的思想则是将上一幅CT图像目标区域分割后得到的轮廓点集投影到下一幅CT图像中，以投影到下一张目标区域中的点集作为新的种子点再继续区域生长都得到新的轮廓[3-5].基于质心的算法的确能够实现一定器官，组织的CT图像序列分割，但不能很好地适用于有拓扑变化的情形，以本文研究的人体主动脉为例，在二维横切CT图像上目标轮廓并不仅仅都是一个连通区域，而是由一个连通区域逐渐变化到2个连通区域，然后再变回到一个连通区域的过程，这时，基于质心的算法就很难将完整的主动脉轮廓提取出来.而基于投影点集的算法则在拓扑变化和分割效果上都有很好的实现，不过，该算法由于利用的是轮廓边缘信息，实现的前提是上一幅目标轮廓的投影点集中必须有包含在下一幅目标区域内的像素点，而这增大了对连续图像间距的依赖性.所以，本文提出了基于区域生长和聚类的CT序列分割新算法.

该算法利用了聚类思想中“合并”和“分裂”的特点，很好地解决了主动脉在二维分割时目标区域数量变化的问题，即根据目标区域的个数确定种子点的个数，为提取完整的主动脉轮廓提供了保障.首先对上一幅CT图像的分割结果进行聚类处理，得到的聚类中心作为下一幅CT图像的新的种子点，再结合基于区域生长思想的算法进行分割.由于主动脉在螺旋CT图像上面积相对较大，使得本文设计的算法减少了对图像间距的依赖，另外在调整好几个重要参数如合并参数的前提下该算法可以适应多种不同的分裂合并情况，提高了分析的智能性.

1算法描述

提出了一种基于区域生长和聚类的序列分割新算法.基于区域生长的算法对主动脉的目标区域进行轮廓提取，然后对得到的目标轮廓在基于isodata的算法上进行聚类处理，将聚类后得到的目标区域的聚类中心作为下一幅CT图像新的种子点再进行区域生长，从而实现仅在选择单幅图像种子点的情况下一次性完成整套CT序列图像的分割.算法流程如图1所示.

1.1　基于区域生长的分割算法

区域生长算法的基本思想是把相同或相似性质的像素聚集成区域.具体是在目标区域内找到合适的种子点，作为生长起点，其次通过已确定好的生长准则，将种子像素周围与其相同或相似的像素并入到种子像素所在的区域内，这些新的像素作为新的种子点继续向四周生长，直到不再有满足条件的像素包含进来，这时目标区域就生成了.区域生长算法要解决的3个核心问题：1)选择合适的种子像素；2)确定能把区域内的像素包含进来的生长准则；3)确定终止条件[6].

种子点的选取直接关系到最终的分割效果，因为本文算法是基于聚类ISODATA算法实现的，所以序列分割时的种子点就是聚类后得到的最终中心坐标.

生长准则的确定：设初始种子点的坐标位置为(x0,y0)，初始种子点的像素值为f(x0,y0)，种子点附近八邻域像素值为

f(xi,yi)=f(x0+Δx,y0+Δy),(Δx,Δy∈[-1,1])，

(1)

则生长准则为

(2)

由上述生长准则可以看出，区域生长中的阈值范围的选取直接决定了最后的分割效果.因为研究的数据是CT图像，各个区域的的灰度差别较大，另外CT图的噪声相对较小，因此本文选择手动设定分割阈值，再结合区域生长便可以比较准确地分割初步目标区域[7].由于省去了阈值计算的过程，也极大地简化了分割过程，缩短了程序运行时间.

关于生长方式，本文选择八邻域生长方式.设定S为出发点，CD为当前搜索方向，它指向要搜索的点方向，CD∈{0，1，2，3，4，5，6，7}.从S点出发，依据搜索设定的方向CD依次寻找下一个边界点.定义下一个点的搜索方向为ND=(CD+7)mod8.在这个搜索算法中，下个边界点能保存上个点的搜索方向，并且逆时针回溯一个方向，再从此方向顺时针搜索.这种搜索方式能优先搜索上个点搜索方向附近的方向，减少了重复搜索，提高了生长速度.

为了避免噪声点的干扰，需要在算法中加入一个约束条件，即在搜索过程中，如果生成了闭合边界，首先计算此闭合边界内部的面积，一旦面积小于一定值，则说明遇到了噪声点，此时抛去噪声点边界，再继续进行.搜索过程一直持续搜索到出发点S，则边界跟踪完成[2].

1.2　基于ISODATA算法的种子点的确定

ISODATA，迭代自组织分析，该算法是机器学习中非监督动态聚类算法中的一种，它是在k-均值算法的基础上，增加了对聚类结果的“合并”和“分裂”2个操作，通过设定初始参数而引入人机对话环节，并使用归并与分裂的机制.当聚类结果某一类中的样本数太少，低于设定值时，或者2个类间的距离太近时(由合并参数调节)，进行合并.而当聚类结果中某一类中样本某个特征类内方差太大，将该类进行分裂.

本文利用ISODATA的算法步骤[8]

Step1: 选取5个对象c1,c2,c3,c4,c5作为聚类中心.

Step2：将每个样本像素(灰度值为255)c1依据公式(3)分派到聚类中心为ci的聚类块Ci中.

(3)

对于任意一个Ci如果满足Ni

Step3：按照下面公式修正聚类中心值

(4)

Step4：计算聚类块Ci中的各个样本像素距离其中心ci的距离

(5)

Step5：计算所有样本像素与其相应聚类中心的平均距离

(6)

情况1：若迭代完成(由参数I决定)，令C=0，不再合并，跳转Step10；

情况3：若迭代次数为偶数，又或者k≥2K错误！未找到引用源，则跳转到Step10，否则跳转到Step8.

Step6：用下面公式求得每一个聚类i对应的标准差σi=(σi1,σi2,…,σin)T

(7)

式中xlm是第l个样本的第m个分量，xl在Ci中，xim是第i个聚类中心的第m个分量；σim是第i个聚类第m个分量的标准差，n是样本的维数.

求出每类具有最大标准偏差的分量σjmax,j=1,2,…,c.

Step7：如果存在σimax>S,i=1,2,…,k，且满足以下2个条件之一

Step8：计算所有的聚类中心两两之间的距离

δij=d(ci,cj),i=1,2,…,k,j=i,i+1,…,k.

(8)

Step9： 比较δij和C，把小于C的δij按大小作升序排列，δi1j1<δi2j2<…<δiLjL，其中L是步骤1给出的每次迭代允许的最大聚类对数.

Step10：从最小的δi1j1开始，对于每个δi1j1合并2个类Cil和Cjl，聚类中心为

(9)

并把聚类数减少，k=k-1.

从第2个δiLjL开始，检查δiLjL是否已在前面合并过程中被合并，如果两者未被合并，则执行合并过程.

Step11：如果迭代完成程序结束，否则根据需要输入参数跳转到步骤1，如输入参数不变，则跳转到Step2.

2实验与分析

2.1　材料

研究对象是人体主动脉CT断层序列图像，该图像是由德国西门子系统医疗集团研制的Sensation 4 多层面螺旋CT 获取，其主要扫描参数为层厚0.6 mm，层距0.6 mm.图2 是利用Sensation 4 多层面螺旋CT 在以上扫描参数下扫描得到的某一病人腹部的断层序列图像，从左至右分别为第176，264，393 张断层图像.从断层序列图像上可以清晰观测到连通区域的主动脉.每个病人样品均采集到700张以上分辨力为512~512的12位断层图像.为了方便在普通显示器上进行处理，本文中的主动脉图像是由12位的DICOM格式转化为8位JPEG格式的灰度图像.

2.2　不同阈值的分割结果

如果图像中含有噪声将会影响区域生长的正常进行，所以在区域生长分割之前，对图像要进行一些边缘保持的平滑处理.为了实现更好的分割效果，本文选择基于阈值的区域生长算法进行分割，所以在算法的最开始需要为CT图像找到合适的分割阈值范围.

序列CT图像中有些区域会在一些特定的状态下与目标区域“粘连”，而它们的灰度值又与目标区域十分相近，对目标轮廓的提取造成了很大的干扰.可以检测出连接区域的灰度值都在[200，250]之间， 又因为主动脉内部像素灰度值没有超过600的，所以设定阈值为[280，600].

现在对特殊位置下的主动脉CT图像进行检测，将不同的阈值添加到单幅区域生长的分割算法中，图3b是阈值[280，600]的结果，图3c是阈值[200，600]的结果.

很明显，不同的阈值范围对目标轮廓的提取影响很大，而检测得到的最佳阈值的确可以很好地提取出这种与干扰区域近似“粘连”状态下的目标区域.这也是本文选择设定阈值的另一个重要原因.

2.3　分割精度对比结果

在ITK的开发平台下，利用面向对象编程语言C++语言实现上述算法，并对实验图像605幅DICOM格式的人体主动脉序列图像进行分割.由于主动脉在医学上主要分为主动脉升部、主动脉弓和主动脉降部，所以本文在这3个部分抽选结果做检测，图4为本文所描述的算法分割后的结果，从上到下分别为主动脉升部、主动脉弓和主动脉降部的结果.从左到右依次为CT原图、医生手动分割、单幅的区域生长分割以及本文的序列分割的结果对比.

需要提出的是，在前2例中，序列分割的结果本是2个连通区域，由于只做主动脉特定部位的比较，所以对序列图像非针对部分做了后期隐去处理.

为进一步定量考察分割算法的性能，本文采用了Jaccard系数[9]来对分割结果进行定量评估，该指标越接近1，说明分割结果越是接近分割标准，也就是说采用的分割算法越好.将A和B看作像素的集合，A是参考图像，B是分割结果，则Jaccard系数的表达式为

对图中单幅CT分割和序列CT图像分割结果进行质量评估，得到表1.分割标准来自有经验的专家对图像进行分割得到的结果.

通过Jaccard系数对比可以看出本文算法在实现快速序列分割的前提下，依然能够保证分割精度，与单幅区域生长分割的效果近似于相同，这也为后期的主动脉三维重构提供了基础和保障.

通过Jaccard系数对比可以看出本文算法在实现快速序列分割的前提下，依然能够保证分割精度，与单幅区域生长分割的效果近似于相同，这也为后期的主动脉三维重构提供了基础和保障.

2.4　多种子点聚类追踪结果

给定初始种子点的前提下，依据图像的连续性利用ISODATA算法的合并与分裂特点对主动脉种子点进行聚类.为了将整个追踪过程更加直观、全面地展现出来，从序列追踪结果中选出了具有代表性的10张结果.

图5为对目标轮廓的聚类结果，圆圈标记处即为聚类中心，为了方便观看，在左上角处将聚类中心的个数及位置坐标显示了出来.

从图5可以看到，首先是对主动脉最顶端，即主动脉弓的单个连通的目标区域进行聚类得到中心，直到目标连通区域的形变超过一定程度，即聚类过程中2个聚类中心的距离超过了设定的合并参数C，这时会分派出2个聚类中心.此时的目标区域主要是主动脉升部和主动脉弓2部分.随着追踪的进行，当目标区域形变成第3幅时，由于左上部似乎又有分裂的趋势，所以这时会增加到3个聚类中心，不过，随着目标区域的恢复回来，并没有再次分裂成3个连通区域，所以上部的2个聚类中心合并，形成第5幅的结果.接着，随着主动脉升部的结束，目标区域只还有1个连通区域，即主动脉降部，当然，此时也只有1个聚类中心.最后当目标区域再次发生分裂时，就可以认为主动脉的聚类过程已经结束.

很明显，基于聚类的多种子点追踪算法可以很好地将主动脉各个阶段的目标区域的聚类中心检测标记出来，无论是针对单个连通的目标区域还是处于分裂状态的多个连通目标区域都能很好地实现出结果，计算出聚类中心，再作为种子点反馈给基于区域生长的分割算法，最终提高完整主动脉序列意义上分割的自动化程序.

3结论

依据主动脉的形状以及在二维CT图像上变化的特点，提出了基于区域生长和聚类的序列分割新算法.在确定好合适的分割阈值的前提下，该算法将聚类中“合并”和“分裂”的2个特点应用到目标连通区域的变化中，得到主动脉在CT图像上各个阶段的聚类中心，再以得到的聚类中心为种子点区域生长得到目标轮廓，从而实现主动脉的序列分割.

结果很好地证实了该算法在保证分割精度的前提下能够实现快速序列分割. 在未来的研究中可以考虑加入自适应阈值算法，即为每一张CT图像都能计算出属于它自己的最优阈值，而并非一个整体的阈值，这样处理后得到的分割结果会更为理想，也为后续的三维重构及虚拟介入手术提供了更好的数据和保障.

参考文献:

[1]胡明. 腹主动脉瘤的病因学分析、治疗评估、经济学评价[D].南宁:广西医科大学,2013.

[2]吴龙海,周荷琴,张鹿，等.从CT图像中自动分割出肺部区域的算法研究[J].航天医学与医学工程杂志,2008,21(5):425-429.

WU Longhai,ZHOU Heqin,ZHANG Lu. et al. Research on algorithm for automated lung segmentation in CT images[J].Space Medicine & Medical Engineering, 2008,21(5):425-429.

[3]陈郁淦,周学成,乐凯.根系CT序列图像区域生长分割的新方法[J].计算机工程与应用杂志,2011,47(28):158-161.

CHEN Yugan,ZHOU Xuecheng,LE Kai.New method for region growing segmentation of plant root CT image sequences[J].Computer Engineering and Applications, 2011,47(28):158-161.

[4]周学成,罗锡文.基于阈值分析与区域生长相结合的根系CT序列图像分割算法[Z]. 中国农业工程学会2005年学术年会, 广州：2005:336-339.

[5]彭丰平,鲍苏苏.一种基于区域生长的CT序列图像分割算法[J].计算机与数学工程, 2007,35(5):1-2.

[6]孙浩,陆雪松,刘海华.CT图像中的主动脉分割方法[J].现代科学仪器,2013(2):45-48.

SUN Hao,LU Xuesong,LIU Haihua.A segmentation method of aorta in CT image[J].Modern Scientific Instruments, 2013(2):45-48.

[7]BALAFAR M A.Review of brain MRI image segmentation methods[J].Artificial Intelligence Review, 2010,33(3):261-274.

[8 ]张语涵,孙劲光,苗锡奎.基于改进的ISODATA算法彩色图像分割[J].计算机系统应用, 2010,19(2):41-45.

ZHANG Yuhan,SUN Jinguang,MIAO Xikui.Color image segmentation based on improved ISODATA[J].Computer Systems & Applications, 2010,19(2):41-45.

[9]詹天明,韦志辉,张建伟,等.脑MR图像分割和偏移场矫正的耦合水平集模型[J].中国图像图形学报,2011,16(11):2017-2023.

ZHAN Tianming,WEI Zhihui,ZHANG Jianwei, et al.Coupling level set model for brain MRI CT image segmentation and bias field correction[J].Journal of Image and Graphics, 2011,16(11):2017-2023.

(责任编辑：孟素兰)

A novel sequence aorta CTA image segmentation algorithm
based on region growing and clustering

WANG Guanglei,YUAN Hao,WEI Shuai, SHI Yasong, LIU Xiuling

(Key Laboratory of Digital Medical Engineering of Hebei Province, College of Electronic

and Information Engineering, Hebei University, Baoding 071002, China)

Abstract:The paper proposed a new sequence images segmentation algorithm based on region growing and clustering applying for human aorta CTA sequence images. After determination of proper thresholding range, contour extractions of interested region of arteries were carried out with the help of region growing. Furthermore, isodata algorithm was used for clustering of the interested contour region. Owing to the spatial continuity of arteries, we employed the extracted clustering center of interested region as a new seed point of the next CT image for initialization of region growing. The main contribution of the paper was the implementation automatic lumen segmentation for sequence CTA with only one seed point. The experimental results demonstrated that the new method could extract contour of arteries in CTA sequence images entirely with less calculation and high segmentation accuracy.

Key words:aorta; CTA sequence images; spatial continuity; region growing; clustering

DOI：10.3969/j.issn.1000-1565.2015.06.015

中图分类号：TP399

文献标志码:A

文章编号：1000-1565(2015)06-0644-07

收稿日期：2015-09-18

基金项目：国家自然科学基金资助项目(61473112; 61203160)

通信作者：刘秀玲(1977-)，女，河北沧州人，河北大学教授，博士，主要从事医学图像处理等方面研究.