CN109003284A - 基于层厚统计信息模型的相干光断层扫描图像的层分割方法 - Google Patents
基于层厚统计信息模型的相干光断层扫描图像的层分割方法 Download PDFInfo
- Publication number
- CN109003284A CN109003284A CN201810758417.2A CN201810758417A CN109003284A CN 109003284 A CN109003284 A CN 109003284A CN 201810758417 A CN201810758417 A CN 201810758417A CN 109003284 A CN109003284 A CN 109003284A
- Authority
- CN
- China
- Prior art keywords
- image
- layer
- statistical information
- region
- gradient
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Withdrawn
Links
- 238000000034 method Methods 0.000 title claims abstract description 44
- 230000001427 coherent effect Effects 0.000 title claims abstract description 14
- 210000004204 blood vessel Anatomy 0.000 claims abstract description 25
- 230000011218 segmentation Effects 0.000 claims abstract description 24
- 238000001514 detection method Methods 0.000 claims abstract description 22
- 230000009467 reduction Effects 0.000 claims abstract description 10
- 239000011159 matrix material Substances 0.000 claims abstract description 9
- 238000012545 processing Methods 0.000 claims abstract description 9
- 230000002708 enhancing effect Effects 0.000 claims abstract description 6
- 230000002146 bilateral effect Effects 0.000 claims description 8
- 238000005516 engineering process Methods 0.000 claims description 8
- 238000001914 filtration Methods 0.000 claims description 8
- 238000010606 normalization Methods 0.000 claims description 4
- 238000011946 reduction process Methods 0.000 claims description 4
- 230000009466 transformation Effects 0.000 claims description 4
- 210000001525 retina Anatomy 0.000 description 28
- 230000004256 retinal image Effects 0.000 description 15
- 238000012014 optical coherence tomography Methods 0.000 description 13
- 238000010586 diagram Methods 0.000 description 7
- 238000012360 testing method Methods 0.000 description 6
- 238000003384 imaging method Methods 0.000 description 5
- 201000010099 disease Diseases 0.000 description 4
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 4
- 230000008569 process Effects 0.000 description 4
- 230000008859 change Effects 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 239000012528 membrane Substances 0.000 description 3
- 230000002207 retinal effect Effects 0.000 description 3
- 230000003595 spectral effect Effects 0.000 description 3
- 210000001519 tissue Anatomy 0.000 description 3
- 208000003098 Ganglion Cysts Diseases 0.000 description 2
- 208000005400 Synovial Cyst Diseases 0.000 description 2
- 239000000654 additive Substances 0.000 description 2
- 230000000996 additive effect Effects 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 238000004587 chromatography analysis Methods 0.000 description 2
- 208000030533 eye disease Diseases 0.000 description 2
- 238000009499 grossing Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 210000003583 retinal pigment epithelium Anatomy 0.000 description 2
- 238000003325 tomography Methods 0.000 description 2
- 230000000007 visual effect Effects 0.000 description 2
- 241000283984 Rodentia Species 0.000 description 1
- 241000396836 Trinodes Species 0.000 description 1
- 239000002671 adjuvant Substances 0.000 description 1
- 230000002411 adverse Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 239000008280 blood Substances 0.000 description 1
- 210000004369 blood Anatomy 0.000 description 1
- 210000005252 bulbus oculi Anatomy 0.000 description 1
- 210000004027 cell Anatomy 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 210000004087 cornea Anatomy 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000018109 developmental process Effects 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 239000000835 fiber Substances 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 238000003706 image smoothing Methods 0.000 description 1
- 238000001727 in vivo Methods 0.000 description 1
- 238000012905 input function Methods 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 230000003902 lesion Effects 0.000 description 1
- 239000002052 molecular layer Substances 0.000 description 1
- 230000000877 morphologic effect Effects 0.000 description 1
- 210000005036 nerve Anatomy 0.000 description 1
- 210000004126 nerve fiber Anatomy 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000004796 pathophysiological change Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000002310 reflectometry Methods 0.000 description 1
- 239000000790 retinal pigment Substances 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000000638 solvent extraction Methods 0.000 description 1
- 238000013179 statistical model Methods 0.000 description 1
- 238000012916 structural analysis Methods 0.000 description 1
- 210000003684 theca cell Anatomy 0.000 description 1
- 230000002792 vascular Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/12—Edge-based segmentation
-
- G06T5/70—
-
- G06T5/92—
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10101—Optical tomography; Optical coherence tomography [OCT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20024—Filtering details
- G06T2207/20028—Bilateral filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30041—Eye; Retina; Ophthalmic
Abstract
本发明公开了一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法，包括步骤：S1、预处理：对图像进行降噪、对比度增强和扁平处理；S2、血管检测：将图像中的血管阴影检测出来；S3、构建邻接矩阵：生成垂直方向的梯度图像，对梯度图像进行归一化处理，并计算邻接矩阵对应的图中各节点间的权重；S4、利用层厚统计信息设置搜索区域：根据层厚统计信息在垂直方向上依次为各层限定搜索的区域；S5、边界检测：先将对比度最高的边界检测出来，根据已检测出的边界迭代地限定每次搜索的区域，再检测剩余边界并进行层分割。本方法具有较好的鲁棒性、更高的分割效率和灵活性。
Description
技术领域
本发明涉及相干光断层图像处理技术领域，特别是涉及一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法。
背景技术
准确、客观、高效的医学图像目标检测技术在临床医学眼科结构分析中具有重大的意义。光学相干断层扫描成像技术(Optical Coherence Tomography,OCT)是一种非常有前途的层析成像技术，特别是在活体检测和生物组织成像领域具有十分诱人的应用前景。作为一种新的层析成像技术，由于其非接触性和高分辨率的特点，它被证实为一种十分有效、适应性较广的成像技术手段。近年来，OCT技术在眼部疾病的诊断和治疗上得到了飞速的发展，为临床眼部疾病的诊断和治疗做出了突出的贡献。时域OCT(TD-OCT)是第一代OCT技术，每秒可获得400张A扫描图，轴向分辨率范围为8到10微米。当前商用的OCT技术是谱域OCT(SD-OCT)，成像速度比时域OCT快10到100倍，轴向分辨率范围为3到7微米，可以获取视网膜和角膜的切面信息。由于成像速度快，分辨率高，谱域OCT广泛应用于眼科疾病的诊断和预测。
视网膜中层的厚度变化是临床疾病状态的一种表现，利用视网膜各层的厚度变化信息来辅助临床医学上对眼科疾病进行诊断，在此基础上，确定眼部视网膜各层组织结构显得尤为重要。然而，常规的手工层分割方法是比较耗时的，并且不同的专家分割出的结果会具有一定的差异性。近年来，一些相关研究人员已经开始研究视网膜层结构的自动分割问题，这些方法主要包括：基于梯度信息的分割方法，活动轮廓模型，基于图论的分割方法，形状与上下文先验信息，基于像素分类的分割方法。
Yang et al.等人利用局部和全局双梯度信息，提出了一种对OCT图像中黄斑区域进行层结构分割的算法。Yazdanpanah利用轮廓信息实现了啮齿类动物图像的视网膜层分割算法。Chiu等人利用图像中的双梯度信息构造出一个图，并通过最短路径算法对图像进行逐层分割。为了实现对脉络膜的自动分割，Vedran Kaji提出了一种基于统计模型的分割方法。基于形态学特征和上下文先验知识的分割方法被Qi Song等人提出，以实现最佳的多表面分割。Vermeer KA将视网膜的分割问题转化为像素分类问题，提出了一种基于像素分类的分割方法。在上述提及的分割方法中，Yang et al.等人、Vedran Kaji、Qi Song等人以及Vermeer KA提出的这四种3D的分割方法显示出了更好的实用性，能更好地反映出视网膜的结构特征，但它同时也有着相对较高的时间复杂度，并且对系统的硬件有着更高的性能要求。
视网膜数据集中，每一帧SD-OCT图像也称为B-scan图像，它的大小均为1000×512像素。由于基于相干光断层扫描成像的视网膜的分辨率较高，图像中的内膜层结构清晰可见。视网膜在人体眼球的正后方，它由多层的纤细透明的层状结构组成。一般地，它被认为是介于视网膜神经纤维层(retinal nerve fiber layer,RNFL)和视网膜色素上皮细胞(retinal pigment epithelium,RPE)之间的水平结构。内界膜(inner limitingmembrane,ILM)是视网膜图像中最顶端的白色区域和黑色区域的边界。由于视网膜色素上皮细胞有着最高的反射率，RPE被认为是视网膜图像中最亮的一条边界。位于ILM和RPE之间的层结构分别为：神经节细胞层(ganglion cell layer,GCL)，内丛状层(inner plexiformlayer,IPL)，内核层(inner nuclear layer，INL),外丛状层(outer plexiform layer,OPL),外核层(outer nuclear layer,ONL)，内感光层(inner photosensitive layer,IS)和外感光层(outer photosensitive layer,OS)，如图3所示。
发明内容
本发明为了克服以上技术的不足，提供了一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法。
本发明中的实验数据集采用了杜克大学公开的相干光断层扫描成像的视网膜数据集。
本发明克服其技术问题所采用的技术方案是：
一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法，包括如下步骤：
S1、预处理：
对图像进行降噪、对比度增强和扁平处理；
S2、血管检测：
将图像中的血管阴影检测出来；
S3、构建邻接矩阵：
生成垂直方向的梯度图像，对梯度图像进行归一化处理，并计算邻接矩阵对应的图中各节点间的权重；
S4、利用层厚统计信息设置搜索区域：
根据层厚统计信息在垂直方向上依次为各层限定搜索的区域；
S5、边界检测：
先将对比度最高的边界检测出来，根据已检测出的边界迭代地限定每次搜索的区域，再检测剩余边界并进行层分割。
本发明优选的，所述步骤S1中，采用双边滤波对图像进行降噪处理，采用对数变换和对比度拉伸方法对降噪后的图像进行对比度增强处理。
本发明优选的，所述步骤S2中，血管检测的步骤包括：
S2.1、计算图像中某一区域内每一列的平均灰度，得到图像对应的一个平均灰度列向量；
S2.2、采用多项式平滑技术对平均灰度列向量进行处理，拟合得到一条适用于整个样本点的曲线P；
S2.3、通过样本点和拟合曲线P的误差判断出血管区域。
本发明优选的，所述步骤S3中，构建邻接矩阵包括步骤：
S3.1、生成垂直方向上的两个梯度图像，一个是由明到暗、另一个是由暗到明，并计算图像中每一个点的梯度；
S3.2、对梯度图像进行归一化处理，将梯度图像中各点的像素控制在(0，1)范围内；
S3.3、通过归一化之后的梯度图像，结合梯度、标准差、血管阴影的位置以及图像的层厚统计信息对搜索图的权值进行初始化。
本发明优选的，所述步骤S4中，为了对层厚统计信息约束项进行建模，图像中某点a和中央凹在水平方向上的距离确定a所在的区域，该约束项表示为：
其中，bx为b点沿x轴(垂直方向)对应的坐标，μ和σ分别对应a点对应层所在区域的数学期望和标准差，k为相关系数因子。
本发明优选的，所述步骤S4中，还定义一个向量来辅助限制Dijkstra算法每次遍历时的搜索区域，该向量的定义为：
其中，代表从源点中的节点指向当前节点的向量，为先前的节点指向当前节点的向量。
本发明优选的，所述步骤S5中，结合层厚统计信息和最短路径算法迭代地对剩余各层进行分割。
本发明的有益效果是：
本发明针对相干光断层扫描图像的分割问题，提出了一种新的分割方法，以Dijkstra算法和动态规划算法为基础，在谱域相干光断层扫描(SD-OCT)的视网膜图像构建搜索图，充分利用各层的层厚统计信息来限制搜索区域，并且迭代地采用最短路径算法对视网膜各层进行分割，实现了视网膜相干光断层扫描图像多层结构的全自动分割，其有益效果表现如下：
(1)鲁棒性：本方法利用平均梯度信息确定图中各节点之间的权重，能避免噪声对分割结果的影响；同时利用血管的位置信息限定区域的权重大小，能够避免血管对分割结果的影响，本方法具有较好的鲁棒性。
(2)高效性：本方法结合层厚度统计信息和最短路径方法，初始化节点时仅考虑5条边，具有更高的分割效率。
(3)灵活性:本方法也适用于由轻微的病理生理变化引起的视网膜层状结构的层分割，具有较高的灵活性和适应性。
附图说明
图1为本发明实施例的流程示意图。
图2为视网膜眼底图像。
图3为相干光断层扫描图像视网膜的层结构示意图。
图4为本发明实施例的观测图像。
图5为本发明实施例降噪后的视网膜图像。
图6为本发明实施例对比度增强后的视网膜图像。
图7为本发明实施例扁平处理后的视网膜图像。
图8为本发明实施例血管检测结果示意图。
图9为本发明实施例归一化后的梯度图(由暗到明)。
图10为本发明实施例归一化后的梯度图(由明到暗)。
图11为本发明实施例ILM边界检测结果示意图。
图12为本发明实施例IS-OS的检测结果示意图。
图13为本发明实施例所有边界的检测结果示意图。
具体实施方式
为了便于本领域人员更好的理解本发明，下面结合附图和具体实施例对本发明做进一步详细说明，下述仅是示例性的不限定本发明的保护范围。
本发明公开了一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法，本实施例以视网膜层为例并结合附图1-13说明其具体实施过程。
步骤一、预处理
a)图像平滑
SD-OCT图像不可避免地受到噪声的污染(一般地，该噪声主要为乘性噪声)，会在一定程度上降低分割算法的准确率。为了减少图像中的噪声并提高图像的视觉效果，分别采用降噪和对比度增强技术对图像进行预处理。针对SD-OCT图像中的斑点噪声，观测图像可以表示为：
G(x,y)＝F(x,y)*N(x,y) (1)
其中，G(x,y)为观测图像(指二维的B-Scan视网膜图像)，F(x,y)为视网膜图像对应的原始图像，N(x,y)为噪声成分。通过对数变换，公式(1)可以表示为：
g(x,y)＝f(x,y)+n(x,y) (2)
其中，g(x,y)＝logG(x,y),f(x,y)＝logF(x,y)，n(x,y)＝logN(x,y)，即噪声n(x,y)可以表示为加性噪声进行处理。
高斯平滑和双边滤波被广泛地应用于加性噪声的处理过程。高斯平滑能较好的抑制全局噪声对图像分割造成的不利影响，然而，这会倾向于丢失掉一些相对重要的局部信息，在SD-OCT视网膜图像的分割问题中，可能会使得视网膜边界的检测出现偏差。利用双边滤波对图像进行降噪处理的时候，能较好的处理局部的斑点噪声，并且较好地保留局部细节；在利用双边滤波进行降噪时，不需要通过整幅图像去计算当前像素点的权重，因为当图像中某个点和当前像素点的距离大于某个阈值时，该点对当前像素点的影响非常小而被忽略，这对SD-OCT图像的分割是比较有利的。通过邻域离散化的双边滤波可以描述为：
其中，k(x)是归一化的输入函数，c是基于空间距离的高斯权重值，s为基于灰度相似度的高斯权重值。
利用双边滤波对图像进行了降噪处理，通过双边滤波，原始图像F(x,y)可以由(x,y)经过指数变化得到。图4和图5分别为观测图像和对应的降噪后图像。
b)图像增强与扁平处理
在图像成像、传输和转换的过程中，受某些客观因素的影响，图像的质量会受到严重的影响。采用对比度增强技术能较好地增强图像中的局部信息并增大图像的特征差异，能将图像以更好的视觉形式呈现出来。利用对数变换和对比度拉伸的方法对降噪后的SD-OCT图像进行了图像增强处理。
在分析OCT视网膜图像时，由于相关临床医生习惯对扁平后的医学图像进行分析，医学图像往往需要进行扁平处理。因为图像中高曲率、高梯度的边界可能会在最短路径搜索时出现偏差，在预处理部分也对SD-OCT图像进行了扁平处理。RPE边界是SD-OCT图像中相对较为平坦的结构，同时，它还以较高对比度的形式呈现在图像中的最下方。所以，将RPE作为图像扁平的一个参考标准是相对比较可取的。利用最短路径方法对RPE进行了检测，RPE边界中垂直方向上的最大值Xmax作为参考基准(图像坐标中原点在左上角)，在RPE边界上方的每一列分别向下平移Xmax-Xi单位个像素得到最终扁平后的图像。图6和图7分别显示了对比度增强的图像和扁平处理后的图像。
步骤二、血管检测
视网膜中的血管在SD-OCT图像中呈现为一种模糊不清的阴影，会使得图像中的边缘变得模糊。这会对视网膜各层边界的检测造成不确定的影响。将血管中的阴影先检测出来，再对SD-OCT图像的各层进行分割，会使得检测出来的结果更加的精确。
视网膜中的血管在SD-OCT图像中呈现为一种垂直方向上的阴影，所以可以沿水平方向对其进行检测。图像中的血管阴影区域的灰度相对较低，而其他区域的灰度相对较高。因此，如果将图像中每一列的灰度进行累加，再求得其平均值，则非血管区域的灰度要远大于血管区域的灰度值。
假设ILM和RPE之间的区域作为研究区域(Region of Interest，ROI)，计算ROI区域内每一列的平均灰度，得到图像对应的一个平均灰度列向量，ROI区域垂直方向上的平均灰度的计算方式为：
其中，和分别表示x在垂直方向上对应的RPE和ILM的边界。为了实现对血管区域的精确定位，考虑到血管阴影在不同区域的宽度不一致，利用多项式平滑技术对该平均灰度列向量进行处理。首先，沿水平方向按相同间距取出适量样本点进行多项式平滑，然后取邻域内样本点的平均值作计算取样样本的值。多项式平滑按以下公式进行计算：
yi＝fmedian(I(xi-ks),......,I(xi+ks)),i＝1,2,......,n (5)
其中，xi可以表示为：
xi＝1+ksi (6)
I(xi)为图像第xi列的平均灰度值，计算邻域内各列平均灰度的均值，ks代表水平方向上取样的步长参数，本实施例中取ks＝1。通过多项式平滑可以拟合适用于整个样本点一条曲线P，该拟合得到的曲线的误差可以通过相关经验的阈值T(T＝0.1)进行评估。由于血管区域的灰度相对较低，其对应的列的平均灰度会相对较小。即样本点和拟合曲线P的误差大于该阈值T时，判断该区域为血管区域。通过上述方法，可将SD-OCT图像中的血管区域检测出来。图8为检测结果，其中黑色区域标记的为血管区域。
步骤三、构建邻接矩阵
a)灰度图像归一化
由于SD-OCT视网膜图像的边界整体上呈现为水平的结构，为了简化相关工作，这里仅获取了图像沿垂直方向上的梯度图像。在图像的两侧分别加入了一列，并赋值为0，左侧和右侧的一列分别代表最短路径中的源点(S)和汇点(T)。这里需要取两个梯度图像，一个由明到暗，用于检测RPE，NFL-GCL，OPL-ONL和IPL-INL；另一个由暗到明，用于检测INL-OPL，ILM，IS-OS和ONL-IS。对于图像中的每一个点，它的梯度为相关邻域范围内垂直方向上的梯度平均值，即：
其中，gradient(x,y)为梯度图像中的点，(x,y)为对应的坐标值，N为邻域内取的像素对的数目，本实施例中取N＝3，F(x,y)为SD-OCT图像在(x,y)处的灰度值。相对于利用单个像素对计算梯度，该方法对噪声更加鲁棒。
考虑到梯度图像中某些点的值可能会负值，有必要先对梯度图像进行归一化处理，将梯度图像中各点的像素控制在(0，1)范围内。这里主要按下列公式对图像进行了归一化处理：
其中，gradient(x,y)为梯度图像中点(x,y)处的值，min(gradient)和max(gradient)分别表示梯度图像像素中的最小值和最大值。图9和图10显示了归一化的梯度图像示意图，它们主要用于对图进行初始化。
b)初始化赋值
得到归一化的梯度图像后，便可以用它们来对搜索图的权值进行初始化。图中的两个节点a和b之间的权重和它们之间的梯度有关，a和b节点的梯度又分别受它们的邻域的影响，所以结合它们的邻域信息来为节点之间的边赋权，对最终的分割是有利的。同理，图像中是否存在着血管阴影，对分割也会造成不同的影响，同时考虑血管阴影为图中节点之间的边赋权。
值得一提的是，由于SD-OCT视网膜图像层结构的特性，在利用Dijkstra算法对图进行遍历时，只需要从左往右进行遍历便可实现分割。也就是说，在为图像中某点和它相邻的点之间的边赋权值时，只需要考虑它右侧的点(仅对右上方、正右方和右下方三个相邻点进行赋权)。它们的权重可以表示为：
Wab＝(2-Gab)+λSab+Bs+T(a,b)+Wmin (9)
其中，Wab代表图中的两个节点a和b之间的权重，Gab为梯度影响因子，它表示亮点梯度的和，即：
Gab＝Ga+Gb (10)
其中，Ga和Gb表示a和b两个节点的梯度，可由公式(7)(8)计算得出。λ为标准差的影响因子系数，取λ＝10。Sab为标准差的影响因子，因为相同层的灰度值会比较接近，它们的标准差的值会较小。Sab的计算公式为：
Sab＝||Sa-Sb|| (11)
其中，Sa和Sb分别代表a和b两个节点领域内的标准差。由于不是所有的图像中都存在血管阴影，血管影响因子Bs可以表示为：
其中，|ax-bx|为节点a和节点b在垂直方向上的距离，k为常量(取k＝0.01)，k|ax-bx|代表血管阴影影响因子，f(s)＝0表示节点a和节点b不在血管区域，f(s)＝1表示节点a或节点b在血管区域，f(s)可由血管检测部分得出，T(a,b)为层厚统计信息约束项，Wmin为设置的误差阈值，可取Wmin＝10-5。
步骤四、利用层厚统计信息设置搜索区域
视网膜的层厚信息能反映出眼部组织的结构特征。近年来，视网膜层厚统计研究领域取得了较大的进展，比较有代表性的有Kaori Sayanag，他将有病变的视网膜图像和正常人的视网膜图像进行对比，分别统计两种情况下视网膜各层的厚度信息。除此之外，还有Raheleh Kafieh，他精确地测量了人体视网膜各层的厚度，如黄斑区附近的NFL厚度约为9±8μm，GCL的厚度约为9±5μm，这项研究为相关眼科疾病的诊断和治疗提供了一个辅助判断的基准。本实施例利用Raheleh Kafieh提供的层厚统计信息辅助视网膜各层的分割。
视网膜各层的层厚信息在不同的数据集上可能略有不同，但它们的厚度的分布情况仍然接近于高斯分布。因此，利用Raheleh Kafieh提供的层厚统计信息进行建模，以提高视网膜各层分割的精度。在建模过程中，首先需要确定黄斑区中央凹的位置，利用先前检测出的ILM和RPE来确认中央凹的位置，ILM上距离RPE最小的点的位置被认为是中央凹的位置；然后以中央凹为基准，在两侧分别沿水平方法将视网膜图像分成三个部分(参考Raheleh Kafieh中的数据)，以各区域中各层的层厚统计信息对搜索区域进行限制。
为了对层厚统计信息约束项进行建模，图像中某点a和中央凹在水平方向上的距离确定a所在的区域，然后该约束项可以表示为:
其中，bx为b点沿x轴(垂直方向)对应的坐标，μ和σ分别对应a点对应视网膜层所在区域的数学期望和标准差，k为相关系数因子，可取k＝3。
除此之外，还定义了一个向量来辅助限制Dijkstra算法每次遍历时的搜索区域，这里主要的依据是假设该向量沿垂直方向的投影的模长不能超过某个阈值(该阈值由层厚信息确定)。当该向量沿垂直方向的投影的模长小于该阈值时，它可以访问其右侧相邻的部分节点。如：假设当前节点坐标为(x,y)，它可以访问a(x+1,y-2),b(x+1,y-2),c(x+1,y),d(x+1,y+1),e(x+1,y+2)5个节点(理论上可以访问右侧的全部节点，在实践中设定为5个)；当该向量沿垂直方向上的投影的模长超过该阈值时，它可能只能访问a、b、c三个节点，此时因为认为访问d和e和层厚统计信息不符合而放弃访问该节点。由此一来，视网膜各层的分割精度能得到一定程度上的提高。该向量的定义为：
其中，代表从源点中的节点指向当前节点的向量，为先前的节点指向当前节点的向量。
步骤五、边界检测
由于ILM和IS-OS边界的对比度最高，它们最先被检测出来。图像中左上角的点被设定为源点A，右下角的点为设定为汇点B。然后迭代地利用Dijkstra算法进行视网膜各层的分割。图11显示了ILM在原图上的检测结果。将ILM边界检测出来之后，搜索区域的上边界设定为ILM的边界，下边界仍然以图像的最下方为准。然后再次利用Dijkstra进行搜索，确定IS-OS边界。图12显示了IS-OS的检测结果。
ILM和IS-OS检测出来以后，结合层厚统计信息和最短路径算法迭代地对其余各层进行分割。所有边界的检测顺序为：ILM，IS-OS，RPE，INL-OPL，NFL-GCL，IPL-INL，OPL-ONL,ONL-IS。在每次检测出新的边界后，动态地限定每次搜索的区域。如当IS-OS和OPL-ONL检测出来后，将ONL和IS层作为搜索区域，在此区域中检测ONL-IS边界。图13显示了所有边的检测结果图。
以上仅描述了本发明的基本原理和优选实施方式，本领域人员可以根据上述描述做出许多变化和改进，这些变化和改进应该属于本发明的保护范围。
Claims (7)
1.一种基于层厚统计信息模型的相干光断层扫描图像的层分割方法，其特征在于，包括如下步骤：
S1、预处理：
对图像进行降噪、对比度增强和扁平处理；
S2、血管检测：
将图像中的血管阴影检测出来；
S3、构建邻接矩阵：
生成垂直方向的梯度图像，对梯度图像进行归一化处理，并计算邻接矩阵对应的图中各节点间的权重；
S4、利用层厚统计信息设置搜索区域：
根据层厚统计信息在垂直方向上依次为各层限定搜索的区域；
S5、边界检测：
先将对比度最高的边界检测出来，根据已检测出的边界迭代地限定每次搜索的区域，再检测剩余边界并进行层分割。
2.根据权利要求1所述的方法，其特征在于，所述步骤S1中，采用双边滤波对图像进行降噪处理，采用对数变换和对比度拉伸方法对降噪后的图像进行对比度增强处理。
3.根据权利要求1所述的方法，其特征在于，所述步骤S2中，血管检测的步骤包括：
S2.1、计算图像中某一区域内每一列的平均灰度，得到图像对应的一个平均灰度列向量；
S2.2、采用多项式平滑技术对平均灰度列向量进行处理，拟合得到一条适用于整个样本点的曲线P；
S2.3、通过样本点和拟合曲线P的误差判断出血管区域。
4.根据权利要求1所述的方法，其特征在于，所述步骤S3中，构建邻接矩阵包括步骤：
S3.1、生成垂直方向上的两个梯度图像，一个是由明到暗、另一个是由暗到明，并计算图像中每一个点的梯度；
S3.2、对梯度图像进行归一化处理，将梯度图像中各点的像素控制在(0，1)范围内；
S3.3、通过归一化之后的梯度图像，结合梯度、标准差、血管阴影的位置以及图像的层厚统计信息对搜索图的权值进行初始化。
5.根据权利要求1所述的方法，其特征在于，所述步骤S4中，为了对层厚统计信息约束项进行建模，图像中某点a和中央凹在水平方向上的距离确定a所在的区域，该约束项表示为：
其中，bx为b点沿x轴(垂直方向)对应的坐标，μ和σ分别对应a点对应层所在区域的数学期望和标准差，k为相关系数因子。
6.根据权利要求5所述的方法，其特征在于，所述步骤S4中，还定义一个向量来辅助限制Dijkstra算法每次遍历时的搜索区域，该向量的定义为：
其中，代表从源点中的节点指向当前节点的向量，为先前的节点指向当前节点的向量。
7.根据权利要求1所述的方法，其特征在于，所述步骤S5中，结合层厚统计信息和最短路径算法迭代地对剩余各层进行分割。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810758417.2A CN109003284A (zh) | 2018-07-11 | 2018-07-11 | 基于层厚统计信息模型的相干光断层扫描图像的层分割方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810758417.2A CN109003284A (zh) | 2018-07-11 | 2018-07-11 | 基于层厚统计信息模型的相干光断层扫描图像的层分割方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN109003284A true CN109003284A (zh) | 2018-12-14 |
Family
ID=64599741
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810758417.2A Withdrawn CN109003284A (zh) | 2018-07-11 | 2018-07-11 | 基于层厚统计信息模型的相干光断层扫描图像的层分割方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109003284A (zh) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112508859A (zh) * | 2020-11-19 | 2021-03-16 | 聚融医疗科技(杭州)有限公司 | 一种基于小波变换的子宫内膜厚度自动测量方法及系统 |
| WO2022007352A1 (zh) * | 2020-07-10 | 2022-01-13 | 温州医科大学 | 一种基于光学相干断层扫描系统的脉络膜三维血管成像及定量化分析方法与装置 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110182517A1 (en) * | 2010-01-20 | 2011-07-28 | Duke University | Segmentation and identification of layered structures in images |
| CN103514605A (zh) * | 2013-10-11 | 2014-01-15 | 南京理工大学 | 基于hd-oct视网膜图像的脉络膜层自动分割方法 |
| CN106097378A (zh) * | 2016-07-24 | 2016-11-09 | 江西理工大学 | 一种融合形状先验的水平集视网膜血管图像分割方法 |
-
2018
- 2018-07-11 CN CN201810758417.2A patent/CN109003284A/zh not_active Withdrawn
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110182517A1 (en) * | 2010-01-20 | 2011-07-28 | Duke University | Segmentation and identification of layered structures in images |
| CN103514605A (zh) * | 2013-10-11 | 2014-01-15 | 南京理工大学 | 基于hd-oct视网膜图像的脉络膜层自动分割方法 |
| CN106097378A (zh) * | 2016-07-24 | 2016-11-09 | 江西理工大学 | 一种融合形状先验的水平集视网膜血管图像分割方法 |
Non-Patent Citations (1)
| Title |
|---|
| 王佳: "相干光断层成像OCT视网膜图像处理方法研究", 《万方学位论文》 * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2022007352A1 (zh) * | 2020-07-10 | 2022-01-13 | 温州医科大学 | 一种基于光学相干断层扫描系统的脉络膜三维血管成像及定量化分析方法与装置 |
| CN112508859A (zh) * | 2020-11-19 | 2021-03-16 | 聚融医疗科技(杭州)有限公司 | 一种基于小波变换的子宫内膜厚度自动测量方法及系统 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Yin et al. | User-guided segmentation for volumetric retinal optical coherence tomography images | |
| Kafieh et al. | A review of algorithms for segmentation of optical coherence tomography from retina | |
| Ghorbel et al. | Automated segmentation of macular layers in OCT images and quantitative evaluation of performances | |
| Lee et al. | Segmentation of the optic disc in 3-D OCT scans of the optic nerve head | |
| US20210319556A1 (en) | Method and system for analysing images of a retina | |
| Dufour et al. | Graph-based multi-surface segmentation of OCT data using trained hard and soft constraints | |
| CN108961261B (zh) | 一种基于空间连续性约束的视盘区域oct图像层次分割方法 | |
| Novosel et al. | Loosely coupled level sets for simultaneous 3D retinal layer segmentation in optical coherence tomography | |
| AU2021202217B2 (en) | Methods and systems for ocular imaging, diagnosis and prognosis | |
| JP6842481B2 (ja) | 深層学習を用いた網膜層の３ｄ定量解析 | |
| WO2015028559A1 (en) | Evaluation of optical coherence tomographic data prior to segmentation | |
| CN113574542A (zh) | 宽场自发荧光图像中患有年龄相关性黄斑变性的患者中地理萎缩模式分割与分类 | |
| EP2437649A1 (en) | Image processing apparatus, control method thereof, and computer program | |
| Cha et al. | High-accuracy retinal layer segmentation for optical coherence tomography using tracking kernels based on Gaussian mixture model | |
| CN109325955B (zh) | 一种基于oct图像的视网膜分层方法 | |
| Rabbani et al. | Obtaining thickness maps of corneal layers using the optimal algorithm for intracorneal layer segmentation | |
| Monemian et al. | Analysis of a novel segmentation algorithm for optical coherence tomography images based on pixels intensity correlations | |
| CN109744996B (zh) | Oct图像的bmo位置定位方法 | |
| CN109003284A (zh) | 基于层厚统计信息模型的相干光断层扫描图像的层分割方法 | |
| Zhang et al. | Fast retinal layer segmentation of spectral domain optical coherence tomography images | |
| US11302006B2 (en) | 3D quantitative analysis with deep learning | |
| Alsaih et al. | Retinal fluids segmentation using volumetric deep neural networks on optical coherence tomography scans | |
| Eichel et al. | Automated 3D reconstruction and segmentation from optical coherence tomography | |
| Guimarães et al. | Three-dimensional segmentation and reconstruction of the retinal vasculature from spectral-domain optical coherence tomography | |
| Leopold et al. | Deep learning for ophthalmology using optical coherence tomography |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| WW01 | Invention patent application withdrawn after publication |
Application publication date: 20181214 |
|
| WW01 | Invention patent application withdrawn after publication |