CN104867151A - 一种基于延展性主动轮廓模型的眼底图像血管分割方法 - Google Patents

Abstract

本发明一种基于延展性主动轮廓模型的眼底图像血管分割方法，属于医学图像、计算机视觉和数字图像处理技术领域，本发明利用眼底血管呈线状这一特征，利用Hessian矩阵的特征向量找到线状目标的方向，构建方向场，并沿该方向对血管进行增强；采用延展性主动轮廓分割模型对眼底血管进行分割，方向场使模型可以跨越由于血管亮度不均引起的断裂血管间隙；本发明利用相应方向场引导曲线演化，可以演血管精确分割；本发明对眼底像进行亮度均衡处理，保证模型在亮度不均匀眼底像中实现对血管的准确提取；本发明利用眼球边缘曲率相等这一特征，有效去除血管分割时眼球边缘的影响，对低对比度图像有很好的分割结果。

Description

一种基于延展性主动轮廓模型的眼底图像血管分割方法

技术领域

本发明属于医学图像、计算机视觉和数字图像处理技术领域，具体涉及一种基于延展性主动轮廓模型的眼底图像血管分割方法。

背景技术

医学眼底视网膜图像的血管检测对心脑血管疾病的诊断和治疗在临床上具有重要的指导意义，血管分割结果的好坏直接影响到眼底检测的准确性，同时也影响到后续其他处理和分析的结果；

目前，研究者在医学眼底像血管分割方面提出了大量的不同算法，大致可以分为如下几类：基于窗口处理的方法、基于分类的方法和基于跟踪的方法等；基于窗口处理的方法首先对每一个像素在给定的模型下用包围该像素的窗口区域中像素进行匹配估计，或者把图像划分成很多同样尺寸的子图像，然后在每个子图像中分别计算，这种方法需要对图像中的每个像素处理，在每个像素处需要进行大量的操作，计算量随着图像的增大而迅速增加，而且这类算法受噪声的影响比较大；基于分类的方法首先利用低层图像处理方法进行简单分割，将图像区域分为血管类或非血管类，然后再通过监督学习的方法进一步确认血管，这种算法可用于分割之后的评估，其缺点是血管的不同尺度特性不能用来驱动精确分割；基于跟踪的方法首先确定一个初始的种子点，然后从种子点开始探测局部的图像性质来迭代跟踪血管，其缺点是现有的跟踪算法不能很好地分割分支和交叉处的血管，其可靠性严重依赖于所确定的初始种子点和方向；

受医学图像成像模式多样、图像噪声、血管形态复杂和血管边缘对比度低等因素影响，目前尚无一种通用的医学眼底像血管分割算法，因此，对眼底像分割方法的研究不仅具有深远的理论意义，而且具有广泛的实用意义。

发明内容

针对现有技术的不足，本发明提出一种基于延展性主动轮廓模型的眼底图像血管分割方法，以达到血管精确分割的目的。

一种基于延展性主动轮廓模型的眼底图像血管分割方法，包括以下步骤：

步骤1、根据图像灰度概率密度，对眼底图像进行亮度均衡；

步骤2、根据Hessian矩阵的特征向量，对眼底图像中的血管部分进行增强；

步骤3、确定增强后的眼底图像边缘点的曲率分布，进而在增强后的眼底图像中剔除眼球边缘轮廓；

步骤4、通过建立延展性主动轮廓模型的方式，对增强后的眼底图像进行分割，从而提取眼底图像血管；

具体步骤如下：

步骤4-1、根据响应函数，确定增强后眼底图像响应方向场；

步骤4-2、根据增强后眼底图像响应方向场，建立延展性主动轮廓模型的延展项；

步骤4-3、将所获延展项与测地线模型相结合，建立延展性主动轮廓模型；

步骤4-4、用差分法将延展性主动轮廓模型进行离散化处理，并设置延展性主动轮廓模型中水平集函数的初始值，将上述初始值代入延展性主动轮廓模型中进行迭代；

步骤4-5、判断相邻两次迭代获得的水平集函数之间的差值是否小于设定值，若是，则获得延展性主动轮廓模型的稳定解，该稳定解中的零水平集函数即为眼底血管边缘，否则，继续进行迭代；

步骤4-6、完成眼底图像的分割，并提取眼底图像血管。

步骤4-1所述的响应方向场，该场中每一个点的方向，由该点的Hessian矩阵特征向量决定。

步骤4-2所述的延展项，公式如下：

$F={\mathrm{\ν}}_o\·(1+\mathrm{\λk})\·|{\hat{u}}_{s,1}\·\stackrel{\→}{N}|---\left(1\right)$

其中，F表示延展性主动轮廓模型的延展项，ν_o表示响应函数，λ表示权重系数，k表示眼球边缘轮廓的曲率，表示某一点的Hessian矩阵的其中一个特征向量，s表示尺度因子，表示单位法向量。

步骤4-3所述的延展性主动轮廓模型公式如下：

${\mathrm{\φ}}_t=g\left(I\right)k|\▿\mathrm{\φ}|+\▿g\left(I\right)\·\▿\mathrm{\φ}+F---\left(2\right)$

其中，为梯度算子，φ_t表示水平集函数，t表示时间，即迭代次数；k表示眼球边缘轮廓的曲率，F表示延展性主动轮廓模型的延展项，g(I)函数是关于眼底图像I梯度的递减函数，称为边缘检测函数，公式如下：

$g\left(I\right)=\frac{1}{1+{\left|\right|G_{\mathrm{\σ}}*\▿I\left|\right|}_{L^2}^p}---\left(3\right)$

其中，L²表示L²范数，p为指数，p≥1，G_σ为Gaussian核函数。

本发明优点：

本发明利用延展性主动轮廓模型进行医学眼底像血管分割，创新点包括两个方面：一方面基于Hessian矩阵的眼底血管增强和具有延展性的主动轮廓模型用于眼底血管分割，利用眼底血管呈线状这一特征，利用Hessian矩阵的特征向量找到线状目标的方向，构建方向场，并沿该方向对血管进行增强；另一方面采用延展性主动轮廓分割模型对眼底血管进行分割，方向场使模型可以跨越由于血管亮度不均引起的断裂血管间隙；本发明思想新颖，分割结果好，本发明利用相应方向场引导曲线演化，可以演血管精确分割；本发明对眼底像进行亮度均衡处理，保证模型在亮度不均匀眼底像中实现对血管的准确提取；本发明利用眼球边缘曲率相等这一特征，可以有效去除血管分割时眼球边缘的影响；该发明对低对比度图像有很好的分割结果。

附图说明

图1为本发明一种实施例的基于延展性主动轮廓模型的眼底图像血管分割方法流程图；

图2为本发明一种实施例的原始眼底图像示意图；

图3为本发明一种实施例的井亮度均衡后的眼底图像示意图；

图4为本发明一种实施例的一段模拟血管示意图；

图5为本发明一种实施例的图5中血管的Hessian矩阵特征向量方向示意图；

图6为本发明一种实施例的图6中矩形框区域的局部放大图像示意图；

图7为本发明一种实施例的眼底像实例的增强图像示意图；

图8为本发明一种实施例的眼底像实例的增强血管图像示意图；

图9为本发明一种实施例的眼底像实例的响应方向场示意图；

图10为本发明一种实施例的图10中矩形框的局部放大区域示意图；

图11为本发明一种实施例的眼底像实例的最终分割结果示意图。

具体实施方式

下面结合附图对本发明一种实施例做进一步说明。

一种基于延展性主动轮廓模型的眼底图像血管分割方法，方法流程图如图1所示，包括以下步骤：

步骤1、根据图像灰度概率密度，对眼底图像进行亮度均衡；具体包括如下步骤：

步骤1-1、获得输入图像I的各级灰度概率密度函数：

$p\left(X_i\right)=\frac{n_i}{n}---\left(4\right)$

其中，n_i表示图像中灰度级为X_i的像素总和，n表示图像各灰度级像素总和；

步骤1-2、获得输入图像各灰度级的累计概率密度函数：

$C\left(X_i\right)=\underset{k=0}{\overset{i}{\mathrm{\Σ}}}p\left(X_k\right)---\left(5\right)$

步骤1-3、计算输出图像的灰度级：

f(X_i)＝[X₀+(X_L-1-X₀)·C(X_i)]  (6)

其中，f(X_i)表示灰度级为X_i时，经转换后输出的灰度级，[]代表取整操作；X₀表示灰度级为X₀，X_L-1灰度级为X_L-1，L表示灰度级总数；

亮度均衡化以后图像的灰度级概率分布将趋于均匀分布，图2为眼底图像原始图，图3是将图2进行亮度均衡后的结果图；

步骤2、根据Hessian(海森)矩阵的特征向量，对眼底图像中的血管部分进行增强；具体如下：

对于图像I中的某一点x₀，计算其邻域的二阶Taylor(泰勒)展开：

$I(x_0+\mathrm{\Δ}{x}_0,s)\≈I(x_0,s)+{{\mathrm{\Δx}}_0}^T{\▿}_{o,s}+\mathrm{\Δ}{x_0}^T{H}_{o,s}\mathrm{\Δ}{x}_0---\left(7\right)$

其中，s表示尺度因子，本发明实施例中s＝5，表示以s为尺度，在图像某一点x₀处的梯度向量：H_o，s表示以s为尺度，在图像某一点x₀处的Hessian矩阵：

${\▿}_{o,s}=\left(\begin{array}{c}{I}_x\left(x_0\right)\\ I_y\left(x_0\right)\end{array}\right),H_{o,s}=\left(\begin{array}{cc}{I}_{\mathrm{xx}}\left(x_0\right)& I_{\mathrm{xy}}\left(x_0\right)\\ I_{\mathrm{xy}}\left(x_0\right)& I_{\mathrm{yy}}\left(x_0\right)\end{array}\right)$

其中，I_x表示图像I对x求偏导，I_y为图像I对y的偏导，I_xx，I_xy，I_yx，I_yy均为二阶偏导；

通过计算得到Hessian矩阵的特征值λ₁，λ₂(其中λ₁表示值较小的特征值)和特征向量特征值和特征向量分别描述该点处二阶导数的大小与方向；

如图4所示，以一段血管为例，阐述Hessian矩阵特征向量和特征值的作用，图5中每个点的方向是图4中对应点的Hessian矩阵特征向量方向(图6是图5矩形框的局部放大图像)，附图5还标记了血管上某点的Hessian矩阵特征向量，可见计算的两个特征向量在空间内相互垂直，两个特征值中λ₁比λ₂小很多，且特征向量的方向为血管方向；根据这一特征，定义响应函数ν_o(s)：

其中，β和c是控制与S对响应函数精度的阈值，本发明实施例中，β＝c＝0.5；

利用响应函数对亮度均衡化后的眼底像进行操作，ν_o(s)作用于血管区域时，灰度显著增强，作用于其他区域则灰度基本不变；如图7所示，增强后的血管图像，由于眼球边缘也是线状，因此会同时被增强，因此提取血管时需要剔除边缘轮廓的干扰。

步骤3、确定增强后的眼底图像边缘点的曲率分布，进而在增强后的眼底图像中剔除眼球边缘轮廓；

本发明实施例中，提取强化图像的边缘，计算每个边缘点(x，y)的曲率k；然后统计所有边缘点的曲率分布曲线，眼球边缘呈圆形，边缘上的点曲率相等，根据这一特点计算所有边缘点的曲率分布直方图，选取直方图峰值对应的边缘点的集合，得到眼球边缘轮廓以及其在增强图像中的对应，从增强图像中剔除眼球边缘，保留血管；如图8所示，是剔除眼球边缘后得到的眼底血管的增强图像。

步骤4、通过建立延展性主动轮廓模型的方式，对增强后的眼底图像进行分割，从而提取眼底图像血管；

具体步骤如下：

步骤4-1、根据响应函数，确定增强后眼底图像响应方向场；

本发明实施例中，首先提取强化图像的响应方向场ν_o(如图9所示，图10是图9中矩形框区域的局部放大图像)，该场中的每一个点的方向由该点的Hessian矩阵特征向量决定，建立主动轮廓模型的延展项：

步骤4-2、根据增强后眼底图像响应方向场，建立延展性主动轮廓模型的延展项；

步骤4-2所述的延展项，公式如下：

$F={\mathrm{\ν}}_o\·(1+\mathrm{\λk})\·|{\hat{u}}_{s,1}\·\stackrel{\→}{N}|---\left(1\right)$

其中，F表示延展性主动轮廓模型的延展项，ν_o表示响应函数，λ表示权重系数，本发明实例中λ＝1，k表示眼球边缘轮廓的曲率，表示某一点的Hessian矩阵的其中一个特征向量，s表示尺度因子，表示单位法向量。

在矢量模为零(非血管区域)或曲线法线方向与矢量场方向相互垂直(血管边缘)的情况下，该牵引速度大小为0，在曲线法线方向与矢量场方向平行(血管中间区域)的情况下，该牵引速度达到最大。而且曲率越大牵引力越大，使得曲线在细小目标处演化时，不会因为正则力而停止演化。为了保证曲线在演化过程中的正则性并接近目标边缘，该延展项与测地线模型相结合，得到延展性主动轮廓模型的演化方程，即延展性主动轮廓模型；

步骤4-3、将所获延展项与测地线模型相结合，建立延展性主动轮廓模型；

延展性主动轮廓模型公式如下：

${\mathrm{\φ}}_t=g\left(I\right)k|\▿\mathrm{\φ}|+\▿g\left(I\right)\·\▿\mathrm{\φ}+F---\left(2\right)$

其中，为梯度算子，φ_t表示水平集函数，t表示时间，即迭代次数；k表示眼球边缘轮廓的曲率，F表示延展性主动轮廓模型的延展项，g(I)函数是关于眼底图像I梯度的递减函数，称为边缘检测函数，公式如下：

$g\left(I\right)=\frac{1}{1+{\left|\right|G_{\mathrm{\σ}}*\▿I\left|\right|}_{L^2}^p}---\left(3\right)$

其中，L²表示L²范数，p为指数，p≥1，本发明实施例中p＝2，G_σ为Gaussian核函数。

延展性主动轮廓模型第一项引导曲线向图像边缘移动，寻找血管边缘，第二项在血管边缘形成凹状速度场，曲线演化到血管边界附近时将会锁定在附近，防止边缘泄露，第三项为延展项，作用是分割过程中引导演化曲线快速沿血管方向延伸。

步骤4-4、用差分法将延展性主动轮廓模型进行离散化处理，并设置延展性主动轮廓模型中水平集函数的初始值，将上述初始值代入延展性主动轮廓模型中进行迭代；

用差分法将演化方程离散化，任意给定初始水平集函数φ₀，本发明实施例中，用半径为R的圆的水平集函数作为初始水平集函数：

步骤4-5、判断相邻两次迭代获得的水平集函数之间的差值是否小于设定值，若是，则获得延展性主动轮廓模型的稳定解，该稳定解中的零水平集函数即为眼底血管边缘，否则，继续进行迭代；

本发明实施例中，经过N次迭代后(迭代终止条件是|φ_t+1-φ_t|＜0.01)得到离散化方程的稳定解φ_N，φ_N的零水平集函数即为眼底血管边缘；φ_N中大于零的点集对应于血管上的点，小于零的点集对应于背景点。

步骤4-6、完成眼底图像的分割，并提取眼底图像血管。

如图11所示，是对眼底像中血管的最终分割结果。

Claims (4)

1.一种基于延展性主动轮廓模型的眼底图像血管分割方法，其特征在于，包括以下步骤：

步骤1、根据图像灰度概率密度，对眼底图像进行亮度均衡；

步骤2、根据海森矩阵的特征向量，对眼底图像中的血管部分进行增强；

步骤3、确定增强后的眼底图像边缘点的曲率分布，进而在增强后的眼底图像中剔除眼球边缘轮廓；

步骤4、通过建立延展性主动轮廓模型的方式，对增强后的眼底图像进行分割，从而提取眼底图像血管；

具体步骤如下：

步骤4-1、根据响应函数，确定增强后眼底图像响应方向场；

步骤4-2、根据增强后眼底图像响应方向场，建立延展性主动轮廓模型的延展项；

步骤4-3、将所获延展项与测地线模型相结合，建立延展性主动轮廓模型；

步骤4-4、用差分法将延展性主动轮廓模型进行离散化处理，并设置延展性主动轮廓模型中水平集函数的初始值，将上述初始值代入延展性主动轮廓模型中进行迭代；

步骤4-5、判断相邻两次迭代获得的水平集函数之间的差值是否小于设定值，若是，则获得延展性主动轮廓模型的稳定解，该稳定解中的零水平集函数即为眼底血管边缘，否则，继续进行迭代；

步骤4-6、完成眼底图像的分割，并提取眼底图像血管。

2.根据权利要求1所述的基于延展性主动轮廓模型的眼底图像血管分割方法，其特征在于，步骤4-1所述的响应方向场，该场中每一个点的方向，由该点的海森矩阵特征向量决定。

3.根据权利要求1所述的基于延展性主动轮廓模型的眼底图像血管分割方法，其特征在于，步骤4-2所述的延展项，公式如下：

$F=v_o\·(1+\mathrm{\λk})\·|{\hat{u}}_{s,1}\·\stackrel{\→}{N}|---\left(1\right)$

其中，F表示延展性主动轮廓模型的延展项，v_o表示响应函数，λ表示权重系数，k表示眼球边缘轮廓的曲率，表示某一点的海森矩阵的其中一个特征向量，s表示尺度因子，表示单位法向量。

4.根据权利要求1所述的基于延展性主动轮廓模型的眼底图像血管分割方法，其特征在于，步骤4-3所述的延展性主动轮廓模型公式如下：

${\mathrm{\φ}}_t=g\left(I\right)k|\▿\mathrm{\φ}|+\▿g\left(I\right)\·\▿\mathrm{\φ}+F---\left(2\right)$

其中，为梯度算子，φ_t表示水平集函数，t表示时间，即迭代次数；k表示眼球边缘轮廓的曲率，F表示延展性主动轮廓模型的延展项，g(I)函数是关于眼底图像I梯度的递减函数，称为边缘检测函数，公式如下：

$g\left(I\right)=\frac{1}{1+{\left|\right|G_{\mathrm{\σ}}*\▿I\left|\right|}_{L^2}^p}---\left(3\right)$

其中，L²表示L²范数，p为指数，p≥1，G_σ为Gaussian核函数。