王胜军 王梦然 黄 艳 康 雁,* 刘积仁,
1(东北大学中荷生物医学与信息工程学院,沈阳 110004)
2(东北大学信息科学与工程学院,沈阳 110004)
3(东软医疗系统有限公司,沈阳 110179)
在医学影像分析中,计算机断层造影增强图像(computed tomography angiography,CTA)是临床检查应用中的一项重要项目。在CTA体数据中,去骨的目的是为了方便临床医师直观地观察血管结构,以检测血管狭窄、钙化和动脉瘤,具有重要的临床价值[1]。然而,在CTA中注入了造影剂的血管,虽然增强了与周围组织的对比度,但是其高亮的CT值与部分骨骼的CT值却十分接近,且部分骨骼的形态与血管的形态也十分相似,对血管的分割造成了严重的影响。因此,在CTA体数据中准确的去除骨骼区域,是一项非常有挑战性的任务[2-4]。
早期基于配准的方法[5],是通过血管造影前后的图像进行配准剪影,最后提取血管。Bani-Hashemi等在二维的基础上提出了三维灵活配准的策略[6],即分别旋转、变换和归一化患者血管造影前后的两套CT扫描数据,通过对这两套数据进行三维配准、剪影,从而达到从CTA体数据中分割血管的目的。这种方法较为直观,但是需要患者血管造影前后对应的两套CT扫描数据集,必然会增加对病人的辐射剂量,同时对具有大量的CTA体数据来说计算量很大,不适用于实际的临床诊断。
文献[7-9]提到的半自动分割算法是目前比较流行的去骨技术,前提是需要用户的手动参与,让用户选择被去除的区域或要保留的区域,然后再辅助一些常规的去骨算法,以达到分割血管的目的。这种方法受操作者个人的经验和偏见的影响比较严重,而且螺旋CT在对骨骼进行扫描时层数太多,需要操作者花大量的时间,计算量也比较大。
为了提高去骨算法的效率,减少去骨技术的时间,Shiffman和Rubin提出了一种混合半自动分割的CTA去骨方法[10-11]。这种方法首先对二维图像进行分割,并对各个分割部分进行不同的标记;然后通过用户选择感兴趣的区域,把有相同标记的区域聚集在一起。这种方法涉及到连通性算法,并且分割的区域依附于前一次分割的结果,虽然避免了操作者太多的主观因素,但是却受限于前一部分的分割结果;而且对于大量的CTA体数据来说,算法的效率仍然会比较低,鲁棒性不是很好,所以也不适用于临床诊断。
针对CTA体数据的去骨问题,本研究提出一种基于造影血管特征和层间相关性的自动去骨方法。此方法只需一套CTA体数据,无需用户交互,是一种全自动的快速去骨方法。首先,利用血管和骨骼的CT值特征,对骨骼进行粗分割;然后,根据血管的空间几何形态特征,进一步分离血管和骨骼,最后利用骨骼的三维空间特性和血管的层间相关性,对血管进行断层修补,在实现快速去骨的同时,完成了三维血管的可视化。
CTA是螺旋CT在临床检查应用中的一项重要项目,虽然在CT扫描之前注入了造影剂,增大了血管区域的CT值,使得血管区域和周围组织的对比度得到了增强,但是由于骨骼的高密度特性,其CT值与血管区域的 CT值很接近,严重阻碍了三维血管信息的提取,所以骨骼仍然是血管三维可视化的一个重要障碍。在CTA体数据中去除骨骼的问题,实质上就是分离骨骼和血管的问题。
利用CTA体数据中软组织器官的低CT值特点,可以首先采用阈值分割的方法,把CTA体数据中的软组织(如皮肤和脂肪等)与骨骼和血管区域分离开来,以得到完整的骨骼和血管区域。然后利用骨骼的强密度特性,对每一层的体数据进行阈值判断,采用区域生长的方法,粗提取出骨骼区域。虽然提取的骨骼区域已经去除了全部的软组织,但是也不可避免地去除了骨骼中的松质骨区域,而且在这部分区域也会包含存在钙化斑块的血管区域。这主要是由于两个原因造成的:一是血管中的钙化斑块的存在,使得血管的部分CT值在皮质骨的CT值之间;二是容积效应。螺旋CT在采集过程中出现的容积效应,使得皮质骨的部分CT值在血管的CT值之间。
血管中的钙化斑块是临床诊断中的重点,是不能忽视的,而骨骼的容积区域在CTA体数据中是以渐变的形式融入背景区域或周围组织的,并不是突然截止的,因此采用阈值分割和区域生长的方法会不可避免地出现骨骼提取不纯的现象。
从骨骼的临床解剖学知识上可知,骨骼中除了有皮质骨外,还有大量的松质骨存在,如图1(a)显示的区域。松质骨经螺旋CT扫描后呈现为较低的CT值,使得在区域生长的过程中容易出现欠分割现象,造成所分割出来的骨骼区域内部存在较多的孔洞,不能包括所有的骨骼区域,如图1(a)中黑色的线所圈起来的图形。因此,在经过阈值分割和区域生长后,可以采用数学形态学方法对封闭的皮质骨区域中未分割出来的松质骨区域进行处理,处理之后的结果如图1(b)所示,进而完成对骨骼的粗提取。图1(b)中黑色的线所圈起来的图形即为阈值分割后,再采用形态学方法得到的骨骼区域。
经CT造影成像后的部分血管区域与骨骼区域虽然在CT值上很接近,但是它们的空间几何形态特征却有着明显的不同。例如,在下肢CTA横断面图像中,一些骨骼的面积明显大于血管的面积,同时骨骼的空间分布与血管也存在很大的差别,因此可以根据血管和骨骼的空间几何形态特征在粗提取的骨骼中分离血管。
图1 形态学方法处理松质骨前后对比。(a)形态学方法处理前;(b)形态学方法处理后Fig.1 Contrast chart before and after morphology processed spongy bone.(a)before morphology processing;(b)after morphology processing
1.2.1 区域面积和周长
从图2中可以看到(黑色标示为血管区域),无论是人体的下肢部位还是盆骨、脊柱,大部分的骨骼区域面积都比血管区域的面积大,而周长也比血管区域的周长要长,因此可以通过区域面积和周长的特征来筛选骨骼和血管。
区域面积为
区域周长为
式中,φ(t)和ψ(t)分别表示分割区域的边界曲线方程。
1.2.2 区域均值和方差
图像均值表示的是图像整体的亮度信息,方差表示的是图像整体的对比度信息。由骨骼的解剖学知识可知,骨骼是由骨膜、骨质和骨髓三部分组成,其中骨质和骨髓中含有大量的皮质骨和松质骨。松质骨在经螺旋CT采集后,表现为低密度形态,即其CT值比较低,从而使得骨骼局部区域变得比较暗。骨骼整体呈现的细节信息比较丰富,则其均值比较小、方差比较大;血管由于造影剂的影响,在CTA体数据中整体区域表现得比较亮、均匀和单一,则其均值比较大,方差比较小。因此,可以根据区域的均值和方差大小来筛选骨骼和血管,如图2所示。
1.2.3 区域圆形度
圆形度是用来描述区域形状接近圆形的程度。从血管的解剖结构中可以看出,血管的横截面大多呈现类圆形,而人体中一些部位(如腹部、胸部等)的骨骼所呈现的形状是不规则的多边形,与圆形相差甚远,如图2所示。因此,可以用区域的圆形度来区分人体部分部位的骨骼和血管区域。区域圆形度为
式中,A为区域的面积,L为区域的周长。
当区域是圆形时,C取最大值1;当区域是细条长形或者形状较为复杂时,C的值将比较小。因此,可以把一部分骨骼和血管分离开来。
1.2.4 区域质心
在分割结果中,不规则图形的质心就是单位面积质量恒定的相同图形的面积中心。在区域的圆形度中,人体部分骨骼区域呈现的横截面大多是不规则的几何形状,而血管横截面基本上都是规则的类圆形。因此,可以在计算区域圆形度的基础上,再加上区域的质心计算来进行判断。如果区域的质心点与区域的几何中心点相距甚远,则说明此区域的结构很不规则,一定不是血管区域。区域质心为
1.2.5 区域偏心率
偏心率可以叫做延伸度,它在一定程度上描述了区域的紧凑性,即目标外形关于区域中心的对称程度。目标的偏心率不受平移、旋转和尺度变换的影响,它只是粗略地描述目标的基本形状和中心的情况。偏心率的计算有多种形式,一种常用的简单方法是区域主轴(长轴)长度DA与辅轴(短轴)长度DB的比值,即
由人体解剖学知识可知,人体中部分骨骼的横截面区域形状是不规则的。利用这种特征,在对区域进行圆形度和质心判断后,还可以辅助判断区域的偏心率,防止单一判断造成骨骼区域被错误排除,保证骨骼区域与血管更好地分离。
经过基于特征的骨骼区域提取后,已经可以把大部分的骨骼和血管分离开来,但是还是会把一部分骨骼区域错当成血管区域,或把一部血管区域错认为是骨骼区域,使得在三维血管可视化中出现多余的骨骼区域和断层的血管。因此,还需要对以上分割结果进行层间相关性判断,以还原残余骨骼体素和修补断层血管。
图2 下肢CTA体数据部分横截面示意Fig.2 Schematic of parts CTA Date of lower limb
由于部分骨骼的亮度和几何形态与血管相似,使得在基于特征的骨骼分割提取时,骨骼与血管混在一起,造成骨骼分割错误。然而,这部分的残余骨骼体素毕竟为少数,且大多数集中于体部数据的脊柱部位和膝盖部位,范围有限,同时由于算法之前的处理均在CTA体数据的横断面上处理,因此这些松质骨体素多呈现为断层现象,即在三维坐标系中呈现X和Y轴方向连通、Z轴方向不连通的现象。
如图3(a)所示,血管在人体中呈连通状态,即无论是X、Y轴方向,还是Z轴方向,都是完全连通的,特别是人体下肢部位,血管在Z轴方向的连通层数远高于断层的骨骼区域。因此,可以对CTA体数据进行骨骼薄厚判断,即通过在Z方向检测连通区域的层厚,将层厚小于一定阈值的区域判断为骨骼区域。图3(b)是经过层间骨骼还原后的血管三维截图,可以看出,经层间厚度检测后,已经完全去除了残余的骨骼体素。
图3 骨骼体素去除前后对比。(a)断层骨骼截图;(b)去除骨骼体系体素后Fig.3 Contrast chart before and after bone removal.(a)schematic of fault bone;(b)schematic of bone removal
由于部分血管与骨骼中的骨膜相连,使得螺旋CT采集后,血管与骨骼不能形成单独的个体区域,造成无法分割,因此需要对血管进行断层修补,以恢复血管的整体拓扑结构。本研究通过血管的中心线,采用逐层检测的方式,对缺失的血管进行修复。
在逐层检测缺失血管的过程中,利用已有的血管分割结果中已经得到的血管中心点进行查找,查找的依据是血管的结构特征以及血管层间的相关性,即相邻两层扫描数据间的血管结构不会出现较大的变化。
步骤1:从当前层血管向下层进行查找,若下层血管存在,即在当前层血管结构中标识。不断重复本步骤1的过程,若下层血管不存在则进入步骤2。
步骤2:若当前层血管不存在而它的CT值符合血管的结构特征,则将这一层作为缺失层,采用步骤3对缺失层进行修补;否则认为查找到了血管的终点,进入步骤4。
步骤3:对缺失层进行修补,修补采用区域生长的方法,利用上层血管的区域范围,以限定缺失层血管的区域生长范围,从而使粘连血管能够从骨骼中分离出来。
步骤4:当查找到某条血管的终点时,需查看是否还有其他血管没有被查找,若有则进入步骤1,否则说明完成了全部血管的层间修补。
图4为使用上面介绍的基于层间相关性的血管断层修补的典型二维横断面图像,其中黑色线所圈起来的区域为分割出来的血管区域。比较图(a)和图(b)可以发现,经基于层间血管结构特征的断层血管连接后,可以很好地把与骨骼粘连的血管分离出来。
图5给出了本算法的整体流程图,其中包括基本特征的骨骼区域提取和基于层间相关性的血管析取。
为了验证本算法的有效性和准确性,采用了中国医科大学第二附属医院提供的来自不同厂商、不同排数CT扫描的胸部到下肢的49套体部CTA数据。实验程序运行的硬件环境是 E7500,2.93 GHz双核处理器,2 GB内存。图6是采用本算法进行去骨后血管分割标记效果图和血管三维显示效果图。
图4 下肢CTA体数据断层修补前后对比。(a)断层修补前;(b)断层修补后Fig.4 Contrast chart of CTA Data of lower limb before and after fault repair.(a)before fault repair processing;(b)after fault repair processing
图5 本算法整体流程Fig.5 Flowchart of bone removal
从图6中可以看到,本算法可以很好分离血管和骨骼,血管轮廓清晰、结构完整。由于本算法只用患者的一套CTA体数据进行操作,因此在临床诊断中只需对病人进行一次螺旋CTA采集,从而减少对患者的伤害,降低患者的痛苦。较其他半自动去骨算法,本方法能够实现快速的、可再现性的血管三维可视化。因为对于大量CTA骨骼数据来说,用户的手动参与会消耗大量的时间,同时也会使血管的三维可视化结果不能多次完全重复再现,对于临床诊断来说是很大的障碍。而本算法是一种完全脱离用户操作的全自动去骨算法,在对同一套数据进行操作的时候,完全可以多次重复得到相同的血管可视化图,从而实现临床诊断过程中对血管疾病的多次查诊和复诊。
图6 三维血管提取效果。(a)血管分割效果;(b)三维血管可视化Fig.6 Volume rendered image after the boneremoval.(a)3D visualization of vascular segmentation;(b)3D visualization of blood vessels
在49套CTA体数据中,对临床诊断时感兴趣的血管和骨骼目标段进行量化误差评估。评估重点是观察血管钙化、血管与骨骼存在粘连的区域,对去骨算法的准确程度划分为极好、好、一般、差4个等级。通过评估实验发现,除去2套由于造影剂太少、使血管区域亮度极弱、造成去骨效果不好以外,其余的47套数据中有45套都达到了极好和好的标准,去骨和血管分割的准确度达到了95%以上。由此可见,本方法能够有效地去除骨骼、分割血管,能够满足临床诊断的要求。
在算法实现的速度上,本算法由于不要求用户的任何输入,实现了CTA体数据的去骨和血管分割的全自动化,因此减少了运行时间,提高了效率。同时,由于本算法不涉及太多复杂的计算,计算量相对较小,对于512像素×512像素×548像素的体数据,处理时间约为1 min 20 s,可以被临床医生所接受。
在CTA体数据中,骨骼与血管具有相似的CT值,同时部分骨骼和血管还存在粘连的情况,对血管的分割造成了严重的影响。针对CTA体数据的去骨问题,提出了一种基于特征和层间相关性的CTA体数据去骨算法。这种方法只需利用一套CTA体数据,不需要用户的任何手动参与,在自动去骨的同时,可以得到完整的血管拓扑结构,实现血管的三维可视化。对49套不同CTA体数据进行实验,结果表明本算法稳定性好、准确率高、速度快,同时对患者的伤害小,可以满足临床诊断中血管三维可视化的常规要求。
[1]董建园,郝重阳.基于3D医学图像的血管三维分割[J].中国生物医学工程学报,2008,27(6):848-853.
[2]Mullick R,Avila R,Knoplioch J,et al.Automatic bone removal for abdomen CTA:A clinical review [J].Proc RSNA,2002,23(4),646-650.
[3]Suryanaranayanan,S,Mullick,R,Mallya,Y,et al.Automatic bone removal for head CTA:A preliminary review [J].Proc RSNA,2003,54(15):152-164.
[4]AlyassinAM,Avinash GB. Semi-automatic bone removal technique from CT angiography data [J].Med Imaging,Proc SPIE,2001,43(22):1273-1238.
[5]Van Tran L,Sklansky J.Flexible mask subtraction for digital angiography[J].IEEE Trans Medical Imag,1992,11(3):407-415.
[6]Bani-Hashemi A,Krishnan A,Samaddar S.Warped matching for digital subtraction of CT-Angiography studies[J].SPIE,1996,27(10):428-437.
[7]Dake MD,McDonnel CH,Song SM,et al.CT angiography with spiral CT and maximum intensity projection [J].Radiology,1992,18(5):607-610.
[8]Mullick R,Kamath VP,Mallya Y.Method and apparatus for removing obstructing structures in CT imaging[P].United States Patent:US7123760B2,2006-10-17.
[9]Thomas M.System and method for automatic bone extraction from a medical image [P]. United States Patent:US7616794B2,2009-11-10.
[10]Shiffman S,Rubin G,Napel S.Semiautomated editing of computed tomography sections for visualization of vasculature[J].SPIE,1996,27(7):140-151.
[11]Daphne Y.System and method for automatically segmenting bones in computed tomography angiography data[P].United States Patent:US7593762B2,2009-09-22.