CN109064476B - 一种基于水平集的ct胸片肺组织图像分割方法 - Google Patents

一种基于水平集的ct胸片肺组织图像分割方法 Download PDF

Info

Publication number
CN109064476B
CN109064476B CN201810819746.3A CN201810819746A CN109064476B CN 109064476 B CN109064476 B CN 109064476B CN 201810819746 A CN201810819746 A CN 201810819746A CN 109064476 B CN109064476 B CN 109064476B
Authority
CN
China
Prior art keywords
image
level set
contour
function
outlines
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.)
Active
Application number
CN201810819746.3A
Other languages
English (en)
Other versions
CN109064476A (zh
Inventor
王蓉芳
杨靖
陈佳伟
郝红侠
缑水平
刘红英
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xidian University
Original Assignee
Xidian University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Xidian University filed Critical Xidian University
Priority to CN201810819746.3A priority Critical patent/CN109064476B/zh
Publication of CN109064476A publication Critical patent/CN109064476A/zh
Application granted granted Critical
Publication of CN109064476B publication Critical patent/CN109064476B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/12Edge-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/13Edge detection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/60Analysis of geometric attributes
    • G06T7/62Analysis of geometric attributes of area, perimeter, diameter or volume
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10081Computed x-ray tomography [CT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20036Morphological image processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing
    • G06T2207/30061Lung

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geometry (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Image Processing (AREA)

Abstract

本发明公开了一种基于水平集的CT胸片肺组织图像分割方法,解决CT图像中肺区域分割边界不精确、需人工干预、结果不稳定的问题。实现过程有:获取CT图像并预处理;构造能量泛函并设定初始零水平集;最小化能量泛函,得零水平集轮廓;从中选出候选肺部区域轮廓;逐个向内候选肺部区域填充轮廓;对填充结果进行形态学开、闭操作,移除小体积连通区域。本发明提取了图像的边缘信息,基于先验知识设计了稳定的轮廓筛选策略,有效的筛选出候选肺部区域轮廓,最后基于轮廓高度信息设计了轮廓填充的优化方案。本发明图像分割结果鲁棒,精度高,是一种全自动的图像分割方法。本发明提取出了CT图像肺部区域,可用于后续对CT图像肺部区域分析。

Description

一种基于水平集的CT胸片肺组织图像分割方法
技术领域
本发明涉及医学影像处理技术领域,尤其是涉及肺部图像分割,具体是一种基于水平集的CT胸片肺组织图像分割方法,用在肺部CT影像处理。
背景技术
近年来,随着计算机断层扫描成像技术的不断提高,可以提供高清晰度、高对比度的CT图像用于人体疾病的诊断。使用肺部CT图像观察肺部结构和功能特征是当今临床针对肺部疾病的重要诊断手段。同时,由于计算机计算能力的增强,深度学习技术在计算机视觉领域取得了突破性的成果。自然的,计算机自动检测肺结节技术也得到了发展。原始的肺部CT图像中包含大量的背景、肌肉、脂肪、病人衣物、骨骼等其他人体组织,为了便于自动检测肺结节,需要先对肺部CT图像处理,提取分割肺组织图像。在提取肺组织图像后,可以大大降低肺结节检测中的假阳性结果,缩小结节搜索空间,提升算法性能。
CT胸片肺组织图像的精确分割是肺结节自动检测与辅助诊断的基础,因此这也是医学图像处理领域的研究热点。目前,CT胸片肺组织图像分割方法有:(1)阈值法,阈值法是最简单的肺组织分割方法,它快速、实现简单,但是不够鲁棒,不能有效去除背景与无关人体组织,易受CT成像时的放射剂量水平影响。(2)区域生长法,区域生长法是一种半自动肺组织分割方法,需要人工参与设置初始种子生长点,人工设置种子点易受主观影响分割结果。(3)基于图像配准的分割方法,该方法是通过与训练模板配准的方法来达到分割的效果,此方法易受分割模板的影响,尤其是难以分割一些具有特异性的肺组织。以上所有方法,主要利用图像的灰度信息来分割,所以在分割结果的边界处会包含有肺部以外的组织,同时容易遗漏位于分割结果边界的肺结节。
现有方法的分割结果不够鲁棒、严重依赖分割模板、分割边界不够精确,分割过程需要人工干预。
发明内容
本发明目的是针对现有技术的不足,提供一种鲁棒、全自动、精确的基于水平集的CT胸片肺组织图像分割方法。
本发明是一种基于水平集的CT胸片肺组织图像分割方法,其特征是,包括有如下步骤:
步骤1获取包含肺组织的胸部CT图像,并进行CT胸片图像预处理;
步骤2在预处理后的胸部CT图像上构造能量泛函方程并进行设定初始零水平集;
步骤3在初始零水平集的基础上最小化能量泛函方程,获取CT胸片图像零水平集轮廓图像;
步骤4挑选出候选肺部区域轮廓:从CT胸片零水平集轮廓图像中进行候选肺部区域轮廓挑选,首先计算出所有零水平集轮廓图像的外接矩形及轮廓的面积,按照轮廓的面积从大到小排序,得到已排序的轮廓列表,计算出轮廓的面积与轮廓的外接矩形的面积之比,然后计算出轮廓的外接矩形的宽度与外接矩形的高度,最后计算出轮廓的外接矩形的长宽比,最终根据轮廓的面积与轮廓的外接矩形的面积之比、轮廓的外接矩形的高度、轮廓的外接矩形的宽度、轮廓的外接矩形的长宽比零水平集轮廓图像中挑选出候选肺部区域轮廓图像;
步骤5根据候选肺部区域轮廓图像的高度逐个向内填充轮廓,获得三维初步分割结果图像;
步骤6对三维初步分割结果图像依次进行三维图像形态学开操作,闭操作,选择使用球形结构元作为基本形态学结构元,结构元半径大小为3,然后对分割结果进行连通域分析,移除体积小于体积阈值w的三维连通区域,得到最终三维分割结果图像。
针对现有方法的缺点,客观上需要提出一种鲁棒、全自动、精确的CT胸片肺组织图像分割方法。
和现有技术方案相比,本发明具有以下优点:
(1)对不同剂量水平的CT图像具有鲁棒性,因为不同的CT图像在成像时具有不同的放射性水平、有不同的成像设备、不同的环境噪声。所以CT图像往往具有不同的灰度值范围,然而传统的CT胸片肺组织图像分割方案往往需要根据CT图像的灰度信息进行分割,所以传统方案极易受到CT图像的灰度值范围的影响,但本方案主要利用了CT图像边缘信息,所以此方案具有很强的鲁棒性;
(2)图像分割过程全自动,传统的阈值图像分割方案、种子生长法图像分割方案都需要人工手动设定图像分割阈值、种子初始生长点。所以图像分割结果极易受到主观经验的影响,也会对图像分割速度造成影响。本方案采用水平集方法,只需设定图像的中心区域作为初始区域,即可通过水平集的自动演化得到所需要的零水平集图像轮廓。
(3)分割结果精度高,传统形态学CT胸片肺组织图像分割方案、阈值肺分割方案分割结果中包含大量的血管,胸腔组织等非肺部区域容易给后续结节检测带来影响。但是此方案利用了肺部区域的边缘信息,通过最终的肺部区域轮廓向内填充来获得最终的图像分割结果,所以此方案分割结果更加精确,基本不包含肺部以外区域,同时对初始图像分割结果进行了三维图像形态学处理,利用了图像的三维信息,使得图像分割边界更精确,进一步减少了肺部以外区域的错误分割结果。
附图说明
图1是本发明实现流程图
图2为本发明实施例2中的待分割肺部CT图像。
图3为本发明实施例3中绘制出的零水平集轮廓图像。
图4为本发明实施例6中经过轮廓处理后的结果。
图5为本发明实施例7中的最终三维分割结果。
图6为本发明实施例9中的仿真实验1结果图。
图7为本发明实施例9中的仿真实验2结果图。
具体实施方式
以下参照附图,对本发明的技术方案和效果做详细描述。
实施例1
现有CT胸片肺组织图像分割方法的结果不够鲁棒、严重依赖分割模板、分割边界不够精确,分割过程需要人工干预。针对现有方法的缺点,需要一种鲁棒、全自动、精确的CT胸片肺组织图像分割方法。为此本发明提出一种基于水平集的CT胸片肺组织图像分割方法,参见图2,包括有如下步骤:
步骤1获取包含肺组织的胸部CT图像,使用MATLAB把图像读取进入计算机,以三维张量在计算机中存储图像,并进行CT胸片图像预处理,图像预处理的方法是给CT图像灰度值加400,这样可以调整原始输入CT图像的灰度直方图分布,使肺部区域在灰度直方图中更加突出,有利于后续最小化能量泛函方程,为了降低算法的复杂度且减少肺部区域以外图像区域的干扰,最后把CT图像的灰度值范围归一化到8位无符号数。
步骤2在预处理后的胸部CT图像上构造能量泛函方程,通过构造能量泛函方程可以有效的把寻找最优水平集函数的过程形式化;在构造能量泛函过程中引入偏置场参数b作为偏置场校正,降低了输入CT图像的多样性对水平集演化过程带来的影响;同时构造能量泛函方程时引入了零水平集弧长函数L(φ)作为正则函数,可以有效避免水平集演化过程中向内收缩的情况。然后开始设定初始零水平集,为了使初始零水平集尽可能准确,也就是说希望零水平集可以位于肺部区域以内,所以根据CT胸片肺部区域先验信息选择将初始零水平集设置在CT图像的中心偏左和偏右的位置。
步骤3在以上能量泛函方程和初始水平集的基础上,当前的主要的主要问题已经转化成为一个优化问题,针对能量泛函方程的特点,选择使用梯度下降算法作为优化算法,逐步向最优的水平集函数演化,在每一次水平集函数更新后都需要计算偏置场参数b和灰度平均值向量参数c,逐步优化直到满足预设的停止优化条件。在优化过程完毕以后可以得到对应于当前图像的最优水平集函数,把它直接应用到原始图像尺度上即可得到零水平集轮廓图像。
步骤4挑选出候选肺部区域轮廓:从CT胸片零水平集轮廓图像中进行候选肺部区域轮廓挑选,首先计算出所有零水平集轮廓图像的外接矩形及轮廓的面积,按照轮廓的面积从大到小排序,得到已排序的轮廓列表,计算出轮廓的面积与轮廓的外接矩形的面积之比,然后计算出轮廓的外接矩形的宽度与外接矩形的高度,最后计算出轮廓的外接矩形的长宽比,最终根据轮廓的面积与轮廓的外接矩形的面积之比、轮廓的外接矩形的高度、轮廓的外接矩形的宽度、轮廓的外接矩形的长宽比零水平集轮廓图像中挑选出候选肺部区域轮廓图像。
步骤5因为候选肺部区域轮廓列表中有可能依然存在没有剔除的肺部区域以外轮廓,但是拍在轮廓列表靠前位置的一定是肺部区域轮廓。肺部区域轮廓以外轮廓一定是存高度差,所以这里选择根据候选肺部区域轮廓图像的高度逐个向内填充轮廓,具体就是如果轮廓高度相差大于30个像素大小就选择丢弃此轮廓,最终获得三维初步分割结果图像。
步骤6三维初步分割结果图像中往往会在肺部区域内包含一些极其细小的空洞,肺部区域以外可能会存在一些小的孤立区域,根据经验选择对三维初步分割结果图像依次进行三维图像形态学开操作,闭操作,针对肺部区域的形态学特征,本方案选择使用球形结构元作为基本形态学结构元,结构元半径大小为3,最后对分割结果进行连通域分析,移除体积小于体积阈值w的三维连通区域,得到最终三维分割结果图像。
本发明采用了基于水平集方法的肺区域边缘轮廓提取技术手段,使得图像分割结果边缘精确;同时在构造能量泛函方程时采用了弧长函数作为正则函数,使得能量泛函方程可以快速有效的收敛到全局最优点,提升了算法的鲁棒性;最后采用了肺区域轮廓的先验信息,从复杂的零水平集轮廓中挑选出了准确的肺部区域轮廓。因为采用了以上技术手段,所以算法得到的图像分割结果精度高、鲁棒性强,相比于传统的形态学CT胸片肺组织图像分割方法、阈值CT胸片肺组织图像分割方法分割结果中包含大量的非肺部区域容易给后续结节检测带来影响,而本发明分割结果精确,基本不包含肺部以外区域。
实施例2
基于水平集的CT胸片肺组织图像分割方法同实施例1-1,参见图2,因为原始输入CT图像的灰度值动态范围是-1024到+1024,且经过统计肺部区域的图像灰度范围在-600到-200,所以步骤1中所述的图像预处理是给CT图像灰度值加400,这样可以将肺部区域的图像灰度范围调整到-200到+200,使肺部区域更突出,方便后续的处理。为了降低算法的复杂度且减少肺部区域以外图像区域的干扰,所以选择将图像灰度值归一化到8位无符号数。
实施例3
基于水平集的CT胸片肺组织图像分割方法同实施例1-2,步骤2中所述的在预处理后的CT胸片肺组织图像上构造能量泛函方程并进行设定初始零水平集,具体包括有如下步骤:
2.1构造能量泛函方程F(φ,c,b)
F(φ,c,b)=ε(φ,c,b)+υL(φ)+μRp(φ)
其中ε(φ,c,b)是能量函数,φ为水平集函数,c为水平集内外的灰度平均值向量,因为CT图像在成像时具有不同的放射性水平、有不同的成像设备、不同的环境噪声,所以在构造能量泛函方程时引入了偏置场校正,偏置场校正可以降低上述因素给图像分割带来的影响,b为偏置场参数,同时在构造能量泛函方程时引入了零水平集弧长函数L(φ)作为正则函数,因为选择的初始水平集都位于图像中心较小的区域,零水平集弧长函数可以有效避免水平集演化过程中向内收缩的情况,尤其是在CT胸片肺组织图像这个领域。CT胸片中肺部区域是独立且内部连续的两个左右对称结构,这就要求最终的水平集也具有类似的结构,也就是说要求能量密度分布较为平均,所以在能量泛函方程中本发明也引入了能量密度惩罚函数Rp(φ),在引入能量密度惩罚函数后可以有效减少图像分割结果中的毛刺现象的出现。υ零水平集弧长函数正则项,μ是能量密度惩罚函数正则项。
能量泛函方程的能量函数具体表示为
Figure GDA0001838614490000061
其中N代表待分割区域数,Mi是成员函数,其中Mi成员函数具体定义为以下形式:
M1(φ(x))=H(φ(x)),M2(φ(x))=1-H(φ(x));
其中,赫维赛德函数H(φ(x))为以下形式:
Figure GDA0001838614490000062
区域函数ei(x)为以下形式:
ei(x)=I21K-2ci(b*K)+ci 2(b2*K);
I其中为待分割图像,ci代表当前水平集区域的灰度平均值,*是卷积算子,K是一个高斯卷积核,指示函数1K(x)为以下形式:
1K(x)=∫K(y-x)dy;
K(y-x)是一个非负的窗函数,满足K(y-x)=1,
Figure GDA0001838614490000063
Oy代表水平集区域y。
零水平集弧长函数L(φ)为如下形式:
Figure GDA0001838614490000064
其中,
Figure GDA0001838614490000065
为梯度算子。
能量密度惩罚函数Rp(φ)为以下形式:
Figure GDA0001838614490000066
其中,势函数p(s)为以下形式:
p(s)=(1/2)(s-1)2
综上,以上即可完成能量泛函方程的构造。
2.2根据CT成像特点以及人体肺部区域位置先验信息,再加上水平集向外演化扩张的特点,可以将初始零水平集设置在CT图像的中心偏左和偏右的位置,且初始零水平集不宜过大,以保证其位于肺部区域以内,所以初始零水平集的设定为:首先计算出CT图像的中心,然后分别对此中心向左和向右偏移m个像素得到两个初始零水平集的中心位置,以此作为两个正方形的中心,正方形的大小为n,最终将这两个正方形作为初始零水平集。
本发明只需设定图像的中心区域作为初始区域,初始区域的计算量小、不需要人工干预,且中心区域的面积较小可保证初始区域中一定包含肺部区域,这样有利于后续的水平集的自动演化最终得到所需要的零水平集图像轮廓。
实施例4
基于水平集的CT胸片肺组织图像分割方法同实施例1-3,步骤3中所述的在初始零水平集的基础上最小化能量泛函方程,梯度下降算法作为一种有效简单的算法,被大规模使用在各类优化问题中,所以最小化能量泛函方程中采用梯度下降算法作为优化算法,包括有如下步骤:
3.1对能量泛函方程中的水平集函数参数φ更新,水平集函数参数φ的更新公式为以下形式:
Figure GDA0001838614490000071
其中,阶跃函数δ(φ)=H′(φ),距离函数
Figure GDA0001838614490000072
3.2对能量泛函方程中的灰度平均值向量参数c更新,参数灰度平均值向量c的更新公式为以下形式:
Figure GDA0001838614490000073
其中N=2,ui(y)=Mi(φ(y)),在计算区域灰度平均值时使用了偏置场校正,减少了输入CT图像的多样性对水平集演化过程的影响,使水平集演化更加稳定,加快演化速度。
3.3对能量泛函方程中的偏置场参数b更新,偏置场参数b的更新公式为以下形式:
Figure GDA0001838614490000074
在其中J(1)与J(2)分别为:
Figure GDA0001838614490000081
Figure GDA0001838614490000082
在最小化能量泛函过程中,选择使用了梯度下降算法作为优化算法,优化过程简洁有效,在具体实现算法时可以降低实现难度;同时在每一次优化过程计算灰度平均值向量时,使用了偏置场校正,使水平集演化过程更加的稳定。
实施例5
基于水平集的CT胸片肺组织图像分割方法同实施例1-4,图3中显示了所有的零水平集轮廓,从图3可见,经过水平集方法提取出的轮廓线中含有大量肺部区域以外的轮廓。例如其中包含有病人身体轮廓、病人衣物轮廓、CT扫描仪器轮廓等。如果直接对所有轮廓直接向内填充来获得分割结果必然只能得到一个不理想的分割结果。所以,如何从这些轮廓当中确定肺部区域轮廓就是需要解决的问题。实际上,在水平集方法提取的轮廓线的基础之上,只需去掉病人身体轮廓、病人衣物轮廓、CT扫描仪器轮廓,然后将所有剩余的轮廓向内填充就可获得基本的肺部区域分割结果,为了论述方便以下统一将病人身体轮廓、病人衣物轮廓命名为外轮廓。分析图3所示的水平集方法提取出的轮廓线可得:a)CT扫描仪轮廓呈狭长形;b)外轮廓的面积绝对大于肺部区域轮廓;c)显然外轮廓的外接矩形宽度必然大于原始输入图像宽度的一半,但是由于肺部区域呈左右对称状,所以可见肺部区域轮廓的外接矩形宽度必然小于原始输入图像宽度的一半。由此本发明设计出轮廓筛选策略。
步骤4中所述的挑选出候选肺部区域轮廓的具体方法是,从已排序的轮廓列表中的每个轮廓依次判断下列条件:
4.1若轮廓的外接矩形宽度小于CT图像宽度的一半则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.2轮廓的面积与轮廓的外接矩形的面积之比小于r则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.3若轮廓的外接矩形的长宽比小于s则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓。
此处面积之比筛选阈值参数r、外接矩形长宽比参数s可以根据实际CT图像的大小,噪声程度来动态设定,以获得更好的图像分割结果。上述筛选策略利用了肺部区域轮廓和外轮廓的先验信息,抓住了肺部区域轮廓与外轮廓的不同点,可以有效的剔除对图像分割结果影响最大的外轮廓。为后续轮廓的填充降低了复杂度。
实施例6
基于水平集的CT胸片肺组织图像分割方法同实施例1-5,由于候选肺部区域轮廓依然存在一些肺部区域轮廓以外的小轮廓,如果对这些轮廓向内填充依然会对图像分割结果带来影响,因为候选肺部区域轮廓是按照轮廓面积大小来进行排序,所以候选肺部区域轮廓列表中的前几个一定是属于肺部区域轮廓。基于这一点,可以先填充的轮廓高度来避免填充后续的小轮廓。步骤5中所述的根据候选肺部区域轮廓的高度逐个向内填充轮廓的具体方法是,确定第一个轮廓的外接矩形右下角点的高度,记为bottom,对候选肺部区域轮廓的每一个轮廓执行以下步骤:
5.1确定此轮廓的右下角点的高度,记为tmp_bottom。
5.2若满足以下条件则丢弃这个轮廓,返回5.1,否则,即不满足下列条件的情况下执行步骤5.3,对轮廓向内填充像素。
(bottom+tmp_bottom)/2-bottom>30or(bottom+tmp_bottom)/2-bottom<-30,本发明实际上以轮廓高度的偏差来选择是否填充此轮廓,此处本发明选择的偏差幅度为30,在实际应用中也可根据实际情况做出调整。
5.3对此轮廓向内填充像素,更新轮廓高度,更新公式如下,最后返回5.1:
bottom=(bottom+tmp_bottom)/2。
按上述方法遍历所有轮廓,获得三维初步分割结果图像。
最终轮廓填充结果参见图4,利用轮廓的高度信息逐个填充轮廓的策略可以有效减少对肺部以外区域轮廓的错误填充,为后续的图像形态学处理降低了复杂度。
实施例7
基于水平集的CT胸片肺组织图像分割方法同实施例1-6,由于最终的填充结果依然会存在一些小的毛刺,肺部区域内部可能会存在一些细小的空洞,肺部区域外部可能会存在一些小的孤立区域,所以这里选择使用三维图像形态学操作来优化最终的结果。特别地,考虑到三维图像分割的特点,图像形态学开操作、闭操作选择使用球形结构元作为基本形态学结构元,结构元半径大小为p;最后为了处理那些图像形态学操作无法处理的孤立区域,选择使用对三维图像分割结果使用连通域分析,因为肺部区域一定是最大的连通区域且体积较大,所以移除体积小于体积阈值w的三维连通区域。由此可得到本发明最终的图像分割结果,参见图5。
通过图像形态学操作可以有效的去肺部区域内部的细小空洞和肺部区域外部的孤立区域;三维连通域分析利用了不同区域的连通信息,就肺部区域图像分割而言,只有一个目标连通区域,所以在连通域分析后移除体积较小的区域是一种合理且有效的优化图像分割结果的方法。
下面给出一个更加详细的例子,对本发明进一步说明
实施例8
基于水平集的CT胸片肺组织图像分割方法同实施例1-7。
步骤1输入原始三维CT图像,从上往下的Z方向取出每一层二维CT图像,对图像进行预处理,给整体图像灰度加上400。然后将CT图像像素灰度值大小归一化到[0,255]。
步骤2基于当前图像的尺度,构建水平集能量泛函方程,并进行设定初始零水平集。能量泛函方程为如下形式:
F(φ,c,b)=ε(φ,c,b)+υL(φ)+μRp(φ)
其中ε(φ,c,b)是能量函数,φ为水平集函数,c为水平集内外的灰度平均值向量,b为偏置场参数,零水平集弧长函数L(φ)为正则函数,Rp(φ)为能量密度惩罚函数,υ零水平集弧长函数正则项,μ是能量密度惩罚函数正则项,能量泛函方程的能量函数具体表示为
Figure GDA0001838614490000101
其中N代表待分割区域数,Mi是成员函数,其中Mi成员函数具体定义为以下形式:
M1(φ(x))=H(φ(x)),M2(φ(x))=1-H(φ(x));
其中,赫维赛德函数H(φ(x))为以下形式:
Figure GDA0001838614490000102
区域函数ei(x)为以下形式:
ei(x)=I21K-2ci(b*K)+ci 2(b2*K);
I其中为待分割图像,ci代表当前水平集区域的灰度平均值,*是卷积算子,K是一个高斯卷积核,指示函数1K(x)为以下形式:
1K(x)=∫K(y-x)dy;
K(y-x)是一个非负的窗函数,满足K(y-x)=1,
Figure GDA0001838614490000111
Oy代表水平集区域y。
零水平集弧长函数L(φ)为如下形式:
Figure GDA0001838614490000112
其中,
Figure GDA0001838614490000113
为梯度算子。
能量密度惩罚函数Rp(φ)为以下形式:
Figure GDA0001838614490000114
其中,势函数p(s)为以下形式:
p(s)=(1/2)(s-1)2
初始零水平集的设定为:首先计算出CT图像的中心,然后分别对此中心向左和向右偏移50个像素得到两个初始零水平集的中心位置,以此作为两个正方形的中心,正方形的大小为50,最终将这两个正方形作为初始零水平集。
步骤3在初始零水平集的基础上最小化能量泛函方程,获取CT胸片图像零水平集轮廓图像。最小化能量泛函方程中采用梯度下降算法作为优化算法,包括有如下步骤:
3.1对能量泛函方程中的水平集函数参数φ更新,水平集函数参数φ的更新公式为以下形式:
Figure GDA0001838614490000115
其中,阶跃函数δ(φ)=H′(φ),距离函数
Figure GDA0001838614490000116
3.2对能量泛函方程中的灰度平均值向量参数c更新,参数灰度平均值向量c的更新公式为以下形式:
Figure GDA0001838614490000121
其中N=2,ui(y)=Mi(φ(y))。
3.3对能量泛函方程中的偏置场参数b更新,偏置场参数b的更新公式为以下形式:
Figure GDA0001838614490000122
在其中J(1)与J(2)分别为:
Figure GDA0001838614490000123
Figure GDA0001838614490000124
步骤4挑选出候选肺部区域轮廓。从图3可见,经过水平集方法提取出的轮廓线中含有大量肺部区域以外的轮廓。为了准确的提取肺部区域轮廓,本发明专门设计出轮廓筛选策略。具体方法是,从已排序的轮廓列表中的每个轮廓依次判断下列条件:
4.1若轮廓的外接矩形宽度小于CT图像宽度的一半则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.2轮廓的面积与轮廓的外接矩形的面积之比小于0.5则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.3若轮廓的外接矩形的长宽比小于4则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓。
步骤5根据候选肺部区域轮廓图像的高度逐个向内填充轮廓,获得三维初步分割结果图像,参见图4。具体方法是,确定第一个轮廓的外接矩形右下角点的高度,记为bottom,对候选肺部区域轮廓的每一个轮廓执行以下步骤:
5.1确定此轮廓的右下角点的高度,记为tmp_bottom;
5.2若满足以下条件则丢弃这个轮廓,返回5.1
(bottom+tmp_bottom)/2-bottom>30or(bottom+tmp_bottom)/2-bottom<-30;
5.3对此轮廓向内填充,更新轮廓高度,更新公式如下,最后返回5.1:
bottom=(bottom+tmp_bottom)/2
最终轮廓填充结果参见图4,利用轮廓的高度信息逐个填充轮廓的策略可以有效减少对肺部以外区域轮廓的错误填充,为后续的图像形态学处理降低了复杂度。
步骤6对三维初步分割结果图像依次进行三维图像形态学开操作,闭操作,选择使用球形结构元作为基本形态学结构元,结构元半径大小为3,然后对分割结果进行连通域分析,移除体积小于体积阈值250000的三维连通区域,得到最终三维分割结果图像,参见图5。
传统形态学CT胸片肺组织图像分割方案、阈值肺分割方案分割结果中包含大量的血管,胸腔组织等非肺部区域容易给后续结节检测带来影响。而本发明利用了肺部区域的边缘信息,通过最终的肺部区域轮廓向内填充来获得最终的图像分割结果,所以分割结果更加精确,基本不包含肺部以外区域,同时对初始图像分割结果进行了三维图像形态学处理,利用了图像的三维信息,使得图像分割边界更精确,进一步减少了肺部以外区域的错误分割结果。
下面通过仿真,对发明的技术效果再做说明
实施例9
仿真实验条件与内容
硬件平台为:处理器为Intel Xeon E5-2609v3@1.9GHz,内存为16.0G,硬盘1TB,操作系统为Ubuntu16.04;
软件平台:MATLAB2015b;
仿真内容与结果
在上述的实验条件下,进行仿真实验。优选地,零水平集弧长函数正则项为62.025,能量密度惩罚函数正则项为4;初始零水平集的设定为:首先计算出CT图像的中心,然后分别对此中心向左和向右偏移50个像素得到两个初始零水平集的中心位置,以此作为两个正方形的中心,正方形的大小为50;初始偏置场与待分割图像大小一致,偏置场所有值为1;高斯核大小为17;梯度下降算法迭代次数为50;面积之比筛选阈值参数r=0.5、外接矩形长宽比参数s=4;结构元半径p=3;体积阈值w=250000,实验结果参见图6。从图6可见本发明的肺组织图像3维分割结果很好地保留肺部区域的边界信息,分割结果精确,可以清晰的把肺部区域和人体其他组织分开。
为了验证本发明的鲁棒性,在以上参数的基础上,改变初始零水平集的位置和大小。具体为,首先计算出CT图像的中心,然后分别对此中心向左和向右偏移100个像素得到两个初始零水平集的中心位置,以此作为两个正方形的中心,正方形的大小为100,实验结果参见图7。对比图6和图7可得,图像分割结果没有差别,由此可以证明本发明的鲁棒性,图像分割结果不受初始区域的影响。
综上所述,本发明公开的基于水平集的CT胸片肺组织图像分割方法,主要解决了现有技术中对CT图像中肺组织区域分割边界不精确、图像分割过程中需要人工干预、分割结果不稳定的问题。其实现过程包括:步骤1获取包含肺组织的胸部CT图像,并进行CT胸片图像预处理;步骤2在预处理后的胸部CT图像上构造能量泛函方程并进行设定初始零水平集;步骤3在初始零水平集的基础上最小化能量泛函方程,获取CT胸片图像零水平集轮廓图像;步骤4挑选出候选肺部区域轮廓:从CT胸片零水平集轮廓图像中进行候选肺部区域轮廓挑选,首先计算出所有零水平集轮廓图像的外接矩形及轮廓的面积,按照轮廓的面积从大到小排序,得到已排序的轮廓列表,计算出轮廓的面积与轮廓的外接矩形的面积之比,然后计算出轮廓的外接矩形的宽度与外接矩形的高度,最后计算出轮廓的外接矩形的长宽比,最终根据轮廓的面积与轮廓的外接矩形的面积之比、轮廓的外接矩形的高度、轮廓的外接矩形的宽度、轮廓的外接矩形的长宽比零水平集轮廓图像中挑选出候选肺部区域轮廓图像;步骤5根据候选肺部区域轮廓图像的高度逐个向内填充轮廓,获得三维初步分割结果图像;步骤6对三维初步分割结果图像依次进行三维图像形态学开操作,闭操作,选择使用球形结构元作为基本形态学结构元,结构元半径大小为3,然后对分割结果进行连通域分析,移除体积小于预设阈值的三维连通区域,得到最终三维分割结果图像。本发明利用水平集方法提取了图像的边缘信息,基于先验知识设计了稳定的轮廓筛选策略,有效的筛选出候选肺部区域轮廓,最后基于轮廓高度信息设计了轮廓填充的优化方案。本发明图像分割结果鲁棒,精度高,是一种全自动的图像分割方法。本发明准确提取出了CT图像肺部区域,可用于后续对CT图像肺部区域分析。

Claims (7)

1.一种基于水平集的CT胸片肺组织图像分割方法,其特征是,包括有如下步骤:
步骤1获取包含肺组织的胸部CT图像,并进行CT胸片图像预处理;
步骤2在预处理后的胸部CT图像上构造能量泛函方程并进行设定初始零水平集;
步骤3在初始零水平集的基础上最小化能量泛函方程,获取CT胸片图像零水平集轮廓图像;
步骤4挑选出候选肺部区域轮廓:从CT胸片零水平集轮廓图像中进行候选肺部区域轮廓挑选,首先计算出所有零水平集轮廓图像的外接矩形及轮廓的面积,按照轮廓的面积从大到小排序,得到已排序的轮廓列表,计算出轮廓的面积与轮廓的外接矩形的面积之比,然后计算出轮廓的外接矩形的宽度与外接矩形的高度,最后计算出轮廓的外接矩形的长宽比,最终根据轮廓的面积与轮廓的外接矩形的面积之比、轮廓的外接矩形的高度、轮廓的外接矩形的宽度、轮廓的外接矩形的长宽比零水平集轮廓图像中挑选出候选肺部区域轮廓图像;
步骤5根据候选肺部区域轮廓图像的高度逐个向内填充轮廓,获得三维初步分割结果图像;
步骤6对三维初步分割结果图像依次进行三维图像形态学开操作,闭操作,选择使用球形结构元作为基本形态学结构元,结构元半径大小为q,然后对分割结果进行连通域分析,移除体积小于体积阈值w的三维连通区域,得到最终三维分割结果图像。
2.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤1中所述的图像预处理是给CT图像灰度值加400,然后将图像灰度值归一化到8位无符号数。
3.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤2中所述的在预处理后的CT胸片肺组织图像上构造能量泛函方程并进行设定初始零水平集,具体包括有如下步骤:
2.1构造能量泛函方程F(φ,c,b)
F(φ,c,b)=ε(φ,c,b)+υL(φ)+μRp(φ)
其中,φ为水平集函数,c为水平集内外的灰度平均值向量,b为偏置场,L(φ)为零水平集弧长函数,Rp(φ)为能量密度惩罚函数,υ零水平集弧长函数正则项,μ是能量密度惩罚函数正则项;
具体地,
Figure FDA0003267843100000021
其中,ε(φ,c,b)是能量函数,ei(x)是区域函数,N代表待分割区域数,Mi是成员函数,为以下形式:
M1(φ(x))=H(φ(x)),M2(φ(x))=1-H(φ(x));
其中,赫维赛德函数H(φ(x))为以下形式:
Figure FDA0003267843100000022
ei(x)为以下形式:
ei(x)=I21K-2ci(b*K)+ci 2(b2*K);
其中I为待分割图像,ci代表当前水平集区域的灰度平均值,*是卷积算子,K是一个高斯卷积核,1K(x)为指示函数,表达为以下形式:
1K(x)=∫K(y-x)dy;
K(y-x)是一个非负的窗函数,满足K(y-x)=1,
Figure FDA0003267843100000023
Oy代表水平集区域y;
零水平集弧长函数L(φ)为如下形式:
Figure FDA0003267843100000024
其中,
Figure FDA0003267843100000025
为梯度算子;
能量密度惩罚函数Rp(φ)为以下形式:
Figure FDA0003267843100000031
其中,势函数p(s)为以下形式:
p(s)=(1/2)(s-1)2
式中,s的值为水平集函数梯度的绝对值;
2.2初始零水平集的设定为:首先计算出CT图像的中心,然后分别对此中心向左和向右偏移m个像素得到两个初始零水平集的中心位置,以此作为两个正方形的中心,正方形的大小为n,最终将这两个正方形作为初始零水平集。
4.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤3中所述的在初始零水平集的基础上最小化能量泛函方程,采用梯度下降算法作为优化算法,包括有如下步骤:
3.1对能量泛函方程中的参数水平集函数φ更新,参数水平集函数φ的更新公式为以下形式:
Figure FDA0003267843100000032
其中,δ(φ)为阶跃函数,δ(φ)=H′(φ),dp(s)是距离函数,
Figure FDA0003267843100000033
3.2对能量泛函方程中的参数灰度平均值向量c更新,参数灰度平均值向量c的更新公式为以下形式:
Figure FDA0003267843100000034
其中N=2,ui(y)=Mi(φ(y));
3.3对能量泛函方程中的参数偏置场b更新,参数偏置场b的更新公式为以下形式:
Figure FDA0003267843100000035
其中,K为非负的窗函数,J(1)与J(2)分别为:
Figure FDA0003267843100000036
Figure FDA0003267843100000041
5.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤4中所述的挑选出候选肺部区域轮廓的具体方法是,从已排序的轮廓列表中的每个轮廓依次判断下列条件:
4.1若轮廓的外接矩形宽度小于CT图像宽度的一半则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.2轮廓的面积与轮廓的外接矩形的面积之比小于0.5则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓;
4.3若轮廓的外接矩形的长宽比小于4则认为它是候选肺部区域轮廓,否则在轮廓列表中去除这个轮廓。
6.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤5中所述的根据候选肺部区域轮廓的高度逐个向内填充的具体方法是,确定第一个轮廓的外接矩形右下角点的高度,记为bottom,对候选肺部区域轮廓的每一个轮廓执行以下步骤:
5.1确定此轮廓的右下角点的高度,记为tmp_bottom;
5.2若满足以下条件则丢弃这个轮廓,返回5.1,否则执行步骤5.3;
(bottom+tmp_bottom)/2-bottom>30or(bottom+tmp_bottom)/2-bottom<-30;
5.3对此轮廓向内填充,更新轮廓高度,更新公式如下:
bottom=(bottom+tmp_bottom)/2
按上述方法遍历所有轮廓,获得三维初步分割结果图像。
7.根据权利要求1所述的一种基于水平集的CT胸片肺组织图像分割方法,其特征在于,步骤6中所述的三维图像形态学开操作、闭操作选择使用球形结构元作为基本形态学结构元,结构元半径大小为3,移除体积小于体积阈值250000的三维连通区域。
CN201810819746.3A 2018-07-24 2018-07-24 一种基于水平集的ct胸片肺组织图像分割方法 Active CN109064476B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810819746.3A CN109064476B (zh) 2018-07-24 2018-07-24 一种基于水平集的ct胸片肺组织图像分割方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810819746.3A CN109064476B (zh) 2018-07-24 2018-07-24 一种基于水平集的ct胸片肺组织图像分割方法

Publications (2)

Publication Number Publication Date
CN109064476A CN109064476A (zh) 2018-12-21
CN109064476B true CN109064476B (zh) 2022-03-04

Family

ID=64835168

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810819746.3A Active CN109064476B (zh) 2018-07-24 2018-07-24 一种基于水平集的ct胸片肺组织图像分割方法

Country Status (1)

Country Link
CN (1) CN109064476B (zh)

Families Citing this family (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109886977B (zh) * 2019-02-19 2020-01-24 闽南师范大学 一种带有邻域约束的图像分割方法、终端设备及存储介质
CN110084824B (zh) * 2019-04-26 2020-03-27 山东财经大学 基于对称水平集的舌体图像分割方法、系统、设备及介质
CN110706236B (zh) * 2019-09-03 2023-01-13 西人马大周(深圳)医疗科技有限公司 血管图像的三维重建方法及装置
CN111429467B (zh) * 2019-10-11 2021-04-06 华中科技大学 改进Lee-Seo模型的水平集三维表面特征分割方法
KR20220074927A (ko) * 2019-10-31 2022-06-03 칼 짜이스 에스엠테 게엠베하 고형상비 구조의 형상 편차를 측정하기 위한 fib-sem 3d 단층 촬영
CN111145185B (zh) * 2019-12-17 2023-12-22 天津市肿瘤医院 一种基于聚类关键帧提取ct图像的肺实质分割方法
CN111627037B (zh) * 2020-05-22 2024-10-29 深圳前海微众银行股份有限公司 图像区域提取方法、装置、设备及存储介质
CN111898600A (zh) * 2020-07-10 2020-11-06 浙江大华技术股份有限公司 字符轮廓提取方法及装置、存储介质、电子装置
CN112102339B (zh) * 2020-09-21 2023-07-07 四川大学 一种基于图集配准的全身骨显像骨骼分割方法
CN112669312A (zh) * 2021-01-12 2021-04-16 中国计量大学 一种基于深度特征对称融合的胸片肺炎检测方法及系统
CN114708277B (zh) * 2022-03-31 2023-08-01 安徽鲲隆康鑫医疗科技有限公司 超声视频图像活动区域自动检索方法和装置
CN115861600B (zh) * 2022-12-20 2023-09-05 西北民族大学 一种spect图像的roi区域识别方法及系统

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102354396A (zh) * 2011-09-23 2012-02-15 清华大学深圳研究生院 基于水平集函数的灰度不均匀图像分割方法
CN104867143A (zh) * 2015-05-15 2015-08-26 东华理工大学 基于局部引导核拟合能量模型的水平集图像分割方法
CN106570867A (zh) * 2016-10-18 2017-04-19 浙江大学 基于灰度形态学能量法的活动轮廓模型图像快速分割方法
CN106997596A (zh) * 2017-04-01 2017-08-01 太原理工大学 一种基于信息熵和联合向量的lbf活动轮廓模型的肺结节分割方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106529555B (zh) * 2016-11-04 2019-12-06 四川大学 一种基于全卷积网络的dr片肺轮廓提取方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102354396A (zh) * 2011-09-23 2012-02-15 清华大学深圳研究生院 基于水平集函数的灰度不均匀图像分割方法
CN104867143A (zh) * 2015-05-15 2015-08-26 东华理工大学 基于局部引导核拟合能量模型的水平集图像分割方法
CN106570867A (zh) * 2016-10-18 2017-04-19 浙江大学 基于灰度形态学能量法的活动轮廓模型图像快速分割方法
CN106997596A (zh) * 2017-04-01 2017-08-01 太原理工大学 一种基于信息熵和联合向量的lbf活动轮廓模型的肺结节分割方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
"Feasibility of automated 3-dimensional magnetic resonance imaging pancreas segmentation";Gou Shuiping 等;《Advances in Radiation Oncology》;20160930;第182-193页 *
"结合全局和局部信息的水平集图像分割方法";刘晨 等;《计算机应用研究》;20171231;第34卷(第12期);第3889-3894页 *
"超声图像乳腺肿瘤分割新方法研究";高梁;《中国优秀博硕士学位论文全文数据库(博士) 医药卫生科技辑》;20140515(第(2014)05期);第E072-131页 *

Also Published As

Publication number Publication date
CN109064476A (zh) 2018-12-21

Similar Documents

Publication Publication Date Title
CN109064476B (zh) 一种基于水平集的ct胸片肺组织图像分割方法
Zhang et al. Fast segmentation of bone in CT images using 3D adaptive thresholding
CN107045721B (zh) 一种从胸部ct图像中提取肺血管的方法及装置
US7315639B2 (en) Method of lung lobe segmentation and computer system
CN107563998B (zh) 医学图像中心脏图像处理方法
CN107527341B (zh) 血管造影图像的处理方法和系统
Pulagam et al. Automated lung segmentation from HRCT scans with diffuse parenchymal lung diseases
KR20150045885A (ko) 초음파 및 ct 영상의 정합에 관한 시스템 및 작동 방법
Chen et al. Pathological lung segmentation in chest CT images based on improved random walker
CN107169975B (zh) 超声图像的分析方法及装置
JP5055115B2 (ja) 識別方法、コンピュータプログラム及びコンピュータプログラム装置
CN116071355A (zh) 一种用于外周血管影像的辅助分割系统及方法
Kaftan et al. Fuzzy pulmonary vessel segmentation in contrast enhanced CT data
CN117237342B (zh) 一种呼吸康复ct影像智能分析方法
CN112712540B (zh) 一种基于ct影像的肺部支气管提取方法
KR101474162B1 (ko) 흉부 컴퓨터 단층 영상의 간유리음영 결절 자동 분할 시스템 및 그 분할 방법
EP1447772B1 (en) A method of lung lobe segmentation and computer system
KR101484051B1 (ko) 능동적 윤곽선 방법을 이용한 cad 시스템을 위한 전처리 장치 및 방법
CN110533667A (zh) 基于图像金字塔融合的肺部肿瘤ct影像3d分割方法
Noviana et al. Axial segmentation of lungs CT scan images using canny method and morphological operation
CN115661020A (zh) 一种dr胸腔图像肋骨检测的方法及装置
JP2020049140A (ja) 肺門部陰影を抑制した肺結節明瞭化法
CN114882060B (zh) 基于水平集的肺野图像分割方法和装置
Wang et al. Liver contour extraction using modified snake with morphological multiscale gradients
Trujillo et al. Segmentation of the Aorta Using Active Contours with Histogram-Based

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
GR01 Patent grant
GR01 Patent grant