基于相位一致和核模糊C均值的眼底血管分割

谭春林，曹 鹏，黄 华

（四川大学 电气信息学院，四川 成都610065）

为避免除数为0，傅里叶幅度和再加上一个小数ε；噪声估计值Tn；多种频率成分在一个点的相位一致性才有意义，故使用一个频率扩展函数W（x）［10］，最后相位一致计算如式 （5）所示

0 引 言

视网膜血管是人体内唯一可进行非创伤性直接观察的深层血管。因为全身性血管病 （高血压等）、血液病 （白血病等）及代谢性疾病 （糖尿病等）严重影响视网膜血管［1］，如视网膜血管形态特征 （如直径，长度等）产生物理变化，故视网膜血管的形态特征对这些疾病的诊断、治疗评价具有重要意义。此外视网膜血管也可用做眼底图像配准［2］、视神经盘定位的参考标准［3］。

常见眼底血管分割方法多是基于图像灰度值，包括基于窗口、分类器及跟踪3 种。Zhang等［4］提出基于窗口的算法，将视网膜血管的横截面近似为高斯形曲线，通过旋转匹配滤波器检测血管；Joes Staal等［5］提出一种基于分类器的算法，在每个像素点处提取多种特征并利用KNN 分类器选取的特征构造从而检测出血管；Sun等［6］提出的基于跟踪的自适应算法，依据血管的中心线、方向、直径等建立血管轮廓模型，再根据灰度值得连续性追踪血管。前述分割方法主要不足包括计算耗时、分割结果严重依赖图像亮度和对比度且易受噪声干扰、二值化时不能采取相同的阈值，因此产生了基于图像频域信息的分割方法。Amin和Yan［7］、陶治江等［8］采用相位一致性算法分割眼底血管。Amin和Yan获得的血管图像有一个亮圆环且含大量分散的白色噪声，这些都会影响分割结果；陶志江等的算法采用最大类间方差法确定分割阈值，导致分割结果易受目标与背景大小比例及噪声的影响。

本文提出一种结合相位一致性算法和基于高斯核函数的模糊C均值聚类算法。相位一致算法利用图像相位信息检测血管边缘特征，有效避免图像亮度和低对比度对分割结果的影响；核模糊C 均值算法自动分类检测结果便于二值化。此外，文中采用的边界模板和面积滤波法提高了分割精度。

1 相位一致性及核模糊C均值

1.1 相位一致性原理

通过对方波信号分析可知，在方波上升沿或者下降沿处，各次谐波分量的相位相等。如果方波占空比为50%，则上升沿和下降沿的相位分别为0°和180°。同理，对三角波做傅立叶分解后，各次谐波在波谷和波峰处的相位分别为90°和270°。方波信号的跳变沿可看作是阶跃边缘［9］。方波和三角波的分解及合成信号分别如图1所示。

图1 方波和三角波前4次谐波及其合成信号

相位一致性认为图像信号相位一致性最大处就是图像特征，从频域来看就是图像信号傅立叶级数分解后各个分量在某点处的相位相同或接近则该点就是一个特征点。Venkatesh和Owens提出相位一致性是傅里叶幅度之和对局部能量归一化得到，如式 （1）所示

式中：An——第n次谐波幅度，E（x）——局部能量。当各次谐波具有相同相位时相位一致性取最大值1，若各次谐波相位不一致则取0。采用傅里叶幅度之和来归一化局部能量使得相位一致算法免受图像亮度、对比度的干扰；图像可采用相同阈值作二值化。

设信号I（x），尺度为n的Log Gabor小波的偶对称和奇对称部分分别为Men和Mon，则尺度n的小波滤波响应、幅度分别如式 （2）、式 （3）所示

所以基于LOG Gabor小波局部能量的计算如式 （4）所示

为避免除数为0，傅里叶幅度和再加上一个小数ε；噪声估计值Tn；多种频率成分在一个点的相位一致性才有意义，故使用一个频率扩展函数W（x）［10］，最后相位一致计算如式 （5）所示

1.2 核模糊C均值原理

聚类分析依据相似性将没分类的目标划分为不同的子类，使类内相似性尽量大。相位一致性检测得到图像像素大致分为3类，即边界区域、目标区域内的背景、目标区域内的血管，因此可在图像二值化之前采用模糊聚类算法提取出血管图像。一般的模糊聚类算法对初始值敏感且结果依赖样本的空间分布，故采用基于核函数的模糊C 均值算法。基于核的聚类算法通过构造新的特征向量将输入数据映射到高维特征空间，将原空间非线性信息转化为线性问题，适合任意形状分布的样本，利用模糊化矩阵来更新核函数的参数就避免了求聚类样本平均值带来的聚类中心的不合理性［11］。

高斯核函数和基于高斯核函数的C 均值模糊聚类算法的目标函数分别如 （6）、式 （7）所示

式 （6）y 为核函数中心，x 为数据点，σ为控制函数径向作用范围的宽度参数；式 （7）中xj（j＝l，2，3…n）表示图像像素的灰度值，c是聚类个数，uij是像素k 属于第i类的隶属度，m 是隶属度加权指数，Vi（i＝1，2，3…c）是当前聚类中心；目标函数的约束条件是1，2…n。样本在映射空间的欧式距离

求目标函数的极值得到隶属度函数uij和聚类中心vi的计算公式如式 （8）、式 （9）所示

综上，基于高斯核函数的模糊C 均值聚类算法的计算步骤如下：

（1）初始化聚类参数：设定聚类数c，迭代阈值T，最大迭代次数iter＿max，初始化迭代次数iter＝0，隶属度矩阵U 和聚类中心V；

2 眼底图像分割

2.1 眼底图像预处理

STARE 图像由TopCon－50眼底照相机拍摄，24 位彩色图像，分辨率605×700，格式PPM，其中20 张含有眼科医生手工分割结果，作为分割评判的参考标准。彩色眼底图像含有红、绿、蓝3个分量。由于红色分量亮度最高，血管和背景对比度低，不易将目标血管和眼底背景区分；蓝色分量的对比度和亮度均较低，且噪声干扰严重，因此也不适合用于眼底血管分割；绿色分量的亮度适中，而且血管与背景对比度较高，能很好反应彩色眼底图像血管分布。图2中 （a）、 （b）分别为彩色眼底图像转换的灰度图和该彩色图像的绿色分量。因此，眼底图像分割通常采用绿色分量作为分割目标。

图2 眼底图像及分割结果

一幅眼底图像大致可分为两个部分，亮度及对比度较高的圆形区域，包含了眼底图像中的全部血管和病变区域，是设计分割算法时的主要考虑对象，文中称为目标区域；另一部分是目标区域的外周区域，从图像来看，该区域不含眼底血管等分割目标称为边界。边界区域内的图像灰度值相近但有一定变化，对目标区域内血管分割结果有较大影响。为去除边界的影响，可经绿色分量提取边界模板，用边界模板对绿色分量做掩码操作可以有效去除该影响。此外，由于经过边界模板掩码操作后的绿色分量图像，目标区域和边界区域的灰度值相差较大，形成一个圆形的边界，在文中的相位一致算法作用下该圆形边界会形成一个量环，因此还需将相位检测结果做掩码操作得到只含有目标区域的眼底图像。边界模板如图2 （c）所示。

2.2 视网膜血管分割

文中视网膜血管分割，主要分为5 个步骤：①边界模板提取，依据边界区域灰度值低于目标区域采用一定阈值即可得到模板；②利用模板去除绿色分量图像的边界；③采用相位一致性算法对去除边界的绿色分量图像做边缘检测，得到眼底图像中血管边缘；④采用高低帽算法增强相位一致检测结果，采用模板去除图像的边界噪声；⑤采用KFCM 聚类算法提取血管图像，图像二值化并进行面积滤波。

图2包括6 幅子图，与分割结果相关的是图2 （d）、（e）、（f），其中图2 （d）为专家手工分割结果，图2 （e）不采用模板处理的分割结果，图2 （f）为采用模板处理的分割结果。

相位一致性算法得到的图像，目标 （血管）灰度值较大故较亮，背景灰度值较低故较暗，可采用高低帽算子分别作用于图像，利用两者差值作为增强图像。高低帽变换（top－hat and bottom－hat transforms）是基于形态学中膨胀和腐蚀的高通滤波算法。高帽算子对图像先腐蚀后膨胀，再用原图像减去高帽算子作用结果，得到原图像的灰度峰值。低帽算子先膨胀后腐蚀，再减去原图像，最后得到原图像的灰度谷值。

算法中用了两次边界模板。因为绿色分量图像边界区域的灰度值并不完全相等，就会影响相位一致性算法检测结果，如图2中 （e）所示的亮圆环外的两条白线及散乱分布的白色团块，故步骤②采用边界模板消除边界区域；又因模板消除边界区域后，其实在边界与目标区域交界处相当于形成了一个阶跃的边界，因此相位检测算法在边界处得到较大的相位一致性，形成一个亮圆环，如图2 （e）所示，故步骤④再次采用边界模板消除边界区域外的图像。图2中 （e）所示，目标区域内还散乱分布着一些小的白色团块，这是由目标区域的噪声引起的，文中采用面积滤波方法去除。边界模板和面积滤波方法可以有效提高分割精度，实验结果如表1和图3所示。

表1 分割结果准确率

图3 分割正确率的柱状

2.3 基于STARE数据库的实验分析

图2 （f）直观的展示了分割结果，并可与参考标准做简略比较。与图2 （e）比较，模板处理后的分割结果图2（f）噪声团块更少，且边界区域也被有效去除；与图2 （d）比较，算法得到较好分割结果。眼底血管分割一般采用分割正确率来评价算法性能。分割正确率等于正确分割像素数与图像像素总数之比。正确分割像素数包含两部分，一部分是正确分割的血管像素数即分割结果和参考标准都含有血管像素，另一部分是正确分割的非血管像素数即分割结果和参考标准都含有的非血管像素。采用提出的算法对STARE图库中的20幅含有专家分割标准的眼底图像分割，分割正确率示于表1，表中ACC1表示不采用模板的分割正确率，ACC2表示采用模板后的分割正确率，Im0001等表示图像在STARE图像库中的编号。

分析表1知，算法对不同图像的分割正确率略有差异，因为不同眼底图像包含有不同大小和方向的血管；比较ACC1和ACC2的数据知，没用模板的分割准确率最低为88% （Im0139），最高达到93.66% （Im0291），均值为90.36%；采用模板后的分割正确率最低为89.45%（Im0044），最高可达95.97% （Im0291），均值为92.65%；边界模板对不同图像分割正确率贡献不同，如分割正确率提高最小为1.68% （Im0002），分割正确率提高最大是2.31% （Im0291），平均提高2.29%；

表1的数据用柱状图示于图3，可以观察出采用模板处理的分割结果要明显优于未采用模板处理的结果；若以分割准确率平均提高约为2.29%计算，则采用模板后，正确分割的像素数增加了9698个。

3 结束语

图2的分割结果表明基于PC 和KFCM 的算法能完成眼底血管视网膜分割，表1中的分割正确率从统计数据的角度证实所提出算法的有效性。PC 算法很好检测出图像边缘特征，且免受图像亮度影响并且检测结果不依赖图像的对比度；利用KFCM 对PC 检测结果做聚类时不必考虑图像二值化阈值从而简化了分割算法。提出的算法能够得到较好分割结果，但是相位一致性算法还需进一步优化，使得边缘定位更准确才能分割出更细小的毛细血管；另外可探索其它方法来消除血管空洞。

［1］Ravishankar S，Jain A，Mittal A.Automated feature extraction for early detection of diabetic retinopathy in fundus images［C］//IEEE Conference on Computer Vision and Pattern Recognition.CVPR.IEEE，2009：210－217.

［2］Deng K，Tian J，Zheng J，et al.Retinal fundus image registration via vascular structure graph matching ［J］.Journal of Biomedical Imaging，2010：14.

［3］Abdel Razik Youssif AH，Ghalwash A Z.Optic disc detection from normalized digital fundus images by means of a vessels’direction matched filter［J］.IEEE Transactions on Medical Imaging，2008，27 （1）：11－18.

［4］Zhang B，Zhang L，Zhang L，et al.Retinal vessel extraction by matched filter with first－order derivative of Gaussian ［J］.Computers in Biology and Medicine，2010，40 （4）：438－445.

［5］Fraz M M，Remagnino P，Hoppe A，et al.Blood vessel segmentation methodologies in retinal images－A survey ［J］.Computer Methods and Programs in Biomedicine，2012，108 （1）：407－433.

［6］Sun K，Chen Z，Jiang S，et al.Morphological multiscale enhancement，fuzzy filter and watershed for vascular tree extraction in angiogram ［J］.Journal of Medical Systems，2011，35（5）：811－824.

［7］Ashraful Amin M，Yan Hong.High speed detection of retinal blood vessels in fundus image using phase congruency ［J］.Soft Couputing，2011 （15）：1217－1230.

［8］TAO Zhijiang，GAO Yuan，HUANG Hua.Retinal vasculature’s segmentation in fundus image based on phase congruency［J］.Computer Engineering and Design，2012，33 （8）：3139－3143 （in Chinese）. ［陶治江，高原，黄华.基于相位一致性的眼底图像视网膜血管分割 ［J］.计算机工程与设计，2012，33 （8）：3189－3143.］

［9］TIAN Weiqing，ZHU Weigang，YANG Shunlong，et al.SAR image edge detection based on phase congruency［J］.Journal of Sichuan Ordnance，2012，33 （9）：89－91 （in Chinese）. ［田微晴，朱卫岗，杨顺隆，等.基于相位一致性的SAR 图像边缘检测［J］.四川兵工学报，2012，33 （9）：89－91.］

［10］Henriksson L，Hyvrinen A，Vanni S.Representation of cross－frequency spatial phase relationships in human visual cortex ［J］.The Journal of Neuroscience，2009，29 （45）：14342－14351.

［11］Yang Miin－Shen，Tsai Hsu－Shen.A Gaussian kernel－based fuzzy c－means algorithm with a spatial bias correction ［J］.Pattern Recognition Letters，2008，29 （12）：1713－1725.