基于模糊先驗(yàn)的懲罰最大似然重建算法

張 麗

(山西大同大學(xué)物理與電子科學(xué)學(xué)院,山西大同037009)

0 引 言

近年來(lái),由PET投影數(shù)據(jù)重建圖像一直是國(guó)內(nèi)外研究的熱點(diǎn)。

目前重建算法主要分為解析法和迭代法。解析法的代表為濾波反投影重建(FBP),其在CT重建中有著廣泛的應(yīng)用,具有形式簡(jiǎn)單和成像速度快的優(yōu)點(diǎn),但是無(wú)法有效的抑制噪聲[1]。迭代法的代表算法為最大似然估計(jì)(Maximum Likelihood,M LEM)重建算法[2]和最小二乘(LS)重建算法[3]等。但是他們由于缺乏有效的先驗(yàn)約束,在迭代重建過(guò)程中會(huì)出現(xiàn)噪聲積累,從而所造成圖像質(zhì)量退化。為了解決迭代中圖像質(zhì)量惡化的問(wèn)題,一些學(xué)者提出在圖像迭代更新過(guò)程中引入合適的圖像先驗(yàn)信息,這種先驗(yàn)具有抑制噪聲、平滑圖像的作用[4],常用的先驗(yàn)為二次型先驗(yàn)。雖然二次先驗(yàn)結(jié)構(gòu)簡(jiǎn)單,但是容易過(guò)度平滑圖像,丟失豐富的邊緣信息,而非二次先驗(yàn)具有保持邊緣的作用。可是先驗(yàn)函數(shù)中的閾值選擇及其麻煩。鑒于此,本文在Hubber先驗(yàn)中運(yùn)用模糊數(shù)學(xué)知識(shí)簡(jiǎn)化閾值的選擇,形成了一個(gè)新的重建算法--基于模糊先驗(yàn)的懲罰最大似然重建算法。

1 懲罰最大似然算法

懲罰最大似然算法可以表示為:

其中

通過(guò)最優(yōu)化方法進(jìn)行求解,同時(shí)運(yùn)用OSL算法[5],可得迭代公式:

本文選用具有保持邊緣的Hubber先驗(yàn)[6],該先驗(yàn)屬于非二次先驗(yàn)。Hubber先驗(yàn)的函數(shù)表達(dá)式為:

其中s為像素j的一個(gè)鄰域,δ為先驗(yàn)的閾值。

2 基于模糊先驗(yàn)的懲罰最大似然重建算法

Hubber先驗(yàn)的閾值選取極其重要,現(xiàn)在常用的方法都是實(shí)驗(yàn)反復(fù)驗(yàn)證來(lái)選取合適的閾值。閾值的選取靠的是經(jīng)驗(yàn),這也增加了算法調(diào)試的難度。而本文將模糊知識(shí)與H ubber先驗(yàn)相結(jié)合來(lái)簡(jiǎn)化先驗(yàn)閾值的選取難度。

當(dāng)圖像受到噪聲污染后,圖像的不確定性在一定程度上表現(xiàn)為模糊性[7],因此,我們引入模糊熵來(lái)表示圖像信息,計(jì)算Hubber先驗(yàn)的閾值。

模糊熵的數(shù)學(xué)表達(dá)式為:

其中J為圖像的像素?cái)?shù),K值為像素?cái)?shù)的倒數(shù),μj表示像素j的隸屬度。本文選取的隸屬度函數(shù)為半升嶺形函數(shù):

其中 ‖ f‖j為像素j處的梯度模,‖ f‖m in和 ‖ f‖max分別為像素梯度模的最小值和最大值。A =‖ f‖min+‖ f‖max,B=‖ f‖max-‖ f‖min。我們經(jīng)過(guò)實(shí)驗(yàn)反復(fù)驗(yàn)證,將Hubber先驗(yàn)中的閾值定義為:,Deep為圖像的像素深度。

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

本實(shí)驗(yàn)利用shepp-logan模型產(chǎn)生所需的實(shí)驗(yàn)數(shù)據(jù)。模型大小為:128×128。投影參數(shù):128個(gè)投影方向(均勻分布在0°～180°內(nèi))并且每個(gè)方向上有128個(gè)探測(cè)器對(duì)。利用公式Y(jié)*=AX*產(chǎn)生無(wú)噪聲的觀測(cè)數(shù)據(jù),用該投影數(shù)據(jù)作為泊松變量的均值生成實(shí)際帶噪聲的投影數(shù)據(jù)[8],文中算法的初始值都為2。

3.1 評(píng)價(jià)準(zhǔn)則

本文采用兩個(gè)評(píng)估準(zhǔn)則對(duì)重建圖像進(jìn)行評(píng)估。

1)信噪比(SNR)

圖1 shepp-logan

其中μ表示重建圖像的均值。

1)歸一化均方誤差(RM SE)

3.2 重建圖像分析

本實(shí)驗(yàn)主要比較M LEM,PM L-QM,PM L-Hubber,PM L-FuzzyHubber四種重建算法的優(yōu)缺點(diǎn)。為了便于對(duì)重建結(jié)果進(jìn)行分析,三種方法都使用相同的均勻初始值2。PM L-QM,PM L-Hubber,PM LFuzzyHubber的β值分別取0.000 055,0.002,0.002 9,其中Hubber先驗(yàn)的閾值取0.002。

由圖2可知,相對(duì)于M LEM算法,PM L-QM算法重建的圖像噪聲少,但是圖像被過(guò)度平滑,丟失了邊緣信息。PM L-FuzzyH ubber和PM L-H ubber算法重建的圖像較平滑,同時(shí)保持了圖像的邊緣,有明顯的細(xì)節(jié)特征,重建圖像的質(zhì)量明顯好于M LEM和PM L_QM算法,而本文提出的PM L-FuzzyHubber算法噪聲更少,圖像質(zhì)量好于傳統(tǒng)的Hubber先驗(yàn),且參數(shù)選取簡(jiǎn)單方便。這也說(shuō)明了PM L-FuzzyHubber算法對(duì)噪聲有很好的抑制作用,且具有保持邊緣的作用,且算法簡(jiǎn)單,容易調(diào)試。

圖3 體模圖像及各算法的重建圖像截面圖

由圖3可知,PM L-FuzzyHubber算法能夠產(chǎn)生更好的結(jié)果,該算法重建的圖像更接近原圖。具有消除噪聲和保持邊緣的作用。

3.3 重建精度的比較

利用3.1中給出的兩個(gè)評(píng)估準(zhǔn)則對(duì)不同算法的重建圖像進(jìn)行定量分析,結(jié)果如表1所示。

表1 不同算法重建圖像的定量分析

圖4、5分別給出各個(gè)算法的RMSE和SNR的變化曲線。

由表1,圖4、5可知,本文提出的算法的質(zhì)量參數(shù)相比于其他算法優(yōu)勢(shì)非常明顯,這也說(shuō)明了算法的有效性。

4 結(jié) 論

基于模糊先驗(yàn)的懲罰最大似然估計(jì)重建算法(PM L-FuzzyHubber)能很好的抑制噪聲,且具有保持圖像邊緣的作用,重建質(zhì)量?jī)?yōu)于PM L-Hubber和PM L_QM算法。相對(duì)于傳統(tǒng)的Hubber先驗(yàn),FuzzyHubber先驗(yàn)無(wú)須手動(dòng)設(shè)定閾值,這也簡(jiǎn)化了算法,降低了算法的調(diào)試難度。仿真結(jié)果也證明該算法達(dá)到了預(yù)期的目標(biāo),能自動(dòng)的選擇閾值,且重建圖像要優(yōu)于傳統(tǒng)Hubber先驗(yàn),這也說(shuō)明了算法的有效性。

[1] Nam Yong Lee,Yong Choi.A modified OSEM algorithm for PET reconstru ction using wavelet processing[J].Com puter Methods and Prog rams in Biomedicine,2005,80:236～245.

[2] L.A.Shepp,Y.Vardi.M axim um likelihood reconstruction fo r em ission tomography[J].IEEE Trans Medical Imaging..1982,1 (2):113～122.

[3] Daube M E,Muehllehner G.An iterative image space reconstru ction algorithm suitab le for volume PET[J].IEEE Trans Medical Imaging..1986,5(5):61～66.

[4] H ongqing Zhu,Huazhong Shu,jian Zhou.Bayesian algorithm s for PET image reconstruction w ithmean cu rvature and Gauss cu rvature diffusion regu larizations[J].Com puters in Biology and Medicine,2007,37:793～804.

[5] Green P J.Bayesian reconstructions from em ission tomography data using amodified EM algorithm[J].IEEE Trans Medical Imaging.1990,9(1):84～93.

[6] Sakeri A lendius,U lla Ruotsalainen.Bayesian reconstru ction for em ission tomog raphy based on m edian root prior[J].Eu ropean Journal of Nuclear Medicine and Molecu lar Image.1997,3(24):258～265.

[7] A.Rosenfeld.The Fuzzy Geometry of Image Subsets[J].Pattern Recognition Letter,1984,2:311～317.

[8] 周健.正電子發(fā)射斷層的圖像重建方法研究[D].東南大學(xué)博士學(xué)位論文.2006,10.