Article Info

1.北京工业大学 a.生命科学与生物工程学院；b.信息学部，北京 100124

[摘要]目的本研究旨在探寻一类性能优异的血管增强算法，并结合阈值水平集分割算法进行肝脏血管系统的三维自动分割。方法首先对原始三维增强CT数据进行S型非线性灰度映射；随后对不同的血管增强算法进行对比分析；最后使用阈值水平集分割算法分割出肝血管系统。选用3Dircadb公开数据集中的20例腹部增强CT数据定量评估了两类经典的血管增强算法，包括血管特征提取算法及扩散滤波算法。结果血管特征提取算法运行效率平均优于扩散滤波算法。血管特征提取算法结果的对比度平均高于扩散滤波算法2 dB以上，导致扩散滤波算法后续的计算复杂度高，准确性降低。阈值水平集分割算法的结果与区域生长算法、形态检测水平集算法和测地线活动轮廓水平集算法相比，准确性达77%以上，高于其余分割算法。结论血管特征提取算法与扩散滤波算法相比，更适合依赖灰度值的血管分割。阈值水平集算法能缓解单纯依赖阈值或依赖血管边界的血管欠分割问题，结合血管增强算法后能更准确的分割出肝脏血管。

[关键词]肝血管分割；血管增强；扩散滤波；血管特征提取；阈值水平集分割

引言

基于腹部增强CT图像的肝血管分割技术是计算机辅助诊断及治疗肝脏相关疾病的核心技术之一，广泛应用于肝肿瘤消融导航及监测、肝区切除或活体移植等肝脏外科手术中，但准确分割肝脏内部血管系统是一项艰巨的任务[1]。除图像噪声及复杂不均匀的肝脏组织影响外，还由于肝脏血管自身结构复杂，以及造影剂的影响[2]，给肝血管分割过程带来了诸多不利因素。

理想情况下，肝脏血管被假设为均匀半径的管状交叉结构[3]，基于这一假设衍生了许多血管增强算法，以血管特征提取和扩散滤波为典型[2]。此外，还有基于梯度向量流保存血管边界的滤波算法[4]，基于最优通量管状检测算法[5]和基于改进的高斯滤波检测邻域曲线结构算法[6]等。

就血管分割算法而言，依据算法原理进行如下分类：① 区域生长算法[7-8]；② 水平集或活动轮廓分割算法[9-10]；③ Graph cuts 算法[11-12]；④ 其他形态学或高斯混合模型的分割算法[13]。但这些算法普遍存在如下问题：① 低信噪比图像复杂血管结构处分割不完整；② 模糊的血管边界处欠分割；③ 对图像噪声敏感；④ 肝区边缘误检。

针对上述问题，本研究希望找到一类血管增强算法，在完整检测并增强血管灰度值的同时，抑制血管周围肝脏组织及噪声。在此基础上，使用阈值水平集算法自动分割肝血管，改良单纯依赖血管灰度值或边界的传统分割方法。本研究取得了较好的结果。

1 材料与方法

1.1 材料

为准确评估现有各种血管增强和分割算法性能，我们选取3Dircadb公开数据集进行算法测试，该数据集来自法国消化系统癌症治疗研究所。其中包含20例对比度增强CT图像（男女各10例），各例数据像素体积在0.56～0.86 mm，层厚为1～4 mm，切片数在64～502片间，单层分辨率为512×512像素。数据集中包含一名放射学专家逐层手工勾勒的肝血管轮廓，以此作为本实验肝脏血管分割的金标准。

1.2 方法

预处理使用Sigmoid非线性映射算法，增强图像目标前景与背景的对比度。这一过程可以有效的降低肝脏边界被血管特征提取算法误检为血管结构的概率，见图1。

其中I表示输入图像，I＇是输出图像，Min、Max是输出图像的灰度极值，α定义了待突出的血管灰度范围，β定义了待突出血管的中心灰度值。本研究α取值25，β取值135。

血管特征提取算法基于三维CT图像的Hessian矩阵特征值分解，分解采用Kopp[14]提出的计算方法。为了增强不同尺度的血管，采用在多个高斯尺度空间响应中求取极大值的多尺度策略。本研究实现了Frangi[15]、Sato[16]、Erdt[17]的血管特征提取算法，见图2。

图1 Sigmoid非线性映射算法结果

注：a.原始图像；b.Sigmoid滤波结果；c.不经预处理的血管增强；d.经预处理的血管增强。

图2 血管特征提取结果

注：a～c分别为Frangi、Sato和Erdt算法结果。d～f分别为Frangi、Sato和Erdt三维重构结果。

Frangi[15]：通过对Hessian矩阵的特征值分进行如下排序

。血管结构的方程为：

其中 α、β、γ是

的敏感系数，分别抑制了滴状、片状结构和噪声。

Sato[16]：通过对Hessian矩阵的特征值分进行如下排序：λ1＞ λ2＞ λ3，令 λc=min(-λ2,-λ3)。其响应函数为：

Erdt[17]：通过对Hessian矩阵的特征值分进行如下排序：λ1＞ λ2＞ λ3，其响应函数为：

扩散滤波需要求解局部偏微分方程ut=div(D×u)，其中Δu是图像梯度，D是控制扩散的扩散张量，如果扩散张量D由标量值g代替，那么这一扩散称为各向同性扩散，而D如果是各向异性张量值，则扩散称为各向异性扩散。扩散滤波结果，见图3。

图3 扩散滤波结果

注：a～c分别为RPM、血管增强滤波算法、混合扩散算法结果；d～f分别为3种算法三维重构结果。

RPM[18]：Perona和Malik[18]引入一种正规化各向同性的非线性扩散，描述方程为：

Δ Δ 旨在将同质区域平滑，而保留高梯度值的边界。Catté等[19]改进了标量扩散张量，提出一种引入高斯微分尺度σ的非线性扩散方式：

血管增强扩散[20]：Manniesing等[20]借助多尺度策略的Hessian矩阵的特征值定义扩散张量Dv。排列Hessian矩阵特征值

是对应的特征向量，则Hessian矩阵可表示为H=Q×Λ×QT。扩散张量定义为Dv=Q×Λ＇×QT，Λ＇阵的对角元素为：

带有平滑开关的混合扩散[21]：因为相干增强异性扩散适用于管状结构的滤波，而边界增强异性扩散适用于平面区域及边界区域的滤波[22]，因此Mendrik等[21]提出了带有平滑开关的混合扩散，主要思想是使用一种投票准则来判断局部区域是否为管状区域。结构分类器的定义为：

其中a=0.001，张量结构

的特征值是 λ1＞ λ2＞ λ3，Kρ是高斯卷积核。扩散张量是 DH=Q×Λ×QT。而对角阵的对角元素λhi混合了边界增强异性扩散和相干增强异性扩散的特征值λei和λci，定义如下：，是对比参数。

对上述6种血管增强滤波器处理后的图像进行血管分割。分割过程主要包括种子点的自动选取及使用阈值水平集算法进行血管分割。

第一步：肝脏掩膜的提取。我们通过数据集中的肝脏掩膜与原数据进行逻辑与运算，得到肝脏感兴趣区域。

第二步：种子点的自动选取。我们将图像分割成有限子集，且每个子集确保有且仅有一条血管通过。考虑到血管的体积及其不对称的血管结构，我们选择5 mm×5 mm×3 mm的子域。之后我们找到每个子域中的灰度极大值点，作为候选种子点。将预处理图像最大灰度极值的0.75倍值作为阈值，当候选种子点灰度大于此阈值时，选为最终的种子点。经试验0.75倍极值能确保种子点分布在各尺度血管当中。

第三步：阈值水平集的血管分割。这是一种基于水平集框架的阈值连通域分割算法，旨在确定一个灰度值范围，

其中

借助水平集方程的传播项在这一灰度值范围内进行传播。通过这种阈值水平集算法，可以将水平集的扩展迭代次数与灰度阈值相结合，防止分割的泄露现象。1.2.4 评价标准

为定量分析血管增强算法，我们首先从3方面将各算法分割的结果与肝血管金标准进行比较：

其中，TP和TN分别是正确分割为血管及背景的体素个数，FP和FN分别是错误分割为血管和背景的体素个数。

为了避免肝脏组织背景像素基数过大（

），会导致血管分割的准确性受TN值左右，因此本研究将血管金标准进行两个像素的膨胀运算并与肝脏掩膜进行逻辑与运算，得到以肝脏血管为主、背景组织为辅的肝血管兴趣域。在此兴趣域中进行血管分割结果的评估，更为客观。

表1 血管增强滤波器最优参数

本文的算法使用Insight Toolkit（ITK）工具包在VS2010平台下实现，阈值水平集算法主要借助itk::Thresh oldSegmentationLevelSetImageFilter类，下限阈值为图像灰度均值与图像灰度均方根值的和，上限阈值为图像灰度极大值。其余参数均为ITK默认值。区域生长算法主要借助itk::ConnectedThresholdImageFilter类，阈值与阈值水平集算法的选择相同。形态检测水平集算法和测地线活动轮廓水平集算法主要借助itk::ShapeDetectionLevelSetImageFilter类及itk::Geo-desicActiveContourLevelSetImageFilter类，参数设置为ITK默认值且与阈值水平集算法一致。

2 结果

血管特征提取算法及扩散滤波算法结果，见图2～3。其评估结果，见表2～4、图4。

表2 血管分割结果评价之准确性（%）

表3 血管分割结果评价之敏感性（%）

表4 血管分割结果评价之特异性（%）

由准确性、敏感性和特异性可以看出血管特征提取算法的准确性略高于扩散滤波，Erdt血管特征提取算法准确性的均值最高，为83.56%。血管特征提取算法的敏感性较扩散滤波算法低，特异性高。除此以外，我们对两种血管增强算法的运行时间及后续运算代价进行比较。随机从20例数据中选出4例数据，进行运行效率的比较，见表5。

表5 血管增强算法评价之运行效率（s）

其中混合扩散算法和血管增强扩散算法的时间代价相对较高，近似其余方法的1.5倍以上。为了分析血管增强算法对于后续血管分割算法运算量的影响，我们首先引入对比度噪声比这一概念，定义如下：

其中μv和μb是血管区域和背景区域的灰度均值，σb是背景区域的灰度标准差，对比度噪声比能定量比较图像前景与背景的对比度信息。同样选择上述4例数据进行分析，对比度噪声比，见表6。血管特征提取算法的对比度噪声比明显高于扩散滤波算法2～3 dB，而血管与背景的对比度越高，后续基于灰度值的分割会越准确，即降低了血管欠分割的几率。此外，按本文的种子点自动选取方式，血管与背景对比度越高，种子点的基数会越少，从而后续分割算法的运算量也会降低。

表6 血管增强算法评价之对比度噪声比（dB）

综上所述，血管特征提取算法在运行时间及对目标与背景对比度增强的性能方面要优于扩散算法。在血管分割算法方面，将本文的阈值水平集分割算法与传统的形态检测水平集分割算法[23]、测地线主动轮廓水平集分割算法[24]和区域生长算法[25]进行对比。首先进行3种水平集分割算法对比，见图5。图5对单层切片选用单一种子点进行血管分割，由上至下分别是形态检测水平集算法、测地线主动轮廓水平集分割算法和阈值水平集算法结果。图5a为S型非线性映射后分割算法的结果对比，此时血管边界对比度高。形态检测水平集算法、测地线主动轮廓水平集分割算法使用单一种子点对血管分割不完整，而阈值水平集分割算法的扩散项在指定灰度区域内进行迭代，能够分割出更完整的血管结构。图5b为直接对原始数据进行分割的结果对比，此时血管边界对比度低。形态检测水平集算法、测地线主动轮廓水平集分割算法会有明显欠分割现象产生，而阈值水平集分割算法在阈值及水平集迭代次数的双重约束下，欠分割现象较轻。同样，在区域生长算法中，单一的灰度阈值面临难以准确分割血管的难题，阈值较高，则血管分割不完整，阈值较低，导致粘连周围组织甚至无法分割。我们对之前选出的4例数据经过Erdt血管特征提取后进行血管分割算法的准确性评估，见表7。

表7 血管分割算法评价之准确性（%）

本研究所有血管增强算法结合阈值水平集算法的血管分割结果，见图6。

3 讨论

由本文结果可以得出：① 结合S型非线性映射的血管特征提取算法效果更佳；② 血管特征提取算法的运行效率平均高于扩散滤波算法，混合扩散算法和血管增强扩散算法的时间代价明显较高；③ 血管特征提取算法结果的对比度高于扩散滤波算法结果的对比度，使血管特征提取算法的种子点基数少于扩散滤波算法，经血管特征提取算法后的分割过程运算复杂度低。且较低的对比度导致扩散滤波算法结合阈值水平集分割算法的分割结果中存在欠分割现象。

本文对比分析了6种主流血管增强算法，并结合阈值水平集分割算法进行血管分割。主要贡献如下：① 结合S型非线性映射，解决了血管特征提取算法将肝区边界误检为血管的问题；② 经过试验，发现血管特征提取算法结果比扩散滤波算法结果的对比度强、效率高，更适合基于灰度值的血管分割工作；③ 综合血管增强算法和阈值水平集算法的肝血管分割方法，避免了依赖单一阈值及在模糊血管边界出现欠分割的现象，取得更高的血管分割精度。

本文研究工作还存在如下不足：① 对于血管增强算法而言，目前实现的6种算法均基于高斯尺度空间，而高斯尺度空间内的高斯变换是一种同向扩散，模糊了边界，并会将临近的血管进行同化，导致欠分割现象。若要解决这一问题，需要寻找合适的算法替代高斯尺度空间算法，结合梯度向量流算法可能会有更佳的效果；② 本文分割的血管结构存在一定的断裂现象，这些断裂主要出现在血管的结点或交点，还在一些对比度较低的大曲率血管结构处。造成这一现象最主要的原因是血管特征检测算法对于滴状（或球状）结构的抑制。因此，可以简化现有的血管特征提取算法中对于滴状结构的抑制项，并将RPM算法作为预处理过程，防止对细小的噪声结构进行增强；③ 在部分容积效应及造影剂影响下，血管增强后内部灰度不均匀，也会出现过分割现象。今后的工作中可以尝试优化血管特征提取算法的响应函数，使血管增强后的血管结构灰度均匀。还可以考虑在分割算法中融入梯度信息，可能会得到更好的血管分割性能。

图6 肝血管金标准及本文血管增强算法结合阈值水平集算法分割肝血管结果

注：a为手工分割的肝脏血管金标准；b～d分别是Frangi、Sato、Erdt提出的血管特征提取算法结合阈值水平集算法的血管分割结果；e～g分别是RPM、血管增强扩散、混合扩散结合阈值水平集算法的血管分割结果。

4 结论

肝血管分割是肝疾病临床手术至关重要的问题，本文基于增强CT图像，发现血管特征提取算法增强肝血管的性能更为优异，更适合基于灰度值的血管分割。此外，本文结合血管增强算法与阈值水平集分割算法，优化了阈值敏感和血管边界灰度敏感的问题。从结果来看，本文提出的算法较传统的区域生长和水平集血管分割算法具有更高的准确性。

[1] Selle D,Preim B,Schenk A,et al.Analysis of vasculature for liver surgical planning[J].IEEE T on Med Imaging,2002,21(11):1344-1357.

[2] Luu HM,Klink C,Moelker A,et al.Quantitative evaluation of noise reduction and vesselness filters for liver vessel segmentation on abdominal CTA images[J].Phys Med Biol,2015,60(10):3905-3926.

[3] Jerman T,Pernus F,Likar B,et al.Enhancement of vascular structures in3D and 2D angiographic images[J].IEEE T on Med Imaging,2016,35(9):2107-2118.

[4] Bauer C,Bischof H.A novel approach for detection of tubular objects and its application to medical image analysis[A].Joint Pattern Recognition Symposium[C].Berlin:Springer,2008:163-172.

[5] Law MWK,Chung ACS.Three dimensional curvilinear structure detection using optimally oriented flux[A].European Conference on Computer Vision[C].Berlin:Springer,2008:368-382.

[6] Xiao C,Staring M,Wang Y,et al.Multiscale bi-Gaussian filter for adjacent curvilinear structures detection with application to vasculature images[J].IEEE T on Med Imaging,2013,22(1):174-188.

[7] Oliveira DAB,Feitosa RQ,Correia MM.Segmentation of liver,its vessels and lesions from CT images for surgical planning[J].Biomed Eng Online,2011,10(1):30.

[8] Wang Y,Fang B,Pi J,et al.Automatic multi-scale segmentation of intrahepatic vessel in CT images for liver surgery planning[J].Int J Pattern Recogn,2013,27(1):1357001.

[9] Freiman M,Joskowicz L,Sosna J.A variational method for vessels segmentation: algorithm and application to liver vessels visualization[A].SPIE Medical Imaging. International Society for Optics and Photonics[C].Washington:SPIE,2009:72610H.

[10] Shang Y,Deklerck R,Nyssen E,et al.Vascular active contour for vessel tree segmentation[J].IEEE T on Med Imaging,2011,58(4):1023-1032.

[11] Esneault S,Lafon C,Dillenseger JL.Liver vessels segmentation using a hybrid geometrical moments/graph cuts method[J].IEEE T on Med Imaging,2010,57(2):276-283.

[12] Pamulapati V,Venkatesan A,Wood BJ,et al.Liver segmental anatomy and analysis from vessel and tumor segmentation via optimized graph cuts[A].International MICCAI Workshop on Computational and Clinical Challenges in Abdominal Imaging[C].Berlin:Springer,2011:189-197.

[13] Wang G,Zhang S,Li F,et al.A new segmentation framework based on sparse shape composition in liver surgery planning system[J].Med Phys,2013,40(5):051913.

[14] Kopp J.Efficient numerical diagonalization of hermitian 3×3 matrices[J].Int J Mod Phys C,2008,19(3):523-548.

[15] Frangi AF,Niessen WJ,Vincken KL,et al.Multiscale vessel enhancement filtering[A].International Conference on Medical Image Computing and Computer-Assisted Intervention[C].Berlin:Springer1998:130-137.

[16] Sato Y,Westin CF,Bhalerao A,et al.Tissue classification based on 3D local intensity structures for volume rendering[J].IEEE T Visual Comput Graphic,2000,6(2):160-180.

[17] Erdt M,Raspe M,Suehling M.Automatic hepatic vessel segmentation using graphics hardware[A].International Workshop on Medical Imaging and Virtual Reality[C].Berlin:Springer,2008:403-412.

[18] Perona P,Malik J.Scale-space and edge detection using anisotropic diffusion[J].IEEE T Pattern Anal,1990,12(7):629-639.

[19] Catté F,Lions P L,Morel J M,et al.Image selective smoothing and edge detection by nonlinear diffusion[J].SIAM Journal on Numerical Analysis,1992,29(1):182-193.

[20] Manniesing R,Viergever MA,Niessen WJ.Vessel enhancing diffusion: A scale space representation of vessel structures[J].Med Image Anal,2006,10(6):815-825.

[21] Mendrik AM,Vonken EJ,Rutten A,et al.Noise reduction in computed tomography scans using 3-D anisotropic hybrid diffusion with continuous switch[J].IEEE T on Med Imaging,2009,28(10):1585-1594.

[22] Weickert J.Anisotropic Diffusion in Image Processing[M].Stuttgart:Teubner,1998.

[23] Koenderink JJ,van Doorn AJ.Local features of smooth shapes:Ridges and courses[J].Proc Spie,1993:2-13.

[24] Caselles V,Kimmel R,Sapiro G.Geodesic active contours[J].Int J of Comput Vision,1997,22(1):61-79.

[25] Conversano F,Franchini R,Demitri C,et al.Hepatic vessel segmentation for 3D planning of liver surgery: Experimental evaluation of a new fully automatic algorithm[J].Acad Radiol,2011,18(4):461-470.

ZHANG Ruia, WU Shuicaia, ZHOU Zhuhuanga, WANG Yuea, WU Weiweib
a.College of Life Science and Bioengineering; b.Faculty of Information Technology, Beijing University of Technology, Beijing 100124, China

Abstract:ObjectiveTo investigate a set of state-of-the-art vessel enhancement algorithms and propose an automatic 3D hepatic vascular segmentation strategy by combining the vessel enhancement and the threshold level set segmentation algorithms.MethodsFirstly, the sigmoid filter was conducted to the original 3D contrast-enhanced CT images. Then, different vessel enhancement algorithms were analyzed and compared on the filtered images. Finally, the liver vessel was segmented by combining the threshold level set segmentation method with the optimal vessel enhancement algorithm. Two groups of vessel enhancement algorithms were quantitatively assessed: vesselness filters and diffusion filters. Experimental data involved 20 cases of abdominal contrast-enhanced CT in the 3Dircadb public data sets.ResultsExperimental results showed that the vesselness filters were more effective than the diffusion filters. The contrast-to-noise ratio of the vesselness filters was higher than that of the diffusion filters by 2 dB. Thus, the diffusion filters yielded higher computational complexity and lower accuracy. The segmentation accuracy of the proposed method was over 77%, which was higher than that of traditional segmentation algorithms including region growing, shape detection level set and geodesic active contours.ConclusionCompared with the diffusion filters, the vesselness filters are more suitable for grey level based on 3D hepatic vascular segmentation. The threshold level set method could alleviate the under-segmentation problem when using only one threshold or vessel boundaries, and could segment hepatic vessels more accurately with the combination of vessel enhancement algorithms.

Key words:hepatic vessel segmentation; vessel enhancement; diffusion filter; vesselness filter; threshold level set

基金项目：国家自然科学基金项目（71661167001）；北京市属高校青年拔尖人才培养计划项目（CIT amp; TCD201404053）；中国博士后科学基金项目（ 2017M620566）。

通讯作者：吴水才，教授，主要研究方向为生物医学电子与信息处理。