冠状动脉血管造影图像三维分割方法

刘昊，王瑜，王怡宁，徐橙

1.北京工商大学人工智能学院，北京100048；2.中国医学科学院北京协和医院，北京100730

前言

作为世界人口主要的死亡原因之一，心血管疾病早期体征不明显且易导致患者猝死。在所有心血管疾病中，多数是冠状动脉问题，主要由左或右动脉狭窄引起［1］。因此，冠状动脉疾病的早期诊断、评估风险和治疗规划至关重要。

近年来，冠状动脉计算机断层血管成像检查（Coronary Computed Tomography Angiography,CCTA）由于其无创性、多视角等特点，已成为心脏病诊断和治疗的重要工具。但是由于CCTA 数据量大，复杂性高，因此即使专家进行图像诊断也通常需要花费大量时间。此外，由于CCTA 图像的噪声影响，冠脉结构复杂且被其他组织环绕，迅速准确的冠状动脉分割难以实现［2］。冠状动脉序列图像分割方法主要包括区域生长算法［3］、水平集［4］、图割方法［5］、分水岭算法［6］和基于深度学习的方法［7］等。其中，分水岭算法在医学图像分割中应用广泛，具有运算速度快、计算简单、边缘定位精准等优点，但易受到微弱边缘的影响，很容易因噪声和对比度低等原因造成过分割和易丢失重要轮廓等问题。为了克服分水岭算法的缺陷，可以结合有效预处理方法先行去除噪声，或对分割后的图像进行区域融合。

为此，本文提出一种新颖的冠状动脉分割方法，基本思路是，利用分水岭算法将冠脉图像分割为子区域，然后交互式分割出种子区域，接着根据形态学和同一组织在断层图像的连续性，以分水岭子区域为单位进行区域生长，最终实现冠状动脉的三维分割。本算法可解决传统分水岭算法造成的过分割及低对比度图像易丢失重要轮廓的问题，且运算速度快、分割精准，同时利用形态学思想解决了少数相邻层冠状动脉区域不重叠的现象。

1 方法

本文提出的方法主要包括3 个步骤：（1）数据预处理；（2）利用分水岭算法将CCTA 图像分割为子区域；（3）利用以子区域为单位的区域生长方法，得到最后的冠状动脉分割结果。

1.1 数据预处理

有效的图像预处理是冠状动脉准确分割的基础，针对CCTA 图像结构的复杂性，本文使用阈值分割方法［8］对原始图像进行预处理，以灰度直方图中冠状动脉最高灰度值为阈值，低于阈值的区域灰度值保持不变，其余区域灰度值设为0。具体公式如下所示：

其中，f(x,y,z)表示原图像(x,y,z)点的灰度值，T表示冠状动脉最高灰度值。

本文根据经验和调试，设置阈值为1 500（数据格式为DICOM）。阈值分割方法能很大程度剔除冠状动脉周围高密度组织的干扰。此外，针对CCTA 图像包含噪声的问题，使用空间域高斯平滑方法［9］对图像进行平滑处理，为后续分水岭分割提供高质量数据。

1.2 基于改进分水岭算法的CCTA图像分割

分水岭算法［10］是一种流行且简单的图像分割方法，其本质源于地理。假设景观被雨水淹没，雨水自然沿着最陡峭的路径流下，并最终落入许多盆地区域。分水岭是这些区域的分界线，是景观的自然分割轮廓。将图像看作是景观，从而利用分水岭算法将图像分割成许多细小子区域。但是，传统分水岭算法对图像中存在的噪声、物体表面细微的灰度变换较为敏感，通常会出现过分割现象，导致分割结果不理想。为此本文根据CCTA 图像中细节丰富的特点，提出一种改进的分水岭算法。首先对于预处理后图像利用Sobel 算子［11］得到梯度幅值图像，然后增强图像中极值的对比度，同时保持图像相对平滑区域，最后利用形态学修正梯度图像，并提取图像的标记修正梯度图像进行分水岭变换，得到分割结果。

1.2.1 梯度图像计算预处理后的冠状动脉CCTA 图像，虽然大量非冠状动脉组织被剔除，且图像噪声得到抑制，但是如果仅仅通过传统分水岭算法的梯度计算方法得到图像的边缘，容易出现信息丢失现象，影响算法后续的分割效果。因此，本文采用sobel 算子计算梯度图像，计算公式如下：

其中，f(x,y)表示图像(x,y)点的灰度值，Gx表示横向梯度检测的灰度值，Gy表示纵向梯度检测的灰度值。该算子较传统梯度计算方法有明显优势，能够找出图像中的细小边缘，并且可以生成噪声水平较低的图像，更适用于CCTA图像。

1.2.2 前景与背景标记虽然利用sobel 算子代替传统的梯度计算方法一定程度上解决了边缘信息丢失以及图像噪声现象，但是图像中仍然存在与目标无关的极小值，导致过分分割的现象。为了解决这个问题，本文使用标记控制的分水岭分割。如果可以提前识别并标记前景对象和背景位置，则分水岭算法的分割效果会更好。本文使用形态学运算和最大类间方差算法［12］来提取图像的标记，具体步骤如下：（1）计算图像的局部极大值；（2）利用公式f=(fΘB)⊕B对局部极大值图像进行关操作，并作为前景标记；（3）利用公式f=f⊕B对图像进行膨胀操作；（4）利用最大类间方差算法将图像二值化；（5）计算二值化图像距离变换，并利用分水岭变换求取其相邻区域间的分界线，以此作为背景标记。

1.3 基于区域生长算法的冠状动脉三维分割

传统区域生长算法［13］以像素单位作为生长的种子点，因此此类算法不仅计算时间长，且对噪声和灰度不均较为敏感，适合一些分割灰度均匀的平滑区域，无法准确高效地提取医学影像中的血管组织。针对上述问题，本文以分水岭子区域为单元，交互式分割种子区域，根据断层图像上下层间的相似性，以及冠状动脉的连通性，利用相邻层图像中对应区域的面积重叠关系以及密度相似性，自上而下进行区域生长，最终实现冠状动脉的三维分割。

1.3.1 种子区域的选择种子区域的选择直接影响最终的分割结果，尤其是在断层图像分割中。CCTA 图像对于同一组织在断层图像的相邻层之间普遍具有相似性。因此，可以根据冠状动脉断层间的相似性来选取种子区域。在本文研究的CCTA 图像中，除了断层图像中部分冠状动脉目标区域像素在前一张断层图像的目标区域中没有重叠的像素点外，其余所有冠状动脉目标区域均能够在前一张断层图像中周围找到重叠的区域像素。因此，选择无法通过上层投影而分割的目标区域作为种子区域，进行序列图像的三维区域生长分割，通过对图像的种子区域集合设置向下投影，从而实现序列图像的三维分割，如此就不会错过部分断层图像中的冠状动脉区域。

1.3.2 区域生长规则的选择在断层图像中，相邻层之间具有良好的连续性，从垂直于断层平面的方向来看，相邻层图像中同一组织普遍存在面积重叠和密度相似关系。但是，在个别情况下相邻层冠状动脉断面相互不重叠。尽管在这种情况下，二者并不重叠，但其相隔距离并不远［14］。针对这种不重叠现象，本文运用公式f=f⊕B进行形态学膨胀，对预分割的两相邻层中的前一层目标区域进行膨胀处理后再利用面积重叠关系合并下层目标区域，其中B为结构元素。本文选用圆形结构元素，半径为3 个像素，如图1所示。其中，圆形虚线是冠状动脉截面进行膨胀后的边界，圆形虚线所包围的区域与相邻下层中的冠状动脉重叠，可以看出形态学膨胀可以克服同一组织不重叠的问题。

图1 冠状动脉断层膨胀处理Fig.1 Coronary artery cross-section dilation

本研究选择面积重叠关系和密度相似性作为区域生长的生长准则，具体如下式所示：

其中，RAnew为新的分水岭子区域，RAbefore‐dilate为上一层已分割目标区域进行膨胀后的区域，RAnew⋂RAbefore‐dilate为RAnew和RAbefore‐dilate相交的像素点个数；RMnew为新的分水岭子区域的灰度均值，RMbefore为上一层已分割目标区域的灰度均值。本文通过上述生长准则，将符合条件的分水岭子区域纳入冠状动脉区域，并计算该层目标区域，进行膨胀处理后作为新的种子区域，进而搜索下一层符合条件的分水岭子区域，直到没有新的分水岭子区域纳入冠状动脉区域为止。

提出的利用交互方法分割出冠脉的步骤如下：（1）将种子区域所在层作为当前层，即层1。（2）对层1目标区域进行膨胀，在层2 中搜索与层1 中膨胀后目标区域满足面积重叠和密度相似关系的分水岭子区域，记录层2 中满足条件的子区域并进行区域合并，其结果作为层2 中目标区域的分割结果。（3）将层2作为当前层，并重复步骤（2）合并下一层中的冠状动脉区域，直到在下一层中没有满足条件的子区域时，停止搜索。依照上述步骤进行循环生长，最终实现冠状动脉各分支的三维分割。

2 实验结果与分析

实验数据集来自中国医学科学院北京协和医院的CCTA 图像，包含231 幅切片，并由专家标注数据，图像为DICOM 格式，大小为512×512 像素。为了验证提出方法的有效性，本文精心设计了一系列实验，对比了现有算法和提出方法的结果，并与专家标注结果进行了比较。实验是在Intel Core (TM)i5‐7400CPU 3.00 GHz 的PC 机上基于MATLAB R2019a软件平台实现的。经过多次优化，实验具体的参数设定如下，冠状动脉灰度值上限是T= 1 510，形态学运算结构元素是半径为3的圆形结构，面积重叠阈值Thδ1= 0.48，相似性密度范围Thδ2= 380，其中选择初始层种子区域的平均密度和膨胀后的区域面积作为实验的初始模式点进行计算，种子区域选择情况如图2所示。

图2 交互式分割种子区域的选择Fig.2 Selection of seed region in interactive segmentation

本文对比了传统分水岭方法［10］、三维区域生长法［13］和本文提出方法的冠状动脉分割结果，3 种方法将在相同预处理和种子点选取的情况下进行实验，具体结果如图3和表1所示。使用有经验的放射科医生手动分割结果作为标准，分别利用Jaccard 指数和Dice 相似系数得分进行客观评测，两种评测标准的范围都为［0,1］，其中1 表示重叠较高，0 表示无重叠，具体公式如下：

表1 冠状动脉分割的评估Tab.1 Evaluation of coronary artery segmentation results

图3 分割结果及医生标注Fig.3 Segmentation results obtained by different methods and manual segmentation results

其中，Y是医生标注结果，Yp是算法分割结果。

从图3中可以看出，传统分水岭算法由于过分割现象严重，分水岭子区域过于粉碎，很难选取合适的面积重叠阈值，导致错分割或欠分割现象严重。而三维区域生长算法的结果同样极易发生欠分割现象，这是因为传统区域生长算法只能从体素灰度特征上设定相似性规则，而对于冠状动脉和周围组织的模糊边界，则很难选取合适的阈值进行有效分割。而本文的方法以分水岭子区域为单位，利用面积重叠关系和密度相似性进行区域生长，可以有效抑制图像中其他组织对于冠状动脉的影响，虽然在细节上和手动分割结果有些误差，但是冠状动脉的主枝区域是十分准确的。

从表1可以看出，在两个定量评价指标中，本文算法的Jaccard 指数为0.936 8，Dice 相似系数得分为0.967 3。两者的指标都接近1，高于传统分水岭算法和三维区域生长算法的定量评价指标，说明本文的分割结果更为精确，即本文算法可以更有效地分割冠状动脉。

3 结语

本文结合分水岭算法与区域生长算法，提出了一种新颖的冠状动脉三维分割方法，首先利用阈值法和高斯滤波对冠状动脉CCTA 图像进行预处理，然后利用改进的分水岭算法将断层图像分割为若干数个子区域，并利用交互式方法分割出种子区域，最后利用形态学方法和断层图像相邻关系，实现冠状动脉的三维分割。实验结果表明，本文方法能有效地分割出冠状动脉，但在血管细小分支部分有待提高，后续工作将进一步研究有效提高血管细小分支部分精准分割的方法。