基于超声超高帧率时空滤波图像的颈动脉血管壁运动追踪*
2022-01-26马祎欣杨雅博徐胜军贠卫国
马祎欣,杨雅博,徐胜军,贠卫国△
(1.西安建筑科技大学信息与控制工程学院,西安 710311;2.西安交通大学生命科学与技术学院 生物医学信息工程教育部重点实验室,西安 710049)
1 引 言
据统计,心脑血管疾病占据非传染性疾病死亡的主要位置[1]。在中国,心血管疾病死亡率居于首位,每5例死亡中有2例死于心血管病[3-4]。颈动脉作为向大脑供血的血管,一旦发生病变,将会直接导致大脑的血液供给不足,危害人们生命健康。其中,尤以斑块和硬化最为常见。颈动脉中斑块和硬化的发生,导致颈动脉血管壁的弹性发生改变[5-6]。因此,通过分析血管壁弹性可以确定斑块和硬化的程度[7]。对于心脑血管疾病的诊断和筛查具有非常重大的临床意义。
图1 块匹配追踪计算过程示意图
动脉血管随着心动周期有序地搏动,随着心脏泵血和静脉回流,动脉血管发生有规律地跳动,正常状况下,其跳动频率与心脏搏动频率一致[8]。动脉血管弹性与血管壁发生搏动时运动的位移相关,对于血管弹性良好的血管,血管壁搏动的位移较大,心脏泵出的血量有效地输送到四周。然而,当血管弹性发生改变时,一次泵血血管壁发生的位移将变小,直接导致输送到四周的血液不足,影响肢体功能。基于此,可以通过颈动脉搏动位移状况反映颈动脉弹性。
基于超声成像实时性的特点,可通过超声图像对颈动脉血管壁进行运动追踪[9-11]。在颈动脉的运动追踪中,块匹配方法最常用[11-13]。但是,目前对块匹配方法的使用,存在以下两个问题:通过常规B超得到的数据进行计算,这类数据帧率较低,往往只有30~40 Hz[14],通过追踪之后得到的结果并不能很好地反映精细的运动信息;块匹配方法是基于图像块之间的互相关运算,实际是计算两幅图像对应的图像块中像素信号之间的互相关值,通过互相关值的大小进行图像块的匹配。但是图像中往往存在着很多随机噪声的干扰,这将给互相关的计算带来很大的误差[15]。
针对以上问题,本研究采用超高帧率平面波发射进行超声成像[16],提高运动追踪的时间分辨率。同时,采用基于奇异值分解(singular value decomposition,SVD)的时空滤波方法[17],利用组织信号的高相关性和随机噪声的低相关性的特征,进行降噪处理,从而提高块匹配方法的准确性,减少追踪结果中的不稳定搏动。
2 超高帧率超声图像块匹配追踪
通过发射平面波进行成像,不需要多次聚焦,从而提高成像的帧率[16,18]。采集到的超声图像fi,i=1,2,3,4,……,N。
块匹配追踪方法总共分为三步:初始化,选择候选中心点,互相关追踪。整个追踪过程见图1。
2.1 初始化
2.2 选择候选中心点
从第二幅图像开始fi,i=2,3,4,……,N。
基于假设:在血管搏动的过程中,相同的时间内血管壁运动产生的位移的绝对值相等。
(1)
由理论中心点产生候选中心点(xcp,ycp)。
(2)
当图像帧率为40 Hz时,相邻帧时间间隔为0.025 s,选择变化范围为-10~10,涵盖颈动脉运动的最大范围,同时可以避免全局计算的时间复杂度较大的问题[9]。
2.3 互相关追踪
(3)
(4)
对应的,fi的匹配块为cwi。
根据上述三个步骤,首先进行初始化,以第一帧图像为基础,选择血管壁上一点,并以此点为中心点建立参考块。从第二帧开始,对每一帧图像进行计算候选中心点和互相关追踪的操作,选择出对应的中心点和匹配块,实现对N帧图像的颈动脉血管壁的运动追踪,中心点变化反映血管壁运动。
3 时空滤波方法滤除随机噪声
对于用来追踪的图像fi,i=1,2,3,4,……,N。直接进行块匹配追踪时,因为图像中存在的随机噪声会对图像的互相关计算带来很大的误差,因此需要进行滤除。
一簇高帧率超声图像中存在的信号可以分为以下两类[20]:组织的回波信号,因为组织的运动频率较低,在帧与帧之间信号具有更高的空间相关性。因此,组织信号在时间上可以看作是做低速运动,空间上具有高相关性;随机噪声信号,这部分信号是随机分布的,相邻两帧图像的相关性很低。因此,时间上可以看作在做高速运动,空间上具有低相关性。一维时间尺度结合三维空间尺度,从而实现四维的时空滤波[21]。
采用奇异值分解的方法[17,22-23],根据组织回波信号与随机噪声信号的不同,在奇异值分解之后,组织信号主要包含于较大奇异值对应的奇异向量中,随机噪声信号主要包含于较小奇异值对应的奇异向量中。将较小奇异值置为0,即可滤除随机噪声成分[24]。
时间序列图像fi(nx,nz),首先转化为二维的时空矩阵F(nx×nz,N),并对F进行奇异值分解[25-26]:
F=UΔV*
(5)
式中,Δ(N,N)为一个对角矩阵,对角线上的元素为奇异值,U(nx×nz,N),V(nx×nz,N)为标准正交矩阵,其中包含各奇异值对应的奇异向量,*表示共轭变换。
对奇异值分解后的F进行滤波处理:
Ff=FVIfV*=UΔfV*
(6)
式中,Ff为随机噪声滤除的结果,If为滤波器,是一个对角矩阵,右侧部分对角元素为0。
经过时空滤波之后,图像中的随机噪声被滤除,通过互相关计算的块匹配追踪可以得到更加准确的结果。
4 颈动脉超高帧率数据采集
采用数字化可编程超声成像实验系统(Vantage 256,美国Verasonics Inc 公司)进行高帧率超声成像。在受试者同意的前提下,对受试者颈动脉进行径向扫描成像。本研究的受试者中包含两位男性、两位女性,分别编号为:受试者1、受试者2、受试者3、受试者4,其中受试者1、2为女性,受试者3、4为男性。
超声图像采集记录颈动脉的运动,因此在数据采集过程中,要尽量避免其他运动成分的干扰。采集成像中会带来运动误差,因此,把握超声探头的手需要保持稳定;受试者自身的运动也是一个重要的干扰源,受试者保持冷静状态,避免因为紧张产生颤抖;受试者在图像采集过程中,保持屏息状态,防止呼吸运动的干扰。
平面波发射进行成像,以每秒400帧图像的成像帧率进行采集,得到一系列超高帧率图像。为了保证对比的一致性,避免两次采集时位置和受试者状态的不同带来的影响,将采集到的超高帧率图像进行十分之一降采样,得到低帧率图像。
通过对颈动脉血管壁的运动追踪可以检测到颈动脉运动的周期性,以此周期作为心动周期计算心率。为了验证改良方法的性能,与使用临床超声设备(DC-40,迈瑞)的PW模式测量的受试者心率进行对比。受试者在两次超声成像的过程中,实验室播放舒缓的纯音乐。且两次超声成像连续,数字化可编程超声成像实验系统采集结束之后,立刻使用临床超声设备进行采集。舒缓的环境和较小时间间隔使得受试者状态保持平稳且一致,排除因为受试者状态不同而对心率产生的影响。
5 不同条件下颈动脉运动追踪结果
通过对低帧率、高帧率、高帧率图像经过时空滤波方法降噪之后再进行追踪(高帧率&时空滤波)三组数组的对比,说明时间分辨率的增加可以反映颈动脉血管壁运动的细节信息,时空滤波方法进行降噪,可减少颈动脉血管壁运动追踪结果中的不平稳波动。
对受试者1颈动脉进行径向扫描成像,颈动脉发生有规律的搏动,选择颈动脉血管壁上一点为参考块中心点,设置参考块和匹配块的尺寸均为3.30 mm×2.46 mm(横向×纵向)进行块匹配追踪。见图2,图中左侧为进行追踪的超声图像,最上方为低帧率图像,对应图2 (a)、图2(d)分别为低帧率图像的追踪结果中,中心点纵、横坐标的变化状况。中间为高帧率图像,对应图2(b)、图2(e)分别为高帧率图像的追踪结果中,中心点纵、横坐标的变化状况。下方为高帧率&时空滤波的结果,在图像中可以发现大部分的随机噪声信号都被滤除,血管下管壁信号得以保留。对上管壁的成像不明显,原因在于上管壁和周围的组织信号声阻抗差异不大。对应图2(c)、图2(f)分别为高帧率&时空滤波的追踪结果中,中心点纵、横坐标的变化状况。
图2 受试者1低帧率、高帧率、高帧率&时空滤波三组追踪结果对比
对比图2(a)和图2(d)、图2(b)和图2(e)、图2(c)和图2(f),运动追踪结果中心点的横坐标数据中无周期性的特点,纵坐标的变化呈现出周期性。这是因为对于径向扫描的血管,颈动脉血管的搏动沿着径向,所以,图像的纵坐标的运动呈现周期性。横坐标为横向运动,无周期性,发生的运动变化为不规则地抖动。
对比图2(a)和图2(b),高帧率成像的时间分辨率增加,反映出的颈动脉血管壁运动情况更加细节。图2(a)中对于颈动脉运动追踪的两点之间的运动成分缺失,两次追踪得到的中心点之间的时间间隔为0.025 s,难以衡量在追踪得到的两个中心点中间某些时刻的运动细节状况。而图2(b)反映出更多的运动细节成分,追踪得到的两个中心点之间的时间间隔为0.0025 s,时间分辨率更高。对于图2(a)无法衡量运动情况的部分,图2(b)结果可以很好地反映发生的运动细节。由追踪结果可以得到颈动脉血管的运动周期。由此分析心动周期,计算心率。
见图2(a)、图2(b)、图2(c),在采集时间内,共有三个心动周期。为避免误差,计算三次心动周期的平均值,并以此进行心率的计算。低帧率时计算得到的心率值为97.3,高帧率图像计算得到的心率值为96.8。高帧率时的心率结果更接近于用临床超声设备得到的实际测量值96,见图3。
对比图2(b)、图2(c)中的虚线框部分,图2(b)中存在不稳定波动信号,这部分信号是因为相邻两帧图像之间的相关性计算被随机噪声信号所干扰,导致相邻两帧之间的追踪存在误差,从而产生波动。通过时空滤波方法进行滤除之后,对于同一位置,不稳定的波动有效减少。同时,经过时空滤波之后,运动追踪得到的中心点的坐标变化曲线的走向与时空滤波之前保持一致,计算得到的心率值也一致,表明时空滤波操作并未对原始信号造成损失。
同样的方法,对受试者2、3、4颈动脉进行径向扫描成像,并分别对低帧率图像、高帧率图像、高帧率图像进行时空滤波之后的结果进行块匹配追踪。绘制纵坐标追踪结果图并计算心率,并用临床超声设备(DC-40,迈瑞)测量心率进行对比。
图4(a)、图4(b)、图4(c)分别为受试者2低帧率图像、高帧率图像、高帧率且时空滤波的追踪结果中,中心点纵坐标的变化情况。图4(d)、图4(e)、图4(f)和图4(g)、图4(h)、图4(i)分别为受试者3、受试者4的低帧率图像、高帧率图像、高帧率且时空滤波的追踪结果中,中心点纵坐标的变化情况。对比图4(b)和图4(c),图4(e)和图4(f),图4(h)和图4(i)中的虚线框部分,可以明显地看出经过时空滤波处理之后,不稳定的波动信号减少。
图3 受试者1临床超声设备测量心率结果
图4 受试者2、3、4低帧率、高帧率、高帧率&时空滤波三组追踪结果对比
从心率值的计算结果来看,受试者2在低帧率时计算得到的心率值为75.5,高帧率时计算得到的心率值为75.0,高帧率并且时空滤波处理之后计算得到的心率值为74.5;受试者3在低帧率时计算得到的心率值为76.9,高帧率时计算得到的心率值为75.9,高帧率并且时空滤波处理之后计算得到的心率值为75.9;受试者4在低帧率时计算得到的心率值为70.6,高帧率时计算得到的心率值为64.5,高帧率并且时空滤波处理之后计算得到的心率值为64.5。见图5,受试者2、3、4通过临床超声设备测量得到的心率值分别为69、75、66。可以看出,在高帧率并且时空滤波处理情况下的测量值更接近用临床超声设备得到的实际测量值。
图5 受试者2、3、4临床超声设备测量心率结果
6 结论
针对颈动脉血管弹性的预测,本研究提出了一种基于超高帧率超声图像,通过块匹配方法追踪颈动脉血管壁的运动,利用追踪块的中心位置变化,反应血管壁运动状况。进一步通过一种基于奇异值分解的时空滤波方法进行降噪处理,滤除相邻两幅图像中相关性较低的随机噪声信号,降低这些信号对相邻图像的互相关计算造成的误差,较好地实现了对颈动脉的运动追踪。
在实验数据采集过程中,采集者把握探头的手保持稳定,避免因为探头抖动造成运动误差。受试者在采集数据的时间段内进行屏息,有效地减少了呼吸运动带来的干扰,高帧率采集数据时间较短,并未造成不适,增强了临床使用的普适性。
通过对低帧率、高帧率、高帧率且时空滤波三种情况的对比,有效地说明了高帧率采集可以通过增加时间分辨率反映出更加细节的运动过程。时空滤波降噪的处理,抑制了随机噪声信号的干扰,图像中可以看出随机噪声成分大幅度减少,使得块匹配过程更加准确,减少了匹配结果中的不稳定波动。在时空滤波之后,追踪结果中,中心点坐标的变化走势未发生改变,计算得到的心率值与未经过时空滤波的心率值相近或一致,说明了时空滤波没有对信号造成损失,进一步说明了时空滤波方法对于改进运动追踪的有效性。所有受试者中,高帧率且时空滤波的追踪结果计算得到的心率值,接近临床超声设备在受试者同样的状态下获得的值,说明了追踪算法的准确性。