張愛桃，陳小茜，肖雨，郭東敏，周旭，李連捷

河北醫(yī)科大學(xué) 醫(yī)學(xué)影像學(xué)院，河北 石家莊 050017

引言

目前，隨著醫(yī)療科技日趨發(fā)達(dá)，醫(yī)學(xué)影像技術(shù)已經(jīng)廣泛應(yīng)用于醫(yī)學(xué)診斷和臨床治療，其中CT影像檢查技術(shù)[1-2]具有圖像分辨率高、操作便利、無創(chuàng)等諸多優(yōu)點(diǎn)，適用于多種疾病的臨床診斷和病情評估。腦部疾病由于不便取樣本探查，腦部CT[3]成像成為腦部疾病臨床診斷最重要的輔助工具之一。但是由于CT成像的機(jī)制，圖像在獲取和傳輸?shù)倪^程中不可避免地要受到噪聲的污染，噪聲導(dǎo)致圖像質(zhì)量下降、組織邊界模糊、細(xì)微結(jié)構(gòu)難以辨認(rèn)，不僅影響醫(yī)生對病情的診斷和治療，而且嚴(yán)重阻礙圖像的后續(xù)處理。因此，抑制CT圖像噪聲[4-5]成為醫(yī)學(xué)圖像處理領(lǐng)域的基礎(chǔ)和熱點(diǎn)之一[6]。對于依賴信號強(qiáng)度的乘性噪聲主要采用基于投影域的方法，近似于白噪聲的加性噪聲主要采用基于圖像域的方法。在日常生活中，最常見的CT圖像噪聲是高斯白噪聲[7]。目前，針對高斯白噪聲的處理方法，主要有空間域中的以平滑為基本思想的均值濾波、高斯濾波[8]、局部濾波[9]等，此外還有頻率域中的維納濾波[10]和小波閾值收縮[11-14]等。局部強(qiáng)度統(tǒng)計特征是衡量區(qū)域內(nèi)像素間的平均相似性，但這一特征難以準(zhǔn)確辨別邊緣與其鄰近點(diǎn)之間的差異，導(dǎo)致了濾波結(jié)果中邊緣信息的模糊。小波閾值收縮算法雖能夠很好地估計信號的噪聲，但是閾值和閾值函數(shù)的選擇上存在不確定性，會造成信號的失真。2005年Buades等[15]提出了去除圖像加性噪聲的非局部均值濾波算法（Non-Local Means，NLM），該算法是利用圖像非局部結(jié)構(gòu)的相似性來去除噪聲，恢復(fù)圖像的主要幾何結(jié)構(gòu)。腦CT 圖像具有很強(qiáng)的對稱性和非局部區(qū)域相似的特性，因此本文根據(jù)腦CT圖像的特性提出一種新的結(jié)合小波變換的非局部均值腦CT圖像濾波算法的改進(jìn)方案。

1 結(jié)合小波變換的自適應(yīng)NLM算法

1.1 非局部均值濾波

非局部均值濾波[15]是考慮到圖像的非局部統(tǒng)計自相似性質(zhì)，利用圖像包含的大量重復(fù)結(jié)構(gòu)進(jìn)行去噪。通過計算非局部像素之間的相似度確定權(quán)值，然后加權(quán)平均來恢復(fù)圖像。NLM的主要特點(diǎn)：將相似像素定義為具有相同鄰域模式的像素，對像素周圍固定大小的窗口內(nèi)的信息進(jìn)行比較，而不只是比較圖像單個像素的信息，因此得到的相似性更加可靠。其基本原理如下：對于一幅含加性噪聲圖像，見公式(1)。

其中X(i)為未受噪聲污染的原始圖像；N(i)為均值為0，方差為σ2的高斯白噪聲，I為圖像域?qū)τ谌我庀袼攸c(diǎn)i，它的非局部均值評估值NL[Y](i)可以通過在圖像中所有的像素加權(quán)平均得到，見公式(2)。

其中，權(quán)值ω(i,j)由像素i,j的相似度決定，見公式(3)。

其中，Ni,Nj是分別以像素i,j為中心點(diǎn)的圖像塊，α＞0是高斯核標(biāo)準(zhǔn)偏差，為高斯加權(quán)的歐幾里得距離，高斯加權(quán)意味著距離中心像素越遠(yuǎn)，權(quán)值就越小。h為濾波參數(shù)，控制著指數(shù)函數(shù)的衰減，也就是控制著歐氏距離對權(quán)重的均質(zhì)化程度。Z(i)為權(quán)值的歸一化常量，見公式(4)。

從定義中可以看出ω(i,j)的值在區(qū)間[0,1]之間，且∑ω(i,j)=1。當(dāng)兩個像素鄰域之間的相似度越高，加權(quán)歐式距離越小，權(quán)值越大。由于這種相似度考慮了像素鄰域的信息而不是用單個像素的值進(jìn)行比較，因此保持了去噪算法的穩(wěn)定性，同時高斯加權(quán)又增加了中心像素點(diǎn)的權(quán)重，使得重構(gòu)信號較好地保留了細(xì)節(jié)特征。

濾波參數(shù)h控制著權(quán)重的均質(zhì)化程度,當(dāng)h較大時，ω(i,j)均質(zhì)化程度高，平滑效果比較好。反之h較小時，均質(zhì)化程度低，圖像細(xì)節(jié)保留程度較高。因此，當(dāng)圖像噪聲水平較高時，則需要相對較大的h平滑噪聲[16]；反之當(dāng)圖像的噪聲水平較低時，就需要相對小的h去保留圖像細(xì)節(jié)。本文中濾波參數(shù)的取值由噪聲方差確定，見公式(5)。

其中，k為調(diào)節(jié)系數(shù)，σ為噪聲標(biāo)準(zhǔn)差。

1.2 仿真數(shù)據(jù)處理

為了定量分析NLM算法的濾波參數(shù)的設(shè)置及去噪效果，本文首先進(jìn)行了仿真實(shí)驗(yàn)。在仿真實(shí)驗(yàn)中，不同程度的高斯白噪聲以式(1)的形式添加到Shepp-Logan頭部模型仿真CT圖像（圖1）中，然后分別采用原始非局部均值濾波、本文算法進(jìn)行去噪處理（鄰域窗口Ni采用77，搜索窗口采用1313），比較去噪效果。為了定量地評估算法的去噪性能，采用去噪前后圖像的峰值信噪比（Peak Signal to Noise Ratio，PSNR）做為度量[17]。

圖1 Shepp-Logan頭模型仿真CT圖像

1.3 真實(shí)腦CT圖像處理

為了驗(yàn)證仿真實(shí)驗(yàn)研究的可行性，本文采集了三幅真實(shí)人腦的CT圖像，分別采用原始非局部均值濾波（h2=100σ2）[18]和本文提出的濾波算法進(jìn)行去噪處理，比較去噪效果。仿真實(shí)驗(yàn)所用圖像噪聲方差是已知的，但是在臨床實(shí)踐中腦CT圖像噪聲方差是未知的，所以還需估計腦CT圖像的噪聲方差，本文采用由Donoho和Johnstone提出的小波域噪聲方差估計[19]，并對其進(jìn)行了驗(yàn)證，證明了其可靠性。

本實(shí)驗(yàn)應(yīng)用MATLAB 2014a軟件完成算法，應(yīng)用SPSS 22.0軟件對數(shù)據(jù)結(jié)果進(jìn)行統(tǒng)計學(xué)分析。采用配對t檢驗(yàn)，對濾波前后的圖像信噪比數(shù)據(jù)進(jìn)行統(tǒng)計學(xué)分析，以P＜0.05為差異有統(tǒng)計學(xué)意義。

2 結(jié)果分析

2.1 仿真實(shí)驗(yàn)的處理結(jié)果

本文采用15組不同噪聲水平的高斯白噪聲樣本仿真圖像（W1～W15）驗(yàn)證算法的有效性。圖2是三組不同峰值信噪比圖像[W6 (PSNR=15.4658)，W10 (PSNR=7.6972)，W15 (PSNR=5.9949)]濾波前后的結(jié)果對比，圖2a、2d、2g分別為W6、W10、W15去噪之前的圖像，將這三幅圖像分別運(yùn)用原始非局部均值濾波后的效果如圖2b、2e、2h，可以看到濾波后圖像邊緣及細(xì)節(jié)出現(xiàn)模糊，且隨著信噪比的降低，濾波后圖像細(xì)節(jié)損失更加嚴(yán)重。通過本文提出的濾波算法處理后圖像（圖2c、2f、2i）邊緣清晰，隨著信噪比的降低細(xì)節(jié)模糊程度較輕，該算法對不同信噪比的圖像具有自適應(yīng)性，且PSNR明顯高于原始非局部均值濾波后的（表1），15組含有不同強(qiáng)度噪聲的圖像通過本文算法濾波后的PSNR遠(yuǎn)遠(yuǎn)高于未處理圖像及原始非局部均值濾波后圖像的PSNR（圖3）。對未處理圖像與經(jīng)過本文算法去噪后的圖像的PSNR（表2）采用配對t檢驗(yàn)進(jìn)行比較分析，差異有統(tǒng)計學(xué)意義（P＜0.001）。

表1 不同噪聲水平的圖像去噪前后PSNR比較

表2 仿真圖像去噪前后PSNR

圖2 不同噪聲水平的圖像去噪效果圖

圖3 仿真圖像去噪前后PSNR統(tǒng)計

通過仿真實(shí)驗(yàn)發(fā)現(xiàn)，本文濾波算法中的調(diào)節(jié)系數(shù)k的最優(yōu)值（即k取該值時處理后的PSNR最高）與去噪前圖像的PSNR具有線性關(guān)系，見圖4。

圖4 最優(yōu)k值與圖像PSNR的關(guān)系

根據(jù)圖像擬合出的函數(shù)關(guān)系式如公式(6)所示。

由此，調(diào)節(jié)系數(shù)k的值可以根據(jù)所需處理圖像的PSNR自適應(yīng)選取。

2.2 真實(shí)數(shù)據(jù)處理結(jié)果

本文選用了三幅PSNR水平不同的真實(shí)人類腦CT圖像來驗(yàn)證算法的有效性，并且將該算法與原始非局部均值算法濾波方法進(jìn)行了對比實(shí)驗(yàn)。

在估計腦CT圖像的噪聲方差時，由Donoho和Johnstone提出的小波域噪聲方差估計[19]，見公式(7)。

其中，MAD是HH子帶（HH子帶表示水平高頻和垂直高頻信息）小波系數(shù)幅度的中值，σ為圖像噪聲標(biāo)準(zhǔn)差。由于在最小細(xì)化對角子帶，即HH子帶中信號能量很小，所以認(rèn)為含噪圖像的HH子帶主要由噪聲組成，可以在該子帶進(jìn)行噪聲方差估計[20]。

本文先采用Shepp-Logan頭模型仿真CT圖像（圖1）對該噪聲估計方法進(jìn)行了驗(yàn)證。給圖1分別加入不同程度的噪聲，利用式(6)對估計的含噪圖像的噪聲標(biāo)準(zhǔn)差和真實(shí)噪聲標(biāo)準(zhǔn)差進(jìn)行對比，結(jié)果如圖5所示，兩者非常接近，差異無統(tǒng)計學(xué)意義（P=0.321＞0.05），所以本文采取了這種方法。

圖5 含噪圖像方差與小波域估計方差對比

在對三幅真實(shí)腦CT圖像進(jìn)行處理過程中發(fā)現(xiàn)，仿真實(shí)驗(yàn)中的調(diào)節(jié)系數(shù)k的取值規(guī)律，同樣適合該真實(shí)數(shù)據(jù)處理。這三幅含噪圖像的PSNR及k的取值如表3所示，去噪結(jié)果如圖6所示，可以看出，原始的非局部均值濾波不僅需要手動調(diào)節(jié)參數(shù)h，還容易導(dǎo)致圖像過于平滑，PSNR降低（表4），而本文提出的濾波算法，可以根據(jù)輸入圖像的PSNR自適應(yīng)選取濾波參數(shù)，在抑制噪聲的同時很好地保留了圖像的細(xì)節(jié)及邊緣，提高PSNR（表4）。所以無論從視覺效果還是圖像PSNR的提高上，本文提出的算法都具有更好的效果。

圖6 真實(shí)數(shù)據(jù)的處理結(jié)果

表3 三幅含噪圖像的PSNR及k的取值

表4 真實(shí)腦部CT圖像去噪前后PSNR的比較

3 討論

由于腦部CT圖像具有很強(qiáng)的對稱性和非局部區(qū)域相似的特性，本文提出了結(jié)合小波的自適應(yīng)非局部均值濾波算法對其存在的高斯白噪聲進(jìn)行去噪處理。仿真實(shí)驗(yàn)表明該算法能夠有效地抑制白噪聲，并在真實(shí)的腦部圖像處理中也驗(yàn)證了其可靠性。

該算法中濾波參數(shù)h和濾波窗口的設(shè)置非常關(guān)鍵。濾波參數(shù)h控制噪聲的抑制程度，當(dāng)h較大時，冪函數(shù)的衰減較慢，平滑效果比較好；反之h較小時，圖像細(xì)節(jié)保留程度較高，因此，h的取值當(dāng)由圖像噪聲水平?jīng)Q定。原始的非局部均值濾波需要手動調(diào)節(jié)參數(shù)h，從而容易導(dǎo)致圖像過于平滑，本文提出濾波參數(shù)的取值由噪聲方差確定，即h2=k×σ2，通過仿真實(shí)驗(yàn)我們得出調(diào)節(jié)系數(shù)k的最優(yōu)值與圖像的PSNR具有線性關(guān)系，并且在真實(shí)的腦部CT圖像中也得到了驗(yàn)證。處理真實(shí)腦部CT圖像時，由于噪聲方差未知，本文采用小波域噪聲方差估計法對其進(jìn)行計算，并且驗(yàn)證了其有效性。因此結(jié)合小波的自適應(yīng)非局部均值濾波算法能夠根據(jù)輸入圖像的PSNR，自適應(yīng)地設(shè)定濾波參數(shù)，在抑制噪聲的同時很好地保留了圖像的細(xì)節(jié)及邊緣。

由于腦部圖像非局部的相似性非常高，搜索窗口應(yīng)該盡量大，因此理想情況下應(yīng)該定義為整幅圖像，但是考慮到計算效率，本文中圖像的搜索窗口統(tǒng)一定為13×13。鄰域窗口的大小應(yīng)該取決于圖像的噪聲水平，信噪比高時應(yīng)減小窗口，反之應(yīng)增大窗口，但是窗口過大又會造成圖像細(xì)節(jié)的損失。本文中鄰域窗口根據(jù)經(jīng)驗(yàn)設(shè)定為7×7，無法根據(jù)噪聲水平自適應(yīng)地改變，因此鄰域窗口的大小與圖像噪聲水平的關(guān)系有待于進(jìn)行進(jìn)一步深入實(shí)驗(yàn)研究。

4 結(jié)論

本文提出的結(jié)合小波的自適應(yīng)非局部均值濾波算法能夠自適應(yīng)處理不同噪聲水平的腦部CT圖像，在有效抑制噪聲的同時保留圖像的細(xì)節(jié)及邊緣信息，提高腦部CT圖像的PSNR，是一種很有前景的腦部CT圖像去噪處理方法，值得做進(jìn)一步的深入研究。