CN109816661B

CN109816661B - 一种基于深度学习的牙齿ct图像分割方法

Info

Publication number: CN109816661B
Application number: CN201910219774.6A
Authority: CN
Inventors: 饶云波; 苟苗; 王艺霖; 张孟涵; 郭毅; 程奕茗
Original assignee: University of Electronic Science and Technology of China; CERNET Corp
Current assignee: University of Electronic Science and Technology of China; CERNET Corp
Priority date: 2019-03-22
Filing date: 2019-03-22
Publication date: 2022-07-01
Anticipated expiration: 2039-03-22
Also published as: CN109816661A

Abstract

本发明属于医学CT(Computed Tomography，计算机断层扫描)图像分割技术领域，涉及一种基于深度学习的牙齿CT图像分割方法。本发明的技术方案，将传统的Level Set算法与U‑net网络模型相结合，利用Level Set算法来解决神经网络所需要训练集问题，使得神经网络可以使用未标记的标签进行训练，同时利用神经网络模型来完成图像的自动分割问题，避开曲线演化的不收敛问题，实现在医学图像训练集不充足的情况下，能获得足够精确的分割效果。

Description

一种基于深度学习的牙齿CT图像分割方法

技术领域

本发明属于医学CT(Computed Tomography，计算机断层扫描)图像分割技术领域，涉及一种基于深度学习的牙齿CT图像分割方法。

背景技术

医学CT图像分割问题由来已久，很多研究人员都在该领域提出各自算法用以实现医学图像的分割，并取得了一定程度上的成果。由于医学图像的诸多限制，目前基于统计的分割、图的分割、活动轮廓分割、深度学习的分割等方法仍存在很多问题，对比分析如下。

基于统计的分割方法包括大津阈值法、多阈值控制分水岭方法、自适应阈值方法、迭代阈值技术等。该方法计算结构简单，容易理解，适合分割比较小且简单的结构。缺点是若区域情况比较复杂，则难以获得良好的分割效果，并且容易出现空洞现象和不准确的分割，如图1所示。

基于图的分割方法如Graph Cut等算法，将像素点转化为无向图，之后构建节点之间的权重函数，利用最大流最小割理论进行求解，如图2所示。该方法要求分割的图形具有宏观上的一致结构，在保证分割精度的前提下，对权重函数的构造和求解便成为影响分割效果好坏的关键步骤。缺点是该方法需要人工参与且对人员技术和经验有要求，需要足够充分的先验模型。

基于活动轮廓Active Contour Model(ACM)的分割方法如snake算法，初始化一条曲线，按照曲线演化理论进行演化，利用能量泛函的极值理论使得曲线收敛于目标边界。该方法具有良好的抗噪性，在目前应用较广。缺点是对初始曲线要求严格，初始曲线极小的差异可能会导致完全不同的演化结果。除此之外，每次演化迭代都会占用大量计算资源，曲线演化也可能停止在不准确的边界，如图3所示。

基于深度学习的分割方法。是当今最为火热的技术，通过采集图像正负样本的相关特征，对数据浅层特征、深层特征进行分析，实现图像分割。该方法拥有极好的数据学习能力，在医学图像分割上也能取得良好的结果，但由于医学图像难以满足深度学习所需要的训练图片数量，过拟合严重。因此深度学习方法在医学图像上的应用并不广泛。

基于以上分析，对基于统计的分割、基于图的分割、基于活动轮廓分割、基于深度学习的分割等方法进行了分析对比，表1是对比分析情况：

表1传统CT图像分割算法对比

发明内容

随着计算机技术的飞速发展、GPU的更新迭代，当今计算速度越来越快，深度学习应用越加广泛，但是深度学习训练模型时所使用的训练集图片数量庞大，医学图像难以提供足够有效的CT图片，且人工标记图片也存在过于严重耗时的问题，容易造成数据过拟合问题导致模型无法应用于训练集外的图片，这也是目前医学图像领域深度学习发展远不如自然图像处理的主要原因。针对这一问题，本发明首先使用Level Set算法思想，结合曲线演化理论，使得初始化曲线收敛于图像边界并输出作为图像标记，用以减低人工负担。在此之上使用U-Net网络进行模型训练，用以减轻医学图像训练图片不足引起的过拟合问题。本发明提出的方法主要包括三个步骤：(1)牙齿CT图片预处理，(2)牙齿CT图像标记，(3)模型训练及优化。

本发明的技术方案是：

步骤A：牙齿CT图片预处理

数据集准备与预处理是训练模型的基础，本研究的原始CT图像数据由四川大学华西医院提供，共5组牙齿上下颌的CT图像数据，每组为一个人的完整口腔扫描结果，共401张原始牙齿CT图像，图片文件为DCM医用格式。DCM文件是遵循DICOM(Digital Imaging andCommunications in Medicine,医疗数字影像和通信协议)标准的一种文件，它定义了能满足临床需求又可用于数据通信的图像格式。目前主流的心电图、核磁共振成像、血管镜、超声心动图、计算机断层扫描均采用该格式进行存储和传输，但因为其文件格式中加入了用于通信协议的相关报头，DCM文件的浏览与修改都需要专业软件来进行，不适用于深度学习的训练，因此需要对DCM文件中的图片数据进行提取与整理。

经研究中，除去DCM格式文件多余的用于通信的数据后，需获取DCM格式文件数据中窗宽和窗位信息。窗宽是指显示图像时所选用的CT值范围，其大小影响到医学图片的对比度。窄窗宽对比度较高，适合用于脑检查等分辨密度较接近的组织。宽窗宽对比度较低，适合分辨密度差异较大的组织，如肺、骨。窗位是指窗宽上下限的CT值平均数，窗位的高低影响图像的亮度。本研究中将DCM文件数据灰度级别映射到0-255，使用最小值进行加窗，其具体实现过程如下：

对原始图片进行预处理。将原有的DCM图片输入到MATLAB软件中转化为一般的PNG格式图片。具体方法是：过滤掉DCM文件中用于通信传输的metadate报头信息，保留存储图像信息的部分，得到401*401的16位int数据，之后根据医学图像相关协议计算图像数据部分的窗宽与窗位，以CT值跨度作为窗宽，最小CT值作为窗位，对数据进行加窗，将数据映射到标准图片的0-255像素级之间，转化为8位的uint类型数据并保存为PNG图像文件。本发明使用了5组原图，转化后有效图片数量为200张左右，如图5所示。

步骤B：牙齿CT图像的标记

DCM文件转化为容易修改编辑的PNG格式后，本研究提出一种新的基于Level Set算法对牙齿CT图片进行自动分割标记的算法。在本发明的算法中，主要步骤如下：

B1首先对原始曲线泛函进行水平集化，显式曲线公式如公式(1)表示：

其中F表示控制曲线演化的速度函数，N表示曲线向内的法向量，s表示曲线演化时任意的参数化变量，t表示时刻。将C(s,t)公式改写成t时刻的一个零水平集φ(x,y,t)，使用φ来表示进行曲线演化，其中(x，y)表示曲线上的点在图像的坐标。φ初始化如公式(2)表示：

其中d为图中点到曲线的最短距离，用高维的曲面能较好的演化曲线的分裂和融合。

B2使用距离约束作为内部能量项、区域像素差作为外部能量项来构建新型水平集能量函数如公式(3)所示：

公式(3)第一项为距离约束项，

表示图像演化边界的梯度变化，用以控制在演化过程中水平集函数和符号距离函数保持一致，第二项和第三项为外部能量项，表示演化曲线内外分别对各自区域灰度均值的距离，当两个外部能量项均接近0时，能量函数E取最小值，演化曲线即为所求分割曲线。第二项接近0而第三项远大于0时，目标边界在演化曲线外，曲线应向外扩张，反之曲线应向内收缩，函数用外部能量项的数值控制曲线演化。其中Ω表示整个图像域，μ₀(x,y)表示点(x,y)处的灰度值，c₁表示演化曲线内部的灰度均值，c₂表示演化曲线外部的灰度均值，μ、λ₁、λ₂均为常数。

B3在本发明中，使用Heavside函数保证能量泛函的连续与平滑，如公式(4)表示：

其中ε是一个趋于0的正数，使用公式(4)带入能量函数公式(3)各项，得到新的能量泛函如下：

∫_φ>0|μ₀(x,y)-c₁|²dxdy＝∫_Ω|μ₀(x,y)-c₁|²H(φ(x,y))dxdy (5)

∫_φ<0|μ₀(x,y)-c₂|²dxdy＝∫_Ω|μ₀(x,y)-c₂|²(1-H(φ(x,y)))dxdy (6)

对公式(7)进行隐式拆分迭代求解可得到函数的数值求解，得出图片边界的分割结果如图6(b)所示，之后将曲线内部作为目标，曲线外部作为背景进行二值化，如图6(c)所示。在本发明中，初始化公式(7)中常数λ1、λ2值为1，μ值为0.02，演化步长η值为0.5，迭代次数取为500。

牙齿CT图像标记的另一部分是由人工使用LabelMe软件制作完成。该软件在Window 10操作系统下，配置python3.6、pyqt5、pillow4.0等运行，在Anaconda的虚拟环境下安装运行LabelMe，软件图形界面如图7(a)所示，牙齿CT图片分割结果如图7(b)所示。

为了保证训练精度，可以通过人工完成部分图片标记的方式，人工标记图片数量与自动标记图片数量比例为1:1，本发明中使用LabelMe进行手工标记制作。

LabelMe所生成的标签图片灰度值范围为0-1，如图7(b)所示，而在深度学习中接收的图片灰度值范围为0-255，LabelMe直接生成的图片在网络中均视为黑色，不具有标记效果，因此在本发明的方案中，需要将其结果映射到0-255的范围。

步骤C：模型训练及优化

在模型训练与优化阶段，本发明提出了基于U-Net的新的深度学习模型，下采样阶段提取图像特征实现分类，上采样阶段使用逆卷积方式恢复特征尺寸、位置，补充高层特征信息以完成定位。U-Net网络模型下采样阶段中由5组节点组成，下采样阶段每经过一个节点，特征图谱通道数就增加一倍，每组节点中包含2个2*2的卷积层和1个2*2的最大池化层，卷积层使用Relu激活函数，在第4个节点中，2层卷积层过后加入一个Dropout层，用以减轻过拟合问题，之后再接2*2的最大处化层，第5个节点中移除最大池化层，只保留Dropout层和卷积层，作为下采样阶段的结果。下采样阶段中第5个节点的输出作为上采样的输入，对其进行逆卷积，卷积核大小、步长也下采样中最大池化层保持一致，逆卷积后的结果再与下采样中对应节点进行裁剪拼接融合，用以补充特征图谱高层位置信息。上采样阶段节点结构与下采样阶段保持一致，每个节点由2个3*3的卷积层和1个2*2的卷积层组成。上采样结束后，在添加一层1*1的卷积层，使用Sigmoid激活函数，将模型输出的64个特征图谱通道转化为需要的结果。在本发明的方案中深度学习模型共30层，其具体结构如图8所示。

本发明主要应用于二分类，在模型训练中使用交叉熵损失函数，其中y为模型输出概率，

为期望概率，如公式(8)、(9)所示：

基于原理以上分析，本发明中所使用的网络模型取消了全连接层，输入图片大小不再受限制，可以使用任意大小进行训练与检测。本研究在Keras框架下训练模型，Keras后台使用TensorFlow，GPU为NVIDIA GeForce GTX 1070，使用GPU加速来减小训练时间。具体实施过程如下：

C1模型训练。将对应的原图和标签以相同的方式命名，并放置于image与label文件夹目录下。使用偏移、旋转、镜像、反转等方式对图片和标签进行数量的倍增，图片和标签使用同种方法以保证扩充后的图片与标签仍能一一对应。模型训练时使用交叉熵损失函数，如公式(8)、(9)所示，TensorFlow中优化模型过程中模型学习率设置为0.01，即模型数值更新在上一次迭代数值的0.01范围内，避免数值更新步幅过大越过最优值，模型采用择优记录的方式每一次迭代只记录最优的结果，模型保存为tooth.hdf5，训练结果如图9所示。

C2模型验证：得到最优训练模型后，使用12张牙齿CT图片进行模型的验证，图片分别按顺序命名为0-11，图片格式为PNG格式，放于test文件夹目录下。使用得到的tooth,hdf5模型进行牙齿的分割，分割结果存储在test文件夹同根目录下的result文件夹下。分割结果如图10所示。

本发明的技术方案，将传统的Level Set算法与U-net网络模型相结合，利用LevelSet算法来解决神经网络所需要训练集问题，使得神经网络可以使用未标记的标签进行训练，同时利用神经网络模型来完成图像的自动分割问题，避开曲线演化的不收敛问题，实现在医学图像训练集不充足的情况下，能获得足够精确的分割效果。

基于Level Set算法的自动标记生成：本发明利用区域信息来构建Level Set能量函数，初始曲线将图像划分为内外两部分，以划分区域内的平均灰度值与各个像素灰度值的差值作为外部能量约束来进行曲线的扩张或收缩演化，弱化初始边界对收敛结果的影响，使得对初始边界的抉择不再艰难。添加水平集梯度相关约束项，保证水平集函数始终和符号距离函数保持一致，尤其是在零水平集附近位置保持一致，用这种方法来规避曲线演化过程中，符号距离函数渐渐脱离水平集函数所引起的重新初始化问题，用以提高曲线收敛速度与精度。

牙齿CT图像特征点提取：本发明基于语义分割中U-Net网络模型，利用其高效性，应用并修改使其适用于牙齿CT图片的特征点提取。使用由华西医院提供的5组完整的DCM图像文件，建立了与其对应的标签图片用以训练网络模型，在训练迭代过程中保留下最优的权重结果，最终得到完成的网络模型，使其可直接用于之后的牙齿分割。

本发明的有益效果为，将传统Level Set算法与U-net网络方法融合，结合双方优点，解决标签图片操作困难问题，减轻人力，使用图像增强方法为神经网络模型训练扩充数据库，使其能避免过拟合问题，使用语义分割方式实现像素级的图像分类，经过迭代获取到适用于牙齿分割的最优模型。

附图说明

图1为大津阈值法分割算法的结果，(a)为牙齿的CT原图，(b)为大津阈值法分割结果；

图2为Graph Cut分割算法，(a)为CT原图，(b)为设置graph cut算法前景、背景，(c)分割结果；

图3为Snake分割算法，(a)为图像原图及初始曲线，(b)为曲线演化停止曲线结果；

图4为本发明提出方法的框架图；

图5为转换后的PNG图示例；

图6为Level set自动标记结果，(a)输入牙齿CT原图，(b)迭代500次时曲线演化结果，(c)经Level Set标记结果；

图7为人工制造训练标签，(a)为LabelMe读入牙齿CT图片，(b)使用LabelMe输出的标记图；

图8为深度学习模型结构；

图9为模型训练结果示意图；

图10为模型分割结果，前两行为测试牙齿CT图像原图，后两行模型预测结果。

具体实施方式

在发明内容部分已经对本发明的技术方案进行了详细描述，在此不再赘述。

Claims

1.一种基于深度学习的牙齿CT图像分割方法，其特征在于，包括以下步骤：

S1、预处理：将原始CT图像转化为位图图像；

S2、采用Level Set算法，对步骤S1获得的图像进行自动分割标记，将获得的标记图片与相对应的原始CT图像作为训练集；Level Set算法的具体过程为：

S21、对原始曲线泛函进行水平集化，显式曲线公式如公式(1)：

其中F表示控制曲线演化的速度函数，N表示曲线向内的法向量；将C(s,t)公式改写成t时刻的一个零水平集φ(x,y,t)，使用φ来表示进行曲线演化，φ初始化如公式(2)：

其中d为图中点到曲线的最短距离；

S22、使用距离约束作为内部能量项、区域像素差作为外部能量项来构建水平集能量函数如公式(3)：

公式(3)第一项为距离约束项，

表示图像演化边界的梯度变化，用以控制在演化过程中水平集函数和符号距离函数保持一致，第二项和第三项为外部能量项，表示演化曲线内外分别对各自区域灰度均值的距离，当两个外部能量项均接近0时，能量函数E取最小值，演化曲线即为所求分割曲线；第二项接近0而第三项远大于0时，目标边界在演化曲线外，曲线应向外扩张，反之曲线应向内收缩，函数用外部能量项的数值控制曲线演化，其中Ω表示整个图像域，c1表示演化曲线内部的灰度均值，c2表示演化曲线外部的灰度均值，μ、λ1、λ2均为常数，μ₀(x,y)表示点(x,y)处的灰度值；

S23、使用Heavside函数保证能量泛函的连续与平滑，如公式(4)：

∫_φ＞0|μ₀(x，y)-c₁|²dxdy＝∫_Ω|μ₀(x，y)-c₁|²H(φ(x，y))dxdy (5)

∫_φ＜0|μ₀(x，y)-c₂|²dxdy＝∫_Ω|μ₀(x，y)-c₂|²(1-H(φ(x，y)))dxdy (6)

对公式(7)进行隐式拆分迭代求解可得到函数的数值求解，得出图片边界的分割结果，之后将曲线内部作为目标，曲线外部作为背景进行二值化；

S3、采用具有编码器-解码器结构的U-Net网络来提取图片的特征图谱，构建深度学习模型，通过迭代实现对模型权重的更新，获取最优模型；

S4、根据获取的最优模型对经过步骤S1处理的CT图像进行分割。

2.根据权利要求1所述的一种基于深度学习的牙齿CT图像分割方法，其特征在于，所述步骤S2中，采用Level Set算法标记一半步骤S1获得的图像，另一半通过人工标记的方式进行标记。

3.根据权利要求1所述的一种基于深度学习的牙齿CT图像分割方法，其特征在于，所述步骤S3的具体方法为：

构建基于U-Net的深度学习模型，下采样阶段提取图像特征实现分类，上采样阶段使用逆卷积方式恢复特征尺寸、位置，补充高层特征信息以完成定位；

U-Net网络模型下采样阶段中由5组节点组成，下采样阶段每经过一个节点，特征图谱通道数就增加一倍，每组节点中包含2个2*2的卷积层和1个2*2的最大池化层，卷积层使用Relu激活函数，在第4个节点中，2层卷积层过后加入一个Dropout层，用以减轻过拟合问题，之后再接2*2的最大池化层，第5个节点中移除最大池化层，只保留Dropout层和卷积层，作为下采样阶段的结果；下采样阶段中第5个节点的输出作为上采样的输入，对其进行逆卷积，卷积核大小、步长与下采样中最大池化层保持一致，逆卷积后的结果再与下采样中对应节点进行裁剪拼接融合，用以补充特征图谱高层位置信息；上采样阶段节点结构与下采样阶段保持一致，每个节点由2个3*3的卷积层和1个2*2的卷积层组成；上采样结束后，在添加一层1*1的卷积层，使用Sigmoid激活函数，将模型输出的64个特征图谱通道转化为需要的结果；

模型训练中使用交叉熵损失函数：

其中y为模型输出概率，

为期望概率。