CN108492285B - 一种基于灰度变化的径向线视盘萎缩弧分割方法 - Google Patents

Abstract

本发明公开了一种基于灰度变化的径向线视盘萎缩弧分割方法，属于图像分割技术领域。该方法能够自动分割眼底图中的视盘萎缩弧，之后简称PPA；并且其分割结果为PPA和面积为主的指标计算提供了便利，该算法可以有效提高医生诊断的效率和准确率。具体步骤如下：首先通过灰度预处理将彩色眼底图像转化为灰度图；其次，采用径向线对视盘及其外侧区域的像素进行采样，选取灰度变化最大的两个点作为PPA和视盘的候选边界点；然后剔除受血管或噪声影响产生的异常点；最后对候选边界点进行椭圆拟合，拟合出的两个椭圆构成的最大闭合月牙形区域即为PPA的分割结果。

Description

一种基于灰度变化的径向线视盘萎缩弧分割方法

技术领域

本发明涉及一种基于灰度变化的径向线视盘萎缩弧分割方法，属于图像分割技术领域。

背景技术

视盘萎缩弧(以下简称PPA)的产生与近视和青光眼相关。已有研究表明，PPA的大小与病理性近视的程度成正相关，即病理性近视越严重，PPA面积越大。但并非仅有病理性近视患者的眼底中存在PPA，轻度和中度近视患者的眼底也存在PPA；另一方面，儿童近视的发病率逐年升高，如何准确地预测和防控儿童近视成为当下亟待解决的问题。所以PPA在预防近视，诊断相应眼科疾病等方面具有越来越重要的意义。

目前，PPA尺寸等指标的检测大多由专业的眼科医生通过肉眼观察实现，这一过程耗时耗力，且无法定量描述PPA面积的变化过程。所以可以完成PPA的自动分割并计算其面积的数字图像处理算法具有极大的应用前景。目前已有的算法通过主动轮廓模型实现PPA的分割，但是这种方法具有较高的灵敏度，不能精确获取PPA的边缘。

发明内容

本发明的目的在于克服现有PPA分割方法灵敏度不高以及不能精确获取PPA边缘的技术缺陷，提供一种基于灰度变化的径向线视盘萎缩弧分割方法，能够精确地实现PPA的分割，并计算得到PPA的面积。

一种基于灰度变化的径向线视盘萎缩弧分割方法，包括如下步骤：

步骤一：灰度图像的预处理，具体为：处理灰度图像，所以应将给定的彩色眼底图像转变为灰度图；为了减少眼底图中血管对PPA分割的影响，采用血管表现较不明显的红色通道作为灰度图像进行处理；

步骤二：采用等角度间隔径向线提取边界候选点：首先在视盘中心构建等角度径向线，使得所有的等角度径向线相交于视盘中心，并且都穿过视盘的边界；此时会有一部分等角度径向线穿过PPA区域；沿着每条等角度径向线对图像进行采样，可得到该等角度径向线对应的灰度曲线，即可得到灰度变化曲线；因为PPA与视盘和视网膜背景之间的灰度差异，通过求解灰度变化曲线中最大的两个点即可确定PPA边界候选点和视盘边界候选点；

步骤二、具体包括如下步骤：

步骤2.1：构建等角度径向线；具体为：以视盘中心为交点，构造N条长度为P的等角度径向线，每条等角度径向线的方向以θ₀＝2π/N为间隔均匀分布；

步骤2.2：沿着每条等角度径向线在图像上计算采样点；对于每条等角度径向线，以p₀为步长在图像中提取若干采样点，采样点坐标(x,y)的计算公式如下(1)：

其中，p_i＝p_i-1+p₀，i＝1,2,…,P/p₀表示采样点(x,y)距视盘中心的距离，其中，P表示等角度径向线的长度；θ＝n×θ₀表示当前等角度径向线n的倾角；x₀，y₀分别表示视盘中心的横坐标和纵坐标；sinθ和cosθ分别表示倾角θ的正弦值和余弦值；

步骤2.3：对采样点进行插值采样；采用双线性插值计算采样点(x,y)对应的采样灰度，记为

其计算公式如(2)：

其中，x₁，y₁，x₂，y₂为点(x,y)的邻近点坐标，u＝x-x₁，s＝y-y₁为两个中间变量，I_r表示图像；

因此，采样点(x,y)处的灰度变化，

可以表示为公式(3)：

其中，

与

分别为第i-1和第i+1个采样点上的采样灰度；

步骤2.4：采样平均；单一等角度径向线上的灰度变化对噪声非常敏感，因此，取邻近的M条等角度径向线L_n-(M-1)/2，…，L_n-2，L_n-1，L_n，L_n+1，L_n+2，…，L_n+(M-1)/2中对应采样点的像素值均值作为等角度径向线L_n上该采样点的像素值，以此减少灰度变化曲线中的波形抖动；

步骤2.5：选取边界候选点；查找得到灰度变化曲线的前两个最大值，这两个最大值对应的位置为当前等角度径向线上的边界候选点；

步骤三：剔除异常点：通过检测灰度变化曲线中的正-负峰值剔除等角度径向线穿过血管时产生的异常点；通过边界点的位置约束剔除噪声引起的异常点；

步骤三、具体包括：

步骤3.1：在灰度变化曲线图中选取所有的局部极小值点，并将变化率的绝对值大于某一阈值的局部极小值点成为显著局部最小值；此时计算每个显著局部极小值点与边界候选点的距离；如果该间隔小于某一阈值，则认为其为一对正负峰值对；该正负峰值对表明当前边界候选点是因血管干扰而产生的异常点，应被剔除；如果该间隔大于或等于这一阈值就保留该边界候选点；

步骤3.2：若当前边界点与视盘中心的距离满足以下公式(4)的条件，则将其判定为异常点：

|D-D_mean|>YH_d (4)

其中，D是当前边界点与视盘中心的距离，D_mean为已检测到的所有边界点到视盘中心的平均距离，TH_d为限定剔除范围的阈值；

步骤3.3：在视盘的鼻侧通常不存在PPA，即鼻侧检测到的PPA边界和视盘边界重合；因此，若检测到的PPA边界点灰度变化值小于一定的阈值，则认为PPA和视盘边界重合，直接选取视盘边界点作为PPA边界候选点；如果检测到的PPA边界点灰度变化值大于或等于这一阈值，则认为PPA和视盘边界不重合，此时不做任何操作，仅跳至步骤四；

步骤四：提取PPA的闭合边界曲线：对提取的所有边界候选点进行拟合，获得PPA的闭合边界曲线；其中，所有边界候选点包括PPA边界候选点和视盘边界候选点；

步骤四、具体包括：

步骤4.1：采用最小二乘法对提取到的视盘边界候选点和PPA边界候选点进行椭圆拟合；

步骤4.2：提取两个椭圆围成的封闭月牙形区域；面积最大的区域即为PPA的分割结果，再根据分割结果计算一些临床可供参考的相关参数；

其中，相关参数主要为PPA的面积、PPA面积与视盘面积的比值；

至此，从步骤一到步骤四，完成了一种基于灰度变化的径向线PPA分割方法。

有益效果

1、本发明提供了一种基于灰度变化的径向线PPA分割方法，可以实现PPA的自动分割；

2、本发明所述方法采用了径向线图像采样方法，对采样结果平均降噪，提高了方法抗噪性；

3、本发明所述方法具有一定的应用价值和商业价值，它可以嵌入于眼底照相机中，实时动态提示眼底PPA的位置与各项指标；

4、本发明所述方法也可应用于临床科研与临床诊断中，自动分割眼底图中的PPA并给出各项指标，从而降低医生的工作量，提高诊断的速度、效率和精度。

附图说明

图1为本发明一种基于灰度变化的径向线PPA分割方法及实施例中的流程示意图。

具体实施方式

下面结合附图并举实施例，对本发明进行详细描述。

实施例1

图1是本发明具体实施方式中PPA分割算法的流程图，具体包括如下步骤：

步骤A：灰度的预处理：选取由北京同仁医院眼科研究所提供的100幅彩色眼底图作为样本，提取该100幅彩色眼底图的红色通道作为待处理灰度图像。因为在红色通道中，眼底图中的血管相对于背景较不明晰。此时血管对PPA分割的影响较小。

步骤B：为了提取大量的PPA边界候选点，首先在视盘中心构建等角度径向线，然后沿着每一条等角度径向线对图像进行采样。此时可以得到每一条等角度径向线对应的灰度变化曲线。因为PPA与视盘和视网膜背景之间的灰度差异，通过求解灰度变化曲线中最大的两个点即可确定PPA边界候选点和视盘边界候选点。

步骤B.1：首先构造等角度径向线。以视盘中心为原点，构造N＝180条长度为P＝200的等角度径向线，每条等角度径向线的方向以θ₀＝2π/N为间隔均匀分布。

步骤B.2：对于每条等角度径向线，以p₀为步长在图像中提取若干采样点，采样点坐标(x,y)的计算如下(5)公式：

其中，p_i＝p_i-1+p₀，(i＝1,2,…,P/p₀)表示采样点(x,y)距原点的距离，P为等角度径向线的长度；θ＝n×θ₀表示当前等角度径向线n的倾角，sinθ和cosθ分别表示当前等角度径向线倾角θ的正弦值和余弦值。注意采样步长p₀不宜过大也不宜过小，若p₀过小，则相邻两采样点间的灰度变化不明显；若p₀过大，边缘点将不能准确定位；因此，在这里设定p₀＝5。

步骤B.3：采用双线性插值计算采样点(x,y)对应的采样灰度。其计算公式如下：

其中，x₁，y₁，x₂，y₂为点(x,y)的邻近点坐标，u＝x-x₁，s＝y-y₁为两个中间变量，I_r表示图像。因此，采样点(x,y)处的灰度变化可以表示为：

其中

与

分别为第i-1和第i+1个采样点上的采样灰度。

步骤B.4：采样平均。单一等角度径向线上的灰度变化对噪声非常敏感，因此，取邻近的5条等角度径向线L_n-2，L_n-1，L_n，L_n+1，L_n+2中对应采样点的像素值均值作为等角度径向线L_n上该采样点的像素值，以此减少灰度变化曲线中的波形抖动。

步骤B.5：选取边界候选点。查找得到的灰度变化曲线的前两个最大值，这两个最大值对应的位置为当前等角度径向线上的边界候选点。

此时从B.2到B.5步遍历所有的等角度径向线，可以找到等角度径向线所有方向上的边界候选点。

步骤C：剔除异常点：当径向线穿过血管时，由血管造成的灰度曲线的变化有可能使算法得到异常的候选边界点。此类型的异常点会在灰度变化曲线中引入一对正-负峰，故通过检测灰度变化曲线中的正-负峰值剔除此类异常点；另外对于噪声引起的异常点，可通过边界点的位置约束剔除。

步骤C.1：对于等角度径向线的每个方向，在其灰度变化曲线图中选取所有的局部极小值点，并称变化率的绝对值超过设定阈值的局部极小值点为显著局部极小值点。此时计算每个显著局部极小值点和候选边界点之间的间隔，如果某一个显著局部极小值和候选边界点的间隔小于某个设定的阈值，则认为其为一对正-负峰值。此时可以剔除该位置处检测到的候选点。

步骤C.2：若当前边界点与视盘中心的距离满足以下条件，则将其判定为异常点：

|D-D_mean|>TH_d (8)

其中，D是当前边界点与视盘中心的距离，D_mean为已检测到的所有边界点到视盘中心的平均距离，TH_d为限定剔除范围的阈值。

步骤C.3：在视盘的鼻侧通常不存在PPA，即鼻侧检测到的PPA边界和视盘边界重合。因此，若检测到的PPA边界点灰度变化值过小，则认为此处PPA和视盘边界重合，直接选取视盘边界点作为PPA边界候选点。

步骤D：提取PPA的闭合边界曲线：对提取的所有边界候选点(包括PPA的边界候选点和视盘的边界候选点)进行拟合，最终获得PPA的闭合边界曲线。

步骤D.1：因为提取到的视盘边界候选点和PPA+视盘区域的外围边界候选点近似椭圆，故对其进行椭圆拟合，以得到光滑连续的边界线。这里采用最小二乘法对其进行椭圆拟合。

步骤D.2：若不存在误差，第一步得到两个椭圆围成的一个封闭月牙形区域即为PPA的分割区域；但是实际上这两个椭圆往往会围成多个封闭月牙形区域。一般而言，面积最大的区域即为PPA的分割区域。

自此，就实现了PPA分割算法的全部过程。实验验证，该方法可以有效实现彩色眼底图中PPA病变区域的分割，并且对于失焦的模糊图像也有较好的分割效果。测试结果表明依据分割结果计算出的PPA面积等指标与医生标定值高度相关，因此能很好地辅助临床诊断。

综上所述，以上仅为本发明的较佳实施例而已，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (5)

1.一种基于灰度变化的径向线视盘萎缩弧分割方法，其特征在于：能够精确地实现PPA的分割，并计算得到PPA的面积；包括如下步骤：

步骤一：灰度图像的预处理，具体为：处理灰度图像，所以应将给定的彩色眼底图像转变为灰度图；

步骤二：采用等角度间隔径向线提取边界候选点：首先在视盘中心构建等角度径向线，使得所有的等角度径向线相交于视盘中心，并且都穿过视盘的边界；此时会有一部分等角度径向线穿过PPA区域；沿着每条等角度径向线对图像进行采样，可得到该等角度径向线对应的灰度曲线，进而可得到灰度变化曲线；因为PPA与视盘和视网膜背景之间的灰度差异，通过求解灰度变化曲线中最大的两个点即可确定PPA边界候选点和视盘边界候选点；

步骤三：剔除异常点：通过检测灰度变化曲线中的正-负峰值剔除等角度径向线穿过血管时产生的异常点；通过边界点的位置约束剔除噪声引起的异常点；

步骤三、具体包括：

步骤3.1：在灰度变化曲线图中选取所有的局部极小值点，并将变化率的绝对值大于某一阈值的局部极小值点成为显著局部最小值；此时计算每个显著局部极小值点与边界候选点的距离；如果该间隔小于某一阈值，则认为其为一对正负峰值对；该正负峰值对表明当前边界候选点是因血管干扰而产生的异常点，应被剔除；如果该间隔大于或等于这一阈值就保留该边界候选点；

步骤3.2：若当前边界点与视盘中心的距离满足以下公式(4)的条件，则将其判定为异常点：

|D-D_mean|＞TH_d (4)

其中，D是当前边界点与视盘中心的距离，D_mean为已检测到的所有边界点到视盘中心的平均距离，TH_d为限定剔除范围的阈值；

步骤3.3：在视盘的鼻侧通常不存在PPA，即鼻侧检测到的PPA边界和视盘边界重合；因此，若检测到的PPA边界候选点灰度变化值小于一定的阈值，则认为PPA和视盘边界重合，直接选取视盘边界候选点作为PPA边界候选点；如果检测到的PPA边界候选点灰度变化值大于或等于这一阈值，则认为PPA和视盘边界不重合，此时不做任何操作，仅跳至步骤四；

步骤四：提取PPA的闭合边界曲线：对提取的所有边界候选点进行拟合，获得PPA的闭合边界曲线；其中，所有边界候选点包括PPA边界候选点和视盘边界候选点。

2.根据权利要求1所述的一种基于灰度变化的径向线视盘萎缩弧分割方法，其特征在于：步骤一为了减少眼底图中血管对PPA分割的影响，采用血管表现较不明显的红色通道作为灰度图像进行处理。

3.根据权利要求1所述的一种基于灰度变化的径向线视盘萎缩弧分割方法，其特征在于：步骤二、具体包括如下步骤：

步骤2.1：构建等角度径向线；具体为：以视盘中心为交点，构造N条长度为P的等角度径向线，每条等角度径向线的方向以θ₀＝2π/N为间隔均匀分布；

步骤2.2：沿着每条等角度径向线在图像上计算采样点；对于每条等角度径向线，以p₀为步长在图像中提取若干采样点，采样点坐标(x，y)的计算公式如下(1)：

其中，p_i＝p_i-1+p₀，i＝1，2，...，P/p₀表示采样点(x，y)距视盘中心的距离，其中，P表示等角度径向线的长度；θ＝n×θ₀表示当前等角度径向线n的倾角；x₀，y₀分别表示视盘中心的横坐标和纵坐标；sinθ和cosθ分别表示倾角θ的正弦值和余弦值；

步骤2.3：对采样点进行插值采样；采用双线性插值计算采样点(x，y)对应的采样灰度，记为

其计算公式如(2)：

其中，x₁，y₁，x₂，y₂为点(x，y)的邻近点坐标，u＝x-x₁，s＝y-y₁为两个中间变量，I_r表示图像；

因此，采样点(x，y)处的灰度变化，

可以表示为公式(3)：

其中，

与

分别为第i-1和第i+1个采样点上的采样灰度；

步骤2.4：采样平均；单一等角度径向线上的灰度变化对噪声非常敏感，因此，取邻近的M条等角度径向线L_n-(M-1)/2，…，L_n-2，L_n-1，L_n，L_n+1，L_n+2，…，L_n+(M-1)/2中对应采样点的像素值均值作为等角度径向线L_n上该采样点的像素值，以此减少灰度变化曲线中的波形抖动；

步骤2.5：选取边界候选点；查找得到灰度变化曲线的前两个最大值，这两个最大值对应的位置为当前等角度径向线上的边界候选点。

4.根据权利要求1所述的一种基于灰度变化的径向线视盘萎缩弧分割方法，其特征在于：步骤四、具体包括：

步骤4.1：采用最小二乘法对提取到的视盘边界候选点和PPA边界候选点进行椭圆拟合；

步骤4.2：提取两个椭圆围成的封闭月牙形区域；面积最大的区域即为PPA的分割结果，再根据分割结果计算一些临床可供参考的相关参数。

5.根据权利要求4所述的一种基于灰度变化的径向线视盘萎缩弧分割方法，其特征在于：步骤4.2中，相关参数主要为PPA的面积、PPA面积与视盘面积的比值。

Images