超聲彈性圖像去噪方法研究

邵黨國，劉 帆，相 艷，馬 磊，易三莉，賀建峰

（昆明理工大學 信息工程與自動化學院和云南省計算機應用重點實驗室，云南 昆明 650500）

超聲彈性圖像去噪方法研究

邵黨國，劉 帆，相 艷，馬 磊，易三莉，賀建峰

（昆明理工大學 信息工程與自動化學院和云南省計算機應用重點實驗室，云南 昆明 650500）

超聲彈性成像能直觀展現各組織的彈性信息，已經成為腫瘤檢測的重要工具。由于成像系統中的固有特性，使得圖像遍布著大量偽影噪聲，降低了圖像的可辨讀性。為了獲取偽影噪聲較少且組織邊緣結構清晰的彈性圖像，本文提出一種分數階偏微分圖像去噪方法，并通過各向異性算法有效平衡了梯度閾值k與微分階數v的關系。人體體模彈性圖像和仿真彈性圖像的測試結果表明，該算法能有效改善彈性圖像質量，增強邊緣結構信息，且去噪后圖像的信噪比及噪聲比例差得到顯著提升。

超聲彈性圖像；分數階微分；各向異性擴散；去噪

0 引言

Ophir等人提出的超聲彈性成像技術在近些年得到迅速發展，已經成為現代醫學成像的新模式[1]。該技術主要是利用生物組織與相關病變組織間的彈性差異進行成像，原理如下：首先給待檢測組織施加一個內部或外部的動態/靜態壓力，受到同等壓力的軟組織塊相比硬組織塊產生的形變程度更大；分別采集組織塊形變前后的回波信號，對其進行分析，計算出各組織的位移分布，最終得到待檢測組織的應變分布及其彈性系數分布并進行成像[2]。

由于超聲成像系統的固有特性，使得彈性成像和其他超聲成像模式產生的圖像一樣，必然會受到噪聲影響，進而降低圖像的可辨讀性。彈性圖像的主要噪聲來源有兩個：解相關誤差和幅度調制誤差。解相關誤差是由壓縮、變形、平滑移動或其他運動而導致運動跟蹤不符合常規的時間延遲估計模型；幅度調制誤差是由信號的隨機波動，使得位移估計的位置不在分析窗口內。各種圖像去噪技術，如直方圖匹配算法[3]、雙邊濾波算法[4]、小波去噪算法[5]以及復合方法[6]等已經被運用到超聲彈性圖像的去噪研究中，并取得了一定的研究成果。噪聲抑制對圖像分析、邊緣檢測、圖像分割、圖像重建等都具有重要意義[7-8]。

分數階微積分算法是近幾年圖像去噪領域的重要工具，其原理是利用分數階微積分在提高信號高頻成分的同時能夠非線性保留信號中低頻成分的特殊性質，將其應用于數字圖像領域，達到抑制噪聲并保留邊緣的目的[9-10]。盡管分數階微分算法在數字圖像中得到了較好的應用，但在超聲圖像中的研究仍較少，特別是彈性圖像。為了獲取邊緣結構清晰且噪聲區域平滑的超聲彈性圖像，本文提出了一種基于分數階微分的各向異性擴散算法（FAD）。

1 分數階微分

分數階微分至今為止并沒有統一的定義，從不同角度推導可得到不同的定義表達式，目前較常用的是 Grünwld-Letnikov定義，Riemann-Liouville定義和 Caputo 定義[11-13]。由于 Grünwld-Letnikov 定義是將連續可導的微分階數由整數階推導至分數階，通過對整數階微分的差分近似遞推式求極限推衍而來，因而該定義更適合應用于圖像處理領域[14]，其推導過程如下：

可微函數f(x)一階導數定義：

其二階導數定義：

推導至函數f(x)的n階導數定義：

根據上述推導，引入 Gamma函數概念，可以得到Grünwld-Letnikov定義表達式如下：

式中，一元信號f(x)的持續期t為[a, t]，h=(t-a)/n為步長，v為微分階數。若將一元信號f(x)的持續期按單位h=1進行等分，則n=[(t-a)/h]h=1=(t-a)，可以得到一元信號f(x)的差分表達式如下：

1.1 濾波模板的構建

針對各向同性擴散模型去噪時出現的邊緣模糊問題[15]，Perona和Malik于1990年提出了一種非線性偏微分方程模型（簡稱 P-M 模型），有效的保護了圖像邊緣，該模型定義如下[16]：

式中，div為散度算子，?為梯度算子，?I表示圖像I的梯度，為擴散系數方程，I0（x,y）為原始圖像，t為表示降噪過程與擴散時間有關的時間算子。Perona等人給出了擴散系數方程，如下所示：

式中，k表示擴散閾值，可阻止方程進一步擴散，圖像梯度?I可視為邊緣檢測器。若?I遠遠大于擴散門限k值，趨于0，擴散得到抑制；反之，若?I遠遠小于擴散門限k值，）趨于1，擴散得到加強，此時可視為高斯濾波。

由一元信號的差分表達式（5）可知，表達式左側前n項的系數總和為一個非零數值，這是分數階微分與整數階微分顯著特性區別之一。圖像相鄰像素間存在著高度的自相關性，為了使濾波模板具有抗旋轉性，本文構建了具有 x軸、y軸正負方向及兩對角線這8個方向的濾波模板；同時對濾波模板進行了部分修改以得到更好的去噪效果，修改后的模板如圖1所示：

圖1 微分掩模算子Fig.1 The fractional order differentiation mask

其中v為微分階數，用上述濾波模板對圖像進行卷積運算，結果發現：當k值不變，v逐漸增加時，圖像的保邊能力逐漸增強；當v值不變，k逐漸增加時，圖像的平滑能力逐漸增強。可見擴散閾值 k和微分階數v是圖像平滑與圖像保邊的重要參數。因此，如何平衡k、v之間的關系是本文重點解決的問題。

1.2 平衡k、v關系

用 FAD算法對體模彈性圖像和數字仿真彈性圖像進行處理，分別求取 k為[1,35],v為[0.1,5.5]時

各圖像所對應的SNRe及CNRe[17-18]，結果如圖2、圖3所示。其中，橫坐標為v值（范圍為[0.1,5.5]，間隔0.01），縱坐標為k值（范圍為[1,35]，間隔1）。其中每個（a）為SNRe的結果，每個（b）為CNRe的結果，圖中顏色越亮意味著SNRe值或CNRe值越大；反之則越小。

圖2 人體體模彈性圖像去噪后的SNRe、CNRe值，橫坐標為v（范圍為0.1-5.5）縱坐標為k（范圍為1-35）（a）、（b）分別為SNRe和CNRe值Fig.2 The de-noising values of SNRe and CNRe for the phantom ultrasonic elastography.The horizontal axis represents v, and k values on the vertical. (a) and (b) are the SNRe and CNRe values respectively

圖3 數字仿真彈性圖像去噪后的SNRe、CNRe值，橫坐標為v（范圍為0.1-5.5）縱坐標為k（范圍為1-35）(a)、(b)分別為SNRe和CNRe值Fig.3 The de-noising values of SNRe and CNRe for the digital simulated ultrasonic elastography.The horizontal axis represents v, and k values on the vertical. (a) and (b) are the SNRe and CNRe values respectively.

由圖2、圖3可以看出k、v之間存在一定的關系，根據實驗所得的SNRe、CNRe值，選取處理效果較好的k、v值進行曲線擬合，結果如圖4所示，得到二者間的關系表達式如下：

式中，v為微分階數，k為擴散門限值。表達式均方根誤差為0.09721，誤差平方和為0.567，擬合系數為0.9943。

2 實驗結果

圖4 彈性圖像k、v擬合關系圖Fig.4 The curve fitting relation between k and v for elastography

為了體現FAD算法的有效性，將FAD算法與傳統PM算法進行結果對比。在PM算法中，除了濾波模板和微分階數 v有所不同，其他參數均與FAD算法保持一致（即積分常數設置為 0.2，傳導系數函數為1，迭代次數設置為100次）。圖5-6中，每個（a）為原始超聲彈性圖像；每個（b）為 PM算法k=6的結果；每個（c）為FAD算法k=6,v=0.56的結果；每個（d）為PM算法k=9的結果； 每個（e）為FAD算法k=9,v=0.85的結果。

圖5 人體體模超聲彈性圖像去噪效果對比Fig.5 The contrast of denoising in phantom ultrasonic elastography

圖6 數字仿真超聲彈性圖像去噪效果對比Fig,6 The contrast of denoising in digital simulated ultrasonic elastography

由圖5-6，可以明顯看到兩種算法均能夠平滑同質區域的偽影噪聲，但FAD算法相比PM算法在邊緣信息保留效果上占有明顯優勢。SNRe及 CNRe值可以定量分析算法的平滑能力。為了減少分數階微分估計誤差，k、v值應滿足關系式：v = 0.001509*k2+0.07144*k + 0.05085，k [6, 16]任意整數。體模彈性圖像中取一個 25*25的區域測試 SNRe及 CNRe值；數字仿真彈性圖像中取一個60*60的區域測試SNRe及CNRe值。

表1和表2為圖5、圖6各圖像的SNRe、CNRe及運算時間。通過實驗分析可知：當k值小于10時，如圖5、圖6的（b）（d），PM算法出現過度保邊現象，進而降低了SNRe、CNRe；當k值大于10時，P-M算法由于平滑過度而丟失了部分邊緣信息，顯然這并不符合圖像去噪的目的。FAD算法平衡了抑制噪聲與邊緣保留的關系，可得到較理想的去噪結果。且當k[6, 16]時，FAD算法在抑制噪聲及邊緣信息保留兩個方面具有明顯效果。若 k值小于 6，FAD算法的平滑噪聲區域能力下降，出現過度保邊的現象；若k值大于16，FAD算法的邊緣保邊能力下降，出現過度平滑，圖像部分邊緣細節信息被丟失。綜上所述，可知梯度閾值k和微分階數v之間的平衡是圖像平滑與保邊的兩個重要參數。大量實驗數據表明，基于FAD算法所提出的k、v關系式，去噪后彈性圖像噪聲抑制效果明顯，且邊緣信息也得到了較好的保留。

表1 SNRe、CNRe以及運算時間Table 1 SNRe、CNRe and the operation times

表2 SNRe、CNRe以及運算時間Table 2 SNRe、CNRe and the operation times

3 結論

超聲彈性圖像彌補了其他醫學成像技術無可比擬的優點，即能提供組織的硬度信息，成為臨床檢測的重要工具。然而超聲成像系統固有的噪聲極大地降低了圖像的質量。為了平滑偽影噪聲并保留特征區域邊緣信息，進而改善圖像質量，本文提出一種基于分數階微分的圖像去噪方法，并通過各向異性擴散算法平衡擴散閾值k與微分階數v之間的關系。人體體模彈性圖像和數字仿真彈性圖像的實驗結果顯示，本文采用的FAD算法能夠有效平滑彈性圖像偽影噪聲，提高SNRe及CNRe值；與此同時，特征區域邊緣結構信息也得到了較好的保留。相比于傳統PM算法，FAD算法占有明顯優勢。

[1] 陳蒙. 醫學X光成像中圖像去噪算法研究[J]. 軟件,2013(11): 62-63.

[2] 王春光, 劉文, 劉東權. 一種基于編碼激勵的超聲彈性成像噪聲抑制方法[J]. 計算機應用研究, 2013, 30(5): 1596-1600.

[3] 何穎妮, 邵黨國, 劉東權. 直方圖匹配算法在超聲彈性成像上的應用研究[J]. 計算機應用研究, 2013, 30(4): 1266-1270.

[4] Zhou X, Liu W, Zhang C, et al. Speckle Reduction of Ultrasound Elastography with Bilateral Filter[J]. Advanced Science Letters, 2013, 11(4): 3025-3027.

[5] 顧偉, 杜景琦. 基于小波閾值的圖像去噪[J]. 軟件, 2017,38(7): 107-110.

[6] Cui S and Liu DC. Noise reduction for ultrasonic elastography using transmit-side frequency compounding: a preliminary study[J]. IEEE Transactions on Ultrasonics, Ferroelectrics,and Frequency Control, 2011, 58(3): 509-516.

[7] 鐘彩. 邊緣檢測算法在圖像預處理中的應用[J]. 軟件,2013, 34(1): 158-159.

[8] 徐紅娟, 韓玉合, 劉文洲, 等. 一種基于Ray Casting算法的醫學影像重建系統[J]. 軟件, 2012, 33(9): 118-118.

[9] 黃果, 許黎, 蒲亦非. 分數階微積分在圖像處理中的研究綜述[J]. 計算機應用研究, 2012, 29(2): 414-420.

[10] 王福興, 蒲亦非, 周激流. 分數階微積分的應用研究[J].無線互聯科技, 2011, (8): 12-14.

[11] Scherer R, Kalla S L, Tang Y, et al. The Grünwald-Letnikov method for fractional differential equations[J]. Computers &Mathematics with Applications, 2011, 62(3): 902-917.

[12] 陳慶利, 蒲亦非, 黃果, 等. 數字圖像的0～1階Riemann-Liouville分數階微分增強模板[J]. 電子科技大學學報,2011, 40(5): 772-776.

[13] 陳慶利, 黃果, 張秀瓊, 等. 數字圖像的Caputo分數階微分增強[J]. 計算機輔助設計與圖形學學報, 2013, 25(4):519-525.

[14] 楊柱中, 周激流, 晏祥玉, 等. 基于分數階微分的圖像增強[J]. 計算機輔助設計與圖形學學報, 2008, 20(3): 343-348.

[15] 蘇永利, 張博, 張書玲. 基于各向同性逆擴散方程的指紋圖像增強方法[J]. 西北大學學報自然科學版, 2009, 39(1):9-12.

[16] 陳婷婷, 沈民奮, 楊金耀. 一種復合各向異性擴散的圖像去噪算法[J]. 數據采集與處理, 2014, 29(5): 757-763.

[17] ]Paul Liu and Dong Liu. Filter-Based Compounded Delay Estimation With Application to Strain Imaging[J]. IEEE Transactions on Ultrasonics Ferroelectrics & Frequence Control,2011, 58(10): 2078-2095.

[18] Shao Dangguo, Liu Paul, Liu DongC. Characteristic matchingbased adaptive fast bilateral filter for ultrasound speckle reduction[J]. Pattern Recognition Letters, 2013, 34(5): 463-9.

Reseach of Elastography Denoising Methods

SHAO Dang-guo, LIU Fan, XIANG Yan, MA Lei, YI San-li, HE Jian-feng
(Information Engineering and Automation, Kunming University of Science and Technology and Key Lab of Computer Technologies Application of Yunnan Province, 650500, Kunming, Yunnan)

Ultrasonic elastography which can directly show the elasticity of the organization's information has become an important tool for tumor detection. The inherent characteristics of the imaging system made the elastography distributed with a large number of artifacts. And then reduced the readability of the elastography. In order to obtain a high quality elastography with few artifacts and a clear edge structure, a fractional partial differential image denoising model was constructed, which can effectively balance k (the gradient threshold) and v (the differential order) by anisotropic algorithm. The results show that the algorithm can effectively improve the elastography quality and enhance the edge structure information. The signal to noise ratio and the noise ratio difference of the elastography after denoising are significantly improved.

Ultrasonic elastography; Fractional order differentiation; Anisotropic diffusion; Artifact denoising

TP391

A

10.3969/j.issn.1003-6970.2017.12.015

本文著錄格式：邵黨國，劉帆，相艷. 超聲彈性圖像去噪方法研究[J]. 軟件，2017，38（12）：81-86

中國博士后科學基金：(2016M592894XB)；云南省科技廳面上項目(KKS0201703015)

邵黨國(1979-)，男，博士后，碩士生導師，主要研究方向為醫學圖像處理、機器學習等；劉帆(1995-)，女，碩士研究生，主要研究方向為醫學信號圖像處理；相艷(1979-)，女，講師，博士研究生，主要研究方向為圖像處理、自然語言處理、數據挖掘、機器學習；馬磊，男，博士研究生，主要研究方向為數據挖掘與分析，機器學習；易三莉(1977-)，女，講師，碩士生導師，主要研究方向為醫學圖像處理；賀建峰(1965-)，男，教授，碩士生導師，主要研究方向為醫學圖像處理、數據分析等。