宫颈细胞核质自动分离算法研究∗

于月娜 梁光明 刘任任

（1.湘潭大学信息工程学院 湘潭 411100）（2.国防科技大学电子科学与工程学院 长沙 410000）

1 引言

宫颈癌是一种发病率很高的妇科恶性肿瘤疾病，已经严重影响到了妇女的健康［1］，目前在国内外病理学研究领域中都很重视对于该疾病的研究。宫颈细胞形态学研究是为了准确检测宫颈细胞的数量、形态种类、分布密度与核质比等，为宫颈癌前病变检测提供有力的依据。其中，在宫颈细胞形态学研究中细胞核与细胞质分割效果的好坏在很大程度上影响特征的提取和分类的准确性。

目前已提出许多宫颈细胞核质分离算法，如用支持向量机方法分割［2］，采用K-means 聚类算法［3］，基于分水岭改进的分割［4］、GVF Snake［5］等。这些方法会出现细胞分割不完整、背景误判为目标等缺陷，而且对于粘连细胞和变形较严重癌变细胞核质分离效果差的问题没有得到有效解决，并且效率也不是很高。

基于以上研究，本文主要研究如何提高宫颈细胞核质分离的准确率和分割效率，提出了一种基于Grab Cut算法结合Canny算子边缘检测高效分离宫颈细胞核与细胞质的方法。对获取宫颈细胞图像进行灰度变化预处理求得灰度直方图，根据直方图求解确定Canny 算子双阈值，进而获取宫颈细胞核和细胞质的轮廓；根据获取的轮廓利用改进的Grab Cut 算法进行核质分离，并且针对Grab Cut 算法在分割效率方面的不足也进行相应改进，使其高斯混合模型的参数的确定效率提升。实验结果表明，本文所提出的算法不仅有较好的分割效果而且效率也很高。

2 Canny 算子和形态学获取细胞核轮廓

2.1 Canny边缘检测算法

边缘检测是图像处理中用于精确定位目标图像的一项重要技术。图形边缘作为图像的基本特征是图像分割和特征提取分析及分类识别的重要依据［6］。经典的边缘检测算子有Robert、Prewitt、Log 等，其优点是简单、易于实现，同样存在对噪声敏感、抗干扰性能差，边缘不够精确的缺点［7］。相比这些算子，Canny 算子具有更好的信噪比和检测精度，在图像边缘检测领域中具有更加广泛的应用范围。

2.1.1 传统的Canny边缘检测

传统Canny 算子检测是通过灰度图像信号的极大值获取到边缘信息。Canny 算子检测的提出是基于准确的信噪比、定位精确性和边缘相应唯一性三个准则；其实现包括图像高斯滤波平滑处理、计算求解梯度幅值和方向、非极大值抑制处理和双阈值边缘检测和连接四个步骤。

传统的Canny 边缘检测算子基于边缘检测准则和简单快递处理特点在图像处理领域应用［8～10］很是广泛，但是传统的Canny 边缘检测算子存在一个重大缺点是双阈值参数的选取困难，通常是根据人工经验确定的，而没有通过图像真实特征信息计算确定，这样就导致了其算法适应性稍差，处理的图像产生边缘信息不完整和假边缘存在；而且人工经验确定需要反复验证，很大程度地降低了分割效率。

针对传统的Canny 边缘检测算子这种缺点本文提出一种快速确定Canny 算子双阈值的方法，可以解决双阈值选取困难的问题，而且可以快速准确地确定高低阈值，更多地保留图像的真实边缘信息，从而更好地获得图像边缘。

2.1.2 改进的Canny边缘检测算子

1）灰度直方图图像预处理

在细胞图像的灰度直方图中不同的峰谷值对应不同灰度分布［1］，直方图图像预处理在医学图像处理中应用非常广泛；特别是在宫颈细胞图像（如图1 所示）中主要包括细胞质、细胞核和背景不同的组织，在图像中可以发现其灰度差异比较明显，非常适合灰度直方图预处理，可以获取很多图像的基本信息。

宫颈细胞的边缘检测是得到细胞核的边缘或者得到整个细胞的边缘，宫颈细胞图像（如图1 所示）灰度直方图如图2 所示；直方图中可以明显看出两个波谷的存在，两个波谷分别表示细胞核分割阈值和细胞与背景的分割阈值大致范围。

图1

图2

2）确定双阈值

根据宫颈细胞的灰度直方图分布特点，两个波谷就是大致的分割阈值G1，G2；而阈值处的像素变化最大，梯度也就最大，经研究分析决定采用如下方法来确定其双阈值：

（1）运用一对卷积阵列（分别作用于x 和y 方向）：

（2）使用下列公式计算梯度幅值：

（3）先求阈值G1处的最大梯度作为双阈值中的高阈值Kh，然后根据Kh=3K1求出来K1。根据双阈值利用Canny算子获取宫颈细胞核的边缘。

（4）同理先求阈值G2处的最大梯度作为双阈值中的高阈值Kh，然后根据Kh=3K1求出来K1。根据双阈值利用Canny算子获取宫颈细胞的边缘。

2.2 形态学后处理

形态学处理是较为常用图像分割定位方法，其基本思想是利用称为结构元素的“探针”收集图像信息；形态学［13］处理主要方法有膨胀和腐蚀以及开运算和闭运算，假设A和B都是Z2中的集合，

由于噪声的影响，宫颈细胞图像在Canny 算子处理过后所得到的轮廓边界通常都存在不平滑和不连通区域。本文对Canny 算子处理过后所得到的轮廓边界用4×4 的椭圆结构元素进行先闭运算后开运算去除噪声，平滑边界；得到的区域作为Grab Cut 分割的前景，这样能够得到更好的分割效果。

3 Grab Cut分离细胞核与细胞质

3.1 Grab Cut算法原理

Grab Cut 算法是基于图论的一种交互式的高效分割算法，是对Graph Cuts 算法的改进算法，引入了迭代的思想融入了传统图像分割的边界信息，通过高斯混合模型（Gaussian Mixture Model，GMM）来标记图像的像素分布情况，通过迭代估计和非完全标记法实现能量最小化进而完成图像分割。Grab Cut 算法能够在复杂背景中快速高效提取目标图像，其分割准确率高、交互量少、分割速度快。

在Grab Cut算法图像分割中，一般都采用高斯混合模型聚类算法来表示前景和背景的颜色分布；设整幅图像为y=( y1,y2,…,yN)，其label 集为α=( α1,α2,…,αN)，用高斯混合模型（Gaussian Mix⁃ture Model，GMM）θ 表示图像数据，那么图像分割就可以转化为数学表达式：

其中，E( α,k,θ,y )为图像的Gibbs能量函数：

公式里的R( α,k,θ,y )和B( α,y )分别表示区域项和边界项。

Grab Cut 算法分割过程是根据RGB 颜色空间和三分图信息，通过计算全协方差GMM 对背景和前景进行建模。每个前景与背景区域创建的GMM有K个高斯分量，根据用户标记的未知区域TU和背景区域TB，设，前景区域TU=Ø；Grab Cut算法通过初始化TU和TB的GMM 参数，并迭代TU区域的GMM 参数使其达到最小化来获取图像最小割。

Grab Cut 算法在处理最小化迭代计算获取最小割的过程中，由于高斯模型进行初始化过程中所采用的聚类方法不同会导致产生不同的效果，也会很大程度上影响算法的效率，常用的聚类算法有K-Means算法［11］、K-Medoids算法［12］、Clarans算法［13］、和EM算法［14］，通常采用EM算法。

3.2 Grab Cut分离算法改进

Grab Cut算法在分割过程中通过利用Alpha 值和聚类算法对高斯混合模型进行初始化构建，然后通过迭代估计的方式确定GMM参数在很大程度上限制了图像分割效率；高斯混合模型的参数的确定通常采用EM 算法，其优点是计算稳定、收敛速度快。

式中，k 是高斯混合模型的成分数，aj是满足的混合比例参数，是第j 个高斯成分密度函数其表示形式为

其中，θj=( μj,Σj)是第j 个高斯分布均值和协方差矩，Θk=( θ1,θ2,…,θk,a1,a2,…,ak-1) 为模型的参数集合。

EM 算法就是我们常说的期望最大算法，最早是在1977 年由Dem pster 等提出的应用于求解参数极大似然估计的一种方法。通过迭代EM 算法求解混合模型似然函数准则，估计模型中的参数，从而解决问题。准则函数可以表示为

其中，L( Θk,x )是似然函数，表示为

EM 算法具体流程为

Step 1：在参数样本I中选取一个合适的初始参数Θ0；

Step 2：根据选取的Θ0，设计并计算辅助函数：

同时更新参数：

如上述所示，EM 算法主要是确定符合高斯混合模型分布的数据集中不同的样本的特定高斯分布。

EM 算法在求解GMM 模型参数的整个过程中，由于似然函数是单调递增的，可以根据初始模型和参数求解得出每个成分下各数据点的概率，然后通过不断修改参数值，不断重复直到收敛某个极值点为止。

EM 算法流程中设计与计算Q 辅助函数其实就是求解似然函数的期望，这是EM 算法的核心步骤，就是为了将完整数据集｛X，Y｝中的隐藏变量Y转化为别的统计量，从而可以消除隐藏变量Y 对后续计算不明确的影响；而求解Q 辅助函数的极大值就是保证Q函数是递增的。

在EM 算法中存在两个问题，一个是算法容易陷入局部最优而导致求解结果不是全局最优解，如图3 所示，图中横坐标表示一个未知参数，纵坐标表示其似然函数L，此时选择的初始点分别为A 和B 时，分别收敛到M1 和M2 是似然函数的全局极大值点和局部极大值点，得到对应的解就分别为全局最优解和局部最优解；另一个问题是传统EM 算法需要预先设定混合分量个数，这样就导致拟合效果不稳定，参数确定过程复杂效率较低。

图3 极值点分布图

根据试验结果分析，在Grab Cut 算法分割过程中确定GMM 参数的时间占整个分割时间的80%以上，经过对高斯混合模型和EM 算法［15］基本原理的学习与研究，本文针对Grab Cut算法在分割效率方面的不足，对传统EM 算法作了相应改进，使其高斯混合模型的参数的确定效率提升，这在很大程度上提高了Grab Cut 算法分割效率。

本文针对传统EM 算法存在的问题进行了改善与优化，针对算法容易陷入局部最优的情况对数据点集进行多个划分，具体根据马氏距离划分为K组，然后根据每组数据点选取合适的初始值进行参数估计收敛得到K 组对应的估计值，最后对K 组极值点进行比较得到最大的一个极值点就是对应的最优解。

在上述过程中所有参数更新过程是并行进行的，这样不会影响到算法的效率但是也不会提高很多，针对这点本文对EM 算法进行相应改进，提高其算法效率。传统EM 算法的Step 2 和Step 3 在更新参数迭代过程中每次都进行参数的更新，这样会一定程度上影响算法的效率。本文对其进行改进改为在每更新完一个成分的参数后，不是进行下一个成分的更新，而是利用更新的参数重新更新似然函数函数L( Θk,x )的期望，然后再更新下一步；同时在迭代过程中将混合分量的权重为0 的成分去除，不再进行更新，这样就可以加快算法的收敛速度，提高算法的效率。

本文改进的EM算法流程为

Step 1：在参数样本I中根据马氏距离将其划分为K 组，并在K 组中分别选取并赋予合适的初始值ΘK0；

Step 2：根据赋予的ΘK0，设计并计算辅助函数：

同时更新参数：

则继续更新下面参数：

Step 6：在K 维数组Qmax中求得最大的对应参数，输出结果。

经过改进的EM 算法可以改善传统EM 算法容易陷入局部最优的缺陷，使得拟合结果更加稳定，同时可以使算法收敛速度更快，效率得到很大的提高。

4 实验结果与分析

本文对200 多幅宫颈细胞图像进行了分割实验。图4 列出了本文算法GVF Snake 模型以及K-means 聚类算法对其中一类癌变宫颈细胞图像的核质分离分割效果图。通过效果图可以看出本文提出的算法对于宫颈细胞的核质分离有很好的分割效果。

图4 分割效果图

GVF Snake 模型基本思想是：把梯度场向图像边缘通过迭代的方式进行扩算形成GVF 力场；然后定义动态轮廓能量函数，通过求解能量函数的局部最小值来实现动态轮廓逐步接近图像真实轮廓。但是GVF Snake 模型存在着明显的不足是无法检测到凹型轮廓或具有较高曲率的凸型轮廓和物体内部轮廓；同时，整个运算过程较为复杂，涉及的参数也难以确定。

K-means 聚类算法是一个采用距离作为相似性指标反复迭代的聚类算法［3］能够准确、快速地处理大数据集，并且有可伸缩、高效的优点［16］。但是，该算法对噪声敏感而且K 值及初始聚类中心有不确定性，使其在图像分割中容易造成图像的误分割，并且在实际应用效率也不高。

图5 不同算法分割时间与准确率对比柱状图

通过对比分析，发现本文算法在速度和准确率上均有较大的提高。

5 结语

本文提出了一种基于Grab Cut算法结合Canny算子边缘检测高效分离宫颈细胞核与细胞质的方法。对获取宫颈细胞图像进行灰度变化预处理求得灰度直方图，根据直方图求解确定Canny 算子双阈值，进而获取宫颈细胞核和细胞质的轮廓；根据获取的轮廓利用改进的Grab Cut 算法进行核质分离，并且针对Grab Cut算法在分割效率方面的不足也进行相应改进，使其高斯混合模型的参数的确定效率提升。实验结果表明，本文所提出的算法不仅有较好的分割效果而且效率与分割准确率也很高，非常适合在实际应用中实现宫颈细胞图像的核质分离。