基于廣義Nakagami分布的醫(yī)學(xué)超聲圖像去斑點(diǎn)噪聲算法

侯建華，朱淑琴，陳少波

(中南民族大學(xué)電子信息工程學(xué)院，武漢430074)

在醫(yī)學(xué)臨床診斷中，超聲圖像與其他醫(yī)學(xué)圖像相比具有非侵入式、使用方便、成本低、對(duì)人體損害小等優(yōu)點(diǎn).由于成像機(jī)制的限制，超聲圖像存在固有的斑點(diǎn)噪聲，它極大地降低了圖像的質(zhì)量，使得正常組織和腫瘤病變組織不易被區(qū)分.因此，斑點(diǎn)噪聲的抑制在醫(yī)學(xué)超聲圖像預(yù)處理過程中是具有非常高的研究?jī)r(jià)值和應(yīng)用前景.

常用的超聲圖像斑點(diǎn)抑制的算法主要有:均值濾波、中值濾波［1］、Jain 提出的同態(tài)維納濾波［2］等.基于斑點(diǎn)噪聲的局部統(tǒng)計(jì)特性，一些學(xué)者研究了自適應(yīng)的低通濾波器，如 Lee濾波器［3］，Kuan濾波器［4］，Loupas提出的自適應(yīng)權(quán)系數(shù)中值濾波［5］等.這些算法的共同點(diǎn)是取一個(gè)滑動(dòng)窗對(duì)圖像進(jìn)行濾波處理，因此濾波后圖像細(xì)節(jié)的分辨率很差.

小波變換良好的時(shí)頻特征和多尺度特性成為圖像去噪的一個(gè)強(qiáng)有力的工具.根據(jù)超聲圖像斑點(diǎn)噪聲的乘性特點(diǎn)，通常先對(duì)含斑噪圖像進(jìn)行對(duì)數(shù)變換，將乘性噪聲轉(zhuǎn)變?yōu)榻频募有愿咚拱自肼?，然后利用常?guī)的小波去噪技術(shù)對(duì)系數(shù)進(jìn)行修正.在對(duì)數(shù)變換和冗余小波變換的框架下，本文選用雙邊廣義Nakagami模型［6，7］、高斯模型分別描述斑噪小波系數(shù)和信號(hào)小波系數(shù)，利用最大后驗(yàn)概率準(zhǔn)則對(duì)每個(gè)高頻子帶的小波系數(shù)進(jìn)行斑點(diǎn)抑制.

1 斑噪與信號(hào)的小波系數(shù)統(tǒng)計(jì)模型

首先建立噪聲模型，超聲圖像一般使用乘性模型描述:

其中y，x分別代表含斑信號(hào)和無噪信號(hào)，n為乘性斑點(diǎn)噪聲隨機(jī)變量，且與無噪信號(hào)相互獨(dú)立.

對(duì)(1)式取對(duì)數(shù)變換，將乘性噪聲轉(zhuǎn)變?yōu)榻频募有愿咚拱自肼?

為了簡(jiǎn)便，公式(2)可寫成:

1.1 斑點(diǎn)噪聲的小波系數(shù)統(tǒng)計(jì)模型

很多統(tǒng)計(jì)模型都能描述不同散射信號(hào)情況下斑點(diǎn)噪聲的特性，包括最常用的瑞利分布以及非瑞利分布，如K分布，Weibull分布，Nakagami分布等.本文采用雙邊廣義Nakagami分布能夠更好的描述超聲圖像各層子帶的斑噪小波系數(shù)的統(tǒng)計(jì)特性，其概率密度函數(shù)(pdf)［6，7］如下:

其中:

其中σ2N(σN＞0)是斑噪小波系數(shù)N的方差，s是形狀調(diào)節(jié)參數(shù)，m是Nakagami參數(shù)，噪聲能量Ω決定了函數(shù)衰減的速度.而參數(shù)m和s可以由最大似然估計(jì)法從各個(gè)高頻子帶系數(shù)求得:

函數(shù)ψ(v)是以v為變量的雙伽馬函數(shù)，v可由(8)式求得:

1.2 信號(hào)的小波系數(shù)統(tǒng)計(jì)模型

很多學(xué)者已經(jīng)指出自然圖像的小波系數(shù)直方圖呈現(xiàn)出很長(zhǎng)的重拖尾并且在零點(diǎn)處出現(xiàn)急劇的尖峰.基于這種分布特性，本文采用均值為0的廣義高斯分布(GGD)對(duì)信號(hào)小波系數(shù)建模.GGD的概率密度函數(shù)［8］如下:

其中:

高斯分布是當(dāng)β=2時(shí)GGD模型的特例，其概率密度函數(shù)為:

2 基于Bayes估計(jì)的去斑算法

2.1 基于MAP的Bayes估計(jì)

在貝葉斯估計(jì)理論中，給出觀測(cè)信號(hào)Y，那么信號(hào)X可以由最大后驗(yàn)概率估計(jì)(MAP)求得:

公式(3)代入，可得:

斑噪小波系數(shù)建模為雙邊GND，信號(hào)小波系數(shù)建模為高斯密度函數(shù)，即公式(4)和公式(11)代入公式(13)可得:

則可推出:

其中A=msσ2X y2s－2，B=2s(2s－1)，C=(2ms－1)ΩσX，公式(15) 被稱為 GNDShrink.

2.2 模型參數(shù)的估計(jì)

斑點(diǎn)噪聲方差 σ2N用 Donoho［9，10］提出的魯棒性中值定理求得:

由于X和N都是獨(dú)立隨機(jī)變量且均值為0，故有:

那么可以得到模型邊緣標(biāo)準(zhǔn)差:是含噪子帶的方差估計(jì)，可由(19)式得到［6］:

其中，Rpj是以p為中心以wj×wj為大小的滑動(dòng)窗口，窗口大小wj是根據(jù)圖像分解層數(shù)來定義的，例如最高層的系數(shù)方差很大時(shí)選擇較小的滑動(dòng)窗可以有效的去除噪聲，粗糙層的系數(shù)方差很小時(shí)選擇較大的滑動(dòng)窗可以保存更多的圖像細(xì)節(jié).所以窗口大小可以由圖像分解層數(shù)j計(jì)算得出:

其中，t1是最高層的滑動(dòng)窗口，尺寸大小為3×3.

3 本文算法步驟

由于正交離散小波變換(DWT)不具備平移不變性，容易使去噪后的圖像出現(xiàn)振鈴效應(yīng)和偽Gibbs效應(yīng)，所以本文采用冗余小波變換(RWT)，雖然它的計(jì)算量有所增加，但是可以顯著的提高圖像去噪質(zhì)量.

基于廣義Nakagami分布的去斑點(diǎn)噪聲算法步驟如下:

⑴對(duì)含斑點(diǎn)噪聲圖像作對(duì)數(shù)變換后進(jìn)行冗余小波分解;

⑵根據(jù)公式(16)估計(jì)噪聲方差^σN;

⑶除低頻子帶外，對(duì)每個(gè)高頻子帶:

(i)由公式(5)～(8)計(jì)算出每個(gè)子帶的形狀調(diào)節(jié)參數(shù)s，Nakagami參數(shù)和噪聲能量Ω;

(ii)根據(jù)公式(18)～(20)估計(jì)模型邊緣標(biāo)準(zhǔn)差^σX;

(iii)根據(jù)公式(15)對(duì)當(dāng)前系數(shù)進(jìn)行處理;

⑷進(jìn)行冗余小波反變換(IRWT)后作指數(shù)變換可得到去噪后的圖像.

4 性能評(píng)估

4.1 去斑噪聲算法客觀評(píng)價(jià)指標(biāo)

去斑噪算法的客觀評(píng)價(jià)指標(biāo)，可以從這幾個(gè)方面來衡量:保留細(xì)節(jié)、邊緣以及紋理的能力.本文中選用以下兩種方法［6］.

⑴信噪比(SNR).

其中，f為原始圖像，^f為去斑噪后圖像.SNR越大說明去噪效果越好.

(2)相關(guān)系數(shù)(CoC).

4.2 實(shí)驗(yàn)結(jié)果及分析

實(shí)驗(yàn)圖像選用合成的無斑噪超聲圖像，如圖1(a)所示，人為地加上斑點(diǎn)噪聲后，用Daubechies4小波作4尺度的冗余小波分解，通過以下5種算法進(jìn)行比較:均值濾波，Lee 濾波［3］，同態(tài)維納濾波［2］，BayesShrink軟閾值去噪算法［8］和GNDShrink算法.

從表1中可以看出，對(duì)人工合成無斑噪超聲圖像分別加噪聲標(biāo)準(zhǔn)差為0.6和0.9的斑點(diǎn)噪聲時(shí)，GNDShrink算法在信噪比SNR和相關(guān)系數(shù)CoC的客觀評(píng)價(jià)指標(biāo)均優(yōu)于其他4種算法.均值濾波，Lee濾波［3］，同態(tài)維納濾波均屬于空域?yàn)V波方法，其本質(zhì)是低通濾波;BayesShrink軟閾值去噪算法和GNDShrink算法屬于時(shí)域?yàn)V波方法，利用子帶領(lǐng)域系數(shù)的相關(guān)性，對(duì)小波系數(shù)進(jìn)行建模，去噪效果比較明顯，尤其是本文提出的算法濾波后的信噪比和相關(guān)系數(shù)都得到了很大的提高.

圖1是噪聲標(biāo)準(zhǔn)差為0.6時(shí)以上5種方法對(duì)人工合成斑點(diǎn)噪聲圖像去噪后的對(duì)比.圖2(a)選用的是大小為的真實(shí)帶斑點(diǎn)噪聲超聲圖像，圖2(b)-(g)用以上5種算法對(duì)其進(jìn)行去噪處理后得到的圖像.從主觀視覺上可以看出，均值濾波、Lee濾波、同態(tài)維納濾波和BayesShrink軟閾值去噪算法雖然濾除了斑點(diǎn)噪聲，但造成圖像過多細(xì)節(jié)的丟失;GNDShrink算法與另外4種算法相比，在有效地抑制斑點(diǎn)噪聲的同時(shí)，能更多地保留圖像的紋理和邊緣.

圖1 (σ=0.6)各種去斑算法的主觀視覺效果比較Fig.1 Subjective visual comparison of various despeckling algorithm

圖2 各種去斑算法的主觀視覺效果比較Fig.2 Subjective visual comparison of various despeckling algorithms

5 結(jié)語

本文實(shí)現(xiàn)了一種具有自適應(yīng)性的去斑噪算法.首先，對(duì)含斑噪圖像進(jìn)行對(duì)數(shù)變換和冗余小波分解變換，然后分別用廣義Nakagami模型描述斑噪小波系數(shù)的統(tǒng)計(jì)特性，用高斯模型描述子帶中的信號(hào)小波系數(shù)統(tǒng)計(jì)特性，最后利用貝葉斯最大后驗(yàn)概率估計(jì)準(zhǔn)則推導(dǎo)出該算法.通過人工合成加斑噪圖像以及真實(shí)含斑噪圖像的客觀評(píng)價(jià)指標(biāo)和主觀視覺上的對(duì)比，該算法比另外4種算法都能更多地保留圖像的紋理和邊緣.

［1］Chen M S，Han J.W，Yu P S.Data mining:an overview from a database perspective［J］.IEEE Transactions on Knowledge and Data Engineering，1996，8(6):866-883.

［2］Jain A K.Fundamentals of digital image processing［M］.Englewood Cliffs:Prentice Hall，1989.

［3］Lee JS.Refined filtering of image noise using local statistics［J］.Computer Graphic and Image Processing，1981，15(1):380-389.

［4］Kuan D T，Sawchuk A A，Strand T C，et al.Adaptive noise smoothing filter for images with signal-dependent noise［J］.IEEE Trans Pattern Analysis and Machine Intelligence，1985，7(2):165-177.

［5］Loupas T，McDicken W N，Allen P L.An adaptive weighted median filter for speckle suppressin in medical ultrasonic images ［J］.IEEE Trans Circuits Systems，1989，36(1):129-135.

［6］Gupta S，Kaur L，Chauhan R C，et al.A versatile technique for visual enhancement ofmedical ultrasound images［J］.Digital Signal Processing，2007，17(3):542-560.

［7］Shankar P M.Ultrasonic tissue characterization using a generalized Nakagamimodel［J］.IEEE Trans Ultrason Ferroelectr Frequency Control，2001，48(6):1716-1720.

［8］Chang G，Yu B，Vetterli M.Adaptive wavelet thresholding for image denoising and compression［J］.IEEE Trans Image Proc，2000，9(9):1532-1546.

［9］Donoho D L，Johnstone I M.Ideal spatial adaptive via wavelet shrinkage［J］.Biometrika，1994，81(3):425-455.

［10］Donoho D L，Johnstone I M.Adapting to unknown smoothness viawavelet shrinkage［J］.Journalof the A-merican Statistical Assoc，1995，90(12):1200-1224.