糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法与流程

本发明涉及一种病变相关性分析的方法，特别是一种糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法。

背景技术：

SD OCT是一种快速、非侵入式的频域光学相干断层成像技术，它可以有效地呈现视网膜各组织层的结构及病变的临床病理特征，临床实验表明SD OCT图像能够用于视网膜的厚度测量、病变的检测与识别，特征分析等。糖尿病性视网膜病变是糖尿病最常见和最严重的微血管并发症之一，是糖尿病患者视力丧失的主要原因。

SD OCT糖尿病性视网膜病变图像中视网膜外丛状层和外核层分布着一些点块状高反射信号，这些高反射信号在糖尿病性视网膜病变的机制还不完全清楚。高反射信号主要是脂质渗出物，蛋白质物质或炎症细胞堆积之后的一种形态学显示，随后可能变成硬性渗出的前驱。目前硬性渗出和高反射信号被广泛的进行研究用于视力下降的预防和治疗。有研究表明高反射信号的增加会导致视力受损，而大部分的高反射信号可能被认为是硬性渗出的前驱，意味着硬性渗出与高反射信号可能存在某种相关性。但是现有技术中尚无相关硬性渗出与高反射信号关系的描述。

技术实现要素：

本发明的目的在于提出一种糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法。

实现本发明的目的的技术解决方案为：一种糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法，包括以下步骤：

步骤1、采集SD OCT视网膜图像和CFP眼底图像；

步骤2、采用双边滤波算法对SD OCT视网膜图像进行去噪；

步骤3、利用二维的层分割方法定位视网膜神经纤维层RNFL和光感受器内节/外节IS/OS，从而限制高反射信号所在区域；

步骤4、利用自适应阈值法确定高反射信号种子点集合，并根据基于人类视觉特性的区域生长法分割出高反射信号区域；

步骤5、根据SD OCT视网膜图像生成投影图像，与CFP眼底图像配准，获得裁剪后CFP眼底图像；

步骤6、利用多尺度的显著性检测法结合阈值法分割出裁剪后CFP眼底图像中的硬性渗出；

步骤7、根据裁剪后CFP眼底图像定位SD OCT视网膜图像中的硬性渗出，在SD OCT视网膜图像中进行硬性渗出与高反射信号的相关性分析。分别得出硬性渗出与高反射信号在面积、数目、高度、灰度这些特征上存在的关系。

本发明与现有技术相比，其显著优点：(1)本发明首次给出了一种糖尿病性视网膜病变图像中硬性渗出和高反射信号的分析方法；(2)本发明通过配准CFP眼底图像和SD OCT投影图像来确定SD OCT视网膜图像中的硬性渗出区域；(3)本发明提取了一组特征(面积，数目，灰度分布，高度分布)用于研究SD OCT视网膜图像中硬性渗出与高反射信号间的关系。

下面结合附图对本发明作进一步详细描述。

附图说明

图1是本发明糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法的流程图。

图2是构建高反射信号约束区域的流程图。

图3是SD OCT视网膜图像中高反射信号所在区域示意图。

图4是SD OCT视网膜图像约束区域。

图5是灰度渐变距离图像。

图6是IS/OS边界图像。

图7是高反射信号约束区域图像。

图8是高反射信号种子点集合图像。

图9是高反射信号区域二值结果。

图10是高反射信号区域轮廓。

图11是CFP眼底图像与SD OCT投影图像配准的结果，其中图11(A)为待配准的CFP眼底图像，图11(B)为SD OCT投影图像，图11(C)为配准后叠加的图像，图11(D)为配准后裁剪得到的眼底图像。

图12是裁剪后CFP眼底图像中硬性渗出分割的结果，其中图12(A)为裁剪后的眼底图像，图12(B)为硬性渗出分割结果，图12(C)为剔除假硬性渗出后的结果。

图13是结合SD OCT视网膜图像和CFP眼底图像分析高反射信号和硬性渗出的结果图，其中图13(A)为裁剪得到的眼底图像，图13(B)为眼底图像中硬性渗出区域和SD OCT投影图像中高反射信号区域对比图，图13(C)和(D)分别为SD OCT视网膜图像其中两帧的硬性渗出和高反射信号区域的对比图。

图14是SD OCT视网膜图像中硬性渗出和高反射信号在灰度分布和高度分布上的统计结果，其中图14(A)为灰度分布统计结果，图14(B)为高度分布统计结果。

图15是SD OCT视网膜图像中硬性渗出和高反射信号在面积和数目上的相关性分析结果，其中图15(A)为两者在面积上的相关性分析结果，图15(B)为两者在数目上的相关性分析结果。(图14和图15中，左边图为对于非增殖性糖尿病性视网膜病变NPDR的分析结果，右边图为对于增殖性糖尿病性视网膜病变PDR的分析结果)

具体实施方式

结合图1，本发明糖尿病性视网膜病变图像中硬性渗出与高反射信号的分析方法包括以下步骤：

步骤1、采集SD-OCT视网膜图像和CFP眼底图像，采用现有的OCT成像设备和CFP成像设备对视网膜图像进行采集；

步骤2、采用双边滤波算法对输入图像进行去噪处理，双边滤波算法的公式为：

式中f和h分别为输入和输出图像，函数c(ξ,x)用于测量邻域中心点x和邻域点ξ之间的空间距离，函数s用于测量两点间的灰度相似性，函数c和函数s都是高斯函数，是归一化函数。

步骤3、利用二维的层分割方法定位RNFL和IS/OS，从而限制高反射信号所在区域，具体为：

步骤3-1、采用图论的方法去估计內界膜(ILM)边界，主要利用垂直方向的梯度和像素灰度构造边的权值，具体权值构造如下：

式中I_a和I_a为相邻两个像素点a和b的灰度值，g_a和g_b表示像素点a和b的归一化梯度图像，λ代表相似度系数，取值为2.0。是一个很小的保真项，取值为1×10^-5，*为卷积操作。

然后采用动态规划策略求解最短路径，即ILM边界，具体最小权值函数定义：

且式中w为窗口大小，这里取作2，n_i,j表示图在第i行，第j列的节点，i-w≥1，i+w≤M，M是图一列节点的数目。

步骤3-2、从ILM边界往下一定距离寻找梯度最大处的像素点，即为RNFL下边界。梯度最大处像素点的集合定义为：

其中i,j表示图像中的第i行，第j列，Δ取值为90，g为梯度算子，即[1,1,1,1,0,-1,-1,-1,-1]。

步骤3-3、根据原图像的灰度渐变图像结合曲线拟合的方法分割IS/OS边界，所述曲线拟合方法为四次多项式曲线拟合，具体为：

步骤3-3-1、首先根据下式构造灰度渐变图G:

步骤3-3-2、求出灰度渐变图的梯度图像，采用滤波器[-1；1]对灰度渐变图像进行滤波得到灰度渐变图像的梯度图。找到梯度图像中每一列灰度最大的像素点，同时排除灰度小于最大灰度一半的像素点，然后对剩下的像素点进行曲线拟合得到曲线f₁；

步骤3-3-3、在曲线f₁上方200像素单位距离内寻找灰度最大的像素点，根据最初找到的灰度最大点初次拟合曲线后，排除曲线上方的像素点，再次对剩下的像素点进行拟合，得到曲线f₂；

步骤3-3-4、在曲线f₂上方30像素范围内寻找梯度最大的像素点，并拟合曲线，排除曲线上方的像素点，再次拟合得到最终的曲线，即IS/OS边界。

步骤4、利用自适应阈值法确定高反射信号种子点集合，并根据基于人类视觉特性的区域生长法分割出高反射信号区域的具体步骤为：

步骤4-1、通过迭代最大类间方差法逐步逼近最佳阈值，不断缩小目标区域的方法来分割出最终目标，即高反射信号的种子点集合(此处为现有技术，可参考文献：王茜,彭中,刘莉.一种基于自适应阈值的图像分割算法[J].北京理工大学学报,2013,23(4):521-524.)；

步骤4-2、对阈值分割后的图像进行膨胀处理，膨胀窗口为[1,1,1；1,1,1；1,1,1]。将刚可分辨亮度差作为区域生长的相似性准则，进行区域生长。刚可分辨亮度差定义为：

ΔI是关于图像灰度I的分段函数，α，β，C为待定参数，α取值在[0.01，0.3]之间，中间区间的起始灰度值a和终止灰度值b可以分别取[45,81]和[180,210]之间。于是，参数取值为α＝0.05，a＝50，b＝200；根据在函数间断点取值相等可以得到CC＝αa，从而计算出C和β。

区域生长之后，剔除亮斑灰度与其边缘外5×5范围像素平均灰度差小于d的假亮斑。边缘外5×5范围像素是指以边缘的每一个像素为中心，半径大小为5个像素的窗口所覆盖且不属于亮斑的像素，d取值为20。

步骤5、根据SD OCT视网膜图像生成投影图像，与CFP眼底图像配准，获得裁剪后CFP眼底图像，具体为：

利用半自动的配准方法来配准CFP眼底图像和SD OCT投影图像，主要根据两者图像中血管分支点相互对应的原则生成一个几何变换模型，基于这个变换模型可以匹配得到眼底图像中对应于SD OCT投影图像的区域，并裁剪相匹配的眼底图像用于下一步中硬性渗出的分割。

步骤6、利用多尺度的显著性检测法结合阈值法分割出裁剪后CFP眼底图像中的硬性渗出，具体为：

步骤6-1、采用显著性检测方法来分割裁剪后CFP眼底图像中的硬性渗出，对于给定的一个尺度，像素的显著性定义为局部区域和其邻域像素特征向量的之间的欧式距离：

C_i,j＝||v₁-v₂||

其中，V₁表示图像局部区域的像素特征向量，V₂表示邻域的像素特征向量，使用CIELab彩色空间来定义像素特征向量，V₁＝[L₁,a₁,b₁]^T，V₂＝[L₂,a₂,b₂]^T，最终在所有不同尺度上的显著性总和为：

步骤6-2、得到裁剪后CFP眼底图像的显著图后，使用阈值法分割出硬性渗出所在区域，将分割出的假的硬性渗出，直接剔除。

步骤7、根据裁剪后CFP眼底图像定位SD OCT视网膜图像中的硬性渗出，在SD OCT视网膜图像中进行硬性渗出与高反射信号的相关性分析，具体为：

首先根据匹配后的眼底图像中硬性渗出的位置确定SD OCT视网膜图像中的硬性渗出区域，然后分别提取SD OCT视网膜图像中硬性渗出和高反射信号(包含硬性渗出)的灰度分布，高度分布，面积和数目特征，利用这些特征分析硬性渗出和高反射信号在SD OCT视网膜图像中的相关性，之后通过相关性系数cc，p值和统计分布来进行分析，最终得出硬性渗出与高反射信号在面积、数目、高度、灰度这些特征上存在的关系。

下面结合实施例对本发明做进一步详细的说明：

实施例

本系统发明以SD OCT视网膜体数据和CFP数据作为输入，采用图像处理和统计分析手段对SD OCT视网膜图像中的硬性渗出和高反射信号进行相关性分析。

本实施例的流程如图1所示，通过OCT成像设备采集到的SD-OCT视网膜体数据大小为1024×512×128，图3给出了一帧SD-OCT视网膜图像的感兴趣区域，图中标注了视网膜的几个主要相关组织结构，如RNFL下边界，IS/OS边界，高反射信号区域。为了便于后续的处理，首先对输入图像进行去噪处理，然后采用图论的方法定位ILM边界。利用垂直方向梯度和灰度构造边的权值，采用动态规划策略求解最短路径，并利用定位的ILM边界估计RNFL下边界。由于高反射信号区域往往出现在视网膜中间层，首先利用灰度渐变距离图像估计脉络膜层，依次估计出Bruch膜(BM)边界和IS/OS边界，根据RNFL下边界和IS/OS边界限制高反射信号检测区域。图4-7说明了IS/OS边界图像和高反射信号检测约束区域图像。得到高反射信号检测约束区域图像后，利用自适应阈值法确定高反射信号种子点集合，如图8所示。然后根据基于人类视觉特性的区域生长法分割出高反射信号区域。高反射信号区域分割结果如图10所示。

根据SD OCT投影图像中和眼底图像中血管分支点相对应的原则，建立几何变换模型，进行图像的配准，并裁剪得到配准后的眼底图像。图11展示了待配准的CFP眼底图像(A)，SD OCT投影图像(B)，配准叠加后的图像(C)和配准后裁剪得到的眼底图像(D)。

对于裁剪得到的眼底图像，运用显著性检测结合阈值法分割出图像中的硬性渗出。图12(A)展示了裁剪后的眼底图像(有棉绒斑)，图12(B)为硬性渗出分割结果，图12(C)为剔除假硬性渗出后的结果。

根据眼底图像中的硬性渗出位置确定SD OCT视网膜图像中的硬性渗出区域，SD OCT投影图像中高反射信号区域与眼底图像硬性渗出区域重叠的部分则为SD OCT视网膜图像中的硬性渗出区域，如图13所示，其中(A)为裁剪得到的眼底图像，(B)为眼底图像中硬性渗出区域和SD OCT投影图像中高反射信号区域对比图，(C)和(D)分别为SD OCT视网膜图像其中两帧的硬性渗出和高反射信号区域的对比图。

对于SD OCT视网膜图像中的硬性渗出个高反射信号分别提取一组特征，包括面积，数目，灰度分布和高度分布。然后针对这些特征分别在这两者之间进行相关性分析。图14(A)和(B)分别展示了硬性渗出和高反射信号在灰度分布和高度分布上的统计结果。图15(A)和(B)分别展示了硬性渗出和高反射信号在面积和数目上的相关性分析结果。

实验结果表明，SD OCT糖尿病性视网膜图像中硬性渗出和高反射信号在面积、数目这两个特征存在正相关的关系，硬性渗出的平均高度要大于高反射信号的平均高度，两者的灰度分布则比较相似。