期刊VIP學術指導 符合學術規范和道德
保障品質 保證專業,沒有后顧之憂
摘 要:低劑量輻射是CT成像中的研究熱點和難點,而不完全角度投影又是低劑量輻射中的一種有效手段。因此對不完全角度投影的重建問題也是一個研究熱點。建立在CS理論基礎上的TV-EM算法雖然能夠克服傳統算法的重建精度不高,去噪能力不明顯等特點,但也對邊緣進行了過度平滑。本文在懲罰正則項中加入中值先驗分布,在抑制噪聲的同時更能保留邊緣信息。加入中值先驗,形成新的TV-EM算法,并將其應用于不完全投影下的重建。通過仿真實驗,選擇合適的懲罰因子,驗證了其超強的抗噪能力和穩健的收斂性。
關鍵詞:不完全投影重建;稀疏懲罰;貝葉斯方法
1 引言
隨著科技的發展,CT在臨床診斷中的作用日顯重要。而人們對健康的重視程度也越來越高。因此CT在醫學診斷中的掃描所產生的射線輻射也受到人們的廣泛關注。畢竟,射線輻射對人體健康是有損害的。據報道,接受超過28次CT掃描的患者致癌幾率比平均水平高12%,而兒童在接受CT檢查時所受的輻射影響會更大。因此,人們都希望接受足夠小的輻射,而能得到精確的CT圖像輔助診斷。因此,低劑量輻射也成為醫學影像中亟待解決的熱點和難點。而有限角度和稀疏角度重建就是降低輻射劑量的一種有效手段。
在不完全角度投影數據條件下,無論是有限角度還是稀疏角度,經典FBP算法的重建結果總是不夠理想,雖然重建速度快,但是存在偽影,甚至有一部分信息嚴重缺失,導致重建圖像質量不盡如人意[1]。這就推動了迭代重建算法的研究。2006年,Candes等[2]提出了壓縮感知(CS)理論。隨后做了大量重建算法的研究,如Sidky等人提出全變分最小化算法[3,4],Chen等[5]于2008年提出PICCS算法,而Wang等[6]于2010年提出RRD算法,Zhang等[7]于2013年提出ADTVM算法,都是基于CS理論,選取一種或幾種稀疏變換作為目標函數,再以成像模型AX=P作為約束條件,再求解出重建圖像[8]。在實際成像過程中,由于數據的隨機噪聲和成像模型的非隨機誤差,造成成像模型AX=P總是不相容的。所以用成像模型AX=P作為約束條件并不總是最合適的。本文研究的TV-EM(全變差-期望最大)算法[9],是用貝葉斯方法建立目標函數,以ML-EM算法作為基礎,以圖像X的稀疏變換提取圖像的先驗知識,加入圖像的中值先驗到每一步迭代中,從而改善重建圖像的質量[10,11]。
2 獲取投影
重建物體的離散模型以網格像素(二維圖像)的形式進行刻畫,獲取投影數據模型如圖1所示[12]。探測器檢測結果為射線方向的累積值,并假設射線沒有寬度。探測器檢測的投影值與射線穿過體素的長度相關。實際上投影值是用像素內的線段長度加權的像素值的和。
這里,像素值xj為第j個像素的灰度值。探測器的第i個探測元發出的射線在每個像素內的線段長度記為aij。如果第i條射線與第j個像素不相交,則aij=0。而第i條射線的投影值pi是探測到的光子數。而投影方程可表示為:
pi=aijxj
3 貝葉斯圖像重建原理
對于一個成像問題來說,假定每一個圖像像素發出的光子數是一個獨立的泊松隨機變量。而實際上第i個探測元理論上探測到像素j發出的光子數的均值是aijxj。而所有的探測元探測到的pi組成一個樣本。那么可建立似然函數[17]:
prob=∏i(∑jaijxj) (2-1)
取對數似然
ln(prob)=∑i[-∑jaijxj+piln(∑jaijxj)-lnpi!] (2-2)
在圖像重建中,ML-EM方法就是尋找一個解,使得探測到的數據最有可能發生,即最大化概率Prob,其中pi為探測到的投影數據,而xj是待求圖像。可用如下優化目標來表示:
=argmaxx≥0ln(prob) (2-3)
為求解目標(2-3),Shepp和Vardi提出ML-EM算法,迭代公式[16]為:
j=∑i (2-4)
其中pi中是投影數據,xj是圖像當前估計值。
4 TV-EM算法
TV約束是基于Donho,Candes和Tao等提出的壓縮感知理論,假定圖像具有分段光滑性質。選擇圖像具有稀疏性的TV范數作為正則項,加入(2-3)作為正則約束。
根據Candes和Donoho提出的理論,圖像重建目標變為:
=argmaxx≥0ln(prob)+?琢TV(x) (3-1)
?琢——正則化參數
TV(x)——圖像的TV范數
TV(全變差)范數是圖像的梯度的l1范數,但由于使用起來不方便,習慣用l2范數近似代替l1范數,這樣就可以對TV范數求導了。但是l2范數在抑制噪聲的同時,把邊緣信息也過度平滑了。引入TV范數為稀疏懲罰項的TV-EM算法迭代公式[9]為:
j=∑i (3-2)
其中,TV范數求偏導是離散化后的公式[12]:
=
+
- (3-3)
參數?著是一個很小的數,為了避免分母為0。
5 中值先驗
中值先驗在降低噪聲的同時,也保持了圖像的邊緣,它的分布如下[18]:
R(x)=aexp-∑j
a——常數;
?茁——正則化參數;
Mj——領域中像素的中間值。
推薦閱讀:低溫建筑技術發表論文要求