基于配准分割的G raph Cuts自动分割算法在肝脏图像中的研究

王建军　谢勤岚

（中南民族大学生物医学工程学院　武汉　430074）

1　引言

肝癌是最常见的6种癌症之一，同时也是全世界范围内导致人类死亡的主要原因。在2012年，全世界大约有78万人被诊断为肝癌，约有74万人死于肝癌［1］。因此，较早的诊断和治疗肝病是非常重要的。肝脏医学影像分割非常重要，它可以辅助医生对肝脏疾病患者的肝脏进行诊断、功能评定和治疗决策。但是在医学图像处理中，肝脏的分割一直都是一项巨大挑战的任务。如图1所示的肝脏解剖结构，其周围有主动脉、骨头、肾等灰度相似的器官，还有肿瘤等疾病，肝脏对比度低、边界模糊，很难将其完全分割出来。

图1　肝脏的解剖结构

肝脏分割方法［2］包括半自动分割和全自动分割两类。目前已经有许多分割算法能对肝脏进行有效的分割。周向荣等使用概率模型估计初始空间的位置并计算肝脏分布的概率［3］，能够实现对肝脏的全自动分割，但分割精度不够理想。Rikxoort等通过使用K最近邻域方法计算肝脏中每个像素的概率和多尺度Atlas配准方法［4］，对肝脏进行有效的全自动分割。Foruzan等提出基于灰度分析和结构信息的方法［5］，使用期望最大值算法，计算肝脏灰度范围的统计参数，结合阈值方法，能够将肝脏与周围的组织有效的分割出来。R Pohle，KD Toennies提出的模型自适应区域增长算法［6］，对于对比度不高的肝脏CT影像，都能获得一个较好的分割结果。但是他们的方法对于初始种子点的选取非常敏感，鲁棒性不够好。

本文利用基于配准分割的Graph Cuts自动分割算法，对人体肝脏CT图像进行自动分割。该方法的分割流程如图2所示，主要的思想是利用Rikxoort提出的基于Atlas的分割方法对肝体进行有效的分割，利用权重投票［7］将所有的Atlas的分割结果融合成一幅预分割图像，然后利用水平集分割中的符号距离函数［8］自动标记预分割图像目标和背景的像素，最后将标记的像素标签送入Graph Cuts图割算法［9］中，完成自动分割。

图2　自动分割流程图

2　基于A tlas的分割

Atlas分割算法的主要思想是利用已经分割的Atlas图像对待分割图像进行非刚性配准，通过配准得到的形变场参数应用到Atlas的二值图像上，该二值图像通过应用形变场参数得到的图像即为待分割图像的分割结果。也就是说，Atlas分割算法是将分割问题转化为了配准的问题。图3是基于Atlas的分割框架图。

图3　Atlas分割框架

Atlas分割算法对肝脏图像分割结果的优劣，主要依靠配准的质量。对于配准，主要分为四个步骤：配准框架（Registration framework）、度量方法（Metrics），插值（Interpolators）和优化（Optimisers）。

本文的配准框架是使用刚性的全局仿射粗配和非刚性的局部B样条精配组合的框架，其框架模型的公式如式（1）。Igabol(x，y，z) 是配准的全局仿射粗配模型，Ilocal(x，y，z)是局部精配B样条自由形变模型（Free-form deformation，FFD）。

度量方法是配准中较为重要的部分，它是衡量两幅图像配准优劣的标准。目前度量的方法有：均方差（Mean squared difference，MSD），归一化相关系数（Normalised correlation coefficient，NCC）［10］，归一化互信息（Normalized Mutual Information，NMI），Kappa统计（Kappa statistic，KS）［11］以及本文所采用的互信息（Mutual information，MI）［12］。MI的公式定义如式（2）。

其中，IF是固定图像，IM是浮动图像，μ是形变参数，LF、LM分别是固定图像和浮动图像的灰度级，pF和 pM是边缘概率，p( f ，m;μ )是两幅图像的联合概率密度。

图像插值的作用主要是在经过变换的像素位置可能不再是整数后获得整数位置的像素值，优化是是变换参数的相似性测度达到最优。

3　权重投票

在Atlas分割实验中，对于待分割图像有多个Atlas图像，就有多少个分割结果。然而，这些分割结果的精度差别较大，并不是所有的分割结果都适合做Graph Cuts实验的预分割。为了得到适合Graph Cuts分割实验的预分割图像，可以采用联合所有的Atlas图像形变分割融合成一幅分割结果图像。权重投票融合公式［13］如下其中，c是分类数，本算法只有两类即目标和背景，故 c为2；lc(x)是Atlas分割结果中，图像的像素x属于目标或者背景的概率；δ[⋅]是克罗内克函数，当c和(Li◦Ti)(x)相等时，其值为1，当 c和(Li◦Ti)(x)不相等时，其值为0；wi是标量权重因子，它是在Atlas试验中，从配准的文件中提取的NMI的值。

权重投票融合，可以将多个Atlas分割结果融合为一幅图像，将融合后的图像作为Graph Cuts分割实验的预分割图像，可以增加实验的鲁棒性。

4　水平集符号距离函数

水平集函数是一种跟踪轮廓和曲面演化的数值方法，而不是直接计算曲面的轮廓。如图4所示，是零水平集函数示意图。图中黑色粗线轮廓部分，是零水平集函数代表的轮廓，它是由更高维度的水平集函数的等于0时对应的曲线轮廓。

图4　零水平集函数示意图

水平集符号距离函数表示像素点到零水平集函数的欧式距离。在本文的实验中，主要是像素点到Atlas分割结果图像中目标的边缘的距离，距离边缘越远的像素，距离越大，距离边缘越近的像素，距离越小。本文符号距离函数的完成是基于作者C.R.Maurer提出的线性时间算法［14］。图5是利用线性时间算法求Atlas分割结果图像的符号距离函数图像的过程，（a）是Atlas分割完成的二值图像，（b）是线性时间算法计算的符号距离函数图像，从图中可以看出，距离二值图像边缘越远的像素点灰度值越亮，表明该到边缘的欧式距离越远。

图5　Atlas分割结果和符号距离函数图像

图6是利用符号距离函数标记目标和背景像素的过程。图（a）是权重投票的结果，将所有的Atlas分割结果融合为一幅预分割图像，图（b）是使用三种不同大小的符号距离标记融合图像的目标和背景，其中的内符号距离3的大小决定着白色目标像素的大小，外符号距离1和外符号距离2决定着背景像素大小。

图6　利用符号距离函数标记前背景像素

5　实验结果

本文所有的实验都是基于IKT（InsightSegmentation and Registration Toolkit）［15］框架医学图像处理工具。其中Atlas分割模块是采用开源配准软件包Elastix，权重投票和水平集符号距离函数等模块是采用ITK的内部类，Graph Cuts算法的基础代码部分是由 Tustison N.和 Yushkevich P.［16］等提供。本实验的数据是来自临床的CT图像数据集3Dircadb。3Dircadb数据集是由20组CT医学影像组成，图像切片均为512×512像素，切片的像素间隔为0.56~0.87，切片层数范围是91~260。

在Atlas分割实验中，用3Dircadb数据集其中一个病人图像作为训练集，另外19个病人图像作为测试集，其结果如表1所示。Patient_ID表示病人CT图像的序号，Overlap是Atlas分割结果图像融合后的分割精度。图7是对3Dircadb数据中第一个病人图像使用Atlas分割算法后将其利用权重融合的示意图。

表1　融合Atlas分割重叠率结果

图7　Atlas分割结果和权重投票融合的例子

将Atlas分割结果融合后，需要对其进行像素标记。图8是利用符号距离函数标记目标和背景像素标签的过程，图8（a）是待分割图像，图8（b）是融合Atlas分割结果图像，图8（c）利用符号距离函数标记像素的前背景像素图像。在图8（c）中，蓝色部分是融合Atlas分割结果，白色部分是标记目标像素，蓝色外的灰白色部分是标记的背景像素。经过大量的实验发现，基于Atlas分割结果图像，当内符号距离为大于5的像素全部标记为目标，外符号距离为0～15的像素标记为背景时Graph Cuts自动分割效果最好，其结果要优于融合Atlas分割结果图像的分割精度。

图8　利用符号距离函数标记像素

图9是对数据集3Dircadb中某一个病人的采用Atlas分割算法和Graph Cuts自动分割算法的对比图，图10是两种算法对数据集20个病人分割结果的线箱图。从图中可以看出，Graph Cuts自动分割算法不尽在割精度上优于Atlas分割算法，而且鲁棒更好。

图9　两种分割结果对比图

图10　两种算法分割结果重叠率Box-plot图

6　结语

本文针对肝脏CT医学影像的分割，采用Atlas配准分割为预分割结果，通过权重投票融合预分割结果图像，最后使用水平集符号距离函数标记前背景像素，实现GraphCut自动分割。最终的实验结果表明，本文提出的基于配准分割的Graph Cuts自动分割算法，相较于传统的Atlas配准分割有较大精度的提升，在操作性方面，优于传统的Graph Cuts分割算法，实现自动分割目标。

［1］Ferlay J，Soerjomataram I，Dikshit R，et al.Cancer incidence andmortality worldwide：sources，methods andmajor patterns in GLOBOCAN 2012［J］.International JournalofCancer，2015，136（5）：359-363.

［2］Mharib A M，Ram li A R，Mashohor S，et al.Survey on liver CT image segmentationmethods［J］.Artificial Intelligence Review，2012，37（2）：83-95.

［3］Zhou X，Kitagawa T，Hara T，etal.Constructing a Probabilistic Model for Automated Liver Region Segmentation Using Non-contrast X-Ray Torso CT images［J］.Med Image ComputComputAssist Interv，2006，9（2）：856-863.

［4］Rikxoort EMV，Arzhaeva Y，Ginneken BV.Automatic segmentation of the liver in computed tomograpy scans with voxel classification and atlasmatching［J］.3d Segmentation in the Clinic A Grand Challenge，2007，147（6）：132-136.

［5］Foruzan A H，Zoroofi R A，HoriM，et al.Liver segmentation by intensity analysis and anatomical information in multi-slice CT images［J］.International Journal of Computer Assisted Radiology and Surgery，2009，4（3）：287-297.

［6］ Pohle R，Toennies K D.A New Approach for Model-Based Adaptive Region Growing in Medical Image Analysis［M］.Computer Analysis of Images and Patterns.Springer Berlin Heidelberg，2001：238-246.

［7］Rohlfing T，Jr C R M.Multi-classifier framework for atlas-based image segmentation［J］.Pattern Recognition Letters，2005，26（13）：2070-2079.

［8］J.A.Sethian.Level SetMethods and FastMarching Methods［M］.Cambridge University，Second edition，1999，7（2）：81-85.

［9］Boykov Y Y，Jolly MP.Interactive Graph Cuts for Optimal Boundary&Region Segmentation of Objects in N-D Images［J］.Iccv，2001，1：105-112.

［10］Knops Z F，Maintz JB A，Viergever MA，et al.Normalized Mutual Information Based PET-MR Registration Using K-Means Clustering and Shading Correction［J］.Medical Image Analysis，2006，10（3）：432-439.

［11］Fitzpatrick JM，Hill D LG，Shyr Y，etal.Visualassessment of the accuracy of retrospective registration of MR and CT images of the brain［J］.IEEE Transactions on Medical Imaging，1998，17（4）：571-85.

［12］Stefan Klein，Marius Staring.Elastix themanual［M］.Image Sciences Institute，2015，9（4）：5-6.

［13］Klein S，Van-Der-Heide U，Lips I，et al.Automatic segmentation of the prostate in 3D MR images by atlas matching using localizedmutual information［J］.Medical Physics，2008，35（4）：1407-1417.

［14］Jr CR M，QiR，Raghavan V.A Linear Time Algorithm for Computing Exact Euclidean Distance Transforms of Binary Images in Arbitrary Dimensions［J］.Pattern Analysis&Machine Intelligence IEEE Transactions on，2003，25（2）：265-270.

［15］Ibanez L，Schroeder W，Ng L，et al.The ITK Software Guide［J］.Computational Statistics&Data Analysis，2005，21：231-256.

［16］N.J TUSTISON，P.Yushkevich，Z.Song.Graph Cuts，Caveat Utilitor，and Euler's Bridges of Konigsberg［J］.Insight Journal，2008，11（14）：128-138.