高玮玮 沈建新* 王玉亮 梁 春 左 晶 2
(南京航空航天大学机电学院,南京 2 10016)2(江苏省中医院眼科,南京 2 10029)
糖尿病是严重影响人类健康的内分泌疾病,致残致死率仅次于心脑血管疾病及癌症。据预测,到2030年全球糖尿病患者将达到总人口数的4.4%[1]。该疾病不仅给人们带来巨大痛苦,且带来众多并发症,其中,以糖尿病视网膜病变(diabetic retinopathy,DR,简称“糖网”)发生率最高,对视力影响也最大,已成为目前,特别是发达国家20~74岁成人致盲的首要原因[2]。由于每个糖尿病患者都有发展为DR的可能,而DR具有进行性、不可逆性,因此,如何准确筛查无明显视力损伤的糖尿病患者是否存在DR,不仅为早期诊断、早期治疗从而挽救患者视功能提供先机,还可节约大量社会医疗资源。目前国内外已有多种筛查方法[3-4],其中免散瞳眼底照相检查,由于具有简单易行、价格低廉、图像易得、直观、易于保存和记录,且与眼底荧光造影诊断结果,具有显著一致性和较高的灵敏度、特异性等优点,被眼科医生普遍认为最适用于DR的筛查和随诊[5]。但现阶段该方法基本依靠眼科医生对眼底图像的肉眼观察,这对实施大规模DR筛查,存在明显不足与局限。因此,若能借助计算机快速、可靠地自动识别出眼底图像中的DR病灶,不仅可将医生从繁重的人工阅片工作中解脱出来,更为DR筛查的大规模实施提供必备的基础条件。临床上通常将DR分为非增殖期、增殖前期以及增殖期[6]。其中在非增殖期,患者视网膜会出现微动脉瘤(microaneurysm,MA)、视网膜内出血(haemorrhage,H)以及毛细血管渗漏引发的硬性渗出(hard exudates,EXs)。MA虽也可在其他疾病中见到,但它仍是DR的特征性表现,是该病最早的可靠特征,它的出现与否以及数目可作为DR等级分类的依据[7]。因此,该病灶的自动检测对于建立基于眼底图像的DR自动筛查系统尤为关键。
针对眼底图像中MA的自动检测,国外许多学者对此进行研究并提出了相关算法,而国内目前没有相关研究的报道。Spencer等[8]、Cree等[9]、以及Frame A等[10]都提出了针对眼底荧光血管造影图像的MA自动检测算法;Hipwell等将此类算法进行了优化,并应用到了高分辨率无赤光眼底图像MA的自动检测中[11];Sinthanayothin等利用循环区域生长法以及人工神经网络实现了散瞳眼底图像中MA的自动检测[12];Larsen等提出了商业化的DR红色病灶(MA和H)自动检测系统,但没有公开相关算法[13];Niemeijer等在对免散瞳眼底图像中每个像素分类从而获取DR红色病灶候选区域的基础上,利用K-近邻分类器获得真DR红色病灶[14];Martins等利用多层感知器(MLP)神经网络,实现免散瞳眼底图像中MA的自动检测[15];Hatanaka等利用double-ring滤波器以及特征分析,实现了免散瞳眼底图像中MA的自动检测[16];Narasimhan等分别通过支持向量机以及贝叶斯分类器,实现眼底图像中MA、H以及EXs的自动识别[17]。
为实施大规模的DR筛查,必须建立以免散瞳眼底图像为对象的DR自动筛查系统。而目前已有MA自动检测算法,多是针对眼底荧光血管造影图像、无赤光眼底图像以及散瞳眼底图像。与它们相比,免散瞳眼底图像质量相对较差,MA自动检测难度更大,故文献[8-12]中的算法并不适用于免散瞳眼底图像中MA的自动检测。文献[13-14]将MA与H作为一个整体进行检测,这在一定程度上忽略了MA作为DR最早期病症对该疾病的预示价值。因此,为构建基于免散瞳眼底图像的DR自动筛查系统,在已实现免散瞳眼底图像中EXs自动检测[18]的基础上,提出了一种新的适用于免散瞳眼底图像中MA自动检测的算法。该算法充分利用免散瞳眼底图像中MA的形态特征,采用扩展极小值(extended-minima,EMIN)变换[19]获取MA候选区域;然后从中去掉EXs和血管以消除由这二者造成的假阳;最后根据MA的尺寸信息获取真正的MA。
源图像由江苏省中医院眼科提供,为Canon CF-60DSI眼底照相机拍摄的100幅分辨率为3 504像素×2 336像素的JPG彩色免散瞳眼底图像(目前,国际上还没有用于评价DR病灶自动检测算法性能的公共标准图像数据库[20])。这些图像具有不同的颜色、亮度和质量。实际使用时,对原图进行了适当裁剪,去掉背景区域并进行了压缩,即实际处理图像的分辨率为750像素×500像素。据眼科医师判断,除45幅健康眼底图像外,其余55幅眼底图像中均出现了MA(在彩色免散瞳眼底图像中MA表现为边界清楚的红或暗红色小圆点,其直径常在15~60 μm,即视盘边缘视网膜静脉直径的1/8~1/2,也偶有较大者,可至视网膜静脉直径,但一般不超过 1 25 μm[7])。MA以及 E Xs在彩色免散瞳眼底图像中的具体形态见图1。对于出现MA的55幅眼底图像,由该眼科医师对其中MA的相关位置进行手工标注,所得标注结果将用于评价算法对该病灶的自动检测性能。
图1 彩色免散瞳眼底图像及MA和EXs细节Fig.1 Color non-dilated fundus image and detail of MA and EXs
图2 算法流程图Fig.2 Flowchart of the proposed algorithm
实验系统配置为:Intel(R)Core(TM)Duo E7500 CPU,6.00 GB RAM的计算机,Matlab R2009a的软件环境。
据临床观察,MA以直径小于125 μm的红或暗红色孤立圆点存在于视网膜上,在免散瞳眼底图像G通道(记为fg,见图3(a))中表现为区域内部具有连续灰度值但其外部边缘像素灰度值确定更高的孤立区域,因此,可利用EMIN变换将其从免散瞳眼底图像G通道中检测出来;但EXs以及血管中也存在满足此特征的区域,故还需将这二者从EMIN变换后的结果中去除;最后在得到的MA候选区域中利用MA的尺寸信息获取真正的MA。因此,文中算法对免散瞳眼底图像中MA的自动检测主要包括4步骤:免散瞳眼底图像G通道预处理;EXs自动检测与去除;血管自动检测及去除;根据尺寸信息实现MA的自动检测。具体算法流程见图2。
1.2.1 免散瞳眼底图像G通道预处理
免散瞳眼底图像在获取、传输和存储的过程中,不可避免地受到各种噪声源的干扰。这些噪声的存在使得所观测到的图像会出现诸如光照不均、对比度不足等缺陷,严重影响了图像的应用效果。因此,为获得理想的分割效果,在进行相关目标分割前,有必要对图像进行预处理。
首先对fg1进行中值滤波以消除噪声,结果记为fg。另外,为增强显示效果,提高图像对比度,需进行灰度均衡化。灰度均衡化可以在整张图像内进行,也可以分区域进行,但都有一定的局限与不足,因此,此处选用对比度受限自适应直方图均衡(contrastlimited adaptive histogram equalization,CLAHE)方法来增强眼底图像的对比度[21]。该方法结合自适应直方图均衡与对比度受限两种方法,因而能够在最大限度提高图像对比度的同时,对图像的视觉表现产生较小影响,从而使得后续相关目标检测的准确性得到显著提高。fg1的CLAHE处理结果记为fg2(如图3(b)所示),预处理前后 M A局部细节分别见图3(c)、(d)。
1.2.2 EXs自动检测
由于EXs常数个或数十个呈簇状堆积,有时相互融合成片,有时排列成环状,因此,它们之间以及内部会有孔洞存在,这些孔洞会被EMIN变换检测出来,故在MA自动检测前必须进行EXs的自动检测。
对于免散瞳眼底图像中EXs的自动检测采用文献[18]提出的算法:首先利用大津法(Otsu)阈值分割以及数学形态学快速提取出视盘,在此基础上利用方差滤波、形态学腐蚀重建等获取EXs候选区域;然后利用形态学膨胀重建等获取精确的病灶轮廓,从而获取真正的EXs。该算法复杂度低、运行速度快、鲁棒性强、分割效果好,图3(a)的具体检测结果,如图4。
图3 免散瞳眼底图像G通道预处理。(a)免散瞳眼底图像G通道;(b)预处理结果;(c)预处理前MA局部细节;(d)预处理后MA局部细节Fig.3 Preprocessing of green channel.(a)Green channel oforiginalnon-dilatedfundusimage;(b)Result of preprocessing;(c)Closeup of MA before preprocessing;(d)Closeup of MA after preprocessing
图4 EXs自动检测结果。(a)检测出的EXs;(b)检测结果叠加在原彩图上的局部细节图Fig.4 AutomateddetectionresultofEXs.(a)Detection result of EXs;(b)Detail of detection result superimposed over the original image
1.2.3 血管自动检测
由MA病理特征可知,其不会出现在血管上。但由于光照不均等因素,免散瞳眼底图像中的血管上会出现孔洞、血管片断等,此类区域亦满足EMIN变换所能检测区域的特性,因此,血管也是必须去除的假阳。对于血管的自动检测,通过对两幅图像的差值进行Otsu阈值分割实现,具体算法描述如下:
(1)对fg2进行灰度形态学闭运算以消除血管以及相关细节信息(结果如图5(a)所示),该形态学操作所用结构元素尺寸要大于图像中血管最大宽度,即
式中,φ表示灰度形态学闭运算,sB表示大小为s的形态结构元素B。
(2)利用灰度形态学腐蚀重建来填充fg2中存在的孔洞(结果如图5(b)所示),例如红或暗红色小圆点,以及血管上的小孔洞,即
表示灰度形态学腐蚀重建,fm为重建运算中的marker,为 m ask。
在图5(c)中存在一些孤立的白色小区域,利用二值形态学开运算将其中小于等于11个像素的区域去除(据统计,在分辨率为750像素×500像素的免散瞳眼底图像中,直径小于125 μm的MA面积不超过11个像素),即可得血管fvessel,如图5(d)所示。
1.2.4 MA自动检测
MA的自动检测具体描述如下:
(1)对 fg2进行 E MIN变换,即
式中,∧表示逐点求取最小值,-表示取反,fEXs为EXs检测结果。
(3)由于在分辨率为750像素×500像素的眼底图像中MA面积不超过11个像素,因此,利用二值形态学开运算获取fg7中面积小于等于11个像素的区域即得MA,结果如图5(g)所示,将其叠加在原彩图上的效果如图5(h)。检测结果细节见图6,其(a)中箭头所指向的MA较不明显,但采用本算法仍可将其检测出来,如(b)中菱形所示。
对MA自动检测算法进行评价时应结合临床要求,故关注的是算法能否将MA快速、可靠地从眼底图像中检测出来,从而准确筛选出DR病人。因此,文中对算法性能的相关评价指标是分别基于图像水平以及病灶区域水平定义的,具体见式(6)~(9)。
图5 MA检测过程。(a)预处理图像的闭运算结果;(b)孔洞填充结果;(c)血管候选区域;(d)血管分割结果;(e)扩展极小值变换结果;(f)MA候选区域;(g)MA检测结果;(h)检测结果叠加在原彩图效果Fig.5 The detection of MA.(a)Result of closing on preprocessed image;(b)Result of filling in holes;(c)Candidate regions of vessel;(d)Segmentation result of vessel;(e)Result of extended-minima transform on preprocessed image;(f)Candidate regions of MA;(g)Detection result of MA;(h)Detection result superimposed over the original color image图6 MA检测结果细节。(a)原彩图;(b)检测结果Fig.6 Detail of detection result of MA.(a)Original color image;(b)Detection result图7 t不同取值对MA检测结果的影响。(a)t不同取值对图像水平检测结果的影响;(b)t不同取值对病灶区域水平检测结果的影响Fig.7 Influence of t on detection result of MA.(a)Influence of t on detection result based on images;(b)Influence of t on detection result based on lesions
式中,TP为真阳性,FP为假阳性,TN为真阴性,FN为假阴性。具体解释为:对于算法检测为包含 MA的眼底图像,根据其与眼科医师人工判断结果的一致、不一致,而分别称为真阳和假阳;对于算法检测为正常的眼底图像,根据其与眼科医师人工判断结果相同、不相同,而分别称为真阴和假阴。将算法对MA区域的检测结果与眼科医生的人工判断结果相比时,由于算法最终给出的都是其判断为MA的区域,因此只存在真阳、假阳,以及漏检的假阴,此时不存在真阴。因此,分别选用图像水平的灵敏度(sensitivity)、特异性(specificity)、准确率(accuracy)以及病灶区域水平的灵敏度、阳性预测值(predict value)评价算法性能。另外,结合处理效率,即处理一幅图像(分辨率为750像素×500像素)的平均时间来评价算法性能。
文中算法利用的是免散瞳眼底图像中MA的形态特征,故除EMIN变换阈值t外,其他相关参数主要与处理对象的分辨率相关,因此,在实际情况下根据处理对象的分辨率选取相关参数即可。对于EMIN变换阈值t的选取,由多幅眼底图像检测结果(见图7)发现,当其取0.05时MA检测结果最佳。
将100幅免散瞳眼底图像依图像质量分为两组,其中60幅被认定为质量较好,命名为 A组;其余40幅质量相对较差,命名为B组。利用文中算法以及文献[15]提出的MLP神经网络方法对这两组图像进行 MA自动检测,具体评价指标见表1。由表1发现,文中算法对两组免散瞳眼底图像的检测结果均优于 MLP神经网络方法,特别是处理效率,文中算法处理一幅图像仅需要9.7s(且包括EXs自动检测过程),明显优于MLP神经网络方法的12.8s。
另外,根据表1数据分别计算算法对两组图像检测结果的相对误差,从而检验算法稳定性。相对误差根据以下公式计算:
式中,cA表示算法对A组图像检测结果的某项指标,cB则为同种算法对B组图像检测结果的对应指标,由此可分别获得算法对两组图像检测结果相应指标间的相对误差,具体见表2。由表2发现,文中算法对两组图像检测结果各项评价指标的相对误差均低于4%,且明显优于MLP神经网络方法。因此,文中算法具有较好的稳定性。
表1 检测结果Tab.1 Detection result
表2 两组图像检测结果的相对误差Tab.2 Relative error of detection result bet ween two groups of images
由前面的算法描述可知,利用EMIN变换获取的MA候选区域中,存在由EXs以及血管造成的假阳。由于前期已实现的EXs自动检测算法能够充分检测出图像中的EXs,因此在最终的检测结果中并未出现由EXs造成的假阳,但出现了一些位于血管上的假阳,这主要是由无法检测的微弱血管片断以及对比度差的血管末梢造成的。因此,文中算法的检测结果在一定程度上依赖血管检测精度。虽然目前已有很多成熟的血管分割算法,但文中血管分割部分并未采用,而是提出了新的简单、快速的血管分割算法,这主要为兼顾效率指标。EXs的自动检测已占用了一定的时间(平均每幅图像需要9.1s),如果血管分割部分为了追求精度继续花费大量的时间,则即使最后的MA检测结果达到100%的准确性,对临床来说也已无甚意义。且依据英国糖尿病协会(British Diabetic Association)于1997年提出的DR病灶自动筛查算法最低标准[22]:图像水平灵敏度80%、特异性95%而言,文中算法已达到该标准的要求。另外,Javitt等[23]在其所建立的健康准则模型中对DR自动筛查技术提出的基本准则是:60%的灵敏度或是在此基础上更大化的成本效率。即对DR病灶的自动检测在灵敏度等相关指标达到一定水平的基础上,也就是说在保证可靠性的基础上,更关注效率指标。因为对于糖尿病患者而言,其视网膜的病变情况会受糖尿病的治疗、控制以及病程等诸多因素的影响,因此,只有保证较高的检查频率,定期检查眼底才能及时、有效地发现视网膜的相关病变。而较高的检查频率需要高效的DR病灶自动检测算法作为保障,因此,文中采用的血管分割算法是合理且有效的。当然,若能进一步研究出更加简单、高效的血管分割算法以及更加快速、有效的EXs自动检测算法,对整个DR自动筛查系统而言都是意义重大的。这也是该系统最终能否应用于临床必须解决的关键问题。
MA是DR的特征性表现,是DR最早的可靠特征。因此免散瞳眼底图像中MA的自动检测,是建立基于免散瞳眼底图像的DR自动筛查系统关键环节之一,直接决定了系统筛查结果的准确性。只有简单、高效的MA自动检测算法,才能满足临床中对筛查速度与精度的要求。文中算法充分利用免散瞳眼底图像中MA的形态特征,采用EMIN变换获取其候选区域,并结合简单有效的EXs以及血管自动检测算法,从而实现免散瞳眼底图像中MA的高效自动检测。该算法在保证较高检测精度的基础上,很好地兼顾了效率指标,且算法稳定可靠,具有很高的实用价值。该研究工作为构建基于免散瞳眼底图像的DR自动筛查系统提供了算法支持,为整个筛查系统奠定了良好的基础。对整个筛查系统而言,在实现了DR非增殖期的EXs、MA自动检测基础上,视网膜内出血H的自动检测将是下一步的研究目标。
[1]Wild S,Roglic G,Green A,et al.Global prevalence of diabetes:estimates for the year 2000 and projections for 2030[J].Diabetes Care,2004,27:1047-1053.
[2]Watkins PJ.ABC of diabetic retinopathy[J].British Medical Journal,2003,7(2):105-107.
[3]Villarroel M,Ciudin A,Hernndez C,et al.Neurodegeneration:An early event of diabetic retinopathy[J].World J Diabetes,2010,1(2):57-64.
[4]彭金娟,邹海东,王伟伟,等.上海市北新泾社区糖尿病视网膜病变远程筛查系统的应用研究[J].中华眼科杂志,2010,46(3):258-262.
[5]蒋升,张莉,克麦尔·艾则孜.免散瞳眼底照相技术在糖尿病视网膜病变筛查中的应用研究[J].国际眼科杂志,2008,8(10):2037-2039.
[6]陈喆,张士胜,朱惠敏.糖尿病视网膜病变的国际临床分类分析[J].国际眼科杂志,2011,11(8):1394-1401.
[7]Dupas B,Walter T,Erginay A,et al.Evaluation of automated fundus photograph analysis algorithms for detecting microaneurysms,haemorrhages and exudates,and of a computerassisted diagnostic system for grading diabetic retinopathy[J].Diabetes& Metabolism,2010,36(3):213-220.
[8]Spencer T,Olson JA,McHardy KC,et al.An image-processing strategy for the segmentation and quantification of microaneurysms in fluorescein angiograms of the ocular fundus[J].Comp Biomed Res,1996,29:284-302.
[9]Cree MJ,Olson JA,McHardy KC,et al.A fully automated comparative microaneurysm digital detection system[J].Eye,1997,11:622-628.
[10]Frame A,Undrill P,Cree M,et al.A comparison of computer based classification methods applied to the detection of microaneurysmsin ophthalmic fluorescein angiograms[J].Comput Biol Med,1998,28:225-238.
[11]Hipwell JH,Strachan F,Olson JA,et al.Automated detection of microaneurysms in digital red-free photographs:a diabetic retinopathy screening tool[J].Diabetic Medicine,2000,17:588-594.
[12]Sinthanayothin C,Boyce JF,Williamson TH,et al.Automated detection of diabetic retinopathy on digital fundus image[J].Diabetic Medicine,2002,19(2):105-112.
[13]Larsen M,Godt J,Larsen N,et al.Automated detection of fundus photographic red lesions in diabetic retinopathy[J].Investigat Ophthalmol Vis Sci,2003,44(2):761-766.
[14]Niemeijer M,Ginneken B,Staal J,et al.Automatic detection of red lesions in digital color fundus photographs[J].IEEE Trans Med Imaging,2005,24(5):584-592.
[15]Martins CIO,Veras RMS,Ramalho GLB,et al.Automatic microaneurysm detection and characterization through digital color fundus images[C]//Proceedings of the II Workshop on Computational Intelligence.Salvador:e Scholarship,2010:45-50.
[16]HatanakaY,Inoue T,Okumura S,etal.Automated microaneurysm detection method based on double-ring filter and feature analysis in retinalfundus images[C]// 25th International Symposium on Computer-Based Medical Systems.Rome:IEEE,2012:1-4.
[17]Narasimhanm K,Neha VC,Vijayarekha K.An efficient automated system for detection of diabetic retinopathy from fundus images using support vector machine and bayesian classifiers[C]//InternationalConference on Computing,Electronicsand Electrical Technologies.Kumaracoil:IEEE,2012:964-969.
[18]高玮玮,沈建新,王玉亮.基于数学形态学的快速糖尿病视网膜病变自动检测算法[J].光谱学与光谱分析,2012,32(4):760-764.
[19]Soille P.Morphological Image Analysis:Principles and Applications[M]//Berlin:Springer-Verlag,1999:170-171.
[20]Sáanchez CI,Hornero R,López MI,et al.A novel automated image processing algorithm for detection of hard exudates based on retinal images analysis[J].Med Eng Phys,2007,30(3):350-357.
[21]Rosenman J,Roe CA,Cromartie R,etal.Portalfilm enhancement:Technique and clinical utility[J].Int J Rdaiat Oncol Biol Physics,1993,25(2):333-338.
[22]British Diabetic Association.Retinal Photography Screening for Diabetic Eye Disease[R].London:British Diabetic Association,1997.
[23]Javitt JC,Canner JK,Frank RG,et al.Detecting and treating retinopathy in patients with typeⅠdiabetes mellitus.A healthy policy model[J].Ophthalmology,1990,97(4):483-494.