左心室超声图像分割及特征分析
2012-02-28刘艳丽
刘艳丽
(承德医学院,河北 承德 067000)
心功能的定量测定是临床对超声心动图检查的一个基本要求。但是由于超声心动图的图像质量较差,噪声高、分辨率低、对比度低、边界缺失和边界不连续等阻碍了心功能超声定量测定的发展。要对左心室心功能进行定量分析,准确地勾画左心室的心内膜边界则是首要条件[1],而较差的超声图像质量使心内膜面的勾画变得很困难。KRIPS R[2]等提出了一种离散小波变换结合time-space算法来实现心肌组织与血液的鉴别,这样不需要事先确定ROI区域,虽实现了鉴别但并没有确切地提取出内膜的边界;LACERDA S G[3]等提出了一个结合图像处理技术的径向搜索算法对左心室进行分割;姚等提出了一个新的全自动勾画心内膜算法[4],该算法应用中心起源的心形搜索,用改进的Sobel算子进行放射形边界探测,并且依据相邻离心距差进行边界点修正。通过对12幅质量不同的超声图像进行实验,准确率达到了56%,已经是比较好的效果。
心脏的运动是永不停息的,所以仅对某一时刻的静态图像分析是远远不够的,本文采用高解析度的M型超声成像,应用类似腹部超声显像中所使用的标准的B扫描声像仪,用心电图作为显波显像的控制闸门,记录心动周期中每一个舒张或收缩状态,集合许多心动周期成一个断面图像[5]。本文以此作为原始材料进行研究,在进行自适应阈值的简单预处理后,应用改进的CV水平集方法进行实验,力求从另外的方向来精确地提取左心室的内、外膜边缘,进而进行参数测量等后续工作。
1 改进的CV水平集算法
活动轮廓模型是基于边缘像素信息的,传统的CV模型是基于区域信息的,但是前提是区域是同质的,即假设同一区域内像素强度均一[6],但实际应用中,由于技术限制或采集时的人为原因,多数医学图像都很难达到区域同质的要求,比如本文讨论的心室心肌纤维层。Li Chunming[7]等提出了一种基于区域扩展性拟合能量的水平集算法,提出了一个RSF拟合区域扩展性能量函数,以一个边界轮廓C和在近似边界两端图像强度的两个拟合函数 f1(x)和 f2(x)来定义,区域扩展的能量依赖于特定参数的核函数,利用可控范围内区域强度信息,从小的领域扩展到整个的图像域,然后将这个能量并入到变水平集公式中,最终曲线经演化获得能量函数的极小化值,应用可控范围的局部区域强度信息来计算这两个拟合函数,引导曲线扩展向目标的边界。
这里引入了一个正性的高斯核函数。
RSF局部强度拟合在给定一个中心点x处的能量函数 εFx(加权平均误差)定义如下:
其中,I(y)为图像亮度,λ1和 λ2为正常数,C 为图像区域中一个闭合边界,Ω1为边界外部,Ω2为边界内部。文中将拟合函数定位在点x附近,是因为任一点y的强度由K(x-y)计算,根据核函数的局域化性质,当点 y远离中心点x时,强度 I(y)对拟合能量的贡献会逐渐消失,高斯核 K(x-y)会同时彻底减少到 0,因此说能量函数是被邻近x的y指针的强度控制的。一般来说,当x-y>3σ,高斯核为0,因此,只有邻近点的强度才对能量函数有支配地位。当对εFx的积分极小化处理后就能找到一个边界C能对所有x的拟合能量最小化,从而获得整个目标边界。为了实现拓扑,则将积分公式转化为规则化的水平集公式进行处理。
其中,M1(φ)=H(φ),M2(φ)=1-H(φ),固定 φ,极小化能量泛函,得到 f1(x)、f2(x),它们的值由尺度 σ 控制,一般来说σ大,则图像评定的范围要大些。
实验过程中首先采用隐性活动轮廓将初始轮廓标记图上,基于此进行水平集演化。此外,参数的设定对图像分割结果影响较大。本文在参数λ1和λ2的选择上进行了部分研究,由于能量函数值是由轮廓内部和外部两部分决定,如果取 λ1=λ2,则拓扑同时发生,朝向整个图像。而为了剔除非心肌组织的干扰,需要拓扑向轮廓限制区域内,即设定了 2λ1=λ2。
2 实验及仿真结果
图1所示为一张包含4个心动周期的心室肌舒张和收缩过程的经食道胃底短轴切面超声图像。从图中可看到,中空部分的上下两侧是由很多层心肌纤维层组成的左心室心肌,这样的层结构导致左心室心肌剖面图像的灰度是明暗相间的,也就是说心肌区域是灰度不均匀的,这在图像分割上是一个很大的难题。
图1 4个心动周期的左心室剖面超声图像
首先是图像ROI区域的选择,实验先选取了一个心动周期的剖面图像,应用传统的CV模型进行仿真实验,发现传统的基于区域同质的分割方法不适用于左心室剖面的图像分割。因此,本文针对ROI图像的特点,应用了基于区域扩展性拟合能量的水平集模型分割,希望能得到好的结果。
应用改进的CV水平集能量拟合方法进行图像分割的结果如图2所示。
图2 基于区域扩展性拟合能量的水平集模型分割结果
以上仿真中σ取3.0,比较结果可得,本文应用基于隐性区域扩展性拟合能量的水平集模型分割取得了预期的结果,分割效果比较满意。实验中参数的设定包括迭代次数、λ1和 λ2,比较图 2(a)和(f),若 λ1=λ2=1,则拓扑向整个图像,这样的多目标寻找与多相位方法优缺参半,对于图 2(e)和(f),若设定 2λ1=λ2=2 和 2λ1=λ2=6的两种情况,心肌边缘清晰可见,完全滤掉了腱索和肌小梁的干扰。实验对原图像4个周期的ROI进行分析,设定 λ1=3,λ2=6,iter=400, 图 3为另外 3个心动周期的仿真结果,除了第4个周期的图像效果较差外,其他的图像都能准确得到心室上边界的内外膜。M型超声随时间变化而高速采样,时间越长,图像的效果越差,造成分割的困难,而对于左心室下半部分,从4个周期可见,边界都比较模糊,暂不做分析。基于隐性的区域扩展性拟合能量的水平集模型分割实现了准确提取上半部分的内外膜目的。这样的分割的方法同时应用于心室边界。分割还需要结合医生的经验,来识别正确的边缘,以免出现误分割。
图3 另外3个心动周期的左心室内外膜边界提取
3 左心室超声图像的特征参数分析
首先依据得到的边界信息来分析内外膜收缩和舒张的速度变化,如图4所示。利用样条插值拟合方法得到了曲线上对应的各点梯度变化,实验证明这种拟合方法更贴近左心室内外膜边界的变化。最终得到内外膜收缩期和舒张期速度情况如表1所示。
图4 心室后壁内外膜与梯度变化
表1 内外膜速度变化参数表
由表1可见,心室外膜的速度变化小于内膜的速度变化,对于心室内膜还可以得出心肌收缩阶段的平均速度变化要小于舒张阶段,舒张过程在加速,但在收缩末期附近速度变化比较均匀,而在舒张的末期会出现速度变化比较大的情况。由4个周期分割后的边界轮廓还可以得到,内膜收缩时间和舒张时间比约为 9:10,这与正常心动周期内舒张过程比收缩过程长是一致的,这样既有利于血液流回心脏,又能使心肌有充分的时间休息。
接着进行了舒张期与收缩期纹理比较。首先采用共生矩阵方法,如表2所示,从4个方向进行纹理灰度共生矩阵统计后各个值。为了统计各参数的影响,选取了4个周期的状态最后各参数取平均值。
表2 4个心动周期特征参数均值
由表2可以分析得出,除了IDM相差不大外,其他参数(如对比度、角二阶矩和相关等)均相差较大,其中在水平方向上,收缩期与舒张期的对比度相差不大,但在其他方向上有COND>CONS,这表明在取样区域内舒张期出现了大的灰度差异频移,而收缩期内的纹理相对细腻柔和,与人眼的感知相符。同时,对于角二阶矩ASMS>ASMD,也反映收缩过程中的纹理是规则变化的。两个过程图像的熵和纹理方差都有一定差异,这都反映出了心肌紧缩阶段的纹理变化比较均匀。分析引起该状况的原因主要是:收缩期心肌增厚,小梁相互交叉,心肌分层达最大值,层与层之间的力学特征基本一致,比较舒张过程中心肌之间的回声一致性下降的结果。
另一方面,由于分形维数直观上与物体表面的租糙程度相吻合,而自然界中的不同物体租糙度有很大差别,可用分形维数作为区分不同类别物体的有效参数。因此如何准确地估计计算分形维数成为研究的重点。本文应用差分盒维数法,采集到的2个心动周期的左心室舒张收缩过程AVI帧进行分析,找到整个过程的分形值的规律,然后再比较舒张过程和收缩过程心室壁的分形维信息,找到它们之间的联系。
在提取了二维B超心室心肌收缩和舒张过程中AVI帧后,选取了一个完整的心动周期的左心室图像,图 5(a)、(b)显示的是舒张末和收缩末左心室的形态及纹理。选取了一个心动周期共55帧B超图像来计算左心室的分形维数,得到的差分分形维数FD值曲线(图5(c)),可以看到分形维的变化:在心肌收缩末点时最小,舒张末取得最大值,有FDmax=1.715 7、FDmax=1.556 7。
从曲线走势来看,在收缩过程中的值是下降的,舒张过程中是上升的,也可以得出舒张过程的时间要长于收缩时间。这与之前分割出来的心室内膜边界的走势是极为相似的,与参考文献[8]中是利用声学造影剂后得出的结论基本一致。
图5 分形纹理的提取
另外,比较如图6(b)(图片来源于维基百科left ventricular pressure)心动周期各时相心室内压所示,心室在心动周期中的压力变化在收缩期和舒张期也体现的比较明显,具体为:心室肌的收缩使室内压突然增加,等容收缩期内心室内压上升的速率加快,当心室收缩结束,立即开始心室舒张,室内压会下降,在等容舒张期内左心室内压急剧下降。这一过程与图5中分析的分形值的反向变化(图 6(a)所示)趋势相关,推测可以进一步应用分形维在估计心室压力曲线的变化值。
为了比较左心室超声心动图心肌纤维层的纹理变化,本文进行了统计,计算心室收缩过程的25帧图像与舒张过程中的30帧图像各自的平均分形维数,结果如下:mean(FDS)=1.644 1 M型超声心动图最早应用于心脏结构的探测,随着现代科学技术的发展,尤其是解剖型、全方向M型超声的发展,使得M型超声图像的每秒扫描取样线大大增加,因而更能准确掌握室壁在心动周期不同时相中的活动规律,了解心壁各层在收缩期跨壁速度梯度的差异。由于图像纹理灰度的差异不同,因而医学分割方法不尽相同。本文针对结构复杂的左心室剖面,在分析前首先找出心肌组织,应用改进的CV水平集方法,基于区域扩展性拟合能量的水平集模型,通过更改参数准确得到了心室的内外膜的边界,从而为后面的特征分析奠定了基础。最后,重点对心肌收缩和舒张过程的纹理变化作比较分析,得出了心肌运动的过程中收缩期和舒张期的时间、共生矩阵各特性参数、分形维等均有差别,可为医生的临床诊断提供参考。 [1]FOLLAND E D, PARISI A F, MOYNIHAM P F, et al.Assessment of left ventricular ejection fraction and volumes by real time two dimensional ec-hocardiography:a comparison of cineangio-graphic and radionuclide techniques[J].Circ,1979, 60 (4): 760-766. [2]KRIPS R,ADAM D R.Left ventricular wall segmentation in ultrasound cross-sectional images computer in cardiology in 1996[C].IEEE Conference, 1999: 141-144. [3]LACERDA S G, DA ROCHA A F, VASCONCELOS D F,et al.Left ventricle segmentation in echocardiography using a radial-search based image processing algorithm[C].30th An-nual International Conference of IEEE EMBS 2008,2008:222-225. [4]姚炜,田建明,赵宝珍,等.关于一个全自动勾画超声左心室内膜算法的初步研究 [J].中国超声医学杂志,2002,18(10). [5]徐智章.左心室超声心动图检查[J].国际心血管病杂志:1975(5). [6]VESE L A,CHAN T F.A multiphase level set framework for image segmentation using the Mumford and Shah model[J].International Journal of Computer Vision, 2002,50(3):271-293. [7] Li Chunming, Kao Chiuyen, GORE J C, et al.Minimization of Region-Scalable Fitting Energy for Image Segmentation[J].IEEE Transactions on Image Processing 2008, 17:1940-1949. [8]SARKAR N,CHAUDHUTIB B.An efficient approach to estimate fractal dimension of textural images[J].Pattern Recognition, 1992, 25(9):1035-1041.