穆 珺,晏峻峰,彭清华
(湖南中医药大学,湖南 长沙 410208)
中医现代目诊以传统目诊为基础,通过对眼科仪器设备检查结果的分析,联系相关临床经验从而做出诊断。而眼底图像分析包括眼底血管分析,属于中医现代目诊的一个重要研究内容[1]。 不仅如此,眼底血管在诸如高血压、糖尿病、青光眼等疾病的检测和诊断中含有丰富的信息[2]。 但对眼底图像血管进行人工分割的方式,通常费时费力,并且需要具有专家知识[3]。 因此,眼底图像血管分割算法的研究一直受到医学图像处理领域的广泛关注[4-6]。
SINGH 等[2]使用基于Gumbel 概率分布函数的滤波器和基于熵的阈值方法对眼底图像血管进行分析。LI 等[3]通过BP 算法训练深度神经网络来完成血管分割。 HU 等[7]提出了一个多尺度卷积神经网络来获取眼底图像血管分割的概率图,并使用条件随机场获取最终分割结果。 DASGUPTA 等[8]提出了一个全卷积神经网络架构来进行眼底图像血管分割。BARKANA 等[9]提出了一个眼底图像血管分割模型,使用多种统计特征,并融合多种分类器来提升和评估分割性能。 RONNEBERGER 等[10]提出了一个名为U-Net 的端到端神经网络架构,并通过U-Net 神经网络进行医学图像分割。 此后的研究者提出了许多基于U-Net 神经网络进行眼底图像血管分割的算法,并显示出优秀的性能[6]。 JIN 等[11]提出了一个可变形的U-Net 神经网络架构来分割眼底血管。 ORUJOV等[12]提出了一种模糊边缘检测的算法来分割眼底图像血管。CHENG 等[13]提出了一个密集连接U-Net 神经网络架构来分割眼底血管。 FENG 等[14]提出了一个跨连接神经网络来分割眼底血管。YU 等[15]使用零相位分量分析和全局对比度规范化进行预处理,然后训练卷积神经网络获取眼底血管分割结果。 尽管研究者们付出了大量的努力,也取得了不少进展,但眼底图像血管分割仍然具有值得进一步改进的空间[4-6]。
本文所提出的眼底图像血管分割算法的总体流程包括:数据增广与预处理、小波系数图像获取、UNet 神经网络训练与预测、阈值分割与后处理。
1.1.1 数据来源 本文所提出的眼底图像血管分割算法是在DRIVE 数据集上进行训练和测试。DRIVE数据集是广泛使用于眼底图像血管分割的公开数据集,其中包含40 张眼底图像,每张图像的尺寸为565 像素宽、584 像素高。 该数据集的40 张图像中,20 张被划分为训练图像,20 张被划分为测试图像。每一张眼底图像都提供了对应的专家手动分割图像。本文采用U-Net 神经网络,将训练图像和对应的分割图像作为训练数据进行训练;将算法对测试图像的分割结果与测试图像的专家分割结果进行对比,以度量算法分割性能。
1.1.2 数据增广 在诸如眼底图像血管分割等医学图像处理问题中,要获取大量的专家手动分割信息以构成训练数据通常是非常困难的。 而训练数据数量的缺乏,很可能造成训练中的模型过拟合问题,从而降低以神经网络等算法为代表的有监督的分割算法的准确性和泛化能力。为解决该问题,本文采用数据增广操作,将原始数据集中的40 张图像增广为160 张图像。数据增广操作对于每一张图像,计算得到其水平镜像图像、竖直镜像图像、水平与竖直镜像图像。 单张原始图像的增广结果如图1 所示。
图1 数据增广
1.1.3 图像预处理 在输入U-Net 神经网络进行训练和测试之前,本文将对原始图像进行一系列预处理操作。预处理操作包括图像灰度化、图像数据标准化、图像对比度受限自适应直方图均衡化[16]、伽玛变换。 图像灰度化是将彩色图像转换为灰度图,其计算公式如下:
GrayScale=0.299×R+0.587×G+0.114×B (1)
其中,R、G、B 分别为每个像素点的红、绿、蓝3个颜色通道的亮度值,通过加权求和,得到GrayScale,即该像素点的灰度值。
图像数据标准化操作的作用是将数据转化为具有均值0 和标准差1 的新数据,以避免具有较大取值范围的变量对计算结果产生剧烈影响,标准化的公式如下:
其中,x 是每个像素点的灰度值,x 是所有像素点的平均灰度值,Sx是所有像素点灰度值的标准差,x′是标准化之后的灰度值。
图像对比度受限自适应直方图均衡化算法的作用是利用局部图像块的直方图对图像的灰度值进行映射,以提高图像对比度,实现图像的增强。伽玛变换的作用也是对图像进行增强,其公式如下:
其中,x 是每个像素点的灰度值,gamma(x)是变换后的灰度值,参数γ 取值为0.833。
原始图像经过预处理后的结果如图2 所示。
图2 图像预处理
1.2.1 哈尔小波变换 小波变换是信号处理与图像处理中一个强有力的数学工具。 小波变换能同时描述图像的空间域和频率域的信息。 哈尔小波是一种应用极为广泛的小波变换[17]。 哈尔小波的计算公式如下:
通过哈尔小波变换,可以得到图像的两种小波系数,分别为近似系数与细节系数。 近似系数可以视为原始图像在低分辨率下的近似;而细节系数可以描述原始图像在各个局部的细微变化。
对二维灰度图像进行一阶哈尔小波变换,可以视为在行和列两个方向分别进行一维哈尔小波变换,从而可以得到以下小波系数矩阵:图像的近似系数矩阵、反映图像水平方向变化细节的水平细节系数矩阵、反映图像垂直方向变化细节的垂直细节系数矩阵、反映图像对角方向变化细节的对角细节系数矩阵。
1.2.2 小波系数矩阵上采样 对灰度图像进行一阶哈尔小波变换,可以得到水平细节系数矩阵、垂直细节系数矩阵、对角细节系数矩阵,这3 个细节系数矩阵可以综合反映出图像各个局部在不同方向上的细微的变化情况。
因此,将这3 个细节系数矩阵作为富含信息的输入,和预处理得到的灰度图像共同组成神经网络训练和测试的输入。
由于小波变换得到的细节系数矩阵尺寸为原始图像的一半,因此,采用双线性插值算法对小波系数矩阵进行上采样,得到尺寸和原始图像一致的小波系数图像。
小波变换的输入图像和上采样得到的3 个小波系数图像如图3 所示。
图3 小波变换的输入图像和上采样得到的3 个小波系数图像
1.3.1 U-Net 神经网络架构 U-Net 是RONNEBERGER 等[10]于2015 年提出的端到端的卷积神经网络。 U-Net 网络由收缩路径与扩张路径两个对称的部分构成。 U-Net 神经网络对输入的图像块进行处理,最终得到同样尺寸的特征图作为分割预测结果输出。
在收缩路径中,原始输入图像块经过两次3×3卷积与ReLU 操作,再进行2×2 最大池化,如此操作重复4 次;得到的特征图像块再经过两次3×3 卷积与ReLU 操作,此后进入扩张路径;在扩张路径中,将收缩路径最末端得到的特征图像块进行上采样,并与收缩路径中处于同层的特征图像块进行拼接,拼接后进行两次3×3 卷积与ReLU 操作,将得到的像块进行上采样,如此重复4 次。 在扩张路径的最末端,通过2 个1×1 卷积将输出特征图映射为每个像素点分别属于血管和背景的预测值。
本文在大体借鉴文献[10]中的U-Net 神经网络架构基础上,对网络架构进行了一些调整和改善。
首先,本文将收缩路径中的池化操作设置为重复2次,相应的扩张路径中的上采样操作也设置为重复2 次,以适应于本文所处理的眼底图像及眼底图像块尺寸有限的特点。 其次,本文在每个3×3 卷积与ReLU 操作之后,加入了批标准化与丢弃操作,以加速神经网络收敛。
本文所采用的U-Net 神经网络架构如图4 所示。
图4 U-Net 神经网络架构
1.3.2 U-Net 神经网络训练 (1)训练图像块采样尽管已经通过数据增广将训练图像的数量从20 增加到80,但是训练图像的数量仍然非常有限。 为避免U-Net 神经网络训练中出现过拟合,本文采用随机采样的方式,在每张原始训练图像中,随机采样9500 个尺寸为48×48 的图像块,从而构造总数量为760 000 的图像块作为训练数据。
(2)损失函数本文对于U-Net 神经网络的最终输出层采用softmax 函数作为激活函数。 softmax 函数的公式如下:
其中,zi(i∈{1,2})表示U-Net 神经网络的输出,即图像块中每个像素点分别属于血管和背景的预测值,并且预测值被排列成2 条向量。 p(zi)就是将预测值归一化后得到的每个像素点分别属于血管和背景的概率值,分别用p(z1)和p(z2)表示,若z1>z2,则p(z1)>p(z2),反之亦然;每个像素点将被算法预测为p(zi)取最大值所对应的类别。
对于每个像素x,损失函数可以衡量算法对像素预测值与真实值之间的差异,本文采用交叉熵损失函数。 交叉熵损失函数的公式如下:
其中,y(x)与p(x)分别是像素x 的真实类别标签和预测值。 y1(x)是表示像素x 是否属于背景点的类别标签,而y2(x)是表示像素x 是否属于血管点的类别标签。 p1(x)是算法预测像素x 属于背景点的概率值,而p2(x)是算法预测像素x 属于血管点的概率值,考虑所有采样的像素点x(x∈Ω),总体的损失函数如下所示:
其中,m 是所有属于训练样本集的像素总数,而Loss(y(x),p(x))是单个像素x(x∈Ω)的损失函数。
1.3.3 U-Net 神经网络预测 在完成训练后, 本文使用训练得到的U-Net 神经网络对于测试图像的所有像素点属于血管的概率进行预测。 预测值越接近1 表示像素点属于血管的概率越大,预测值越接近0 表示像素点属于背景的概率越大。
由于训练中采用了尺寸为48×48 的图像块,在测试中,也需要将测试图像划分为48×48 的图像块。 对每一张原始测试图像,本文使用步幅为2 像素的间隔进行图像块采样。 对于原始测试图像中的同一个像素点,将多次出现在不同采样块的不同位置;每个采样图像块中分别得到的同一个像素点的预测值,再将同一个像素点的多个预测值求平均作为该像素点的最终预测值,以此减轻采样位置的差异对于预测结果的影响。
原始测试图像和U-Net 神经网络预测得到的概率图像如图5 所示。
图5 原始测试图像与U-Net 的预测结果
1.4.1 阈值分割 通过U-Net 神经网络预测,可以得到测试图像所有像素点属于血管的概率预测值,该预测值取值为0 到1 之间,预测值越接近1 表示像素点属于血管的概率越大。
本文采用阈值分割方法将每个像素点属于血管的概率预测值转化为像素值是否属于血管的分类标签。 具体计算公式如下:
其中,PredictionMap(x)表示每个像素点x 属于血管的概率预测值,T 是分割的阈值,本文中设置T=0.5。而seg(x)表示阈值分割的结果,seg(x)=1 表示像素x 为血管点;而seg(x)=0 表示像素x 为背景点。
阈值分割的结果如图6(a)所示。
图6 阈值分割与后处理结果
1.4.2 后处理 在眼底图像中,许多局部区域的对比度很低,细小的血管非常不清晰,因此,即使借助于U-Net 神经网络,许多血管尤其是细小血管的分割也是极其困难的。 表现在阈值分割结果上,就是分割出的血管有很多断裂的情况。 为解决这一问题,本文在阈值分割的基础上,引入闭运算作为后处理,从而将血管分割结果中出现的断裂处进行融合。
后处理的结果如图6(b)所示。 直接对比阈值分割结果与后处理结果,两者差异不明显;但在“2.4”的实验分析中,通过对一些测试图像的阈值分割结果与后处理结果进行局部放大,可明显观察到后处理结果相对于阈值分割结果有所改进。
本文所有实验均运行于笔记本计算机,该计算机配置了型号为英特尔酷睿i7-9750H、主频为2.60 GHz的CPU,并配置了16 GB RAM 和6GB型号为NVIDIA GeForce GTX 1660Ti 的GPU。 该计算机所使用的操作系统为64 位Win10。本文算法实现所使用的编程语言为python,用到的第三方库主要包括TensorFlow、Keras、sklearn、numpy、opencv、pywt 等。本 文实验所用到的眼底图像数据集为DRIVE 数据集,该数据集的介绍在“1.1.1”。
本文用来对算法性能进行评价与度量的主要指 标 包 括:准 确 率(Acc)、特 指 度(Spe)、灵 敏 度(Sen)、精度(Pre)、F1 值与ROC 曲线的AUC 指标。上述评价指标的定义如下:
其中,TP 表示被模型正确预测的正样本数;FP表示被模型错误预测为正样本的负样本数;FN 表示被模型错误预测为负样本的正样本数;TN 表示被模型正确预测的负样本数。 在本文实验中,正样本是指对应血管的像素点;负样本是指对应背景的像素点。
Acc、Spe、Sen、Pre、F1 值,是从不同角度衡量分类模型的性能,上述指标的取值范围都在0 到1 之间,取值越接近1,说明模型性能越好。 ROC 曲线是接受者操作特征曲线,ROC 曲线下方的面积,简称为AUC 指标,可以综合地评价模型的平均性能。 AUC指标的取值范围也在0 到1 之间,取值越接近1,说明模型性能越好。
为分析本文所使用的小波系数图像对于分割算法性能的影响,本文在其他实验设置完全相同的情况下,对比分析了不加入小波系数和加入小波系数时,血管分割结果的性能差异,如表1 所示。 实验结果显示,从Acc、Spe、Pre、F1 值、AUC 指标等5 个方面比较,加入小波系数的结果均优于不加入小波系数的结果;唯一在Sen 指标上,加入小波系数的结果略差于不加入小波系数的结果。综合而言,本文所使用的小波系数对分割结果具有一定的提升作用。
表1 不加入小波系数和加入小波系数的性能比较
分析对比阈值分割结果与后处理结果,从图像全局观察来看,两者不明显。 为此,本文将一些测试图像的阈值分割结果与后处理结果进行了局部放大和对比分析。 阈值分割结果与后处理结果的局部放大图像如图7 所示。 其中,列(a)为阈值分割结果局部放大图,列(b)为后处理结果局部放大图。通过对于细小血管所处局部图像进行放大,可以明显观察到,引入闭运算作为后处理所得到的分割结果中,原先断裂的细小血管实现了融合。
图7 阈值分割与后处理局部放大图
根据现有眼底图像血管分割算法所普遍采用的评价指标,从Acc、Spe、Sen、ROC 曲线的AUC 指标4 个方面对本文算法和现有算法进行性能比较,如表2 所示。 本文算法分割结果的ROC 曲线如图8 所示。实验结果显示,本文算法均超过或者不弱于现有的眼底图像分割算法。
表2 本文算法与现有算法性能比较
图8 本文算法分割结果的ROC 曲线
需要说明的是,本文算法为有监督的分割方法,需要将数据集中的原始图像划分为训练图像和测试图像。本文采用DRIVE 数据集的发布者所给出的划分,将DRIVE 数据集中的40 张原始图像划分为20张训练图像和20 张测试图像。 此外,文献[3]、文献[8]、文献[9]、文献[7]、文献[11]、文献[13]、文献[14]、文献[15]的算法也都是有监督的分割方法,上述论文中均是采用相同的划分方法,将DRIVE 数据集中的40 张原始图像划分为20 张训练图像和20 张测试图像。
本文算法所得到的眼底图像血管分割结果如图9 所示,其中列(a)为原始眼底图像,列(b)为本文算法得到的最终分割结果。
针对中医现代目诊的一个重要研究内容,即眼底血管分析,本文提出了一种基于小波变换和UNet 的眼底图像血管分割算法。 该算法主要通过数据增广与预处理、小波系数图像获取、U-Net 神经网络训练与预测、阈值分割与后处理等步骤完成眼底血管的分割。本文从不同角度对本文算法进行了性能分析,实验结果显示,本文所使用的小波系数对分割结果具有一定的提升作用,并且本文算法超过或者不弱于现有的眼底图像分割算法。
综上所述,本文提出的眼底图像血管分割算法能得到满意的分割结果,在中医现代目诊中具有潜在的应用价值。