奇异值分解滤波器在超声造影成像的应用及性能分析∗

章希睿 桑茂栋 杜宜纲 林穆清 朱磊

(深圳迈瑞生物医疗电子股份有限公司 深圳 518057)

超声造影成像(Contrast-enhanced ultrasound,CEUS)能够完成针对扫查目标的高增强、动态实时显影，可呈现病灶与正常组织在血流灌注模式及形态结构上的差异，以达到辨别病灶良恶性的诊断目的[1−3]。超声造影在肝脏肿瘤尤其是肝细胞癌(Hepatocellular carcinoma,HCC)诊断方面带来了革命性突破，诊断准确率可达90%，已广泛应用于临床，并纳入多个国家和地区的肝癌诊疗规范[4]。

从临床应用的角度来看，穿透力(灵敏度)和本底残留是医生评价增强显影能力的两个标准，对其做进一步的改善，必将利于造影临床价值的提升。在技术层面，穿透力(灵敏度)可抽象为信噪比(Signal-to-noise ratio,SNR)，SNR越高说明穿透力越好、对造影剂的检测灵敏度越高；本底残留则由造影-组织残留比(Contrast-to-tissue ratio,CTR)来衡量，该值越大说明本底残留越低。理论上，在未注射造影剂的造影图像中只存在背景噪声，但实际情况并非如此，例如图1。首先，组织自身也会产生与造影剂类似的非线性成分，易与造影剂的非线性回波一同被检测到；其次，前端电路在发射脉冲序列时难免出现幅度、相位不一致的情况，导致组织成分无法被完全抵消；此外，血管壁、脏器包膜和骨头等强反射面极易导致信号饱和，在图像上同样体现为组织残留。

图1 造影图像中的组织残留Fig.1 Residual tissue in a CEUS image

目前为止，造影图像SNR和CTR的提升主要依托于新型非线性成像方法及其配套的频域滤波处理[5−14]。超声造影剂是一种由惰性气体和脂质外壳构成的微气泡(Microbubble,MB)，其在声场激励下会产生收缩及膨胀行为，可形成频率成分较为丰富的非线性散射回波，包括二次谐波[6]、非线性基波(频率位于基波段的非线性成分)[7]、次谐波(频率为1/2倍基频的非线性成分)[8−12]和超谐波[13−14]等。较早应用于临床的是二次谐波方法，系统噪声多为白噪声，且组织残留和MB均可产生二次谐波成分导致频带重叠，配套的频域滤波对SNR和CTR的提升极为有限[6]。非线性基波和次谐波方法的投入使用[7−12]，一定程度改善了频带重叠问题，性能较二次谐波有所提升，但仍无法消除由前端电路不对称性和信号饱和导致的宽带残留信号。近年来，一类基于回波联合空时域矩阵奇异值分解(Singular value decomposition,SVD)的滤波器被广泛应用于血流成像及超分辨率微血管成像，被统称为SVD滤波器[15−19]。该类滤波器为非频域方法，不存在上述频带重叠的问题；且为自适应滤波器，处理效果与信号内容(不同的切面、运动速度等)相关，因而较传统壁滤波方法(高通滤波，且截止频率固定)具有更高的血流检测灵敏度[20−22]。另外，基于造影剂微泡定位的超分辨率微血管成像也将SVD滤波器作为微泡信号提取的主要方法，替换经典的造影脉冲序列技术以获取更高的成像帧率，而且已被证明在高频应用场合具有更优的性能[23]。图2提供了一个简单的仿真实例(采用B-mode下的平面波发射，成像帧率高达242 Hz；共采集2185帧波束合成后数据，使用SVD滤波器对其进行处理，并呈现了其中1帧数据在滤波前后的效果)，图2(a)为含有造影剂的二维B-mode图像数据，图2(b)为图2(a)的数据经过SVD滤波后的结果。相比图2(a)，图2(b)中血管的位置及形状变得清晰可辨。所以，SVD滤波能同时达到提取微泡信号及抑制组织信号的目的。

图2 对超高帧率B图像数据做SVD滤波前后的对比Fig.2 Before and after SVD filtering with ultrafast B-mode data

鉴于此，本文考虑将SVD滤波器应用于常规造影成像，以能否进一步提升造影图像的SNR和CTR作为评价标准，明确在造影脉冲序列处理后增加SVD滤波器的可行性及临床价值。

1 方法及其仿真实现

1.1 联合空时SVD滤波原理

本文采用的是最常见的联合空时SVD滤波器。其中，空域指代的是由轴向深度及横向线构成的二维空间数据，时域则体现为不同时刻采集得到的多帧二维空间数据。因此，该滤波器为多帧批处理方法，其原理详述如下。

假设共有N帧待处理数据，每次批处理其中的M帧(M≤N)，记为sk,l(t1),sk,l(t2),···,sk,l(tM)，索引k∈[1,K]和l∈[1,L]分别表示纵向第k个深度采样点及横向第l根接收线。如图3所示，为了方便进行矩阵SVD，第一步需要将尺寸为K×L×M的三维数组转换成KL×M维Casorati矩阵，记为S并可表征为

由式(1)可知，矩阵S的行维对应于数据的空间位置，列维对应于数据的采集时刻。随后，对矩阵S进行SVD，可得

其中，运算符“(·)H”表示矢量/矩阵的共轭转置；λ1>λ2>···>λP为降序排列的P个奇异值，P为矩阵S的秩；{u1,u2,···,uP}和{v1,v2,···,vP}分别是上述奇异值对应的空间奇异矢量和时间奇异矢量，维度分别为KL×1和M×1，且相互正交[24]。需要强调的是，造影成像在不同时刻的回波通常是非相干信号，此时S为列满秩，满足P=M。

图3 联合空时SVD滤波器的原理示意图Fig.3 Schematic diagram for a spatiotemporal SVD filter

从物理角度来看，SVD将一个大的信号空间表征为若干空时域子空间的线性组合，奇异值为线性加权因子，数值大小对应各子空间的强度权重与相干程度。具体而言，子空间的奇异值越大，所对应成分的强度越大、运动速度越慢(速度越慢，经多帧采样后的时域相干度越高)。SVD滤波性能的优劣取决于奇异值的合理选取。经过脉冲序列技术处理后的造影数据中，主要包括组织残留、微泡信号和噪声，本文采用了“三段式”优选策略：

(1)组织残留：多为静止或低速运动成分，且强度一般高于微泡信号和噪声。设置组织残留阈值α，当λp>α时被认为是组织残留对应的奇异值，应舍弃；

(2)系统噪声：接近白噪声，各帧数据间的噪声趋近统计独立。设置噪声阈值β(β<α)，当λp<β时被认为是噪声对应的奇异值，同样应舍弃；

(3)微泡信号：强度及运动速度介于上述二者之间。当α≥λp≥β时被认为是微泡对应的奇异值，应予以保留。被保留的奇异值及其奇异矢量分别记为、，对其进行线性加权重建，进而获取SVD滤波后的Casorati矩阵：

如图3所示，对矩阵¯S的列矢量做矩阵化处理，即可得到SVD滤波后的M帧K×L维数据：

1.2 仿真实现

考虑到联合空时SVD滤波器为批处理方法，为了后续更好地分析批处理帧数M对图像性能的影响，本文使用迈瑞公司Resona 7平台的线阵高帧率造影功能(L11-3U探头，造影帧率为48 Hz)采集了比格犬肝脏在注射造影剂后(比格犬体重为12 kg，注射SonoVue，剂量为0.05 mL/kg)的波束合成基带数据，共计103帧。图4呈现了第1帧数据在未做SVD滤波的仿真结果，图中3个不同颜色的虚线框从上到下分别标出了组织残留(绿色框)、微泡信号(黄色框)和噪声(红色框)。

图4 未做SVD滤波的犬肝脏双实时造影图像(基于波束合成基带数据的仿真结果)Fig.4 A dual-live CEUSimage of dog liver without SVD filtering(simulation results using baseband beamforming data)

首先，将SVD滤波器置于波束合成之后、包络检测之前，直接对复数信号进行处理，其余处理环节与不做SVD滤波的造影成像保持一致。在本文的仿真实例中，矩阵的SVD由Matlab函数“svd(·)”提供数值解[24]。另外，可供调整的SVD滤波系数包括批处理帧数M、组织残留阈值索引Iα以及噪声阈值索引Iβ。三者的详细定义分述如下：

(1)批处理帧数M的取值范围为[2,N]，其中N是待处理的总帧数，在本仿真实例中为103。需要指出的是，本文采用了非重叠的批处理子集划分方式(只针对M

(2)本文默认奇异值采用降序排列，索引小于组织残留阈值索引Iα∈[1,Iβ)的奇异值，应舍弃。

(3)类似地，索引大于噪声阈值索引Iβ∈(Iα,M]的奇异值，应舍弃。

综上可知，最终得以保留的是索引为Iα,Iα+1,···,Iβ的奇异值。

图5比较了第60帧数据在进行SVD滤波(M=103、Iα=2、Iβ=36)前后的图像。可见，处理后犬肝脏造影图像的近场组织残留得到了显著的抑制，且在背景噪声不变的前提下达到了增强微泡信号的效果。综合以上两点，该结果从图像层面初步验证了SVD滤波器能够同时提升造影图像的SNR和CTR。第2节将针对3个参数展开性能分析，力求得到更加明确的结论。

图5 犬肝脏造影数据在进行SVD滤波前后的造影图像对比Fig.5 Comparison between CEUSimages before and after being SVD filtered with dog liver data

2 性能分析与讨论

本节以SNR和CTR作为评价指标，分析SVD滤波器3个系数对造影成像性能的影响，并且确定其优化策略。在上述仿真实例中，共采集到103帧波束合成基带数据，其在SVD滤波前后的SNR和CTR可分别由式(5)和式(6)计算得到。

其中，ROISignal、ROINoise和ROITissue分别表示微泡、噪声和组织残留对应的信号区域，选取方法参见图4(a)中的三色虚线框。

造影成像为动态模式，其性能应考虑所有帧的结果以作综合评价。因此，本文以所有帧的SNR和CTR均值来分别评价造影成像的SNR及CTR性能。另外，在仿真过程中发现SVD滤波器会导致动态图像出现闪烁现象，因而以所有帧SNR和CTR的方差值来量化该现象的程度，并以此作为滤波器系数优化的参考标准之一。

2.1 批处理帧数

在分析批处理帧数M对性能的影响时，令Iα=2、Iβ=36。表1给出了M在不同取值下的性能指标。可见，SVD滤波器能够同时提升造影图像的SNR和CTR，但会导致图像出现闪烁现象(SNR和CTR方差变大)，其严重程度与M呈负相关。此外，批处理帧数对CTR的影响程度明显高于SNR，而且对二者的影响呈相反趋势。因此，在优化该滤波器系数时，需在SNR、CTR以及图像闪烁程度中寻求妥协与平衡。针对本文的应用实例，M=103为最佳选择，理由是SNR最高，CTR虽处于劣势，但仍可超过未作处理时的CTR约达7 dB，且图像闪烁程度最低。

2.2 组织残留阈值

在分析组织阈值索引Iα对性能的影响时，令M=103、Iβ=36。表2给出了Iα在不同取值下的性能指标。当Iα=2时，图像性能出现拐点：CTR显著提升，SNR几无下降。在此之后，随着Iα的变大，SNR和CTR呈同时下降趋势，图像闪烁程度也相应提高。图6绘制了本文应用实例中奇异值随奇异值索引变化的曲线，可以发现奇异值曲线的拐点同样出现在Iα=2处。综上，组织残留阈值的优化准则相对简单，选择奇异值曲线的拐点即可。其中，拐点的判定可采用阈值法：设置一个阈值(例如本文的阈值为10 dB)，当某个索引处的奇异值较上一个奇异值的下降量大于该阈值时，则认为是拐点。

表1 M在不同取值下的SNR及CTR指标(Iα=2、Iβ=36)Table 1 SNRs and CTRs under different values of M(Iα=2,Iβ=36)

表2 Iα在不同取值下的SNR及CTR指标(M=103、Iβ=36)Table 2 SNRs and CTRs under d ifferent values of Iα(M=103,Iβ=36)

图6 本文应用实例中奇异值随奇异值索引变化的曲线Fig.6 Singular values versus singular value indices for the application scenario in this paper

2.3 噪声阈值

在分析噪声阈值索引Iβ对性能的影响时，令M=103、Iα=2。表3给出了Iβ在不同取值下的性能指标。可见，造影图像的SNR会随着造影阈值索引的变大相应降低，而CTR受到噪声阈值的影响则十分有限。SNR方差和CTR方差的结果显示，噪声阈值的选取对闪烁程度的影响小于另外两个参数。若仅从SNR和CTR的角度来看，噪声阈值索引越小越好。但在实际选择时，不建议采取过分激进的噪声抑制策略，否则可能导致与噪声强度接近的微泡小信号被抑制。

2.4 小结与讨论

上文以SNR和CTR作为量化指标，逐一分析了SVD滤波系数对造影成像性能的影响，可总结为如下3点：

(1)舍弃最大奇异值(组织残留阈值索引Iα≥2)是造影图像CTR得以提升的基础。在满足此条件后，批处理帧数是决定CTR性能的关键因素，这一结论与文献[17–19]相吻合；且批处理帧数对CTR和SNR的影响呈现相反的趋势，原因是：当采用较小的批处理帧数时，组织残留和微泡信号更容易集中在极少数大奇异值对应的子空间里，在对大奇异值作阈值化处理时虽然强化了组织残留抑制能力(CTR提升)，但更容易出现微泡信号被当成组织残留进行滤除的情况(SNR降低)；在相同的阈值条件下，采用较大的批处理帧数，组织残留和微泡信号更偏向于分布到多个较大奇异值对应的子空间里，降低了微泡信号被错误滤除的概率(SNR提升)，但组织残留的抑制能力被相应削弱(CTR下降)。噪声阈值的选取对CTR性能几乎没有影响。

表3 Iα在不同取值下的SNR及CTR指标(M=103、Iα=2)Table 3 SNRs and CTRs under different values of Iα(M=103,Iα=2)

(2)造影图像的SNR性能由批处理帧数、组织残留阈值和噪声阈值共同决定。具体而言，批处理帧数越多，SNR性能越好；组织残留阈值索引在越过奇异值曲线拐点后，SNR呈现逐渐下降的趋势，采用奇异值曲线拐点作为组织残留阈值是一种普适性的优化策略，该结论在文献[17]也得到了印证；较为激进、严格的噪声阈值选取策略固然能够提升图像SNR，但会相应提升微泡小信号被抑制的风险，导致造影图像丢失细节信息。

(3)SVD滤波器会导致造影图像出现闪烁现象，批处理帧数是决定其程度的关键因素，组织残留阈值和噪声阈值的影响比较有限。因此，可适当增大批处理帧数，将闪烁现象控制在可以接受的范围内，同时还须考虑到批处理帧数对图像SNR和CTR的影响，力求达到三方面性能的平衡。此外，还可使用帧相关等时域平滑类算法以及采用重叠式批处理子集划分进行SVD滤波，进一步改善图像闪烁现象。

另外值得一提的是，作为一种非频域的处理方法，SVD滤波器的性能受发射频率的影响较小[23]。相比之下，脉冲序列技术更容易受到发射频率的影响[23]。以应用最为广泛的SonoVue造影剂为例，脉冲序列技术在低频条件下的造影性能优于高频条件，其应用场合受到限制。本文在脉冲序列处理后增加SVD滤波器的方法，一定程度弱化了常规造影成像性能对发射频率的要求，应用面得以拓宽。

综上所述，在造影脉冲序列处理后增加SVD滤波器，能够同时提升造影图像的SNR和CTR，尤其对后者的提升较为显著，特别适用于乳腺造影、血管造影等组织残留信号较强的临床应用场合。

3 结论

本文提出了一种在脉冲序列处理后增加联合空时SVD滤波器的造影图像增强方法，由于是非频域方法，可有效解决传统频域滤波方法难以克服的频带重叠问题，能够达到同时提升造影图像SNR和CTR的目的。基于犬肝脏造影数据的仿真实例验证了所提出方法的有效性，并且通过详尽的性能分析实验，更加全面、充分地论证了SVD滤波器在提升造影图像性能方面的优越性。此外，本文还总结得到适用于常规造影成像的SVD滤波器系数优化策略，为该方法后续的产品实现和应用打下了扎实的基础。

致谢感谢迈瑞超声Resona 7开发团队在前端数据采集与解析方面提供的技术支持，为本文的方法验证及性能分析奠定了坚实的基础。