基于全连接条件随机场模型的眼底图像血管分割

（沈阳工业大学信息科学与工程学院，沈阳110870）

1 引言

眼底图像又被称为视网膜图像，是人眼内表面的彩色投影图像。大多数的疾病可以通过观察眼底图像中的血管树、视神经盘、中央窝等内部构造的变化来诊断、预防和治疗[1]。使用计算机辅助方式对眼底图像中血管网络的结构特征进行准确检测和定量分析是医生完成眼底图像分析的基础，其准确性和实时性影响着对疾病的诊治效率[2]，因此眼底图像的血管分割一直都是医学图像处理和机器智能领域的研究热点。

2002年Zana和Klein基于数学形态学[3]的思想和线性处理的方法提出了一种视网膜血管分割方法。该方法首先定义了四种不合需要的模式的几何模型，以便能够从各种噪声环境中分离出血管；然后运用数学形态学方法从视网膜图像中提取线状结构，接着使用Laplacian算子计算出曲线的不同属性；最后根据曲线的曲率差异提取血管。这种方法的优点是速度快，对噪声不敏感，缺点是血管的其他有效特征没有得到充分利用且分割结果中血管的连续性较差。

2006年Soares等人[4]利用小波变换理论提出了一种基于2D Gabor变换的快速分割方法。该方法采用了像素灰度值和各尺度下的2D Gabor小波变换响应[5]结果作为特征向量，使用贝叶斯高斯混合分类器对血管特征进行识别。此方法的缺陷是在像素特征选取方面仅考虑了图像的2D Gabor滤波信息，忽略了血管的形状和结构信息，但其分割结果中容易出现血管融合的现象。

2013年Nguyen等人[6]提出了一种从视网膜图像中自动提取血管的有效方法。该方法基于对基本线性探测器的长度的改变，使其可在不同尺度上进行分割，同时又保持了原探测器的强度，消除了各探测器的缺陷，由线性探测器在不同尺度上的响应的线性组合产生了最终分割结果，但其在血管检测中出错的问题比较严重。

2015年王晓红等人[7]以视网膜2D Gabor小波变换为基础，针对视网膜粗血管和细血管分别进行提取，提出了一种基于多尺度2D Gabor小波变换的视网膜血管自动分割方法。该方法只使用了低分辨率图像，并没有扩展到高分辨率眼底图像。

2016年吴奎等人[8]针对单一的2D Gabor算法性能的缺陷，提出了一种基于2D Gabor小波变换和组合线检测算子的视网膜血管自动分割方法，但同样该方法并没有应用于高分辨率眼底图像。

实际上，眼底图像中血管像素之间不是相互独立的，而是具有互相影响的特点，这一特点恰恰可用条件随机场理论中的线性链结构进行描述，故此提出一种基于全连接条件随机场（Fully Connected Conditional Random Field，FC-CRF）模型的眼底图像血管分割方法。该方法不仅考虑了眼底图像的像素特征，同时还考虑了像素之间的相互作用，很好地改善了细小血管处的错误分割、血管融合等问题，提高了血管的连续性以及分割的准确性。研究中采用不同的支持向量机（Support Vector Machine，SVM）实现了对不同像素的分类以及对此方法所涉及参数的自动调节。另外，通过调整不同数据集上的最佳特征参数，解决了在不同分辨率图像中特征提取不成比例的问题。本研究的精华之处在于，将一种用于大面积分割事物的方法成功应用于眼底图像之中，使之不仅可以大面积分割目标图像中的事物，还可以分割例如血管之类的微弱结构，并通过SVM实现像素的分类以及参数的自动调节，最终应用于不同分辨率的眼底图像。本研究将以两个公开数据集中的眼底图像为实验对象，对所提算法进行验证，并与现有方法进行了对比。对比结果表明从分割图片和评估值上均可体现出该算法的高效性和实时性。

2 选用数据及实现方法

2.1 数据集的选用

根据实际情况，选取两个公开的数据集DRIVE和HRF中的图片来进行研究和验证。

DRIVE[9]数据集中的所有图片来源于荷兰的一个糖尿病筛查项目。该项目的筛查对象为25岁至90岁之间的400人。数据集包含40张彩色眼底照片，其中有7张是有轻度糖尿病病症的图像，33张无任何病变；全部眼底图像被平均分为训练集和测试集。这些图像使用Canon CR5无需散瞳的3CCD照相机以45°为拍摄角度进行拍摄，得到像素为565×584、格式为JPEG的彩色眼底图像。

HRF[10]数据集中包含15张健康的彩色眼底图像、15张糖尿病患者视网膜病变图像、15张青光眼患者视网膜病变图像。照相机以60°为拍摄角度进行拍摄，得到像素为3304×2336、格式为JPEG的彩色眼底图像。

本试验中，将每种图像的前五张组成训练集，其余的图片作为测试集。

2.2 实验方法

研究所采用的方法一共包括三个部分：预处理、特征提取、血管分割。整体方法结构如图1所示。

图1 方法结构图

(1)预处理

大多数图像处理方法中的特征提取采用的都是原始彩色眼底图像的绿色通道的灰度图像[11]。本研究考虑到相机边缘光圈会造成边缘的错误检测，使用[4]中所提出的方法对图像进行预处理，对绿色通道图像边界进行扩充，然后取其反转图像。实验选用的原始图像如图2(a)所示，预处理后的图像如图2(b)所示。

(2)特征提取

在特征提取部分，选用了广泛应用于眼底图像血管分割领域的特征提取方法，分别对一元特征和成对特征进行提取，即采用Nguyen等人提出的多尺度线性检测器响应以及二维Gabor小波响应来提取一元特征，采用Zana等人提出的形态学方法来提取成对特征。

Nguyen等人提出的多尺度线性检测器利用血管为细长结构的属性，将一条长为l的直线以不同的角度通过目标像素P来计算像素的平均强度[6]。从所有考虑的方向中选择具有最大平均强度的线Ll(P)，通过做差的方式得到像素的谱线强度Sl(P)=Ll(P)-Ns(P)，是以P为中心，边界长度为s的矩形窗口的平均强度。图2(c)所示为多尺度线性检测器所获得的结果。

二维Gabor小波变换具有检测取向特征的能力，并且可以调整特征频率[4]。由于血管是以不同的大小和方向出现的，所以这种方法能够有效地增强眼底图像中的血管。此处选用Soares等人在[4]中提出的方法，在不同的尺度α计算特征，将图像在不同尺度α的小波反应作为特征，如图2(d)所示。

Zana等人利用血管是线性的、连续的且沿脊线变化是平滑的等特性来增强眼底图像中的血管[3]。首先，使用长度为l的线性结构元素对图像进行开运算操作以减少图像中的噪音。再用同样的结构进行多个顶帽形态操作，并对每个角度相应的响应求和。这个操作可以减少亮度较小的噪点，并且提高线性结构的对比度。然后通过使用窗口大小为7×7像素、标准差为7/4的拉普拉斯算子计算横向曲率来检测曲率线性相干的结构。最后，用连续的形态学开运算、闭运算去除对过亮或过暗的非线性的不规则区域和背景的错误检测（以上三步使用相同的长度为l的线性结构）。由于连续使用形态学操作会降低血管与背景之间的对比度，因此为了提高眼底图像的对比度，在使用此方法前会对眼底图像进行额外处理：从原始彩色图像的绿色通道中减去图像的估计背景。估计背景是通过绿色通道图片进行中值滤波后得到（滤波器的内核的大小要能够确保模糊的图像中不含有像血管一样明显的结构）如图2(e)的结果。图2(f)所示则为形态学方法所获得的结果。

图2 图像预处理及一元特征、成对特征对比

(3)血管分割

CRF模型已被应用于道路检测领域[12]，但是在道路检测中只需要考虑近距离像素之间的相关性（例如四邻域、八邻域）就可以完成道路检测。在医学领域中，FC-CRF模型[13]也已经被用在肝脏和大脑肿瘤的分割中，但是目前FC-CRF模型只能用来定位事物的轮廓或者大面积的分割事物，应用于分割血管之类的微弱结构并不成功。不同于大面积的分割事物，眼底图像中的血管非常细小且分支结构比较明显和复杂，有时还会因为背景中存在不同的病变以及照明问题导致分割错误，因此要想准确分割出连续的、清晰的血管类微弱结构难度会大大增加。因此，考虑到眼底图像在分割过程中像素之间存在不同转移特征和状态特征的特点，在此提出一种基于全连接条件随机场模型的眼底图像血管分割方法。该方法通过计算图像中远距离像素之间的转移特征函数和状态特征函数，不停训练并调整权值向量以及特征图谱，有效地增强了对于微弱结构的连结性的检测能力。不同于[13]的是此算法能够成功应用于分割整个眼底图像，包括对其中分支结构明显的细小血管的分割，而不只是对目标图片中的事物进行大面积分割或检测。

对应全连接条件随机场模型理论，可把眼底图像看作由点和线组成的无向图。无向图中的一个节点表示原眼底图像的一个像素。根据连接规则，将各个节点和它的邻居节点相连，如图3所示。

图3 CRF模型节点连接方式示意图

在该无向图中，任何两个结点间均有边连接的节点子集称为团。若一个团中不能再加进任何一个节点使其成为更大的团，则称其为最大团。根据不同的连接方式在眼底图像中会形成不同的最大团。CRF模型的连接方式可以分为以下两种：

(1)局部邻域CRF(Local neighborhood based CRF，LNB-CRF)：常用的局部邻域连接方法一共有三种：四邻域、八邻域以及D邻域。图3(a)为四邻域的连接示意图。图3中节点为眼底图像的像素点，黑色点是根据不同的连接规则进行能量计算的像素点。

(2)FC-CRF：全连接方法是假设无向图中的任意节点都与无向图中的其余节点相连，如图3(b)所示。

如果无向图模型能够表示成一系列最大团的非负乘积形式，这个无向图模型的概率分布就称为吉布斯分布，于是求最大后验概率被转化为求能量函数的最小值。因此在CRF中可将分割问题看作是能量最小化问题。

为了分割血管像素和非血管像素，首先定义一个标签空间L={-1,1}用于对眼底图像I中所有像素进行标识，1表示该像素是血管，-1表示该像素为非血管，此处用y={yi}表示眼底图像I中所有像素的标签。这样，由对眼底图像I和其中所有像素的标签y组成的条件随机场{I,y}可以用吉布斯分布表示为：

式中，Z(I)是一个归一化常数，G是与I相关的无向图，CG是G中的团，yc为无向图中最大团的像素标签，Ψc为势函数。式(1)表示在给定图像I的情况下，求得标签y的条件概率函数。定义式(1)中吉布斯能量函数为：

通过最小化式(2)中得到的能量，可以得到眼底图像血管的二值分割[14]：

为了表示方便，用Ψc(yc)表示Ψc(yc|I)。

上述的吉布斯能量函数可以分为一元能量函数和成对能量函数两种。一元能量函数是通过纹理、颜色等一元特征计算得到的；但是眼底图像具有复杂性和多样性（例如：光照不均匀、病变等因素），仅仅提取一元特征是不够的，还需要结合像素之间的相互作用（成对特征）进行分析，成对能量函数即是通过成对特征计算得到的。因此，在此将能量函数分别对一元能量函数和成对能量函数来进行分析。假设在无向图G中的标签为y，其总能量是通过对一元能量函数和成对能量函数求和得到的（Ψu和Ψp分别表示一元能量函数和成对能量函数），如下式：

式中xi和fi分别表示从眼底图像中提取的一元特征和成对特征。

因此根据上文所述，能量函数计算流程如图4所示，且有下述三点论证：

图4 能量计算流程图

1)一元能量函数

一元能量函数(Unary Potential,UP)是通过提取眼底图像的一元特征计算得到的。在CRF模型中一般采用纹理、颜色、位置、形状等图像特征来训练分类器，然后在分类器的分类结果上计算CRF模型所需的一元能量函数值。

在LNB-CRF模型和FC-CRF模型上一元能量函数的计算公式是相同的，可根据下式计算得到：

2)成对能量函数

成对能量函数(Pairwise Potential)是通过提取的眼底图像的成对特征计算得到的。传统的CRF模型一般采用邻域连接的结构（实验中采用四邻域连接的方式），但是这很难将较远处的同类像素与之联系起来，而FC-CRF模型则考虑了所有像素点之间的联系，提高了眼底图像中血管的连通性，并且有助于细小血管的分割。

成对能量函数定义为高斯核函数的线性组合，因此，成对能量可以根据下式计算得到：

式中每个k(m)都是一个任意特征f(m)的固定函数，wp(m)是权重的线性组合，滋(yi,yj)是表示标签兼容性的函数，利用Potts模型将两个相似的像素分为不同的标签，避免了不利因素的存在。高斯内核通过f(m)确定连接的像素之间的相似度，LNB-CRF和FC-CRF有不同的内核函数计算公式。

LNB-CRF模型定义在网格图上，以四邻域结构连接为例，可根据下式计算得到核函数：式中，θ(m)是带宽，调节不同图像特征的相关性，f(m)表示提取的眼底图像特征。

FC-CRF模型是假设无向图中每一个节点都与图中其余节点相连接。此方法不仅考虑节点周围的邻域信息，还能够考虑远距离像素间的相互作用，从而提高了分割的精确度，但是也使得计算量和推理过程大大增加，因此实验中选用[13]等人最近提出来一种基于原始CRF的平均近似场理论，可以提高计算速度，在几秒钟内产生准确的血管分割图像。全连接模型的核函数有以下形式：

式中pi和pj是像素i和j的坐标向量。内核宽度θp和θ(m)调节式(8)中像素和特征之间的相关性程度；表示提取的眼底图像特征，例如，随着θp的增加，考虑的更多是远距离之间的相互作用，与之相反，若θp减小，只有局部邻域结构会对分割结果产生影响。

3)像素分类

线性分类器有三大类：感知器准则函数、支持向量机(SVM)[14]、Fisher准则。SVM是很经典的算法，优化目标是最大间隔（margin），又称最大间隔分类器，是一种典型的线性分类器，通过使用核函数可解决非线性问题。因此实验选用SVM进行像素分类。

对于式(5)中涉及到的权值向量憎以及式(6)中涉及到的特征图谱渍(s,y)，在实验中是通过两种不同的SVM训练得到的：

(1)第一种训练方式

用具有边缘调整能力的支持向量机SVM来对全连接模型进行训练分类，并完成参数的自动调节。这种训练方式可以应用于多个领域的高度复杂且精确的模型中，有助于眼底图像的血管分割[15]。具体过程为：

首先，使用[16]中提出的具有边缘调整能力的SVM的计算得到憎。先假设一个训练集S={(s(1),y(1)),...,(s(n),y(n))}，其中，n表示训练的眼底图像的数量，y(i)表示的是训练集中第i张图片所对应的专家人工分割图像。每个集合S(i)={x(i),β,f(i))中包括每张图像的一元特征集合x(i)，偏差常数茁=1，以及成对特征集合f(i)。若特征图谱函数与损失函数之间满足下式那么权值向量憎为：

其次，计算特征图谱：

则眼底图像特征图谱分别通过下式计算得到：

(2)第二种训练方式

选用结构化最小支持向量机（Structural Least Squares Support Vector Machine，SLSSVM）对仅考虑一元特征时的眼底图像进行训练分类，并完成参数的自动调节。

SLSSVM是将SVM中的不等式约束改为等式约束，通过求解方程组得到两类超平面，即若特征图谱函数与损失函数之间满足：

那么权值向量w为：

SLSSVM通过将不同的数据分成不同的簇结构，引入样例均值和协方差矩阵来考虑数据的分布情况，这样也就使分类器的泛化能力得到增强。

2.3 比例模型

比例模型被用来调整在不同分辨率图像中特征不成比例的问题。一般来说，各个方法中所涉及的参数会影响特征的提取结果，这些参数通常与血管属性（例如血管口径）相关，同时也与图像分辨率相关。例如，2D Gabor小波的响应取决于其尺度琢。同样，Nguyen等人的线检测器和Zana等人的形态学方法分别取决于检测器或线性结构元件的长度l。大多数的特征参数最初是根据低分辨率图像（例如DRIVE数据库）设置的。当将这些参数应用于高分辨率图像时，特征提取过程可能会不成比例，那么系统的性能就会显著降低。

因此，当系统应用于不同分辨率图像时，需要对其中的一些参数进行适当的调整，由于所包含的参数及其组合属于高维空间，所以手动调整参数的方法并不可行。因此，针对所选的分割方法提出一种简单的调节方法，通过调整每个数据集上的最佳特征参数来代替调整每个数据集的分辨率，保证在高分辨率图像中也能够成比例的提取特征从而保证了系统性能不会降低。然后，通过乘以修正系数籽=xnew/xtraining来调整这些参数，其中xtraining表示眼底图像中视场(Field of View,FOV)的平均宽度，xnew表示新图像中FOV的平均宽度。由于分辨率的变化与FOV区域内相关联的像素数量的变化有关，这种简单的方法近似于有关缩放特征计算的不变性。

3 实验结果

3.1 评价指标

将分割结果与数据集中提供的标准分割相比较，可得到七种不同的测量结果，这些结果都是根据真阳性 TP(true positive)、真阴性 TN(true negative)、假阳性 FP(false positive)、假阴性 FN(false negative)的数量计算得到，并且只考虑FOV区域内，具体参数定义如表1。

表1 血管分割的性能指标

其中N=TP+TN+FP+FN是图像中像素总数，S=(TP+FN)/N,P=(TP+FP)/N。灵敏度（Se也称为召回率Re）衡量方法正确检测血管的能力，特异性（Sp）则是区分其它非血管结构能力的指标。然而，尽管需要较高的Se值，但是必须与Sp结合进行分析，因为如果将所有像素标记为血管，则可以将Se最大化。精度（Pr）正确识别出血管像素所占的比例。此外还使用马修斯相关系数（MCC），F1-score（F1）以及 G-mean(G)作为整体性能的评估指标，因为相比于准确率（ACC）更适用于评估非平衡数据（非平衡数据集是指阳性数据在整个数据集中所占比例很小）。MCC是人工分割与预期分割之间的相关系数，已经广泛用于评估视网膜血管分割方法。F1-score(F1)是精度和召回率的谐波均值，是一种统计量，F1越高说明实验方法越有效[17]。同样，几何均值G-mean(G)也是非平衡数据集学习常用的评估指标，通过获取和的几何平均数来进行评估，返回一个介于0和1之间的数。最后，生成ROC曲线，通过计算得到曲线下的面积(AUC)。根据ROC曲线的定义，ROC曲线下面积应在1.0和0.5之间。在AUC>0.5的情况下，AUC越接近于1，说明诊断效果越好。AUC在0.5耀0.7时有较低准确性，AUC在0.7～0.9时有一定准确性，AUC在0.9以上时有较高准确性。AUC=0.5时，说明诊断方法完全不起作用，无诊断价值。AUC<0.5不符合真实情况，在现实中很少出现。

3.2 测试结果

表2分别给出了仅用一元能量和使用实验算法所得出来评估结果。训练LNB-CRF模型时，由于SVM给局部能量函数的权重wp一个几乎为0的值，同时又不会更改一元能量函数的权重，导致使用LNB-CRF所获得的分割结果与只使用一元函数获得的结果相同，所以表中不包含使用LNB-CRF所获得的结果以及其评估值。与之相反，当训练FCCRF模型时，具有边缘调整能力的SVM为wp分配了一个非零值，并且还会修改一元特征和偏差项的权重，这意味着成对能量函数会对其他参数造成影响并且有助于预测函数。

表2 在DRIVE集和HRF集上的评估结果

从DRIVE数据集上的评估值中可以观察到，虽然成对能量函数在提高Se值的同时略微的降低了Sp值，但是从F1值、G值、MCC值的改善可以看出，Se的改善程度是要高于Sp的降低程度的。同样的，在HRF数据集中，一些相应的值都得到了改善。

图5所示为实验算法对DRIVE集中的20幅测试图像所得的性能指标。包括各个图像的敏感度、特异性以及分割的准确度。从图中可以看出本文的敏感度和特异性值较高，这说明实验方法的分割更准确，其准确性的均值可以达到0.9486，这高于Gabor分割方法的平均准确性0.9466、形态学分割法的平均准确性0.9320以及多尺度线性检测器方法的平均准确度0.9407。

图5 血管分割的性能结果

图6所示的分割结果中包含原始数据集中的彩色图像、专家手动分割结果、仅考虑一元特征的分割结果以及全连接条件随机场的分割结果。从中可以看出，相比于图6(c)、6(g)仅考虑一元能量的分割结果，图6(d)、6(h)使用FC-CRF模型分割的血管图像中包含更多的细小血管，改善了血管结构的连通性，并且与专家手工分割的血管更为相似。图5中的图像分别取自DRIVE集和HRF集，因此可以证明实验方法无论是应用于低分辨率图像还是高分辨率图像都有很好的分割结果。

图6 分割结果

在图7中，图7(a)为形态学方法的分割图像，图7(b)为Gabor小波变换后的分割图像，图7(c)为多尺度线性检测器的分割图像,图7(d)为使用本研究中提出的实验方法所获得的分割图像，从图中可以看出与形态学方法的分割图像（如图7(a)）相比，实验方法所获得血管的连通性更好。与Gabor小波变换的血管分割图像（如图7(b)）相比，血管融合的情况得到了改善。与线性检测器增强的血管分割图像（如图7(b)）相比，错误的血管检测较少。因此，相比于其他的血管分割方法本文方法更接近专家的手动分割结果。

图7 分割结果

图8是实验方法生成的ROC曲线。曲线下面积的均值可以达到0.9476，这说明此方法有较高的准确性。M.Niemeijer等人[18]的方法获得的曲线下面积为0.9294，M.A.Amin等人[19]的方法获得的曲线下面积为0.9360，T.Chakraboriti等人[20]的方法所获得的曲线下面积为0.9419，均小于此方法获得的曲线下面积值。

图8 DRIVE集上的ROC曲线

表3是将实验方法（属于监督式）与其它现有监督式、非监督式方法进行比较所获得的评估参数值。从表中的评估值可以看出相比于Soares、Wang等人的方法，此法的敏感度有所提升，虽然在特异性方面没有其他方法所获得的数值高，但是从F1、G以及MCC的值可以看出，实验方法敏感度的改善程度要高于特异性的降低程度。尽评估值没有很大程度的提升，但已相当接近于专家手动分割值。

表3 在DRIVE集上评估值对比

上述数据均说明本研究所用的实验方法的分割结构在性能参数上相对优势较高，对血管和非血管的定位更加准确。

4 结束语

结合条件随机场理论与眼底图像血管分割问题，提出一种基于全连接条件随机场模型的眼底图像血管分割方法，具体工作包括：使用条件随机场模型实现眼底图像血管分割;使用不同的支持向量机对分割过程中不同的像素进行分类并对实验涉及参数自动调节。通过调整不同分辨率图像的最佳特征参数避免特征提取不成比例的问题，确保了分割结果的准确性。在两个有着不同分辨率的数据集DRIVE和HRF上进行评估，结果表明，此方法既可以用于低分辨率图像又可以用于高分辨率图像。在评估值方面，不仅使用了敏感度、特异性以及准确率来证明方法的有效性，还使用更具代表性的评估值F1、G-mean以及MCC进行评估。通过与现有方法评估值的对比可以看出，此方法可以较好地提高敏感度、改善特异性。