黄斑中心视网膜OCT图像的分层分析研究 眼科学专业

目录摘要 1Abstract 2第一章绪论 31.1课题背景与研究意义 31.2相关技术介绍 41.2.1OCT技术 41.2.2图像分割技术 41.2.3眼底OCT图像分层 41.3研究内容及章节安排 5第二章相关算法介绍 62.1Canny算子 62.1.1Canny算子原理 62.1.2Canny算子改进 72.2.相关处理 82.2.1中值滤波 82.2.2曲线拟合 82.2.3高斯滤波 92.2.4膨胀算法 10第三章OCT图像分层算法 113.1算法流程 113.2算法分步 123.3错误校正及结果优化 163.3.1错误校正 163.3.2结果优化 19第四章实验结果与分析 224.1误差评价标准 224.2实验结果分析 224.3总结 24参考文献 25致谢 26摘要医学成像技术的成熟发展促使越来越多的电子成像技术应用于临床医学。光学相干断层扫描技术(OCT)因其较高分辨率和不与眼部接触的特性而被广泛应用。然而仅依靠医生手工对OCT眼底图像进行分割，既耗时又耗力,再加上不可忽视的人为因素的影响，因此研究一个OCT视网膜图像自动分层算法是非常必要的。本文首先介绍了OCT技术、图像分割技术的发展现状及眼底OCT图像的分层，然后介绍了研究内容及章节安排。第二章介绍Canny算子的原理及其改进，以及相关的算法原理。在第三章本文给出了整个自动分层算法的流程图，并在此基础上对整个流程进行分步详细阐述，基于不同参数下的Canny算子初步分层后，再对分层过程中的出现的错误比如跳变误差和中央凹进行了校正，最后对各分层结果进行优化。本文描述的方法可以找到视网膜OCT图像的7个边界，以专家手动分割的结果为标准，误差分析显示，该自动分层算法的各层误差控制在2个像素左右，算法整体检测性能良好。关键词：OCT，自动分层，Canny算子，边缘检测

AbstractThematuredevelopmentofmedicalimagingtechnologyhaspromptedmoreandmoreelectronicimagingtechnologytobeappliedtoclinicalmedicine.Opticalcoherencetomography(OCT)iswidelyusedbecauseofitshighresolutionandnocontactwiththeeye.However,itisverynecessarytostudyanOCTretinalimageautomaticlayeringalgorithm,whichreliesondoctors'manualsegmentationofOctfundusimage,whichistime-consumingandconsuming,coupledwiththeinfluenceofhumanfactorswhichcan’tbeneglected.ThispaperintroducesthedevelopmentofOCT,imagesegmentationandthestratificationofOCTimagefirstly,thenintroducestheresearchcontentandchapterarrangement.Thesecondchapterintroducestheprincipleandimprovementofcannyoperator,andsomerelatedalgorithmprinciples.Inthethirdchapter,wegivetheflowchartofthewholeautomaticlayeredalgorithm,onthisbasis,thewholeprocessiselaboratedstep-by-step,basedontheinitiallayeringofthecannyoperatorunderdifferentparameters,theerrorssuchasjumperrorandcentralconcaveinthelayeredprocessarecorrected,andtheresultsareoptimized.Themethodpresentedinthispapercanfind7boundariesofOCTimage,andtheresultofexpertmanualsegmentationisthestandard,theerroranalysisshowsthattheerrorcontroloftheautomaticlayeredalgorithmisaround2pixels,andthealgorithmhasgoodperformanceinoveralldetection.Keywords:OCT,automaticlayersegmentation,Cannyoperator,edgedetection.

第一章绪论1.1课题背景与研究意义眼底是眼球内后部的组织,如图1.1所示，眼底包括脉络膜、视网膜、黄斑、眼底血管、视神经纤维和视神经乳头等部分。图1.1眼底结构图通过检查眼底可以发现很多疾病，毛细血管瘤、出血点和渗出物为早期糖尿病的征兆；通过检测眼底分泌物、血管病变图像或眼底相关部位的病变，可诊断出眼底疾病，如青光眼、白内障、老年性黄斑病变等，因此眼底被眼科医生誉为“检目镜”[1]。获得眼底图像的方法有很多，主要包括超声影像、红外图像和新兴的光学相干断层扫描图像。传统方法是依靠有经验的医生观察病人的医学影像中的眼底特征信息发现病变，主要特征为视神经盘、中央窝、黄斑区域视网膜内层厚度等。但随着人口老龄化的加剧，眼部疾病的发病率也越来越高，而眼底病病种繁多,病因复杂,使得眼底检查这项工作既耗时又耗力,再加之不可忽视的人为因素的影响，因此，如何快速有效地处理医学影像并检测出眼部疾病有着重大意义，也具有非常大的市场潜力[2]。1.2相关技术介绍1.2.1OCT技术OCT全称是光学相干断层成像技术(OpticalCoherencetomography)。在眼科应用治疗中，由于其具有很高的分辨率，可以精确辨别内膜性质和支架结构，能够清晰地显示病变区域，具有非接触性检查等优点，所以OCT技术给我们提供了一个全新的窗口，让我们能够全面观察眼底疾病，尤其是眼底黄斑疾病。经过长时间的改进，OCT已经发展为一种临床实用的检测工具，制成了商品化仪器，并于1995年开始正式用于眼科临床。眼科OCT应用前景十分广阔，未来发展包括能包括高速、三维视网膜成像，功能的成像及使用不同波长范围的光源增强穿透力至脉络膜成像等[3]。1.2.2图像分割技术图像分割是指是在预处理的基础上对信息进行组织和加工的过程，实质就是建立在像素之间的连续性和非连续性上将像素分类。图像分割技术作为数字图像处理中的一个重要组成部分，不但能够按照一定的准则将所选定的区域划分成具有一定功能的区域，也能将图像中特定的边缘、区域等特征提取出来。图像分割技术应用很广泛，几乎出现在图像处理的各个领域，但是还没有一种能够适应各种图片的处理方法，所以还有待于不断地深入研究[4]。1.2.3眼底OCT图像分层视网膜位于眼球壁的内层，是一层透明的薄膜。下图1.2描述了正常人眼的三维OCT视网膜图像分层的结果。由上而下分别为：色素上皮层，外节层、连接纤毛层，内节层，外核层，外丛状层，内颗粒层，内丛状层，神经节细胞层，神经纤维层，这10个层次共形成11个边界。视网膜OCT图像各层之间灰度值不同，存在明显的暗亮变化，由暗到亮或者由亮到暗的位置区域就是我们要寻找的边界。本文所寻找的7个边界，自上而下分别为surface1，surface2，surface4，surface6，surface7，surface10和surface11，为了方便行文，以下叙述将用c1、c2、c3、c4、c5、c6和c7分别对应这7个边界。图1.2正常人眼的OCT视网膜图像分层结果图1.3研究内容及章节安排本文的主要内容是设计一个自动分层算法，在误差较小的前提下，该算法能够检测视网膜OCT图像中的多个边界。主要行文安排如下：第一章首先介绍OCT技术发展，图像分割技术现状以及所要研究的OCT视网膜分层图像。第二章首先介绍Canny算子原理及其改进，再介绍分层算法中的中值滤波,高斯滤波,曲线拟合和膨胀算法。第三章先给出自动分层算法的流程图，再详细阐述分层算法的步骤，最后进行错误校正和结果优化。第四章先对自动分层算法寻找的边界明确一个结果评价标准，再从两个维度进行误差分析和总结。

第二章相关算法介绍2.1Canny算子因为边缘总是出现在相邻的两个灰度值不同的区域连接处，特征表现为灰度值的不连续。而图像边缘信息主要集中在高频段，通常说图像锐化或检测边缘，实质就是高频滤波。Canny算子是一个优良的边缘检测算子，是JohnF.Canny于1986年提出来的。Canny算子具有如下优点：较高的信噪比；定位精度高；具有较强单边响应能力[2]。Canny算子具有较强的抑制噪声的能力，可以较完整地检测边缘，但相比其他传统的边缘微分算子计算复杂，运算速度慢。另外，Canny算子的双阈值是根据全局特征信息来决定的，一方面导致无法消除局部噪声的干扰，另一方面又会丢失灰度值变换缓慢的边缘。2.1.1Canny算子原理利用Canny算子进行边缘检测,包括以下4个步骤：（1用高斯滤波器平滑图像首先用二维高斯函数对图像进行平滑，函数如下[6]。Gx,y,σ公式(2-1)是求权重的公式，σ值是高斯函数的分布参数，它控制着平滑程度，σ值越小，平滑效果越差，但边缘定位精度高；σ值越大，平滑效果越好，但边缘定位精确度差，所以平滑图像时要根据情况选择适当的σ值令g(x,y)为平滑后的图像，用Ggx,y=G(2)用一阶有限差分计算梯度幅值和方向.平滑后的数组g(x,y)的梯度可以用2×2邻域一阶有限差分近似式分别表示x、y方向的偏导数，即Gx≈(g[Gy≈(gx,y-在这个邻域内求有限差分的均值，可得到该点x方向与y方向的偏导数梯度。梯度幅值和梯度方向角表达式如下[7]：G=Gx2+Gy2θ=arctan⁡(GxGy)(3)对梯度幅值进行非极大值抑制梯度幅值越大，对应的图像梯度值也就越大，但这并不足以确定边缘。为了确定边缘，必须保留幅值局部变化最大的点，这一过程称为“非极大值抑制”。非极大值抑制通过比较当前点的梯度强度和正负梯度方向点的梯度强度，如果相比同方向的其他点的梯度强度是最大，则保留其值，否则抑制。(4)用双阈值算法检测和连接边缘Canny算法应用双阈值，τ1和τ2，其中τ1<τ2，可得到阈值边缘两个图像T1[i,j]和T2[i,j]，由于图像T2[i,j]是高阈值得到的，去除了大部分噪声，但在轮廓上会有间断，所以我们把图像T2[i,j]的边缘连接，当检测到断点处，就在图像T12.1.2Canny算子改进因为视网膜OCT图像分层涉及到两种类型的边界：由暗到亮和由亮到暗，所以我们要对Canny算子做一个改进，以满足不同类型的边界区分。Canny算子中的梯度大小如公式2-1所示，因为自动分层算法中检测的边缘是水平方向的，所以我们只需考虑从上向下看的y方向的梯度差，现修改如下：边界由暗到亮时，取y方向正梯度；M1=max⁡(Gy边界由亮到暗时，取y方向负梯度的绝对值。M2=|min⁡(Gy,0)2.2.相关处理2.2.1中值滤波中值滤波是选用一个含有奇数个像素的滑动窗口，该窗口在图像上扫描，把其中所含的像素点按照灰度值的大小进行排序，用中间的灰度值来代替窗口中心点的灰度值。用数学表达为gi=Medfi-式中，Med(…),表示取序列的中值，m表示有m个像素，i表示窗口的中心位置。图2.1可看出中值滤波有较好的平滑效果[8]。图2.1中值滤波效果图2.2.2曲线拟合曲线拟合（curvefitting）是指通过得到的数据来拟合一条曲线方程，并尽可能使最多的离散的点分布在所拟合的曲线方程上。效果图可见图2.2。图2.2曲线拟合效果图在MATLAB中，利用函数polyfit进行多项式拟合。函数polyfit根据观测数据及用户指定的多项式阶数得到光滑曲线的多项式表示，polyfit的一般调用格式为：P=polyfit(x,y,n)。其中x为自变量，y为因变量，n为多项式阶数[9]。2.2.3高斯滤波高斯滤波可以消除高斯噪声，是一种线性平滑滤波。高斯滤波所得到的每一个像素点的值，都由该像素本身和其邻域内的其他的像素值经过加权平均后得到。一维零均值高斯函数为：gx=e-x22其中，高斯分布参数σ决定了高斯函数的宽度。下图2.3可见高斯滤波的效果。图2.3高斯滤波效果图2.2.4膨胀算法膨胀对原图像有扩张的作用,可以较好地填充物体内部的空洞以及连接间距小于结构元素的相邻目标区域[10]。对X中的每一点膨胀定义为：Y=X⊕B={y:B(y)∩X≠∅}（2-用B（y)对X进行膨胀，就是把结构元素B平移后，把B与X的交集非空的点构成集合，过程如图2.4所示。图2.4膨胀运算过程第三章OCT图像分层算法3.1算法流程视网膜Canny算子从亮到暗检测,确定c7结束结果显示误差分析对各层进行修正，处理跳变误差和中心凹在c5和c7之间寻找像素最小值,确定cCanny算子从亮到暗检测,确定c7结束结果显示误差分析对各层进行修正，处理跳变误差和中心凹在c5和c7之间寻找像素最小值,确定c6Canny算子从亮到暗检测,确定c4、c3和c2Canny算子从暗到亮检测，确定c1和c5读图像，输入待检测的帧数开始图3.1分层算法流程图由流程图可看出，我们设计的这个自动分层算法是按照寻找边界的难易顺序决定的。先确定比较明显的边界位置（c1和c5），再根据边界之间的位置关系，缩小下一个边界所检测的区域，就可以找到大致的边界位置了。但是分层算法最重要的是它的普遍性，即对大部分眼底OCT图像都适用。对一个图像适用的参数，其他图像就不一定适用。还有各个边界所处的位置不同，或者所选取的帧不同，结果也会不同，因此有必要针对各个边界的进行算法修正。3.2算法分步（1）显示原始图像（以图像4为例）首先输入待检测的图像，把图像读成一维数组，然后再重组为三维数组，输入所要检测的帧数（以第36帧为例）。图3.2OCT原始图像（2）Canny算子由暗到亮检测将Canny算子参数设置为低阈值=0.1，高阈值=0.5，σ=6，所检测到的两个边界就是c1和c5（图3.3），因为c1和c5是最明显的两条从暗到亮的边界，所以较为容易检测到。我们将自上而下检测到的每列第一个白点位置标记为c1，如此便确定了c1的位置[11]。图3.3Canny算子找到的暗到亮边界确定第一层位置后，就可以进一步确定第五层的位置，首先我们要用形态学的膨胀算法将第一层消除，然后再检测消除过c1的图像，从上而下找的每列第一个白点就组成第五层的位置了。为了完整的将c1消除，我们找到宽度为1的第一层（图3.4a）,将其膨胀后（图3.4b）再取反，然后与原图像进行与运算，即可。(b)图3.4像素宽度为1的第一层膨胀前后（a）膨胀前（b）膨胀后将膨胀结果与图3.3进行逻辑与运算得到图3.5，找到每列的第一个白点即可确定第五层位置。图3.5第五层位置（2）Canny算子由亮到暗检测①将Canny算子参数设置为低阈值=0.1，高阈值=0.5，σ=2，此时阈值较大。结果如图3.6所示。图3.6较大阈值Canny算子找到的亮到暗边界因为此时阈值设置较大，故只能依上而下检测到第二层和第七层。在y方向上，从下往上第一个白点的位置即可确定为第七层位置。但是第七层下方有血管阴影，存在干扰，下图3.7中可看到，我们利用之前找到的第五层，将搜索区域限制在第五层下方30个像素以内，这样很好的限制了第七层下方的干扰[12]。②将Canny算子参数设置为低阈值=0.05，高阈值=0.2，σ=3，此时阈值较小。结果如图3.7所示。图3.7较小阈值Canny算子找到的亮到暗边界因为此时阈值设置较小，自上而下可检测到第二层、第三层、第四层和第七层。上文中我们已经找到了第七层，因此，为了避免第七层对第四层的位置的影响，我们将搜索区域限制在第一步所找的到的第五层上方，这样的话，依据从下而上找的第一个白点即可确定第四层的位置；如此，找到第四层后，用同样的方法将搜索区域限制在第四层以上10个像素（根据以往大量参考资料和经验所得）的位置，即可确定第三层位置；再如此，在第三层10个像素以上可找到第二层[13]。（3）确定第六层位置利用Canny算子很难找到合适的参数检测到第六层，第五层和第七层之间的区域有限，因此我们可以利用像素值大小来寻找，我们将搜索区域限制在第五层下方五个像素和第七层上方两个像素之间，此举是为了排除第五层和第七层自身宽度和其它层的影响，结果如图3.8a所示。但是所找到的像素值并不在一条光滑的曲线上的，我们可以利用MATLAB中的拟合函数，对上面找到的像素值进行一个曲线拟合处理，如下图3.8b所示，如此便确定了第六层的位置。

(b)图3.8拟合前后的第六层（a）拟合前（b）拟合后（4）结果显示我们已经找出来了视网膜OCT图像的七层，整体结果展示如图3.9所示。其中涉及到的各层之间的错误修正及优化，将会在3.3节讲述。图3.9OCT图像自动分层结果3.3错误校正及结果优化3.3.1错误校正OCT自动分层算法检测过程中，由于图像本身的特性，会产生两种错误，一是跳变误差，包括边界缺损和中间缺损两类，二是中心凹，所以需要修改自动算法方法。（1）处理跳变误差跳变误差又可分为两类。一类是边界缺损，主要是原图边界不清晰造成的，如下图3.10所示。处理此类跳变的会选择水平延伸的方法，我们把第6个像素的值直接赋值给前五个像素，把倒数第六个像素的值直接赋值给倒数五个像素，因为视网膜的层边界几乎是水平的，所以位置靠前或者靠后的像素点的高度差也较小。图3.10边界缺损下图3.11a和3.11b体现了修正前和修正后的差异。(a)(b)图3.11边界缺损修补前后（a）修补前（b）修补后另一类是中间缺损，主要原因是血管的存在或边界的不清晰。具体可见下图2.12。一般情况下通常采用线性插值来处理此类缺损，我们在已获取的边界上从左往右检测，设置检测的宽度为2，当检测到前一个像素值和后一个像素值差值大于10，就将前一个点标记为A,同理，再次检测到差值大于10的情况，就将后一个点标记为B。然后连接A、B两点即可。这是因为缺少的像素点数目较少，所以两个端点之间的高度差也较小。图3.12中间缺损下图3.13a和3.13b体现了修正前和修正后的差异。（a）（b）图3.13中间缺损修补前后（a）修补前（b）修补后（2）处理中心凹寻找c2时，因为c2层次非常薄，难以检测，常常会产生中央凹，其边界不是水平的，此时若根据已知边界延伸方向，就会发生错误。因为中央凹一般出现在30--36帧的第二层，所以我们把这部分单独设置，在已经检测到的第一层边界的最低点附近，如果检测到第二层是断掉的，那就把第一层该处的像素值再加一，作为第二层此处的像素值，然后再用均值滤波来处理第二层边界，使其连接处更平滑。问题图如图3.14，校正结果如图3.15a和3.15b。图3.14中心凹（b）图3.15中心凹修补前后（a）修补前（b）修补后3.3.2结果优化如果仅仅是简单地处理跳变误差和中心凹，原有的自动分层算法直接检测的曲线是不够平滑的，例如较小的斜跳变误差，或者有些区域凹凸不平，因此，在自动检测完成每一层后，都会对已获取的边界进行一个中值滤波的处理。对中间的分界不是不太明显的二、三和四层也会再有一个高斯平滑处理，尽管肉眼观察区别不是很大，但是在4.2节误差分析中，可看出误差有所减小[14]。下图3.16（a）--（l）将给出在此自动分层算法下的6张不同的人眼图像在不同帧处的分层效果图。（a）（b）（c）（d）（e）（f）（g）（h）（i）（j）（k）（l）图3.16六张不同眼图不同帧数下OCT自动分层结果（a）（b）第12帧的图1（c）（d）第18帧的图2（e）（f）第24帧的图3（g）（h）第30帧的图4（i）（j）第36帧的图5（k）（l）第48帧的图6第四章实验结果与分析4.1误差评价标准本次实验用了6幅不同人眼的视网膜OCT图像，通过所设计的自动分层算法的结果与手动分割的结果进行比较，得到误差。因为图像的分辨率是480×512×64，计算误差的公式定义为：：E=1512x=1其中，ya是手动分割的数据，ym每张图片有64帧，我们每隔6帧选取一幅图像，计算每层的误差，再将选取的这10帧图像的误差求平均值，得到每张图像各层误差。4.2实验结果分析首先，我们对一张人眼视网膜OCT图像不同帧时进行误差分析，将自动分层结果与手动分割结果相对比，误差如下：表4.1图4不同帧处平均误差帧数C1C2C3C4C5C6C7Y=61.50201.15721.18851.08010.62301.36571.0234Y=121.53911.25201.53611.04100.85940.71910.7422Y=181.35941.05571.37501.84960.68551.20870.7520Y=241.14651.83201.73731.41210.72851.41440.8984Y=301.03132.95121.72070.94730.70902.16700.6230Y=361.05082.68751.71291.31451.00002.36830.6855Y=421.39262.34861.34381.57320.79691.15940.6328Y=480.72851.59860.98051.33790.54300.83430.5449Y=560.97661.43161.21091.28470.99220.97200.9102Y=600.88482.61911.87601.04440.84181.21290.9160误差1.16121.89341.46821.28850.77791.34220.7728图4.1一张图不同帧平均误差表4.1和图4.1是计算出同一张图的不同帧处的平均误差，可看出平均相对误差都比较小，但是不同帧之间误差分布略有差异。然后，我们针对6张不同人眼的视网膜OCT图像求各层的平均误差，将自动分层结果与手动分割结果相对比，误差如下：表4.2图1~6不同层平均误差C1C2C3C4C5C6C7图11.16121.89341.46821.28850.77791.34220.7728图21.24242.38252.02151.74990.89141.33751.0576图31.04652.14491.54731.59650.69571.56751.4531图41.03072.10761.56231.54840.73571.53070.9969图51.16972.43321.77281.93290.71522.39270.9441图61.40192.54212.03971.75571.41882.71881.4604平均1.17542.25061.73531.64530.87251.81491.1142图4.2六张图不同层平均误差表4.2和图4.2是6张不同人眼视网膜OCT图像计算出的平均误差，可看出平均相对误差都比较小，但是有的不同层之间误差分布略有差异。根据以上两种不同维度的误差分析，再分析图4-1和4-2，不难看出，通过自动分层算法找到误差比较大是c2和c6，因为存在中央凹，而c2本身就比较薄，很难找到合理的参数使Canny算子准确找到c2，因此误差较大，c6是通过c5和c7之间像素值最小确定的，然后再曲线拟合，与手动分割相比也存在较大误差。4.3总结尽管OCT视网膜自动分层算法的研究现状已经较为广泛，但是本次设计花了大量的时间和精力却也只找到了11层边界中的7层，精确度也并没有达到理想效果，现总结了以下两点原因：1．Canny算子本身是存在局限性的，对噪声比较敏感，孤立边缘点和伪边缘会被检测出来，造成干扰[15。2．后期的边界修正需要考虑大量的情