医学图像的分割

第六章医学图像分割
医学图像分割是医学图像处理和分析的关键步骤，也是其它高级医学图像分析和解释系统的核心组成部分。

医学图像的分割为目标分离、特征提取和参数的定量测量提供了基础和前提条件，使得更高层的医学图像理解和诊断成为可能。

本章首先对医学图像分割的意义、概念、分类及其研究现状进行了概述,然后分别对基于阈值、基于边缘、基于区域和基于模式识别原理的各种常见医学图像分割方法作了详尽而系统的介绍，接着在对图像分割过程中经常用到的二值图像数学形态学基本运算作了简单叙述之后，较为详细地讨论了医学图像分割效果和分割算法性能的常用评价方法。

第一节医学图像分割的意义、概念、分类和研究现状
医学图像分割在医学研究、临床诊断、病理分析、手术计划、影像信息处理、计算机辅助手术等医学研究与实践领域中有着广泛的应用和研究价值，具体表现为以下几个方面：(1) 用于感兴趣区域提取，便于医学图像的分析和识别。

如不同形式或来源的医学图像配准与融合，解剖结构的定量度量、细胞的识别与计数、器官的运动跟踪及同步等；(2)用于人体器官、组织或病灶的尺寸、体积或容积的测量。

在治疗前后进行相关影像学指标的定量测量和分析，将有助于医生诊断、随访或修订对病人的治疗方案； (3)用于医学图像的三维重建和可视化。

这有助于外科手术方案的制定和仿真、解剖教学参考及放疗计划中的三维定位等；(4)用于在保持关键信息的前提下进行数据压缩和传输。

这在远程医疗中对实现医学图像的高效传输具有重要的价值；(5)用于基于内容的医学图像数据库检索研究。

通过建立医学图像数据库，可对医学图像数据进行语义学意义上的存取和查找。

所谓医学图像分割，就是根据医学图像的某种相似性特征（如亮度、颜色、纹理、面积、形状、位置、局部统计特征或频谱特征等）将医学图像划分为若干个互不相交的“连通”的区域的过程，相关特征在同一区域内表现出一致性或相似性，而在不同区域间表现出明显的不同，也就是说在区域边界上的像素存在某种不连续性。

一般说来，有意义的图像分割结果中至少存在一个包含感兴趣目标的区域。

区域(Region)作为图像分割中像素的连通集合和基本分割单位，可以按照不同的连通性来定义：4连通区域和8连通区域。

区域的连通性是指在一个区域中任意两个像素之间，都存在一条完全由属于这个区域的元素所构成的连通路径。

如果只依据处于四正位(上、下、左、右)或四角位(左上、左下、右上、右下)的相邻像素确定区域的连通性，就称为4连通；如果同时依据处于四正位和四角位相邻的像素确定区域的连通性则称为8连通。

在数学上，医学图像分割可以用集合论模型予以描述：已知一幅医学图像I 和一组相似性约束条件
i C (Λ,2,1=i )，对I 的分割就是求取它的一个划分的过
程，即: N
j 1=Y j R =I , ],1[,,,N k j k j R R k j ∈≠∀=φI （6.1）
其中， j R 为同时满足所有相似性约束条件i C (Λ,2,1=i )的连通像素点的集合，即我们所谓的图像区域；N 为不小于2的正整数,表示分割后区域的个数。

在如上集合论模型描述中，如果保持区域连通性的约束被取消，那么对图像所属像素集的划分就称为分类(Pixel classification)，其中每一个像素集合称为一类(Class)。

在本章后面的讨论中，为了描述上的方便，我们往往不加区分地将经典的区域分割和像素分类统称为图像分割。

通常，医学图像分割方法可以划归为三大类：基于阈值的分割方法、基于边缘的分割方法和基于区域的分割方法。

在理想情况下，医学图像中的每一个区域都是由相应的封闭轮廓线包围着。

原则上，使用边界跟踪算法可以得到区域的边缘（或封闭的轮廓线）；反过来，使用简单的区域填充算法也可以得到边缘所包围的区域。

但在实际的医学图像中，很少能够从区域中得到连续、封闭的边缘，反之亦然。

由于受人体内外环境中种种确定性、不确定性因素的干扰和成像噪声的影响，实际所获得的医学图像不可避免具有模糊、不均匀等缺陷；另外，人体的解剖结构比较复杂而且因个体的病理或生理差异有很大的不确定性，这在医学图像中引入了新的复杂性，同时也给医学图像分割带来了很大的困难；还有，现有医学图像分割的基本方法大多数是针对2D 图像进行的，当推广到3D 乃至4D 医学图像分割应用场合时，在数据结构和算法处理上不可避免导致更大的复杂性，使得医学图像的分割更为困难。

近年来大量学者致力于将新概念、新思想和新方法
应用于复杂二维医学图像和高维医学图像或者图像序列的分割，其中包括数学形态学、模糊理论、神经网络、遗传算法、蚁群算法、粗糙集理论、水平集理论、支持向量机、马尔科夫随机场理论、核函数方法、小波分析和小波变换等，其间有很大一部分属于基于模式识别原理的医学图像分割方法。

各种分割方法或数学工具的有效应用，极大地改善了医学图像的分割效果。

纵观医学图像分割技术的发展历程，实际上是一个从人工分割到半自动分割和全自动分割逐步发展的过程。

人工分割是指由经验丰富的临床医生在原始胶片图像上直接勾画出组织的边界，或者通过图像编辑器用鼠标在计算机显示器上勾画出组织的边界或感兴趣区域。

半自动分割技术是随着计算机科学的发展而产生的，它把计算机强大的数据处理、算法分析和智能计算能力以及自动存贮和记忆功能与医学专家的知识和经验有机地结合起来，通过人机交互的方式完成图像分割。

全自动分割则彻底脱离了人为干预，完全由计算机自动完成图像分割的全过程。

由于全自动分割方法不存在人为因素的影响，为图像中感兴趣区域的自动精确测量奠定了基础。

但是绝大多数自动分割算法实现复杂，运算量较大，在很多情况下，分割结果尚不理想，分割速度和性能也有待提高。

从目前图像分割技术在临床上的应用情况来看，人工分割的精度在所有分割方法中是最高的，被视为金标准，但该方法费时、费力，其分割结果的优劣完全取决于操作者的经验和知识，且分割结果难以重现。

半自动方法与人工分割相比，分割速度有明显提高，但其分割结果在很大程度上仍然依赖于操作者的主观经验和知识，这在一定程度上影响了半自动分割技术在临床上的广泛应用。

研究高效、实用的全自动分割方法并最终取代繁琐的人工分割和主观性依然很强的半自动分割一直是人们追求的目标，也是近年来图像分割方法的研究重点。

由于医学图像经常存在对比度低、组织特征的可变性大以及不同软组织之间或软组织与病灶之间边界模糊、微细血管或神经结构形状分布复杂等特点，加上成像中种种客观因素的制约，使得医学图像分割成为医学图像处理过程中的一个经典难题，到目前为止既不存在一种通用的图像分割方法，也不存在一种准确评价分割成功与否的客观标准[1]。

第二节基于阈值的图像分割方法
阈值法是一种常用的将图像中感兴趣目标与图像背景进行分离的图像分割方法，常用的阈值法基本上都是基于一维灰度直方图统计特征的分割方法(见图
6-1(2))，它简单地用一个或几个阈值将图像灰度直方图分成两段或多段，而把图像中灰度值在同一段内的所有像素归属为同一个物体。

很明显，这类分割方法基于如下前提与假设：对应于特定物体或背景的像素灰度呈现峰状分布特征并且基本集中于不同的灰度区间内。

设原始图像为),(y x I ，阈值法按照一定的规则在),(y x I 中确定若干个门限值,,,,21N T T T Λ其中1≥N ，利用这些门限值将图像分割为几个部分。

分割后的结果图像可表述为
⎪⎪⎪⎩⎪⎪⎪⎨⎧<<≤<≤≥=--10
211
11),(),(),(),(),(T y x I L T y x I T L T y x I T L T y x I L y x R N N N N N 若若若若M
（6.2） 其中，N L L L ,,,10Λ为结果图像的灰阶。

如果2=N ，上述分割步骤也被称为基于阈值法的医学图像二值化过程。

需要注意的是，医学图像的二值化通常是指通过一定的方法使得医学图像上的所有点的灰度或颜色值只有两种取值。

图像二值化的目的主要就是为了简化图像，并使图像具有鲜明的对比性，以便于对图像进行后续处理。

实际上医学图像分割和医学图像的二值化是既有区别又有联系的两个概念，如式（6.2）在2>N 的情况下，图像分割的结果并不一定只具有两种颜色或灰度值，此时该步骤便只能被称为图像的多阈值分割过程，而不能称作图像的二值化。

传统的阈值分割法通常直接利用图像的一维灰度直方图进行分割，或基于图像本身的灰度或颜色分布特征，来确定一个或多个门限N T T T ,,,21Λ，把N 个门限将灰阶所对应数轴划分成1+N 个区间。

对于图像中的任意一个像素，如果它的值处于第i 个区间内，它就属于第i 类目标。

该类方法适用于内容不太复杂且灰度分布较为集中的图像。

理想情况下，从灰度直方图上根据经验直接选取合适的门限即可很好地区分开不同种类的组织（比如CT 图像中皮肤、血管、骨骼等硬、软组织的分割），但在绝大多情况下,简单阈值法并不能得到正确有效的分割。

其原因在于图像的一维直方图一般是图像中各像素灰度值出现频数的统计结果，它只反映出图像中不同灰度值出现的频率，而不能反映某一灰度值的像素所对应的
位置及其邻域特征, 因此它孤立地对每个像素进行运算而仅仅相关于图像中某点像素的信息，却忽略了其空间邻域信息,这使得该方法对于噪声和灰度不均匀性很敏感,此外对于各物体不存在明显灰度差异或各自灰度范围有较大重叠的图像而言，在灰度直方图中，可能没有显著的统计特征，如直方图呈现大量毛刺、没有明显多峰、谷底表现比较平坦等，这种情况下就非常不利于阈值的选取,此时如何确定最佳分割阈值将是一个极为关键和困难的问题。

如果阈值选取不合理，则会把一些本来不是目标的像素也当作目标，造成误识；或者把一些目标漏掉造成漏识。

由于在绝大多数情况下物体和背景的对比度在图像中各处不是完全一样的，很难用一个统一的阈值将物体与背景分开。

在阈值分割技术具体实施过程中，人们往往还需要通过控制阈值选取范围的方法实现局部分割阈值的选择，即将原始图像划分成较小的图像，并对每个小图像根据图像的局部特征分别采用不同的阈值进行分割。

实际处理时，既可按照具体问题将图像分成若干子区域分别选择阈值，也可动态地根据某点邻域内的图像特征选择每点处的阈值而实现图像分割，这往往被称为自适应阈值分割。

当然在图像中背景和目标具有明显对比度差异的情况下，直接对整幅图像采用统一的阈值做分割处理即可，这也是通常所说的全局阈值分割。

阈值法作为一种古老的图像分割方法，因其实现简单、计算量相对较小、性能较稳定已经成为图像分割中最基本和应用最广泛的分割技术之一。

在运用阈值法进行图像分割的过程中由于每个像素的所属类别只依赖于它的数值,能够实现并行的快速实时操作，这使得阈值法常常作为关键的预处理步骤被用在各种图像处理过程之中。

为了提高图像中感兴趣目标和背景的分割精度和效率，目前人们运用信息熵、最优化方法、模糊集合论、数学形态学、小波变换等数学理论或工具发展了各种各样的基于直方图统计特征的阈值选取和分割技术，本节只就p-分位数法、双峰法、迭代法、最大熵法、矩量保持法、大津法等常见的图像阈值分割方法进行详细的介绍。

一、 p-分位数法
p-分位数法(也称p-tile法)是最古老的阈值选取方法之一。

其基本思想是使医学图像中目标所占图像像素数的比例等于其先验概率p来设定阈值。

在很多医
学图像中像素属于目标的先验概率可以根据解剖结构中背景和目标的相对比例基于临床经验预先估计出来。

然后根据先验概率直接在图像直方图上找到合适的阈值，把大于阈值的像素作为目标，小于阈值的像素作为背景，最终实现医学图像的快速分割。

在病理条件下，图像中感兴趣目标和背景之间像素数目可能并不存在相对固定的比例关系，此时就不容易估计出目标的先验概率p ，也就不能用p -分位数法进行阈值分割。

实际上，即使在先验概率能够预先估算出来这一先决条件满足的情况下，也只有当图像背景和目标差别比较显著时候，p -分位数法才可能得到又快又好的分割效果；而对于目标比较多，目标和背景对比不是很明显的图像其分割效果就未必理想，甚至可能完全失效。

二、 双峰法
双峰法原理比较简单，它假设图像是由前景和背景组成，且灰度直方图呈现明显双峰结构：一个与目标相对应，另一个对应于背景。

通过在双峰之间的最低谷处选择阈值即可实现图像分割。

在一般情况下，由于各种噪声和图像细节的干扰，图像的直方图往往有很多毛刺，只可能具有明显的双峰趋势，并不呈现过渡光滑自然的双峰结构，此时需要首先对直方图曲线进行平滑处理，然后再进行阈值选取和分割。

常用的方法是采用特定窗口大小的一维均值滤波器对直方图中每一灰度级处直方图统计值进行修正。

该方法计算量不大，且能有效地平滑直方图曲线，使得直方图双峰特点更加突出和鲜明。

三、 迭代法
迭代法基于最优逼近的思想，通过迭代的过程选择一个最佳阈值，实现图像的分割。

其基本算法如下：
统计图像中各像素灰度的最大值和最小值，分别记为max G 和min G ，置迭代控
制变量0=k ，令阈值2/)(min max G G T k +=，并将其作为初始阈值；
根据阈值
k T 将图像分割为目标和背景区域，分别求出两区域的平均灰度值为0M 和1M ；
令迭代控制变量1+=k k ，求出新的迭代阈值2/)(101M M T k +=+；
若k k T T =+1，则终止迭代，1+k T 即为所得阈值，否则转（2）继续迭代。

迭代法本质上是一种简单的两均值聚类技术，每一次迭代对应的目标和背景区域的均值0M 、1M 分别为相应的聚类中心，而分割阈值k T 与0M 和1M 的距离分别为各自聚类半径，以两区域均值基本不变时的阈值作为迭代法最终分割阈值，完成图像分割。

对于直方图呈现双峰形状且峰谷特征比较明显的图像，迭代方法可以较快收敛到满意结果，此时迭代所得的阈值分割图像结果很好，能较好区分图像的前景和背景的主要区域，但是对于图像直方图双峰特征不明显，或目标和背景比例差异悬殊情况下分割效果可能不理想。

对某些特定图像，迭代过程中微小数据的变化甚至会引起分割结果的巨大变化，导致分割失效，这是由非线性迭代系统对初始条件的敏感性也即俗称的“蝴蝶效应”造成的。

四、 最大熵法
熵是信息论中的一个术语，是所研究对象平均信息量的表征。

其定义为： ⎰+∞∞--=dx x p x p H )(lg )(，其中)(x p 是随机变量x 的概率密度函数。

基于最大熵原则进行阈值选择一直是最重要的阈值分割方法之一。

这种方法的目的在于将图像的灰度直方图分成两个或多个独立的类，使得各类熵的总量最大。

从信息论角度来说就是使这样选择的阈值能获得的信息量最大。

根据最大熵原理进行图像阈值分割，人们通常的做法是选取一个阈值t ，使图像用这个阈值分割出的两部分的一阶灰度统计的信息量最大，即一维熵最大。

一维最大熵阈值图像分割法的基本方法如下：
统计图像中每一个灰度级出现的概率，计算该灰度级的熵，假设以灰度级t 分割图像，图像中低于t 灰度级的像素点构成目标物体0，高于灰度级t 的像素点构成背景B ，那么各个灰度级在本区的分布概率为：
O 区：t i p p t
i ,,2,1,Λ= B 区：1,,2,1,1-++=-L t t i p p t
i Λ （6.3） 上式中∑==t i i
t p p 0，这样对于数字图像中的目标和背景区域的熵分别为
∑=-=i t i t i o t i p p p p H ,
,,2,1),/lg()/(Λ (6.4)
∑-++=---=i t i t i B L t t i p p p p H 1
,,2,1)],1/(lg[)]1/([Λ (6.5)
对图像中的每一个灰度级分别求取B O H H w +=，选取w 最大的灰度级作为
分割图像的阈值。

由于基于熵的阈值选取法受到目标大小的影响很小，所以可以处理目标较小情况下的分割，但一维最大熵法由于涉及对数运算，速度较慢，实时性较差，对于细节较多、噪声较大的图像分割效果也不理想。

五、 矩量保持法
也叫保矩法或矩保持法(Moment-preserving) 。

其基本思想是使阈值分割前后图像的矩保持不变。

设每个像素处灰度}1,,2,1,0{-∈L g Λ ，其中L 为灰度级总数，假定i N （}1,,2,1,0{-∈L i Λ）表示第i 个灰度级在图像中出现的次数，则归一化直方图中第i 个灰度级出现的频率可表示为：
N N p i i /=，其中 N 为图像中像素总数。

图像的第k 阶矩
k m 定义为： Λ,2,1,11
00===∑-=k i p m m L i k i k (6.6)
图像阈值法实现图像二值化分割，意味着保持前3阶矩不变，即存在如下矩量保持方程组：
011000m Z p Z p =+ (6.7) 1111100m Z p Z p =+ (6.8)
2211200m Z p Z p =+ (6.9)
3311300m Z p Z p =+ (6.10)
这里0Z 和1Z ，表示二值化后每个类别对应灰度值，0p 和1p 代表二值化后两
个灰度值的分布概率，经简化后
2/10211
0)4(c c m d p --=， (6.11)
其中
])4[(21,,12/102121232112122310c c c d m m m m m c m m m m m c --=--=--=
然后从灰度直方图中按照-0p 分位数法选取阈值t ，即可实现图像阈值分
割。

矩量保持法算法原理简单，易于实现，可以很容易推广为多阈值分割，广泛应用在各类图像分割过程之中，关于其用于多阈值图像分割的实现方法，具体可参考有关文献。

六、 Ostu 方法
Ostu 法是利用方差来度量图像灰度分布均匀性进行分割的一种阈值化方法。

该方法在一定条件下不受图像对比度与亮度变化的影响，被认为是阈值自动选取的最有效方法之一。

其基本原理是将直方图在某一阈值处分割成两组，一组对应于背景，一组对应于目标。

当被分成的两组物体组间方差最大时，得到最佳分割阈值。

这种确定阈值的方法，也称为最大类间方差阈值分割方法或大津法。

Ostu 法选取出来的阈值非常理想，在各种情况下的表现都较为良好，可以说是很稳定的分割方法，
其基本原理如下：
设图像具有L 级灰度，任一门限值t 将图像中的像素划分为两类：
}1,,2,1{},,,1,0{10-++==L t t C t C ΛΛ (6.12)
若t y x f ≤),(，则0),(C y x ∈；若,),(t y x f >则1),(C y x ∈
对图像的直方图进行归一化便得到灰度级的概率分布：
N n p i i /=（N 为图像中的像素总数） (6.13)
在图像中
0C 和1C 类的出现概率分别为
∑∑=-+=-====t i L t i i i t w p w t w p
w 01110);(1),( (6.14)
C 和1C 类的类均值分别为
,/000∑=⋅=t
i i w p i μ (6.15)
,/)]([/11111w t w p i L t i T i ∑-+=-=⋅=
μμμ (6.16)
式中，∑-=⋅=1
0L i i T p i μ为整个图像的均值。

因此对]1,0[-L 间的任何t 值，有
1100μμμw w T += （6.17）
C 和1C 类的方差分别为 ∑∑-+==⋅-=
⋅-=111212100202
0/)(,/)(L t i i t i i w p i w p i μσμσ (6.18)
这两类的类内方差2w σ，类间方差2B σ和总体方差2T σ分别为：
类内方差 2112002σσσw w w
+= (6.19) 类间方差 210102112002)()()(μμμμμμσ-=-+-=w w w w T T B
(6.20) 总体方差 ∑-=-=10
22
)(L i i T T
p i μσ (6.21) 三者的关系为 222T B w σσσ=+ (6.22)
式中，2w σ和2B σ是门限t 的函数，2T σ与t 无关。

为了使类间方差最大，可以选用统计量22/T B σση=，由于方差是灰度分布均
匀性的一种度量，当η越大时，说明构成图像的两部分差别越大，当部分目标错分为背景或部分背景错分为目标时都会导致两部分差别变小，所以使类间方差最大的分割意味着错分概率最小，此时的分割阈值自然也就是图像分割的最佳阈值。

因此以使η最大时阈值t 作为分割阈值,便可实现图像的理想分割。

大津法本质上属于单阈值的图像分割方法，即只能将图像分为两类，当图像中目标与背景大小比例很小时该方法可能会失效。

在利用该方法求得最佳阈值的过程中，大津法要求用穷举法对每一灰度级进行一次类间方差的计算，这样无疑增加了运行时间。

尽管如此，它仍不失为一种性能优良的自动化阈值分割方法，
在图像分割领域中得到了极为广泛的普及和应用。

图6-1(3)-(6)分别给出了用迭代法、最大熵法、矩量保持法和Ostu方法对肺部CT图像进行阈值分割的结果
（1）原始肺部CT图像（2）一维灰度直方图（3）迭代法分割结果
（4）最大熵法分割结果（5）矩量保持法分割结果（6）Ostu法分割结果
图6-1 图像阈值分割示例
第三节基于边缘检测的图像分割方法
图像最基本的特征是边缘，它是图像局部特性不连续(或突变)的结果。

基于边缘检测的图像分割方法的基本思想是先检测图像中的边缘点，再按照一定的策略连接成轮廓，通过轮廓跟踪完成区域分割。

实际上，在数字图像中，边缘总是以某种图像特征所对应数值发生突变的形式出现的，往往体现为图像局部特性的不连续，比如像素灰度、颜色、纹理等特征的突变。

图像的边缘包含了物体形状和目标结构的重要信息，它常常意味着一个区域的终结和另一个区域的开始。

因此，基于边缘检测的图像分割方法可以通过检测出不同均匀区域之间的边界来实现分割。

作为所有基于边界分割方法的第一步,经典的边缘检测方法是通过构造对图像灰阶变化敏感的差分算子来进行图像分割的。

常见的对于灰度值不连续(或突变)的检测主要借助于空域微分算子进行，通过将各种边缘检测模板与原始图像卷积完成计算。

这类方法大多是基于局部空间信息的处理方法，一般利用图像一
阶导数的极大值点或二阶导数的过零点信息来提供判断边缘点的基本依据。

常见的边缘检测算子有梯度算子、Robert 算子、Sobel 算子、Prewitt 算子、kirsch 算子、Laplacian 算子、Marr 算子等。

一、梯度算子
对于二维图像中的每一个像素),(j i ，计算它的梯度 22]/),([]/),([),(j j i f i j i f j i g ∂∂+∂∂=， (6.23)
其中
)
1,(),(/),(),1(),(/),(--≈∂∂--≈∂∂j i f j i f j j i f j i f j i f i j i f 梯度),(j i g 的大小反映了图像灰度局部变化的强弱，因此可以作为检测边缘点的依据。

梯度算子可以有效地检测出阶跃型边缘，而对屋脊型边缘容易产生双边缘效果。

此外，Robert 提出了另一种类似算子，该算子的形式是:
})1.1()1,1(,)1.1()1,1(max{),(+---+---++=j i f j i f j i f j i f j i R (6.24)
Robert 算子基于可通过任意一对互相垂直方向上的差分来计算梯度的原理，采用两对角线方向相邻像素之差近似梯度幅值来检测边缘。

它检测斜向边缘的效果好于水平和垂直边缘，具有计算简单，定位精度高，对噪声敏感等诸多特点。

二、Sobel 算子
该算子定义为
),(),(),(j i f j i f j i S j i ∆+∆= (6.25) 其中
),
1,1()1,(2)1,1()1,1()1,(2)1,1(),(),1,1(),1(2)1,1()1,1(),1(2)1,1(),(-+------+++-++-=∆+-------+++++-+=∆j i f j i f j i f j i f j i f j i f j i f j i f j i f j i f j i f j i f j i f j i f j i
Sobel 算子比梯度算子和Robert 算子的抗噪声能力要强一些，这是因为 Sobel 算子是一个8-邻域算子，而梯度算子是一个4-邻域算子。

Sobel 算子对噪声具有一定平滑作用，能提供较为精确的边缘方向信息，但它同时也会检测出许多的伪边缘，边缘定位精度不够高。

当对边缘定位精度要求不高时，不失为一种较。