基于Hankel-SVD的非平穩(wěn)超聲血流成像雜波抑制技術研究

王錄濤王 微金 鋼

①(電子科技大學自動化工程學院 成都 611731)

②(國家空氣動力研究與發(fā)展中心 綿陽 621000)

基于Hankel-SVD的非平穩(wěn)超聲血流成像雜波抑制技術研究

王錄濤*①王 微②金 鋼②

①(電子科技大學自動化工程學院 成都 611731)

②(國家空氣動力研究與發(fā)展中心 綿陽 621000)

有效抑制由血管或血管周圍組織時變運動引起的非平穩(wěn)雜波對于提高診斷超聲彩色血流成像中血流動力學參數(shù)描述的準確性有著極其重要的意義。該文基于奇異值濾波技術提出一種改進的非平穩(wěn)雜波自適應抑制方法。該方法逐次利用單個慢時多普勒回波采樣矢量構建Hankel矩陣，然后根據(jù)奇異值分解后得到的正交Hankel主成份所代表的頻域內容，動態(tài)選取高階Hankel主成份重構多普勒血流信號，實現(xiàn)非平穩(wěn)雜波的有效抑制。為驗證算法的有效性，分別對多普勒回波仿真模型合成數(shù)據(jù)與利用彩色超聲設備(Sonix RP)采集的頸動脈血流基帶回波信號進行濾波處理，然后采用滯一自相關估計法計算血流平均速度與功率并進行成像。處理結果表明，相對于傳統(tǒng)IIR濾波方法與多項式回歸濾波技術，利用該文所提算法可對高強度、非平穩(wěn)雜波進行充分抑制，提高血流估計精度，此外，該算法具有空間自適應性，無需人為設定閾值參數(shù)以估計雜波空間維數(shù)，與現(xiàn)有基于特征分解的自適應濾波方法相比，可以有效提高組織空間高強度時變運動時血流與組織的區(qū)分能力。

彩色血流成像；奇異值分解；自適應雜波抑制；非平穩(wěn)雜波

1 引言

在超聲彩色血流成像技術中，雜波指由血管或血管周圍慢動組織反射引起的多普勒回波分量，該分量的強度通常比血流信號強度要高40～100 dB，即使當采樣容積位于血管中心位置處，由于采樣容積的3維特性、合成聲束旁瓣泄露、聲波多重反射以及距離模糊等因素的影響，血流回波信號中也包含大量的雜波信號[1,2]。雜波信號的存在導致血流動力學參數(shù)估計偏差，嚴重降低了彩色血流成像系統(tǒng)的血流鑒別能力，因此充分抑制雜波成份對于提高彩色超聲成像質量有著至關重要的意義[3,4]。

傳統(tǒng)的雜波抑制濾波器是一種高通濾波器，包括有限沖激響應(FIR)濾波器、無限沖激響應(IIR)濾波器與多項式回歸濾波器等。在彩色超聲血流成像系統(tǒng)中，有限的多普勒回波矢量采樣數(shù)據(jù)導致FIR濾波器的幅頻特性無法滿足高強度雜波抑制需求。與FIR濾波器相比，相同階數(shù)的IIR濾波器具有更加陡峭的過渡帶滾降特性與阻帶抑制特性，但當樣本數(shù)據(jù)有限時，IIR濾波器的幅頻響應由暫態(tài)響應主導[5]。多項式回歸型濾波器采用正交多項式對雜波分量進行最小二乘擬合以消除雜波成份的影響，為雜波濾波器的設計提供了新的途徑[6]。當雜波為平穩(wěn)信號時，傳統(tǒng)靜態(tài)雜波濾波器可獲得良好的雜波抑制效果，但當組織存在加速運動時，不充分的雜波抑制將導致血流參數(shù)估計產(chǎn)生較大偏差[7]。

針對雜波信號的非平穩(wěn)特性，文獻[8]提出了基于特征值分解的雜波濾波器(Eigen filter)設計方法。該方法采用對回波多普勒信號自相關矩陣進行特征分解獲得的特征矢量構造雜波濾波器，將特征矢量作為回歸濾波器的正交基矢量可使得期望信號的擬合具有最小均方誤差性質，但特征值分解濾波器假設整個觀察范圍內的多普勒矢量信號是平穩(wěn)的，不具有空間自適應性，此外，雜波空間維數(shù)需要根據(jù)經(jīng)驗值選擇，不恰當?shù)倪x擇將導致血流參數(shù)估計嚴重偏離真實值[9?11]。

作為一種數(shù)據(jù)特征提取的有效方法，奇異值分解技術可以將數(shù)據(jù)主導特征用少數(shù)分解系數(shù)集中體現(xiàn)，具有良好的穩(wěn)定性和不變性，近年來在特征提取、信號去噪等領域的應用引起了廣泛重視[12,13]。本文介紹了一種基于Hankel奇異值分解(Hankel-SVD)的非平穩(wěn)雜波抑制技術。該技術采用單個慢時回波多普勒矢量構建Hankel矩陣，具有空間自適應性；根據(jù)Hankel主成份代表的頻率內容自適應選擇雜波主成份階數(shù)，通過摒棄雜波Hankel主成分，保留血流Hankel主成份并重構血流多普勒信號，實現(xiàn)非平穩(wěn)雜波的有效抑制。最后采用基于多普勒回波仿真模型合成數(shù)據(jù)的血流參數(shù)估計結果與人體頸動脈血流成像實例，對本文所提算法針對高強度、非平穩(wěn)雜波抑制的有效性進行了驗證。

2 Hankel-SVD濾波

超聲血流成像技術利用超聲聲束對待成像區(qū)域進行掃描，在每一個掃描方向形成有限個發(fā)射聲束，通過接收血紅細胞散射的多普勒信號估計血流動力學參數(shù)。Hankel-SVD濾波器采用KL(Karhunen-Loeve)變換將解調后的回波多普勒矢量信號表示為一組正交基矢量的線性疊加，通過保留與血流信號相對應的正交成份并進行信號重構，可有效降低雜波對血流參數(shù)估計的影響[11]。對于包含N個采樣值的慢時回波多普勒矢量x，其Hankel-SVD逼近可表示為

其中P(P＜N)為用于逼近x的正交基矢量的階數(shù)，T表示轉置，γk,?k分別為k階展開系數(shù)與正交基矢量，滿足下述正交條件：

為計算式(1)中的正交成份，將多普勒矢量x分割成部分重疊的數(shù)據(jù)段并重新組織以構造Hankel矩陣：

式(3)的奇異值分解可表示為

其中，Ak表示秩為1的第k階Hankel成份，σk,uk, vk分別為奇異值、左奇異矢量(N維)、右奇異矢量(N?P+1維)。式(4)中的P個奇異值滿足σ1＞σ2＞…＞σP，因此Ak又可稱為第k階Hankel主成份。

超聲血流多普勒回波信號由雜波、血流信號與噪聲3種成份構成，其中雜波成份是由心壁、瓣膜、血管壁及靜止的或慢速運動的肌肉引起的，平均頻率較低但強度較高，在回波信號中占主導地位，而血流信號是由血紅細胞散射的多普勒信號構成，強度較低，但頻率成份高于雜波信號[1]。文獻[6]對超聲多普勒回波自相關矩陣進行特征值分解，得到了超聲多普勒信號特征值譜，并結合超聲回波多普勒信號特性對特征值譜分布與回波信號不同成份的能量、頻率特性構成之間的關系進行了研究，認為與包含低頻成份的特征矢量相對應的較大特征值反映了雜波成份構成，與包含高頻成份的特征矢量相對

其中Kc為雜波Hankel主成份階數(shù)。

分析Hankel矩陣的構造過程可知，在構造Hankel矩陣A時，輸入信號矢量中的每個元素被沿反對角線方向重復配置，雖然經(jīng)奇異值分解，反對角線方向的元素不再相等，但仍可沿Ak反對角線方向對代表同一輸入數(shù)據(jù)的所有元素求平均值來重建血流多普勒回波信號所需的正交基矢量γk?k。應的數(shù)值較小的高階奇異值則集中體現(xiàn)了血流多普勒成份特征，而噪聲信號能量則均勻分布在所有奇異值上。因此，當雜波分量包含的Hankel主成份階數(shù)已知時，可由Ak計算得到的正交基矢量γk?k，進而按照式(1)實現(xiàn)血流多普勒信號的準確重構：

3 雜波Hankel主成份階數(shù)估計

為充分抑制雜波成份，必須準確判斷Hankel主成份Ak是否為雜波成份。基于特征值分解的雜波濾波器設計方法根據(jù)雜波成份與血流成份對奇異值譜的貢獻不同，通過將分解獲得奇異值與預先設置的經(jīng)驗閾值比較，估計雜波Hankel主成份階數(shù)[11]。當雜波成份在回波信號中處于主導時，可實現(xiàn)雜波子空間維數(shù)的準確估計，但對于高頻超聲血流成像，由于血流信號的Rayleigh散射特性，雜波與血流信號的功率比值(CBR)接近0 dB，采用閾值法將導致雜波子空間維數(shù)估計產(chǎn)生較大偏差[15]。此外，組織空間生理特性差異與組織運動的空間、時間不確定性，導致回波信號隨空間、時間變化存在巨大差異，因此很難通過一個特定的閾值準確實現(xiàn)雜波信號與血流多普勒信號的分離。

由超聲多普勒回波成份構成與奇異值分解原理可知，Hankel-SVD的正交基矢量包含了多普勒回波不同成份的頻率信息，雜波成份通常為低頻分量，而血流信號則由高頻多普勒成份構成[6,16]。因此，可根據(jù)?k所反映的頻率信息判定Hankel主成份Ak是否包含于雜波子空間，進而估計雜波Hankel主成份的階數(shù)。首先利用?k計算第k階平均多普勒頻率：

其中，TPRI為發(fā)射脈沖重復間隔。fD,1為一階Hankel主成份，可用于區(qū)分當前多普勒矢量是否由雜波成份所主導，如為雜波成份主導，則可以根據(jù)當前多普勒矢量的瞬時雜波頻率估計雜波帶寬：其中，fi,ct為瞬時雜波頻率，可由當前采樣的復多普勒矢量計算得出：

其中I與Q分別為復多普勒信號采樣值的實部與虛部。根據(jù)fD,1與fth,ct可確定雜波頻率范圍，進而對Hankel主成份Ak的性質進行判決：

當fD,k(k＞1)滿足式(9)時，說明Ak仍為雜波成份，此時，需增加k值進行下一步搜索，即將fD,k+1代入式(9)，驗證fD,k+1是否位于雜波頻率范圍內，雜波Hankel主成份階數(shù)為使式(9)不成立的最小k值，即Kc=k。增加k值使得雜波成份逐次從回波多普勒信號中抽取出來，當k＞Kc時，則雜波成份被有效剝離，Ak以血流成份為主導，根據(jù)式(5)可準確重構血流多普勒信號，而累加大于Kc的奇異值并進行歸一化，則可得到濾波后的血流平均功率：

4 血流參數(shù)估計仿真

為了對本文設計的濾波器的雜波抑制性能進行量化評估，采用文獻[16,17]提出的回波合成模型產(chǎn)生仿真多普勒回波信號，仿真參數(shù)如表1所示，濾波處理后的數(shù)據(jù)經(jīng)滯一自相關估計處理以估計血流平均多普勒頻率與功率。為分析不同流速時的濾波器雜波抑制性能，設計了20組血流信號，其平均多普勒頻率在0～0.5發(fā)射脈沖重復頻率(PRF)范圍內均勻分布。

表1 壁濾波器性能比較仿真參數(shù)取值表

圖1給出了200次獨立仿真實驗得到的血流平均多普勒頻率與功率估計結果。在圖1中，IIR-Prj為六階投影初始化Chebyshev IIR濾波器，通帶歸一化截止頻率為0.15，為保證相近的阻帶抑制性能，Reg-Pol多項式回歸濾波器采用五階Legendre 多項式重構血流信號，Eigen濾波器雜波空間維數(shù)設為5。

圖1(a)給出了估計平均多普勒頻率比較結果，虛線表示理想估計平均多普勒頻率。從圖1(a)中可以看出采用Reg-Pol與Eigen得到的平均多普勒頻率估計精度優(yōu)于IIR-Prj，這是因為采用正交基函數(shù)重構血流信號可有效提升濾波器的雜波抑制性能，有效降低因雜波殘留導致的頻率估計偏差。而對比Reg-Pol與Eigen則可以發(fā)現(xiàn)，Eigen可以獲得更小的逼近誤差，這是由于Eigen將回波信號作為隨機過程處理，根據(jù)雜波的統(tǒng)計特性構造雜波子空間的基矢量，當雜波子空間維數(shù)小于信號空間維數(shù)時，可使得雜波逼近最小均方誤差最小[6]。相對于其它3種被計較濾波算法，當輸入血流平均頻率從100 Hz增加到200 Hz時， Hankel-SVD可以顯著提高血流速度估計精度。

圖1(b)示出了平均功率估計比較結果，其中IIR-Prj濾波處理后的數(shù)據(jù)仍保留較多雜波成份，導致功率估計高于其它3種雜波濾波器。Reg-Pol濾波器具有最強的阻帶抑制性能，但過渡帶較寬，而在更寬的頻率范圍內，Hankel-SVD功率估計值更接近于理想值。

5 頸動脈血流成像比較

為驗證本文算法在活體診斷應用中的有效性，利用Sonix RP (Ultrasonix Inc., Vancouver, Cannada)超聲成像系統(tǒng)采集人體頸動脈基帶多普勒回波數(shù)據(jù)，探頭中心頻率為6.6 MHz，脈沖重復頻率為2.5 kHz，發(fā)射脈沖重復次數(shù)為10，系統(tǒng)采樣頻率為40 MHz，軸向采樣750點，側向掃描50線。對雜波抑制后的采樣數(shù)據(jù)，采用滯一自相關估計法計算每個空間采樣點的血流速度與功率并進行成像。

圖2 給出了一個完整心動周期50幀數(shù)據(jù)中心臟收縮期的一幀數(shù)據(jù)經(jīng)4種雜波濾波器處理后的，未設置速率門限與功率門限的彩色多普勒血流成像結果。 其中，IIR-Prj為二階、歸一化截止頻率為0.05的投影初始化IIR濾波器，Reg-Pol為二階Legendre多項式回歸濾波器，構建Eigen濾波器自相關矩陣的慢時多普勒矢量數(shù)為750，雜波空間維數(shù)固定為2。

在圖2中，不同的顏色代表不同的速度估計結果。在遠離血管區(qū)域雜波成份被充分抑制，噪聲成份使得速度估計趨向于系統(tǒng)抗混疊速度，在圖像中表現(xiàn)為紅、藍交錯分布。綠色區(qū)域表示未能被充分抑制的雜波成份，在血管壁、血管附近組織以及靠近內血管壁的血流區(qū)域均有所分布。對比圖2可知，圖2(d) 中雜波殘留區(qū)域最小，圖2(c)中最大，圖2 (b)優(yōu)于圖2(a)。這是因為Reg-Pol具有更加陡峭的過渡帶，因而雜波性能優(yōu)于IIR-Prj，但在血管壁上緣區(qū)域仍殘留一定數(shù)量的雜波； Eigen濾波器假設雜波沿聲束軸向是平穩(wěn)的，當過多的多普勒矢量用于構建自相關矩陣且濾波器階數(shù)選擇不當時，雜波抑制性能嚴重惡化，大量殘留的雜波導致血流估計速度明顯低于其它3種方法；Hankel-SVD根據(jù)慢時方向單個多普勒回波矢量的雜波成份所反映的頻率特性，動態(tài)調整Hankel雜波主成份階數(shù)，因而能有效抑制由于組織空間高強度時變運動導致的非平穩(wěn)雜波，提高血流與組織區(qū)分能力。

圖3為濾波后功率多普勒成像比較，顯示為歸一化多普勒功率，濾波器參數(shù)設置同圖2。從圖3中可以看出，經(jīng)4種濾波器進行濾波處理后，雜波成份均得到了不同程度的抑制，遠離血管區(qū)域的雜波抑制最為充分，而在血管壁及其附近區(qū)域，血管的非勻速搏動使得在血管壁邊緣仍殘留一定部分的較強功率雜波。圖3(d)為Hankel-SVD濾波處理結果，可以看出雜波功率抑制最為充分，因而血流區(qū)域平均功率值明顯高于其它3種濾波算法。

圖1 濾波器雜波抑制性能比較

圖2 彩色多普勒血流成像

圖3 功率多普勒成像

6 結論

針對超聲彩色血流成像中的非平穩(wěn)雜波，本文提出了一種改進的基于奇異值分解的雜波抑制算法。該方法逐次采用慢時方向單個回波多普勒矢量構建Hankel矩陣，利用Hankel-SVD分解獲得的Hankel主成分作為回歸濾波器的正交基函數(shù)重構血流多普勒信號，可使雜波擬合最小均方誤差最小，因而雜波抑制性能優(yōu)于傳統(tǒng)采用多項式函數(shù)作為正交基函數(shù)的回歸型雜波抑制方法。根據(jù)Hankel主成份包含的頻率內容區(qū)分雜波與血流信號，當雜波成份為主導時，通過估計雜波瞬時帶寬與一階Hankel主成分的平均多普勒頻率確定雜波多普勒頻率范圍，并與高階Hankel主成分平均多普勒頻率進行比較，動態(tài)選擇雜波Hankel主成份階數(shù)，實現(xiàn)非平穩(wěn)雜波的有效抑制。基于多普勒回波模型合成數(shù)據(jù)的仿真實驗與頸動脈血流與功率成像結果表明該算法在有效抑制由血管加速運動導致的非平穩(wěn)、強雜波的同時充分保留血流成份，改善低速血流參數(shù)的估計精度。由于僅采用單個慢時方向的多普勒回波矢量構建Hankel矩陣，本文算法具有空間自適應性，與Eigen雜波濾波器相比，無需人為設定閾值參數(shù)以估計雜波空間維數(shù)，因而可有效提升組織空間高強度時變運動時的血流成像質量。

[1] Oeltze S, Lehmann D J, and Kuhn A. Blood flow clustering and applications in virtual stinting of intracranial aneurysms[J]. IEEE Transactions on Visualization and Computer Graphics, 2014, 20(5): 686-701.

[2] Park G Y, Yeo S M, and Yoon C H. New adaptive clutter rejection based on spectral decomposition and tissue acceleration for ultrasound color Doppler imaging[C]. IEEE International Ultrasonics Symposium, Prague, 2013: 1484-1487.

[3] Masuda K, Bossard A, and Sugano Y. Reconstruction and error detection of blood vessel network from ultrasound volume data[C]. IEEE 26th International Symposium on Computer-Based Medical Systems, Porto, 2013: 497-501.

[4] Yap C H, Thiele K, and Wei Q F. Novel method of measuring vascular regurgitation using three-dimensional nonlinear curve fitting of Doppler signals within the flow convergence zone[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 2013, 60(7): 1295-1311.

[5] Kadi A and Loupas T. On the performance of regression and step-initialized iir clutter filters for color Doppler systems in diagnostic medical ultrasound[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 1995, 42(3): 927-937.

[6] Bjaerum S, Torp H, and Kristoffersen K. Clutter filters design for ultrasound color flow imaging[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 2002, 49(2): 693-704.

[7] Yoo Y and Kim Y. New adaptive clutter rejection for ultrasound color Doppler imaging: in vivo study[J]. Ultrasound in Medicine & Biology, 2010, 36(3): 480-487.

[8] Kruse D E and Ferrara K W. A new high resolution color flow system using an eigen decomposition-based adaptive filter for clutter rejection[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 2002, 49(12): 1384-1399.

[9] Yu A C H. Eigen-based signal processing methods for ultrasound color flow imaging[D]. [Ph.D. dissertation], University of Toronto, 2007.

[10] Song F, Zhang D, and Gong X. Performance evaluation of eigendecomposition-based adaptive clutter filter for color flow imaging[J]. Ultrasonics, 2006, 44(12): 67-71.

[11] Yu A C H and Lovstaken L. Eigen-based clutter filter design for ultrasound color flow imaging: a review[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 2010, 57(5): 1096-1111.

[12] 劉松, 張水蓮, 李世凱, 等. 一種低信噪比損失SVD 濾波部分頻帶干擾抑制算法[J]. 電路與系統(tǒng)學報, 2013, 18(1): 70-75.

Liu Song, Zhang Shui-lian, Li Shi-kai, et al.. A low distortion algorithm based singular value decomposition for partialband interference rejection[J]. Journal of Circuits and Systems, 2013, 18(1): 70-75.

[13] 徐彥凱, 雙凱. 自適應奇異值分解瞬變信號檢測研究[J]. 電子與信息學報, 2014, 36(3): 583-588.

Xu Yan-kai and Shuang Kai. Detection of transient signal based on adaptive singular value decomposition[J]. Journal of Electronics & Information Technology, 2014, 36(3): 583-588.

[14] Yu A C H and Lovstakken L. Eigen-based clutter filters for clolor flow imaging: single-ensemble vs. multi-ensemble approaches[C]. IEEE International Ultrasonics Symposium, New York, 2007: 1101-1104.

[15] Park G Y, Yeo S M, and Yoon C H. New adaptive clutter rejection based on spectral analysis for ultrasound color Doppler imaging: phantom and in vivo abdominal study[J]. IEEE Transactions on Biomedical Engineering, 2014, 61(1): 55-63.

[16] Yu A C H, Johnston K W, and Cobbold R S C. Frequency based signal processing in ultrasound color flow imaging[J]. Canadian Acoustics, 2007, 35(2): 11-23,

[17] Yu A C H and Cobbold R S C. Single ensemble based eigenprocessing methods for color flow imaging part 1 the Hankel-SVD filter[J]. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control, 2008, 55(3): 559-573.

王錄濤： 男，1979年生，博士后，研究方向為陣列信號處理與醫(yī)學成像處理.

王 微： 男，1982年生，博士生，研究方向為雷達信號處理.

金 鋼： 男，1958年生，研究員，博士生導師，從事雷達信號處

理、圖像處理與流體力學方向研究.

Non-stationary Clutter Rejection Based on Hankel-SVD for Ultrasound Color Flow Imaging

Wang Lu-tao①Wang Wei②Jin Gang②

①(School of Automation Engineering, University of Electronic Science and Technology of China, Chengdu 611731, China)

②(China Aerodynamics Research & Development Center, Mianyang 621000, China)

Effective rejection of the time-varying clutter originating from slowly moving vessels and surrounding tissues is very important for depicting hemodynamics in ultrasound color Doppler imaging. In this paper, a new adaptive clutter rejection method based on Hankel Singular Value Decomposition (Hankel-SVD) is presented for suppressing non-stationary clutter. In the proposed method, a Hankel data matrix is created for each slow-time ensemble. Then the orthogonal principle Hankel components can be obtained through the SVD of the Hankel data matrix. It achieves non-stationary clutter suppression by reconstructing the flow signal with only the high order principle Hankel components, which are estimated from the frequency content carried by the principle Hankel components. To assess its efficiency, the proposed Hankel-SVD based method is applied to synthetic slow-time data obtained from a Doppler flow model and carotid arterial complex baseband data acquired by a commercial ultrasound system (Sonix RP). The resulting flow and power images show that the proposed method outperforms the traditional IIR and polynomials regression filter in attenuation of high intense non-stationary clutter signal. It is also adaptive to highly spatially-varying tissue motion and can automatically select the order of the filter, which leads to improved distinguishing between blood and tissue regions compared to other eigen-based filters.

Color flow imaging; Singular Value Decomposition (SVD); Adaptive clutter rejection; Non-stationary clutter

TN911.7; R445.1

： A

：1009-5896(2015)04-0830-06

10.11999/JEIT140893

2014-07-07收到，2014-11-24改回

*通信作者：王錄濤 wltuestc@163.com