基于粒子群优化的Otsu肺组织分割算法
2015-10-22孟亚州高林爽
孟亚州，马 瑜，白 冰，高林爽
（宁夏大学研究生院，宁夏银川750021）
基于粒子群优化的Otsu肺组织分割算法
孟亚州，马 瑜∗，白 冰，高林爽
（宁夏大学研究生院，宁夏银川750021）
为减少肺组织分割算法的运算时间，提出了一种基于粒子群优化的Otsu肺组织改进自动分割算法.针对传统粒子群优化的二维Otsu算法中二维直方图计算量大、粒子搜索容易陷入局部最优解的缺陷，使用灰度级-梯度二维直方图减少二维直方图的计算量，并减小粒子搜索范围，采用基于粒子空间对称分布的改进粒子群获取最佳阈值.算法实现过程中利用孔洞填充算法去除背景，基于形态学操作去除噪声、修补病变区域产生的孔洞.仿真实验结果显示，本文算法对图像尺寸为512像素×512像素CT图像的阈值分割时间约为0.2 s，比基于灰度级-邻域均值二维直方图的粒子群优化的Otsu算法的阈值分割速度提高了约16%.较好地实现了胸腔CT图像的肺组织自动分割，与传统算法相比较，本算法在保证分割精度的基础上分割速度明显提升.
二维Otsu；改进粒子群；孔洞填充；形态学操作；肺组织自动分割
1 引 言
肺组织分割是计算机肺部疾病辅助诊断系统的关键技术之一.随着计算机辅助医疗诊断技术的快速发展，肺组织图像分割算法对分割的精度和速度提出更高要求.
近年来，肺组织分割算法的研究日益得到关注.2009年，孟琭等人提出了基于Snake模型的病变CT肺部图像分割算法，无论肺内有无病变，分割效果均比较理想，但算法运行时间相对较长［1］.2012年，Tetyana Ivanovska等人先初步分离出肺和气管，然后将肺和气管分离开［2］. 2013年，Ying Wei等人利用三维区域连通标记、三维区域生长法对肺组织进行分割，分割精度较高，但算法复杂度相对较高［3］.Shiloah Elizabeth Darmanayagam等人利用多层前馈神经网络判断肺组织分割是否完全，利用阈值法、形态学处理修补肺组织.虽然对病变区域分割效果理想，但是算法处理时间较长［4］.
传统阈值法实现简单、快速，但分割精度较低，对有病变的胸部CT图像分割效果较差.区域生长法属半自动分割方法，而且对种子点及生长合并规则的选取比较敏感.基于活动轮廓模型的分割法，效果较好，但是分割需花费大量时间.基于模式分类的胸部CT分割方法，训练样本需要人为选取，而且由于要提取多个特征，处理时间较长.
针对肺部疾病计算机辅助诊断技术对全自动性、实时性要求较高，因此肺组织分割算法不仅仅要求分割精度，还要求算法必须全自动、快速分割.
综上考虑，本文提出了一种基于改进粒子群优化Otsu的胸腔CT图像快速自动分割算法，使用灰度级-梯度二维直方图，大大减少了直方图的计算量，同时减少了粒子的搜索范围.使用改进粒子群算法来搜索最优分割阈值，不仅加快了传统OTSU算法的运算速度，也避免了过早收敛陷入局部最优解的现象.
2 基于粒子群优化改进的二维Otsu分割算法
2.1二维Otsu分割算法
Otsu算法是日本学者Otsu提出的一种自动阈值分割算法［5］.该算法在图像的灰度直方图基础上，将目标和背景的类间方差最大化作为自动确定分割阈值的依据.一维Otsu算法在信噪比低的图像中分割效果不理想，而二维Otsu利用了图像的局部信息，有效地避免了噪声的影响［6］.考虑到胸腔CT图像也存在噪声，因此本文采用二维Otsu算法.
图1 两种二维直方图的比较Fig.1 Comparison of two different two-dimensional histogram
如图1所示，给定一个分割阈值（s，t），二维直方图分成0、1、2、3四个区域，0、1区域对应图像的背景类、目标类，2、3区域对应图像的边界和噪声.在图（a）中，灰度级-邻域均值二维直方图存在错分的情况:对应图像背景类、目标类的0、1区域的左上角、右下角可能存在边界和噪声，对应边界和噪声的2区域的右下角处、3区域的左上角处可能存在目标或背景点.计算二维直方图需要对整个二维直方图计算；而灰度级-梯度二维直方图如图（b）所示，信息主要分布在直方图的下方，不存在错分情况，计算二维直方图时，只需计算0、1区域，并且粒子群的粒子搜索范围为0、1区域，有效减少了粒子的搜索范围，进而缩短算法运算时间［7］.基于以上考虑，本文采用灰度级-梯度二维直方图.
二维Otsu具体算法概述如下:设一幅大小为M×N的图像有L个灰度级，灰度值为i、其梯度为j的像素点个数为nij，则其概率pij为:
以阈值变量（s，t）将灰度图像的像素分为背景类A、目标类B，则背景类A、目标类B出现的概率PA、PB分别为:
则背景类、目标类的均值向量uA、uB分别为:
图像的总均值向量u为:
离散度矩阵为:
则离散度测度为:
由于离散度测度越大代表两类间的差别越大，因此选取离散度测度达到最大值时对应的分割阈值作为最佳阈值（s∗，t∗）［8］，即
2.2改进粒子群算法
粒子群算法是Kennedy和Eberhart受到飞鸟集群活动的规律性启发提出的一种基于群智能的随机优化算法［9］.传统粒子群算法存在收敛速度慢、容易陷入局部最优解的问题［10］，改进算法则有效解决了这两类问题.本文中粒子群算法选用基于粒子空间分布的粒子群算法，利用胸腔CT图像的灰度级-梯度直方图分布特性，减少粒子的搜索范围，惯性权重、学习系数采用线性变化策略.粒子在运动初期，提高全局搜索能力，更多向群体学习；粒子在运动后期，提高局部搜索能力，加强自我学习能力［11］.
图2 粒子空间位置分布Fig.2 Space distribution of particles
粒子搜索在迭代过程中，会出现粒子位置在局部最优解一侧的情况（如图2（a）所示），这样会导致粒子容易陷入局部最优解.为避免此类情况，迭代过程中需要调整粒子位置提高种群的多样性，使粒子的位置尽量在局部最优解的周围（如图2（b）所示）［11-12］.
迭代过程中出现当前全局最优解左右两侧的粒子个数不均衡时，则需调整粒子位置使粒子在当前全局最优解两侧的粒子个数均衡.如图3所示，在j维上，出现粒子在当前最优解右侧的个数为6，而在左侧的粒子的个数为3的情况时，将右侧的两个粒子，对称调整到左侧.
图3 粒子基于空间对称分布调整示意图Fig.3 The sketch map of adjusting particle's position based on symmetrical distribution
改进粒子群算法概述如下:
Step 1.粒子群初始化.设种群共有N个粒子，其搜索空间为d维空间，则第i个微粒的位置可表示为xi=（xi1，xi2，...，xid），其飞行速度表示为vi=（vi1，vi2，...，vid）.每个粒子所经历过的最好适应值的位置为pi=（pi1，pi2，...，pid），种群所有微粒所经历过的适应值最好的位置为全局最好位置，记pg=（pg1，pg2，...，pgd）.
Step 2.判断是否满足迭代条件，满足则转到Step 3，否则转到Step 7.
Step 3.更新pi、pg，计算适应度函数.
Step 4.速度、位置更新.对粒子群算法的每一次迭代，速度和位置的更新方程为:
其中:
其中，w是惯性权重，T是当前迭代次数，Tmax是最大迭代次数，j是粒子所在的维度，c1、c2为学习系数，cmax、cmin为c1、c2的最大值、最小值，r1、r2为介于［0，1］的两个相互独立的随机数.
Step 6.判断粒子的位置是否搜索过，若搜索过则粒子当前位置变为附近未搜索过的位置.转到Step 2.
Step 7.终止操作，输出全局最优解.
2.3基于改进粒子群优化的二维Otsu算法
二维Otsu算法是通过搜索到能够使背景类和目标类差别最大的阈值，进行阈值分割.传统Otsu算法采用穷举的方式去搜索，计算量大、耗时多.利用粒子群算法优化的Otsu分割算法，能够大大减少传统Otsu算法的计算量，但是传统粒子群算法存在过早收敛陷入局部最优解的缺点.二维直方图选择上，为节省算法运行时间，避免出现使用灰度级-邻域均值二维直方图出现目标或背景与噪声或边界错分的情况，本文算法采用灰度级-梯度二维直方图；粒子群算法选择上，采用基于粒子空间对称分布的改进粒子群算法去优化二维Otsu算法来避免陷入局部最优解.利用图像在灰度级-梯度二维直方图的分布特性，不仅减少了直方图的计算量，也缩小了粒子的搜索范围，在保证了分割精度的同时，又有效提高了分割的速度.
具体算法的基本步骤如下:
Step 1.生成m个粒子，其搜索空间为二维空间，则第i个微粒的位置可表示为xi=（xi1，xi2），考虑到图像在灰度级-梯度二维直方图的分布特性，其搜索范围为图1（b）图的0、1区域，即xi1在［Gmin，Gmax］随机生成，xi2在［Gradmin，Gradmax］随机生成，每个粒子的速度在［-Vmax，Vmax］随机生成，设置最大迭代次数N，惯性权重最大值、最小值，学习系数最大值、最小值.其中，Gmin、Gmax为图像的灰度级最小值、最大值，Gradmin、Gradmax为图像的梯度最小值、最大值.
Step 2.计算图像的灰度值-梯度二维直方图.
Step 3.判断是否满足迭代条件，满足则转到Step 4，否则转到Step 8.
Step 4.根据式（6）计算粒子中每个粒子的适应值.
Step 5.根据式（8）、（9）更新第t+1次迭代后各个粒子的位置、速度.
Step 6.基于粒子空间分布对称，根据式（13）对粒子进行位置调整.
Step 7.判断粒子的位置是否搜索过，若搜索过，则粒子当前位置变为附近未搜索过的位置.转到Step3.
Step 8.将改进粒子群搜索到的全局最优解作为分割阈值，对图像进行阈值分割，输出图像.
图4 基于改进粒子群的二维Otsu算法流程图Fig.4 Flow chart of the two-dimensional Otsu based improved PSO
3 基于粒子群优化Otsu的肺组织改进自动分割算法
针对传统肺组织分割算法普遍存在分割速度慢的缺陷，本文基于胸腔CT图像中心区域近似孔洞特性，利用2.3节中改进粒子群优化的Otsu算法对胸腔CT图像进行肺组织分割，以此来提高分割速度.
算法具体步骤如下:
Step 1.图像预处理:利用胸腔CT图像中心区域近似孔洞特性，进行孔洞填充操作后得到背景区域，与原图像相减，得到去除背景的预处理图像；
Step 2.阈值分割:采用基于改进粒子群优化的Otsu算法，对胸腔CT图像进行肺组织阈值分割；
Step 3.形态学操作:经过阈值分割后，利用开操作去除噪声，利用闭操作填补病变区域的孔洞，得到二值图像，进而得到分割图像.
4 实验结果和分析
为验证基于粒子群优化的二维Otsu改进算法的分割精度和速度，比较本文算法和传统Otsu算法、文献［13］中基于灰度级-邻域均值二维直方图的粒子群优化的Otsu算法的分割精度和速度.采用来自宁夏第四人民医院的肺结核病临床数据.图像尺寸为512像素×512像素，实验硬件环境为Windows8.1系统、主频2.6 GHz、内存4 GB的PC笔记本，软件开发平台为MATLAB R2008a.
图5 不同阈值算法分割结果对比Fig.5 Results of segmentation by different algorithm
表1 不同算法的最佳阈值、分割时间对比Tab.1 The best thresholds and processing time by different algorithm
图5为使用传统二维Otsu算法、使用文献［13］中的基于灰度级-邻域均值二维直方图的粒子群优化Otsu算法、使用本文的阈值分割算法对图5（a）原图像分割的结果对比.
由图5可得:本文算法与传统二维Otsu算法、文献［13］中使用灰度级-邻域均值二维直方图的粒子群优化的Otsu算法的阈值分割结果近似，均不影响后续的肺组织形态学分割.
图6 肺组织分割算法Fig.6 Lung segmentation algorithm
表1为使用图5中3种不同阈值分割算法获取的最佳阈值、分割时间对比.由上可得:文献［13］中，粒子群优化的传统二维Otsu算法处理时间为传统二维Otsu算法处理时间的1/3左右，由于使用的灰度级-邻域均值二维直方图本身存在缺陷，而本文算法使用的是灰度级-梯度二维直方图，无论在直方图的计算，还是粒子群的搜索，都能减少算法的处理时间，因此本文算法处理时间又比文献［13］算法节省16%.
图6描述的是对胸腔CT图像进行肺组织分割的过程图.图6（a）为胸腔CT原图像，图6（b）为对（a）图填充孔洞操作后的图像，图6（c）为图6（b）减图6（a）得到的去除背景的图像，图6（d）对图6（c）利用本文算法进行阈值分割后的图像，图6（e）对图6（d）进行形态学操作去除噪声、填充孔洞后的图像，图6（f）为最终分割结果.
图7是对3幅不同的胸腔CT图像进行肺组织分割的结果对比.图7（a）、7（c）、7（e）为胸腔CT原图像，图7（b）、7（d）、7（f）分别为对图7（a）、7（c）、7（e）的肺组织分割结果，图7（f）、7（g）、7（h）为利用半径分别为13、22、28的结构元素对图7（e）进行形态学操作后得到的分割结果.
由图6、图7可得:本文算法能够有效地对类似图7（a）、7（c）以及正常的胸腔CT图像分割出肺组织，对类似图7（e）进行分割的分割效果较差，通过增大形态学处理所用结构元素的半径，可以修补大部分孔洞（如图6（f）、7（g）、7（h））.
图7 不同CT图像的肺组织分割结果Fig.7 Lung segmentation results of different CT image
表2为使用图5中的3种Otsu算法分别对胸腔CT原图像1、2、3进行阈值分割的处理时间对比.
表2 不同算法的阈值分割处理时间对比Tab.2 Processing time of different threshold segmentation algorithm
由表2可得，本文算法对于不同胸部CT图片的阈值分割处理时间都较少，相比于使用灰度级-梯度直方图的粒子群优化的OTSU算法，本文算法的处理速度平均提高了16%.
5 结 论
提出了一种基于改进粒子群优化Otsu的胸腔CT图像快速自动分割算法.使用灰度级-梯度二维直方图减少直方图的计算量，减小粒子群的搜索范围，达到了节省分割时间的目的.对图像尺寸为512像素×512像素CT图像的阈值分割时间约为0.2 s，相比于文献［13］中使用灰度级-邻域均值直方图的粒子群优化的Otsu算法，本文算法的处理速度平均提高了16%.而且在分割精度上，与传统的二维Otsu算法以及文献［13］粒子群优化的Otsu算法相比，阈值分割精度几乎一致，不影响肺组织分割结果.
［1］ 孟琭，贾同，赵大哲，等.基于Snake模型的病变CT肺部图像分割［J］.系统仿真学报，2009，21（15）:4603-4606，4612.
Meng L，Jia T，Zhao D Z，et al.Pathological lung segmentation for CT images based on Snake model［J］.Journal of System Simulation，2009，21（15）:4603-4606，4612.（in Chinese）
［2］ Ivanovska T，Hegenscheid K，Laqua R，et al.A fast and accurate automatic lung segmentation and volumetry method for MR data used in epidemiological studies［J］.Computerized Medical Imaging and Graphics，2012，36（4）:281-293.
［3］ Wei Y，Shen G，Li J J.A fully automatic method for lung parenchyma segmentation and repairing［J］.Journal of Digital Imaging，2013，26（3）:483-495.
［4］ Darmanayagam S E，Harichandran K N，Cyril S R R，et al.A novel supervised approach for segmentation of lung parenchyma from chest CT for computer-aided diagnosis［J］.Journal of Digital Imaging，2013，26（3）:496-509.［5］ Otsu N.A threshold selection method from gray-level histograms［J］.IEEE Transactions on System，Man and Cy-bernetic，1979，9（1）:62-66.
［6］ 刘健庄，栗文青.灰度图像的二维Otsu自动阈值分割法［J］.自动化学报，1993，19（1）:101-105.
Liu J Z，Li W Q.The automatic thresholding of gray-level pictures via two-dimensional OTSU method［J］.Acta Automatica Sinica，1993，19（1）:101-105.（in Chinese）
［7］ 吴成茂，田小平，谭铁牛.二维Otsu法阈值法的分割快速迭代算法［J］.模式识别与人工智能，2008，21（6）:746-757.
Wu C M，Tian X P，Tan T N.Fast iterative algorithm for two-dimensional Otsu Thresholding method［J］.Pattern Recognition and Artificial Intelligence，2008，21（6）:746-757.（in Chinese）
［8］ 黄岗，李俊，潘金贵.基于粒子群优化方法的2维Otsu快速图像分割算法［J］.中国图象图形学报，2011，16（3）: 377-381.
Huang G，Li J，Pan J G.A fast 2D Otsu image segmentation algorithm based on particle swam optimization algorithm［J］.Journal of Image and Graphics，2011，16（3）:377-381.（in Chinese）
［9］ Kennedy J，Eberthart R C.Particle Swarm Optimization［C］.Proc.IEEE Int.Conf.on Neural Networks，Piscataway:IEEE，1995:1942-1948.
［10］ 王洪涛，李丹.基于改进粒子群算法的图像边缘检测研究［J］.液晶与显示，2014，29（5）:800-804.
Wang H T，Li D.Image edge detection research based on improved particle swarm optimization algorithm［J］. Chinese Journal of Liquid Crystals and Displays，2014，29（5）:800-804.（in Chinese）
［11］ 孙越泓.基于粒子群优化算法的图像分割研究［D］.南京:南京理工大学，2011.
Sun Y H.Image segmentation research based on particle swarm optimization algorithm［D］.Nanjing:Nanjing University of Science&Technology，2011.（in Chinese）
［12］ 孙越泓，魏建香，夏德深.一种基于粒子对称分布多样性的PSO算法［J］.模式识别与人工智能，2010，23（2）: 137-143.
Sun Y H，Wei J X，Xia D S.An improved PSO based on diversity of particle symmetrical distribution［J］.Pattern Recognition and Artificial Intelligence，2010，23（2）:137-143.（in Chinese）
［13］ Helen R，Kamaraj N，Selvi K，et al.Segmentation of pulmonary parenchyma in CT lung images based on 2D Otsu optimized by PSO［C］.Proceedings of 2011 International Conference on Emerging Trends in Electrical and Computer Technology（ICETECT），Tamil Nadu:IEEE，2011:536-541.
Improved lung segmentation algorithm based on 2D Otsu optimized by PSO
MENG Ya-zhou，MA Yu∗，BAI Bing，GAO Lin-shuang
（School of Graduate，Ningxia University，Yinchuan 750021，China）
This paper develops an automatic segmentation method of lung image using two dimensional Otsu method based on particle swarm optimization in order to reduce operation time.Aiming at the shortage of the traditional two dimensional Otsu based on particle swarm optimization，which the calculating quantity is large and standard particle swarm algorithm is easy to fall into local optimum，an improved lung segmentation based on 2D Otsu optimized by PSO is proposed in this paper.Using the grayscale-gradient two-dimensional histogram，not only reduces the amount of histogram's calculation，but also narrows the searching area of particles.The algorithm uses the improved PSO which based on diversity of particle symmetrical distribution to search optimal threshold.In the process of algorithm，the region filling algorithm is used to remove background in order to make the threshold segmentation of lung better，and the morphology operations to remove noise and repair holes which in the target image.The threshold segmenting time in this algorithm is about 0.2 s，increased about 16%than the speed of the traditional Otsu threshold segmentation optimized by PSO，the size of the CT images is 512×512 pixels in this experiment.The segmentation algorithm in this paper can seg-ment the lung in CT image automatically，not only ensures the accuracy of the segmentation，but also improves the speed of the segmentation in comparison with conventional algorithms.
two-dimensional Otsu；improved PSO；region filling；morphological operations；automatic lung segmentation
TP394.1
A doi:10.3788/YJYXS20153006.1000
1007-2780（2015）06-1000-08
孟亚州（1991-），男，山东邹城人，硕士研究生，主要从事图像处理、模式识别方面的研究.E-mail:xingqing_ myz@163.com
马瑜（1974-），男，宁夏西吉人，工学博士，副教授，硕士生导师，主要从事三维医学图像可视化、模式识别等方面的研究.E-mail:mayu95@163.com
白冰（1990-），男，硕士研究生，主要研究方向为图像处理.E-mail:fengfubaibing@163.com
高林爽（1990-），女，硕士研究生，主要研究方向为图像处理.E-mail:1046935271@qq.com
2014-12-04；
2015-01-06.
宁夏回族自治区2012年科技攻关计划项目（No.2012ZYG011）
Supported by Key Science and Technology Program of Ningxia Hui Autonomous Region，China（No. 2012ZYG011）
∗通信联系人，E-mail:xingqing_myz@163.com
