CN110033861B

CN110033861B - 适用于octa图像的血管及黄斑无血管区定量分析方法及系统

Info

Publication number: CN110033861B
Application number: CN201910232818.9A
Authority: CN
Inventors: 许夏瑜; 雷剑琴; 徐峰; 丁文祥
Original assignee: Xian Jiaotong University
Current assignee: Xian Jiaotong University
Priority date: 2019-03-26
Filing date: 2019-03-26
Publication date: 2021-03-23
Anticipated expiration: 2039-03-26
Also published as: CN110033861A

Abstract

本发明公开了一种适用于OCTA图像的血管及黄斑无血管区定量分析方法及系统。本发明通过多尺度线匹配方法实现OCTA图像上的大血管及毛细血管的区分，并采用图论建模方法提升了黄斑无血管区的分割效果，具有稳定性好、精确度高的特点。本发明为OCTA图像的高精度分割提供了有效途径。

Description

适用于OCTA图像的血管及黄斑无血管区定量分析方法及系统

技术领域

本发明属于计算机辅助诊断领域，具体涉及一种适用于OCTA图像的血管及黄斑无血管区定量、自动分析方法及系统。

背景技术

黄斑水肿是糖尿病视网膜病变中一种严重威胁视力的病变。对糖尿病患者的常规眼底检查方法包括眼底彩照、血管造影图像以及相干光断层成像(OCT)。眼底彩照和OCT图像无法对视网膜毛细血管循环系统进行很好的捕捉，而血管造影图像由于其有创性，无法大规模推广使用。

OCTA是一种非侵入性的眼底影像检查技术，可高分辨率识别视网膜及脉络膜的血流运动信息，对活体组织视网膜及脉络膜微血管循环成像。在视网膜及脉络膜血管疾病改变的管理随访和治疗效果检测等方面具有独特优势。尽管国内外已提出多种算法对OCTA图像特征(Characteristics)进行分割及定量(Quantification)，但由于黄斑水肿内部缓慢流动的液体，导致其在OCTA图像上具有与血管类似的信号强度及纹理，使得单纯基于信号强度的分割方法容易产生错误。

对具有黄斑水肿病变的OCTA图像进行正确前处理(Preprocessing)，是利用图像特征有效预测威胁视力的糖尿病视网膜病变的关键技术。目前亟待提出一种针对黄斑水肿病变的OCTA图像分析方法。

发明内容

本发明的目的在于提供一种适用于OCTA图像的血管及黄斑无血管区定量分析方法及系统。

为达到上述目的，本发明采取的技术方案为：

一种适用于OCTA图像的血管及黄斑无血管区定量分析方法，包括以下步骤：

1)通过多尺度线检测对OCTA图像上的大血管进行加强，然后进行阈值分割，得到大血管分割图像，对大血管分割图像进行图像细化得到大血管骨架图像；

2)通过Hessian矩阵对所述OCTA图像上的血管网络进行加强，然后进行血管网络分割，然后利用步骤1)中的大血管分割图像滤除血管网络中的大血管，得到毛细血管分割图像，对毛细血管分割图像进行图像细化，得到毛细血管骨架图像；

3)对步骤1)中的大血管分割图像进行距离变换后建立图割模型，利用图割模型得到黄斑无血管区的分割图像。

优选的，所述步骤1)具体包括以下步骤：

1.1)多尺度线检测算子的设计：选取m个尺度的正方形算子，每个尺度Li选取d个方向，i＝1,2,…,m，得到m×d个线检测算子，d大于等于4；

1.2)利用线检测算子滤波：在每个尺度Li下，利用对应d个方向的线检测算子各自计算得到的灰度平均值中的最大值，对原始图像I的像素灰度进行匹配滤波，得到该尺度下的滤波结果；其中，原始图像I为所述OCTA图像；

1.3)多尺度滤波结果的融合：将各尺度线检测算子滤波结果的平均值与原始图像I的对应像素灰度值求和，得到融合图像J；

1.4)提取大血管图像：利用阈值分割处理融合图像J后进行去噪，得到大血管分割图像(记为V_large)，对V_large进行连续细化，得到大血管骨架图S_large。

优选的，所述步骤1)还包括以下步骤：根据V_large及S_large进行大血管定量分析，得到大血管密度LVD以及大血管骨架密度LSD。

优选的，所述步骤2)具体包括以下步骤：

2.1)基于Hessian矩阵的血管增强：以任意像素点(x,y)为中心对原始图像I进行Hessian滤波，得到Hessian滤波的结果图像H，Hessian滤波的结果为多个尺度下像素点(x,y)的最大灰度值，其中，Hessian滤波的尺度集合为{σ}，尺度σ的取值≥1，σ的步长为1～3；

2.2)去除大血管：对Hessian滤波的结果图像H进行阈值分割得到所有血管分割图像V，在V上减去V_large，得到毛细血管分割图像(记为V_cap)，对V_cap进行连续细化得到毛细血管骨架图S_cap。

优选的，所述步骤2)还包括以下步骤：根据V_cap及S_cap进行毛细血管定量分析，得到毛细血管密度CVD以及毛细血管骨架密度CSD。

所述步骤3)具体包括以下步骤：

3.1)提取前景及背景种子：对V_large进行距离变换得到距离变换图像D，在距离变换图像D上寻找灰度最大值I_max对应的像素点p(x,y)，以p(x,y)为中心在距离变换图像D上建立选区，得到选区内像素点的集合S_f，所述选区的形状选自边长为1～2倍I_max的正方形、该正方形的外接圆或该圆的其他(除正方形外)内接多边形；将距离变换图像D上位于所述选区外的部分与V_large取交集，得到像素点集合S_b；

3.2)黄斑无血管区分割：以像素点的集合S_f为前景种子，以像素点的集合S_b为背景种子，建立图割模型并进行优化，得到黄斑无血管区分割结果。

优选的，所述步骤3)还包括以下步骤：对根据步骤3.2)确定的黄斑无血管区FAZ进行定量分析，得到区域面积(area)、区域周长(perimeter)，以及区域圆度(R₁)。

一种适用于OCTA图像的血管及黄斑无血管区定量分析系统，包括OCTA图像定量分析软件，该分析软件包括分别用于原始图像导入(Input image)、自动分析及结果显示的功能模块；其中，原始图像导入模块具体用于载入一个或一组待分析的OCTA图像，自动分析模块具体包括用于执行以上步骤1)、步骤2)及步骤3)的若干个子模块(例如，图像前处理子模块、血管特征提取子模块、定量分析计算子模块)，结果显示模块具体用于输出以上各步骤的图像及定量分析结果。

上述分析软件还包含账户管理、结果存储和远程传输等功能模块。

本发明相比于现有技术，具有以下有益效果：

本发明通过多尺度线匹配方法首次实现了在黄斑水肿OCTA图像上对大血管的单独提取，不仅为大血管定量分析及毛细血管网络的分割提供了有效的分析方法，而且结合采用图论建模方法解决了现有单纯基于信号强度的分割方法无法有效提取黄斑水肿OCTA图像的黄斑无血管区的问题，实现了对黄斑水肿病变下的黄斑无血管区的有效分割，具有稳定性好、精确度高的特点。

本发明提出的分析系统在自动实现对黄斑水肿病变下的黄斑无血管区的有效分割的同时，还提高了对黄斑水肿引发的视网膜病变的分析方法的可操作性，降低了临床应用的难度，可以实现方法从工作站到用户端的转移，便于用户使用。

附图说明

图1为本发明系统的工作流程图；其中：(a)显示了原始图像I，(b1)显示了不同尺度(m＝3)线检测滤波结果图像，(b2)显示了融合图像J，(b3)显示了距离变换灰度图像D，(b4)显示了Hessian滤波的结果图像H，(b5)显示了所有血管(即血管网络)分割图像V，(c)显示了大血管分割图像V_large，(d)显示了毛细血管分割图像V_cap，(e)显示了黄斑无血管区分割图像。

图2为基于线检测的大血管分割流程图；其中：(a)显示了原始图像I，(b)显示了不同尺度的滤波结果图像，(c)显示了融合图像J，(d)显示了大血管分割图像V_large。

图3为基于图割方法的黄斑无血管区分割流程图；其中：(a)显示了原始图像I，(b)显示了背景种子(顶层)、原始图像I(中间)、前景种子(底层)的叠加，(c)显示了图割模型的建立，其中S表示图割模型中的源点，T表示图割模型中的汇聚点，虚线表示优化后的分割结果，(d)显示了基于图割模型分割的黄斑无血管区，其对病变区进行了精确的分割，(e)为分割结果在原始图像I上的显示，其中线圈指出了传统单纯基于强度的分割结果，其位于(d)的分割结果之内，且所标注的黄斑无血管区的范围明显小于实际范围。

具体实施方式

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

如图1所示，适用于OCTA图像上的血管及黄斑无血管区定量分析方法，包括以下步骤：

步骤一：基于多尺度线检测的大血管分割及定量分析

在步骤一中，本发明提出了一种利用多尺度线检测方法对OCTA图像上的视网膜大血管(即前毛细血管，Precapillary)进行单独提取的方法，包括多尺度线检测算子的设计、利用线检测算子滤波、多尺度滤波结果的融合、阈值分割提取大血管分割图像，及大血管的定量分析，具体流程如下。

针对不同物理尺寸的OCTA图像，尺度个数(m)可以适度调整(即在≥2的范围内选取尺度个数)，但尺度个数过大，则不仅增加计算消耗，也对图像分割精度没有进一步改善作用，设计中选取的最大尺度的线检测算子，其检测窗口应当能覆盖原始图像中最大血管的径向范围，某尺度下线检测方向的选取以沿0～180°均匀分布为准。以原始图像I采用物理尺寸为6×6mm²的黄斑水肿OCTA图像为例，图像大小为1024×1024像素。首先在多尺度线检测算子的设计中，如图2所示，一共设计了三个尺度的正方形算子，具体尺度大小为L1＝27×27像素，L2＝37×37像素，L3＝47×47像素。在每个尺度下，都有八个方向的线检测算子(0°、23°、45°、67°、90°、113°、135°、157°)。因此，共有24个线检测算子。

在每个尺度的检测窗口下，用八个方向的线检测算子分别对原始图像I进行逐像素匹配，并根据八个方向计算中的最大值得到滤波结果，单个像素的匹配滤波的具体操作方法是：在每个检测窗口中，沿不同线检测(Line detection)算子方向分别计算原始图像I上该检测窗口内所有像素的灰度平均值，并确定沿各方向计算的灰度平均值中的最大值

i＝1,2,3，之后执行以下公式，并将J^Li作为窗口中心像素点的滤波结果：

式中，

为尺度L_i窗口下原始图像I的平均灰度值。

在三个尺度下重复上述操作，可得到三幅对原始图像进行不同尺度线检测匹配滤波的结果。通过融合(Combination)各尺度下的滤波结果，得到最终的多尺度线检测(Multiscale line detecion)结果(图像J)：

即最终的多尺度线检测滤波结果为三个尺度线检测结果(J^L1、J^L2、J^L3为对应尺度滤波结果图像上任意一点的灰度值)的平均值与原始图像I上对应点的灰度值求和。

用最大熵阈值(Maximum Entropy)分割处理图像J，之后进行去噪处理(去除面积小于200像素的斑点区域)，即得到大血管分割图像V_large，对V_large进行图像连续细化(连续细化方法参见Gonzalez等人的著作《数字图像处理(第三版)》)，即可得到大血管骨架图S_large。

对大血管的定量分析方法如下：

式中，X(x,y)为原始图像I中的环形区域面积，该环形区域的中心一般位于图像的中心，内径的选取以能够覆盖黄斑无血管区为准，外径与原始图像I的尺寸相应，故内径取为1mm、外径取为6mm；v_lar(x,y)为V_large图像中对应面积下的非零(灰度值)像素数，s_lar(x,y)为S_large图像中对应面积下的非零(灰度值)像素数。

步骤二：毛细血管分割及定量分析，包括基于Hessian矩阵的血管增强滤波、阈值分割获得血管网络图像、去除大血管，及对毛细血管网络进行定量分析，具体流程如下。

Hessian矩阵是一个多元函数的二阶偏导数构成的方阵，描述了函数的局部曲率，Hessian矩阵的表达形式为：

式中，尺度σ的最小值为3，最大值为30，σ的步长为1，I_xx(x,σ)、I_xy(x,σ)、I_yx(x,σ)及I_yy(x,σ)表示对原始图像I在x和y方向上的一阶偏导I_x(x,σ)和I_y(x,σ)在对应方向进一步求偏导而得到的二阶偏导。

在每个σ下，以任意像素点(x,y)为中心对原始图像I进行滤波，并取所有尺度下的滤波最大值(灰度值)作为像素点(x,y)的滤波结果。遍历原始图像I的所有像素点，得到Hessian滤波的结果图像H，对H进行Otsu阈值分割，得到所有血管(All vessels)的分割图像V，在V上减去大血管分割图像V_large，得到毛细血管(Capillary)分割图像V_cap，对V_cap进行图像连续细化，即得到毛细血管骨架图S_cap。

对毛细血管的定量分析方法如下：

式中，X(x,y)为原始图像I中的环形区域面积，该环形区域的内、外径等选取参照大血管的定量分析方法；v_cap(x,y)为V_cap图像中对应面积下的非零(灰度值)像素数，s_cap(x,y)为S_cap图像中对应面积下的非零(灰度值)像素数。

步骤三：对黄斑水肿病变下的黄斑无血管区进行分割及定量分析，包括提取前景及背景种子、建立图割模型、对模型进行最优化，及对黄斑无血管区进行定量分析，具体流程如下。

参见图3，首先确定前景及背景种子点：以大血管分割图像V_large为基础，对V_large进行距离变换，得到距离变换图像D，在D上寻找灰度最大值I_max对应的像素点p(x,y)，以p(x,y)为中心建立边长为1.5I_max的正方形区域。距离变换图像D上该正方形区域内的像素点构成集合S_f，即前景种子；将距离变换图像上该正方形区域外的部分与V_large取交集后，得到像素点集合S_b，即背景种子。

将S_f与S_b输入图割模型并进行最优化(图割模型建立及优化方法参见Boykov等人的“Interactive Graph Cuts for Optimal Boundary&Region Segmentation of Objectsin N-D Images”IEEE Computer Society,2001年)，得到黄斑无血管区的分割结果。

对黄斑无血管区FAZ的圆度R₁的定量分析方法如下：

式中，area为FAZ的区域面积，perimeter为FAZ的区域周长。在本实施例中，area根据黄斑无血管区的非零像素数计算，perimeter根据黄斑无血管区的边缘像素个数计算。

基于Windows的OCTA血管定量分析系统，包括OCTA图像定量分析软件。此软件包含以下几个模块：1)账户管理模块：系统允许多个用户在同一个手机上使用，每个用户都可以通过注册(Register)获得账号和密码，分析前首先进行用户登录(Login)；2)载入OCTA图像模块：用户可以从指定文件夹载入图像或载入整个文件夹的图像；3)自动分析模块：以载入的OCTA图像为原始图像，完成自动的图像分析，分析内容包括上述大血管分割及定量分析、毛细血管分割及定量分析、黄斑无血管区的分割及定量分析；4)结果的显示和存储模块：图像处理的结果可以输出显示和存储，并可以通过电子邮件的方式发送；5)远程传输模块：可选择将分析得到的结果图像及定量分析结果等数据传输至指定的远程工作站或客户端(例如，医生或患者的手机、电脑、平板等智能终端)；6)使用说明模块：此模块中提供了详细的操作说明和对处理结果的解释。

在测试阶段，将利用本发明分析得到的黄斑无血管区的分割结果(简称自动分割结果)与临床医生的手动分割结果进行了比较测试。在一组30张OCTA图像的数据集分析中，由临床医生手动分割并计算了所有测试图像上的黄斑无血管区面积，与此同时通过本发明介绍的自动分析方法对测试图像进行了独立处理。结果用Dice系数(DC)来评测：

式中，A1和A2分别为手动分割结果和自动分割结果。

本发明方法在30幅测试图像上的平均DC值达到0.90，精度与手动分割相当，分析时间较手动分割显著缩短。

总之，本发明针对现有对视网膜血管进行无差别提取和分析的方法进行了改进，首次实现了黄斑水肿OCTA图像上的大血管及毛细血管区别提取，优化了黄斑水肿情况下的黄斑无血管区的分割提取，并实现了基于软件平台的自动分析，为疾病状态下的OCTA图像的高精度分割及定量分析提供了有效、便利的途径。

Claims

1.一种适用于OCTA图像的血管及黄斑无血管区分析方法，其特征在于：包括以下步骤：

1)通过多尺度线检测对OCTA图像上的大血管进行加强，然后进行阈值分割，得到大血管分割图像，所述大血管为前毛细血管；

2)对大血管分割图像进行距离变换后建立图割模型，利用图割模型得到黄斑无血管区的分割图像；

所述步骤2)具体包括以下步骤：

2.1)提取前景及背景种子：对大血管分割图像进行距离变换得到距离变换图像，在距离变换图像上寻找灰度最大值I_max对应的像素点p(x,y)，以p(x,y)为中心在距离变换图像上建立选区，得到选区内像素点的集合S_f，所述选区的形状选自边长为1～2倍I_max的正方形、该正方形的外接圆或该圆的其他内接多边形，将距离变换图像上位于所述选区外的部分与大血管分割图像取交集，得到像素点集合S_b；

2.2)黄斑无血管区分割：以像素点的集合S_f为前景种子，以像素点的集合S_b为背景种子，建立图割模型并进行优化，得到黄斑无血管区分割结果。

2.根据权利要求1所述的分析方法，其特征在于：所述步骤1)具体包括以下步骤：

1.2)利用线检测算子滤波：在每个尺度Li下，利用对应d个方向的线检测算子计算的灰度平均值中的最大值，对所述OCTA图像的像素灰度进行匹配滤波，得到该尺度下的线检测算子滤波结果；

1.3)多尺度滤波结果的融合：将各尺度下的线检测算子滤波结果的平均值与所述OCTA图像的对应像素灰度值求和，得到融合图像J；

1.4)提取大血管图像：利用阈值分割处理融合图像J后进行去噪，得到大血管分割图像。

3.根据权利要求1所述的分析方法，其特征在于：所述步骤2)还包括以下步骤：根据黄斑无血管区的分割图像，对黄斑无血管区进行定量分析，得到黄斑无血管区的区域面积、区域周长，以及区域圆度。

4.一种适用于OCTA图像的血管及黄斑无血管区分析系统，其特征在于：该分析系统包括原始图像导入模块及自动分析模块，所述原始图像导入模块用于载入OCTA图像，自动分析模块包括图像前处理子模块及血管特征提取子模块，图像前处理子模块包括多尺度线检测算子，图像前处理子模块通过多尺度线检测以及对各尺度下的线检测结果进行融合，使OCTA图像上的大血管得到加强，并对大血管分割图像进行距离变换得到用于建立图割模型的距离变换图像，所述血管特征提取子模块用于对图像前处理子模块输出的融合结果进行阈值分割，并生成所述大血管分割图像，以及利用所建立的图割模型得到黄斑无血管区的分割图像，所述大血管为前毛细血管；

所述利用所建立的图割模型得到黄斑无血管区的分割图像具体包括以下步骤：在距离变换图像上寻找灰度最大值I_max对应的像素点p(x,y)，以p(x,y)为中心在距离变换图像上建立选区，得到选区内像素点的集合S_f，所述选区的形状选自边长为1～2倍I_max的正方形、该正方形的外接圆或该圆的其他内接多边形，将距离变换图像上位于所述选区外的部分与大血管分割图像取交集，得到像素点集合S_b；以像素点的集合S_f为前景种子，以像素点的集合S_b为背景种子，建立图割模型并进行优化，得到黄斑无血管区分割结果。

5.根据权利要求4所述的分析系统，其特征在于：所述图像前处理子模块还包括Hessian滤波器，图像前处理子模块通过对所述OCTA图像进行Hessian滤波，使OCTA图像上的血管网络得到加强，并对加强后的血管网络进行分割，所述血管特征提取子模块通过利用大血管分割图像滤除血管网络中的大血管得到毛细血管分割图像；所述分析系统还包括定量分析计算子模块，定量分析计算子模块用于根据黄斑无血管区的分割图像对黄斑无血管区进行区域面积、区域周长以及区域圆度的计算，以及在对大血管分割图像和/或毛细血管分割图像分别进行细化后完成对血管密度VD以及血管骨架密度SD的计算。