谢勤岚,潘先攀
(中南民族大学 生物医学工程学院,武汉 430074)
人体组织器官,如肝脏、肾脏和脾的位置、体积以及形态的变化与多种疾病息息相关,分析这些器官的位置、体积与形态可为多种疾病的临床研究提供相关支持[1].对这些组织器官的精准分割是进行相关定量分析的前提,因此,寻找一种有效的分割算法来实现腹部CT图像的分割显得极其重要.
针对医学图像分割问题,研究人员提出了许多分割方法.赵于前[2]在提取肝脏初始轮廓的基础上,采用5个不同方向的结构元素检测不同方向的边缘,并结合不同尺度的结构元素对梯度图像进行加权平均,实现肝脏边缘的准确检测;何兰[3]等人通过3D卷积神经网络对数据集训练提取肝脏,以此为形状先验结合Graph cuts算法分割出肝脏,但是该算法对训练数据的数量需求比较大,常常会因为缺少足够的训练数据而造成难优化的问题;潘晓佩[4]提出一种基于Snake balloon模型和Canny算子的融合算法,将Canny算子提取的边界结果作为Snake balloon模型在超声肾脏图像分割中的收敛信息,使Snake balloon模型的迭代轮廓停在肾脏图像的弱边界上,有效分割出肾脏;Jayaram[5]等人通过定向主动外观模型提取目标器官的大致形状,将该形状作为形状先验融入图割算法,实现肝脏、肾脏和脾脏的分割,但主动外观模型对目标特征点标定的准确度要求比较高.
本文根据腹部CT图像的成像特点,提出了一种自适应形状约束的Graph cuts分割算法,实现从腹部CT图像中提取肝脏、肾脏以及脾等组织器官的目标,其方法流程如图1.首先使用基于多图谱配准的分割方法处理待分割图像,得到目标区域的初始分割结果,将它作为待分割区域的形状先验;然后由形状先验的轮廓计算符号距离图,将其通过形状比较函数加入到Graph cuts能量函数中;接着根据配准分割过程中产生的目标概率图选择能量函数中每条边的形状约束项系数,实现自适应形状约束;最后使用最大流最小割算法[6]优化Graph cuts能量函数,得到分割结果.实验表明,该算法可以从对比度低,灰度不均,组织边界模糊,背景复杂的CT图像中较为准确高效地分割出目标组织器官.
图1 方法流程图Fig.1 Method flow chart
Graph cuts是Boykov[7]等人提出的一种交互式图像分割算法,该算法将图像的分割问题转化为最小化能量泛函的问题.Graph cuts能量函数如方程(1):
E(A)=λR(A)+B(A),
(1)
其中,R(A)表示区域项能量,B(A)表示边界项能量.A=(A1,A2,…,AP)表示一个二值向量,AP表示给像素点P分配的标记,当P标记为目标,AP值为1,当P标记为背景,则AP值为0.λ是一个非负常数,用于平衡边界项B与区域项R之间的相对重要性关系.方程(1)可以进一步写成:
(2)
其中,N为像素P的领域点的集合.区域项能量Rp(Ap)表示给像素P分配“前景”或“背景” 标记所产生惩罚.
方程(2)中,定义区域项Rp(Ap)可以使用负对数似然性模型[7],如式(3)和(4):
Rp(o)=-lnp(Ap|o),
(3)
Rp(b)=-lnp(Ap|b),
(4)
式中o代表目标标记,b代表背景标记.若P标记为目标,则
(5)
其中,σ表示人工交互过程中标记的所有目标像素的标准差,μ表示标记的所有目标像素的均值,给P标记为背景则同理.
Derive the transfer function of the integrator in the z-domain as:
方程(2)中,边界项能量定义为:
(6)
dist(p,q)表示像素点p和q之间的欧氏距离,σ是一个非负常数,Ip、Iq分别代表p、q点的像素值.该边界项能量表示给像素p和q分配标签所产生的惩罚.
本文采用基于多图谱配准的分割方法[8]获取目标器官的形状先验.基于多图谱配准的分割方法,是指将多个图谱分别与待分割图像进行配准,获取多组形变参数,分别对各图谱对应的标记图像进行形变,再采用合适的图像融合技术,对所有形变后的标记图像进行融合,得到分割的结果.
图2 基于多图谱配准的分割流程图Fig.2 Segmentation flow chart based on multi-atlas registration
为在Graph cuts算法中加入目标器官的形状先验,本文将配准分割结果的轮廓作为零水平集,计算符号距离图,将其通过形状比较函数加入到能量方程(2)中.
符号距离函数(Signed Distance Function,SDF)定义了空间中任意一个点到特定区域边界的距离,该边界对应于水平集函数中的零水平集.符号距离函数同时对距离的符号进行定义:点在区域边界内部为正,外部为负,位于边界上时为0.假设d是空间X的一种度量,那么SDF用数学公式表达为:
(7)
式中(x)为符号距离函数,Ω为空间中的某个区域,∂Ω是区域Ω的边界,d(x, ∂Ω)为点x到区域边界∂Ω的最短距离.
本文符号距离函数的计算使用了C.R.Maurer提出的线性时间算法[9].在算法中我们使用无符号距离函数,故在所得的符号距离图中,对每个点的符号距离取绝对值.如图3所示,(a)为多图谱配准分割结果,(b)为由(a)计算的符号距离图像,(c)为无符号距离图像.
图3 符号距离图的计算Fig.3 Calculation of signed distance map
在Graph cuts算法中加入形状约束后,能量方程(2)改可改写为:
(8)
其中Sp,q表示形状约束项能量,ω是衡量形状约束项能量的在能量函数中比重的常数.
在本文方法中,我们引用Freedman和Zhang[10]的方法,将形状先验以水平集的形式结合到Graph cuts算法中.该方法中,
(9)
方程(8)中形状约束项的权重系数ω可写为:
ωp,q=e-(αp-αq)2,
(10)
其中α为目标概率图,αp、αq分别表示p、q点属于目标的概率,即α中对应p、q点的值.
将式(9)、(10)带入方程(8)可得新的能量函数方程:
(11)
通过使用最大流最小割算法优化能量方程(11),便可得到最终分割结果.
为验证本文提出的算法的可靠性及有效性,我们将分割方法在3Dircadb数据集上进行实验,依次分割出肝脏、肾脏和脾脏,并将分割结果与专家手动分割的金标准比较.3Dircadb数据集包含了20组病人的腹部CT图像,每个病人的影像切片大小为512×512像素,层数为91~260层.本文算法代码的实现基于ITK( Insight Segmentation and Registration Toolkit)医学图像处理软件包.
本文使用国际通用的Dice系数(Dice coefficient)、体积相似度(Volume similarity)和假阳性误差(False Positive Error)[12]三种评价指标对分割结果进行评估,其计算公式分别如下:
(12)
(13)
(14)
其中,S表示分割结果,T表示专家手动分割的金标准,ST为分割结果减去金标准与分割结果重叠的部分,│·│表示体素.三种评价指标的中,Dice系数越大,表示分割结果越接近金标准,即分割结果越准确,相反,体积相似度和假阳性误差越小,说明分割结果越准确.为方便比较不同方法分割结果的VS均值,本文对体积相似度VS取绝对值.
图4展示了部分实验结果.从上至下,第一行为3Dircadb数据集中第3个病人的肝脏分割,第二行为第3个病人的肾脏分割,第三行为第11个病人的脾脏分割.从左到右,(a)列为待分割的原始图像,红色区域为目标区域,(b)列为传统图割算法的分割结果,(c)列为多图谱配准分割的结果, (d)列为本文分割方法的分割结果.
从图4 (b)列的分割结果来看,当目标轮廓不规则或者背景比较复杂时,传统的Graph cuts算法容易出现欠分割或者过分割现象,导致这种结果的原因可以从Graph cuts算法的能量方程中分析:当边界项在能量方程中占比比较高,而两个不连通的区域同属于一个目标时,分割结果可能不同时包含两个区域,这会导致欠分割,如(b)列中的肝脏分割;当区域项在能量方程中占比比较高,而图像中其他非目标区域的灰度与目标区域相似时,分割结果会包含这一非目标区域,这会造成过分割,如(b)列中的肾脏分割和脾脏分割.在(c)列分割结果中,基于多图谱配准的分割方法虽然能大致分割出目标区域,但是仍然不能避免过分割现象.从(d)列的分割结果来看,本文方法通过加入形状先验约束分割结果后,欠分割和过分割现象得到了很大的改善,且相比配准分割结果,本文方法的分割结果在轮廓上保留了更多细节.
为更直观地展示本文方法的分割效果,在原始图像上叠加了金标准的轮廓和本文方法分割结果的轮廓,如图5.从上至下,第1到第3行分别为肝脏、肾脏和脾脏的分割结果,从左到右分别为相应分割结果的横截面、矢状面和冠状面.绿色线表示金标准的轮廓,红色线表示本文方法分割结果的轮廓,红色线越接近绿色线,分割的结果越准确.结果表明,本文方法分割结果的轮廓与金标准的轮廓重叠程度很高,说明本文提出的分割方法是有效的.
图4 实验结果展示1Fig.4 Experimental results 1
图5 实验结果2Fig.5 Experimental results 2
为分析形状约束在Graph cuts算法中的抗噪声干扰能力,实验分别使用传统Graph cuts算法和本文方法对具有低对比度和高噪声特征的合成图像分割,并对实验结果对比分析.实验结果如图6所示,第一行为3Dircadb数据集中第9个病人的肝脏分割,第二行为第6个病人的肾脏分割,第三行为第11个病人的脾脏分割. (a)列为低对比度的待分割图像, (b)列为(a)中加入高斯噪声后的图像,(c)列为传统Graph cuts算法对(b)的分割结果, (d)列为本文分割方法对(b)的分割结果.图中绿色线表示目标区域轮廓,红色线表示算法分割结果的轮廓.从实验结果可以看出,由于图像对比度低,传统Graph cuts算法分割结果的错分现象比较严重,尤其是脾脏的分割,灰度相似的背景区域会被归入目标区域;图像中的噪声导致传统Graph cuts算法分割结果产生许多空洞,而且背景像素去除的不是很干净.相比之下,本文方法能较好地解决以上问题,可以得到相对理想的分割效果.
图6 实验结果3Fig.6 Experimental results 3
实验对比了5种不同分割方法的分割效果,每种方法都对3Dircadb数据集中的20幅CT图像(pat1-pat20)分割,并统计每种分割方法的分割精度均值,统计结果如表1.Graph cuts为传统Graph cuts算法分割,Level Set[13]为水平集方法分割,MARS (Multi Atlas Registration Segmentation)为基于多图谱配准的分割,LS-GC为使用水平集分割结果作为形状先验的Graph cuts分割(注:形状约束项系数由模糊处理后的待分割图像自适应选择),MARS -GC为本文分割方法,即使用多图谱配准分割结果作为形状先验的Graph cuts分割.从统计结果来看,不论是肝脏分割、肾脏分割,还是脾脏分割,相比其他几种分割方法,本文分割方法的分割精度更高,错误率更低,分割效果更好.
表1 不同分割方法的分割精度均值统计Tab.1 Segmentation accuracy mean statistics of different segmentation methods
本文根据腹部CT图像具有对比度低,灰度不均,低信噪比,不同组织之间或软组织与病灶之间边界模糊等特点,提出了一种自适应形状约束的Graph cuts算法,采用基于多图谱配准的分割方法分割原图像得到初始分割结果,将初始分割结果作为形状先验加入到图割算法中,约束分割结果的形状.实验结果表明,该方法能够很好地解决传统Graph cuts算法分割时出现的欠分割和过分割问题,且相比水平集分割和基于多图谱配准分割,本文分割方法能够得到更好的分割结果和更高的分割精度.