王天宇 石征锦* 黄 钲 宋国立 赵忆文
1(沈阳理工大学自动化与电气工程学院 辽宁 沈阳 110159) 2(中国科学院沈阳自动化研究所 辽宁 沈阳 110016) 3(中国科学院机器人与智能制造创新研究院 辽宁 沈阳 110016) 4(中国科学院大学 北京 100049)
病理性近视(PM)是高度近视的一种,典型症状为屈光度进行性加深、眼轴不断增长、视网膜及脉络膜组织进行性损害引起的视功能障碍,多伴有眼后极部的巩膜变薄,并产生多发病理性改变,如后巩膜葡萄肿、漆裂纹、Fuchs斑、色素上皮萎缩、视网膜下出血、脉络膜新生血管等,严重时将导致视网膜脱落而造成失明[1-3]。研究表明,该病的发病率在逐年增加并且趋向于年轻化[4]。其在中国大陆造成完全失明的人数多达670万,患有视力障碍的人数多达710万[5]。一般而言,多数眼部疾病都有一定的治疗窗口期[6],超过治疗窗口期容易造成不可逆性视力障碍,因此针对病理性近视进行及时的检查以及提前预防具有重要的临床意义。
目前针对病理性近视的诊断大多依靠资深眼科医生根据患者眼部的全面体检结果进行人工诊断,该过程不仅耗时耗力,而且在资深眼科医生紧缺、医疗条件不足的发展中国家或者国内的相关贫困地区将很难做到准确诊断,从而因为病情的延误而造成不可逆的视力损失。据统计全球大概有89%的视力障碍患者生活在低收入和中等收入国家[7],因此欠发达国家和相关贫困地区的视力障碍和失明现象依然严峻,所以需要一种高效、自动的机器诊断方法,在不需要大规模的人力参与以及医疗器械介入的情况下可以辅助医生做出及时的诊断决策,为以后远程辅助医疗打下坚实的基础。
近年来,随着人工智能的发展以及眼底图像技术的成熟,国内外学者开发了大量针对眼部疾病的计算机辅助诊断系统,包括对糖尿病性视网膜病变、年龄相关性黄斑病变、青光眼等疾病的诊断。Liu等[8]提出PAMELA(Pathological Myopia Detection Through Peripapillary Atrophy)系统,其自动接收视网膜眼底图像,并执行感兴趣区域(ROI)提取和视神经乳头的分割随后根据毛细血管萎缩(PPA)这一特征在含有80幅眼底图像的数据集上使用支持向量机(SVM)对病理性近视进行自动诊断。Benghai等[9]描述了PAMELA系统中的纹理分析和灰度分析模块。然后使用决策引擎对两个独立的预测结果进行融合,以获得整体分析,最终进一步证明了该系统的优越性能。Zhang等[10]利用MRMR(Minimum Redundancy-Maimum Relevancy)特征选择技术,对候选特征集进行选择并且排列,最后使用SVM分类器进行分类,从而对病理性近视进行诊断。Xu等[11]提出使用特征袋(bag-of-feature)和稀疏学习的方法进行病理性近视的自动诊断。其对数据集中的图像提取SIFT特征,把训练集的所有SIFT特征聚类为K类,构成码本;针对所有数据集中的每一幅图片的SIFT特征,统计视觉词的出现次数,得到相应的直方图;将直方图进行稀疏学习作为样本向量来构建SVM的训练数据和测试数据,从而进行病理性近视的自动诊断。目前来看针对病理性近视的自动诊断方法具有以下3个缺陷:
(1) 分类模型构建过程中所使用的样本数量较少,所以泛化能力相对较差。
(2) 大多采用机器学习的方法,其需要对特征进行手动提取以及筛选,工作量相对较大。
(3) 采用传统的眼底图像进行病理性近视的识别,而传统的眼底图像检测区域较窄。
本文采用深度学习技术在超广角眼底图像上实现病理性近视的高效自动诊断,在一定程度上解决了传统人工诊断的烦琐、低效,以及因医疗设备不足、医生匮乏等因素造成病情延误等问题。相比于先前的研究以及在临床上的应用,本文具有以下3个优点:
(1) 采用较大的数据集,模型的泛化能力较强。
(2) 卷积神经网络强大的特征自动提取能力,省去了特征工程这一烦琐的操作,且深度学习泛化能力较强。
(3) 数据采用超广角眼底图像,超广角图像检测区域更广泛,所以针对超广角眼底图像的自动识别相对于传统眼底图像的识别有着更加重要的意义。
本文方法首先进行图像分割,识别视盘并以其为中心对原始图像进行裁剪,以去除高干扰区域,提取目标区域;然后使用数据增强方法对数据进行扩充,减少过拟合的风险,提高识别准确度;最后利用SRM-SE-DenseNet网络模型对数据集进行自动识别得出诊断结果。
超广角眼底成像技术是一项近年来兴起的眼底图像采集技术,相较于传统的眼底图像,其具有免散瞳、范围广(可达200°范围的成像)、速度快等优势[12],在近几年被大范围应用于临床诊断,有效地提高了诊断水平,如图1所示。
本文所采用的超广角眼底图像由爱尔眼科医院提供,如图2所示。针对每例患者一共有三幅眼底图像,分别为绿激光眼底图像、红激光眼底图像和伪彩色眼底图像。
因为绿激光的穿透能力较弱,所以绿激光图像呈现的是视网膜层面的血管和结构,因为红激光的穿透能力较强,所以红激光图像呈现的是脉络膜层面的血管和结构,而伪彩色图像则是由红激光图像和绿激光图像合并之后得出的伪彩图,既包含脉络膜结构也包含视网膜结构,因此是进行识别分析的最佳选择。
本数据集对病理性近视的标定是根据屈光度和眼轴长度来进行区分,即等效球径小于等于-8D或者眼轴长大于等于26.5 mm就认为是病理性近视[13-14]。
一共收集到5 007例18至57岁之间的男女眼底图像,其中非病理性近视有2 297例,病理性近视有2 710例。在病理性近视中性别比例以及眼别比例相对平衡,而在非病理性近视中女性占比较高于男性。总体而言,眼别数量相对平衡,女性数量稍多于男性,病理性近视数量稍多于非病理性近视,其具体分布如表1所示。
表1 清洗后数据分布
因为患者的伪彩色眼底图像综合了视网膜和脉络膜的结构以及血管所包含的特征,针对病理性近视的识别主要采用患者的伪彩色眼底图像。研究方法主要包括以下几个部分:图像预处理、数据增强、网络结构。
本文主要基于超广角眼底图像进行病理性近视的自动诊断,由于其成像范围较广,所以图像内包含大量对识别产生干扰的区域以及噪声(包括眼睫毛、反光点等),因此需要对原始图像进行预处理。主要分为以下两部分:
(1) 视盘识别。视盘全称视神经盘,在超广角眼底图像中为靠近图像中心的高亮圆盘区域,在灰度图中该区域的灰度值可达到最大灰度级即255。因此可认为距离图像中心点最近的高亮点群即为视盘。所以我们先对图像进行灰度化,然后保留灰度值为255的像素,其他像素的灰度值设置为0。公式如下:
(1)
式中:Iin(x,y)为输入像素点;Iout(x,y)为输出像素点并由此得出像素点的坐标值。随后计算图像中心点坐标,其公式如下:
xm=w/2,ym=h/2
(2)
式中:xm、ym为中心点坐标;w、h为图像的宽度和高度。然后计算每个最大灰度值点与图像中心点的距离,该距离公式如下:
(3)
式中:(xi,yi)为第i个最大灰度值点坐标;(xm,ym)为图像中心点坐标;di为第i点与中心点间距离。为了防止其他高亮点对识别的干扰,我们使用一个距离阈值来对其他干扰点进行筛选。即当两点距离在阈值之内时此点将被保留,否则剔除。该距离阈值设为w/6。随后在被保留的最大值中寻找中位数,以去除孤立噪声点,最后以该中位数坐标作为视盘的位置。对于无法识别视盘的图像给予删除。
(2)
图像裁剪。当识别出视盘坐标点之后,将以其为中心,以边长为400的矩形窗对原始图像进行裁剪。这样既保留了视盘、黄斑区等重要识别区域,也剔除了眼皮、眼睫毛等无关高干扰区域。然后将图像尺寸转换为224×224大小以直接输入给神经网络。
由于总体数据集相对较少,所以需要采用数据增强方法对训练集数据进行扩充,在一定程度上解决过拟合的问题并进一步提高识别准确率。
由于光照条件的不同而导致数据集中图像明暗不均,这将使得训练集具有较大的类内方差,因此,本文通过使用Gamma校正来实现光照不变性,其具体流程为使用不同的γ值生成不同明暗程度的数据对数据进行扩充。其公式如下:
s=crγ
(4)
式中:s为输出;c、γ为正常数,一般情况下c的取值为1;r为输入。γ取值为0.7、0.9、1.1、1.3,增强后图像如图3所示。最终经过数据增强之后的训练集和测试集数据分布如表2所示。
表2 训练集和测试集数据分布
本文采用深度学习方法对超广角眼底图像进行病理性近视的自动诊断。深度卷积神经网络是一种常用的深度学习方法,其相对于传统的机器学习方法省去了特征工程这一烦琐的操作,即利用卷积对图像的特征进行自动提取,从而进行分类。本文所提出网络的基础结构选用DenseNet[15],并结合SE_block(Squeeze-and-Excitation Networks)[16]、SRM_block(Style-based Recalibration Module)[17]对图像进行识别。网络总体结构如图4所示,其中SE_Dense_Block是在原来Dense_Block的基础上加入SE_block后所得结构。如图5所示,即在Dense_Block内部的批量归一化层(Batch Normalization,BN)[18]、ReLU[19]激活层与卷积层的组合(BN_ReLu_Conv)之后添加SE_block,并在过渡层(Transition Layer)之后添加SRM_block。
DenseNet网络是由Huang等[15]提出,其网络出色的Dense Block结构对每一层输出的特征图(feature map)进行充分复用,一定程度上减轻了梯度消失的问题,网络更易于训练并且有一定的正则化效果,相对于ResNet[20]减少了参数量。该网络首先经过卷积层(Convolution)和池化层(Pooling),然后输入给Dense Block。其一共包含四个Dense Block,每一个Block结构中含有一个1×1卷积层和一个3×3卷积层。随后接入过渡层减小特征图尺寸。
SE_block(Squeeze-and-Excitation Networks)由Hu等[16]提出,可以认为是一种通道Attention机制。其利用两个全连接层(FC)训练产生针对特征图中每一个通道的权重,进而对特征进行重新校准,凸显有用特征并抑制无用特征。事实上,SE_block能够嵌入到现有的任意基础网络模型中以进一步提高网络性能。SE_block结构如图6所示,图7(a)为其网络架构。该结构首先对输入特征图进行全局平均池化以把每个特征通道变成具有全局感受野的实数,即Squeeze操作,假设忽略样本个数,其针对每一个通道的特征图计算公式如下:
(5)
式中:g∈Rc,gc为第c个通道的输出;f∈RH×W×C,fc为第c个通道的输入;H为高;W为宽;Fsq(·)为Squeeze操作函数。随后接入excitation操作,以捕获通道之间的相关性,训练生成针对于每一个通道的权重,其计算公式如下:
e=Fex(g,W)=σ(W2δ(W1g))
(6)
Uc=Fscale(fc,ec)=ec·fc
(7)
式中:Fscale(·)为权重和输入特征的融合函数。最终输出特征为U=[U1,U2,…,Uc]。
SRM_block(Style-based Recalibration Module)是由Lee等[17]提出,其Block也可以嵌入在现有的基础网络模型中,与SE结构不同的是其关注的是每一个特征通道的style特征,并利用该信息对通道进行重新校正。其网络结构与SE_block类似,网络架构如图7(b)所示。
其输入特征X∈RN×H×W×C,style特征T∈RN×C×2由平均池化(AvgPool)和标准差池化(StdPool)计算得出,其计算公式分别如式(8)和式(9)所示。
(8)
式中:N为数据个数;C为通道数;μnc为第n个样本中第c个通道所计算得出的均值;xncwh为输入特征图中的单个特征点。
(9)
式中:snc为第n个样本中第c个通道所计算得出的标准差。该通道所得style特征向量如下:
tnc=[μnc,snc]
(10)
该特征向量tnc∈R2,由μnc、snc连接而成。随后利用通道全连接层(channel-wise fully connected CFC)对特征进行编码,公式如下:
znc=Mc·tnc
(11)
式中:Mc是CFC第c个通道的权重。M∈RC×2为CFC的权重参数,Z∈RN×C为CFC层输出,随后经过BN层以及Sigmoid函数进行整合。式(12)-式(14)为BN层计算公式;式(15)为Sigmoid函数计算公式。
(12)
(13)
(14)
(15)
(16)
本文所提出的网络结构是在DenseNet的基础上结合SE_block和SRM_block,DenseNet的基本网络参数保持不变,增长率k=24,网络其他的具体示例以及参数如表3所示。
表3 网络结构参数列表
数据增强后的图像一共20 690幅,其中训练集20 030幅,测试集660幅,具体数据分布如表2所示。本文采用准确率(ACC)、灵敏度(SEN)、特异性(SPE)三个度量指标对网络模型进行全面的衡量。其表达公式如下:
(17)
(18)
(19)
式中:TP为正确分类为正样本的数量;TN为正确分类为负样本的数量;FN为实际为正样本但是预测为负样本的数量;FP为实际为负样本但是预测为正样本的数量。在本数据中,病理性近视为正样本,而非病理性近视为负样本。
一般而言,在基础网络之上不同模块的添加以及组合对分类的性能有着重要的影响,本文在DenseNet的基础网络之上结合SE_block、SRM_block,组合出5种网络结构分别为:
(1) SE_DenseNet(添加的模块均为SE_block)。
(2) SRM_DenseNet(添加的模块均为SRM_block)。
(3) SRM_SE_DenseNet(本文所提出的网络结构)。
(4) SE_SRM_DenseNet(将SRM_SE_DenseNet中的SE_block和SRM_block互换位置,即Dense_Block内部为SRM_block,Transition Layer后接SE_block)。
(5) SS_DenseNet(每个添加的模块均由SE_block和SRM_block组成)。
训练过程的初始学习率设为0.001,L2权重衰减系数为0.000 1,DenseNet的增长率为24,SE_block的衰减率为8,每次迭代传输16幅数据。损失函数为交叉熵损失函数并使用Adam[21]优化方法对参数进行优化。经过一定次数的迭代之后,其测试最高准确率分别为:78.94%、79.85%、80.03%、79.85%、79.11%,而不加模块的DenseNet准确率只有77.59%,可见增加模块确实提高了网络的性能。其ACC、SEN、SPE如图8所示,ROC曲线以及AUC面积如图9所示,图10为混淆矩阵。
根据实验结果可知,SRM_SE_DenseNet网络结构达到了最佳分类准确率,其特异性和灵敏度分别为88.55%、69.54%。DenseNet网络取得了最大的特异性值为96.93%,但是灵敏度却最小仅有53.97%,而SRM_DenseNet却得到最大的灵敏度值为70.86%和最小的特异性值86.59%。从ROC曲线以及混淆矩阵可以看出SE_DenseNet和SS_DenseNet网络结构的AUC面积要大于SRM_SE_DenseNet,但是这两个模型的特异性和灵敏度差值过大,可见对识别具有较大的倾向性,即对正样本识别较差,而对负样本识别较好,相较与SRM_SE_DenseNet来说这种倾向性减小了很多。其次可知SRM_DenseNet和SRM_SE_DenseNet的倾向性都相对较低但是SRM_SE_DenseNet网络的准确率、AUC、特异性均要好于SRM_DenseNet网络。综合考虑,最理想的网络结构即为SRM_SE_DenseNet。
为了充分证明本文提出的网络性能,本文也与几种经典的卷积神经网络进行了实验对比,分别为Alexnet[22]、VGG16[23]、GoogleNet[24]、ResNet50[20]。其实验结果如表4所示。从实验结果可以看出本文所提出的网络结构相比于其他几种网络取得了最高的准确率,其灵敏度仅次于Alexnet,总体性能最佳。
表4 与经典网络的对比实验数据表
与此同时本文也在CPU:Intel Xeon Gold 5118,内存:62 GB,显卡:GeForce RTX 2080Ti的设备上测试了不同组合的网络结构和经典网络结构在该数据集上的识别效率(时间)。分别是整体测试集的识别时间和单幅图像的识别时间,实验结果如表5所示。可以看出,SRM_SE_DenseNet在整体测试集上的识别时间能达到18.131 9 s,针对单幅图像的识别时间能达到0.027 5 s,与其他网络的识别效率相差不多。
表5 识别效率时间表
续表5
识别效率与多种因素有关,比如模型的层数、复杂度、计算量等,各参数不同则识别效率不同,相对来说,模型越简单识别所用时间越短但是识别性能则会相对下降。由于本系统为离线诊断系统,并且图像识别的准确性关乎病人的生命健康,所以在识别时间差距不大的情况下,更注重网络的识别率,因此SRM_SE_DenseNet的综合性能较为突出。
本文提出一种在超广角眼底图像上对病理性近视进行自动诊断的方法,做到辅助医生进行有效的诊断,相对减少医生的工作量,为远程辅助医疗打下了一定的基础,以减少因为资深眼科医生以及医疗设备不足而造成病人不能及时诊疗的事情发生。该方法主要通过图像处理以及深度学习的卷积神经网络实现,最终的实验结果也证明了该方法的有效性。但是由于超广角眼底图像的特殊性质,在裁剪过程中细节特征的丢失是在所难免的,以及对病理性近视的分类方法并不是根据眼底的实质性病变进行判断的,所以这都是造成判断错误以及灵敏度较低的原因。但是总体来说,本文所提出的网络结构表现出来的性能也已经充分证明了深度神经网络在眼底图像上特别是针对病理性近视诊断上有着重要的应用价值。未来的工作将会围绕以上提出的两个问题进行研究,尽量在保留细节特征的前提下,使用根据眼底实质性病变进行病理性近视分类的数据集进行诊断,进一步提高识别的准确率。