基于分區域處理的低劑量CT重建算法

趙 霞，趙金龍，趙榮格，陳 燕，桂志國，劉 祎+

(1.中北大學 山西省生物醫學成像與影像大數據重點實驗室，山西 太原 030051； 2.中北大學 信息與通信工程學院，山西 太原 030051)

0 引 言

計算機斷層掃描(CT)已被廣泛應用于各種臨床應用。雖然CT[1]檢查得到了廣泛的認可，但是當過多的X射線照射劑量停留在人體內時，會給人體健康帶來極大的隱患。因此CT檢查應遵循輻射防護最優化的ALARA(as low as reasonably achievable)原則，盡可能地保護受檢者。然而，X射線劑量的降低會使投影數據中產生大量的噪聲，不能得到較為滿意的重建圖像[2]。因此，用最小的代價獲得最佳CT重建圖像，有著重要的實際應用價值。

圖像域統計迭代重建算法是低劑量CT成像的一種方法，該算法是將所獲得的投影數據進行統計建模，然后用某種迭代方法對目標函數進行求解，其最優解為待重建圖像[3-12]。懲罰加權最小二乘算法[13-17]是最常見的一種統計迭代重建算法，在各種懲罰模型中，MRF是一種廣泛使用的模型，但MRF沒有明確地考慮非均勻性，不是一個令人滿意的模型，因此，Zhang等對整個圖像使用高斯混合MRF來考慮分布復雜度[18]。TV先驗模型也是一種常用的模型，尤其是在投影數據缺失的情況下，He等基于加權方差的加權因子與TV模型相結合提出自適應加權全變分(AWTV)模型[19]，Zhang等先把梯度保真約束項和能夠區分圖像平滑區和細節區的邊緣指示函數應用到TV模型中得到基于梯度保真項的自適應全變分模型[20]。Zhang等提出了基于絕對差值排序(ROAD)和中值先驗(MP)模型的TV低劑量CT重建算法[21]。

TV正則化方法在圖像去噪的同時能較好地保持圖像的邊緣，但是容易使圖像的平坦區域產生階梯效應。高斯MRF正則化方法使用固定的平滑系數導致圖像邊緣較為模糊。根據TV正則化和高斯MRF正則化的優缺點，本文提出了一種基于分區域處理的聯合先驗低劑量CT統計迭代重建算法。該算法首先對低劑量CT圖像進行中值濾波，再利用濾波后圖像的梯度值劃分出圖像的邊緣區域和平坦區域，然后將平坦區域中的圖像塊提取出來用高斯MRF正則項處理，邊緣區域中的圖像塊提取出來用TV正則項處理，最后將這兩種正則項作為聯合先驗應用到懲罰加權最小二乘法中，使用SOR迭代重建算法進行求解。

1 噪聲模型

在低劑量CT重建中，投影數據和噪聲模型是至關重要的一個環節。以往的研究揭示了CT投影數據噪聲產生的兩個主要來源，即X射線量子噪聲和系統電子噪聲。探測器檢測獲得的對數變化前的CT投影數據可以用統計獨立的泊松分布加上統計獨立的高斯或正態分布來描述

(1)

文獻[17]中指出，當投影數據中有電子噪聲時，噪聲方差可以定義為

(2)

2 基于分區域處理的低劑量CT統計迭代重建算法

2.1 懲罰加權最小二乘算法

參見文獻[19]中的加權最小二乘(weighted least square，WLS)算法，其在圖像域的代價函數為

Φ(u)=(y-Gu)′∑-1(y-Gu)

(3)

Φ(u)=(y-Gu)′∑-1(y-Gu)+βR(u)

(4)

式中：第一項為式(3)中的WLS，第二項R(u)表示懲罰項，β是平滑參數。一般的懲罰項有高斯MRF模型、TV模型等。

2.2 改進算法

在基于高斯MRF的懲罰加權最小二乘法中，由于目標函數中正則項的權重是固定的，沒有考慮平坦區域和邊緣區域具有不同的特點，使得圖像邊緣較為模糊。TV正則化方法雖然在去噪的同時能較好地保持圖像的邊緣，但是容易使圖像的平坦區域產生階梯效應。因此，本文提出了一種基于分區域處理的聯合先驗低劑量CT統計迭代重建算法，該算法首先對低劑量CT圖像進行區域劃分，再對不同區域的圖像塊進行不同的處理，然后將改進的正則項應用到懲罰加權最小二乘法中，最后使用SOR迭代重建算法進行求解。

2.2.1 區域劃分并分區域處理

本文算法首先對低劑量CT圖像進行平坦區域和邊緣區域的劃分，由于低劑量CT圖像信噪比低，進行準確劃分比較困難，所以先對其進行中值濾波，將濾波后的圖像進行梯度變換，然后根據圖像的像素值確定一個閾值，將大于閾值的像素都等于255，小于等于閾值的像素都等于0，即可劃分出圖像的平坦區域和邊緣區域。

對于二維CT圖像ui,j來說，梯度變換公式為

(5)

圖1 Shepp-Logan模型區域劃分

圖1(a)是Shepp-Logan模型，圖1(b)是對Shepp-Logan模型進行FBP算法重建后的圖像，圖1(c)是劃分區域后的圖像。

根據劃分區域后的圖像，將低劑量CT圖像中所對應的邊緣區域的圖像塊提取出來進行TV正則化處理，將平坦區域的圖像塊提取出來進行高斯MRF正則化處理。高斯MRF正則項的表達式如下

(6)

對于二維CT圖像ui,j來說，它的全變分表達式為

(7)

將圖像值非負作為約束條件，且可通過如下最優化問題來得到圖像u

(8)

采用梯度下降法來求解該最優化問題，其中全變分的計算公式為

(9)

式中:ε是一個非常小的正參數，防止分母為0。

2.2.2 基于分區域處理的懲罰加權最小二乘算法

本文算法在圖像域的目標函數為

Φ(u)=(y-Gu)′∑-1(y-Gu)+β1R(u1)+β2TV(u2)

(10)

在重建實驗中，兩個懲罰系數β1和β2是可以獨立設置的，首先固定β1，并搜索β2使得RMSE最小，然后固定β2再搜索最優的β1來確定β1和β2的取值，u1和u2分別是從圖像平坦區域和邊緣區域提取的圖像塊，R(u1)表示根據式(6)對平坦區域的圖像塊進行高斯MRF正則化處理，TV(u2)表示根據式(9)對邊緣區域的圖像塊進行TV正則化處理。

(11)

(12)

其具體實現步驟如下：

步驟1 使用FBP算法對低劑量CT投影數據進行重建后的圖像作為初始化圖像，記為u0。

步驟2 將u0先進行中值濾波，然后利用濾波后圖像的梯度值將其分為平坦區域和邊緣區域。

步驟3 分別將u0所對應的平坦區域和邊緣區域的圖像塊提取出來。

步驟4 根據式(6)對平坦區域的圖像塊進行高斯MRF正則化處理，根據式(9)對邊緣區域的圖像塊進行TV正則化處理。

步驟5 利用式(12)計算出迭代重建圖像u。

重復步驟2～步驟5，經過一定次數目標函數的解達到收斂后的重建結果作為最終的重建圖像。

3 實驗結果與分析

為了驗證本文所提算法的有效性和可行性，本文用模擬數據和臨床數據來進行實驗，然后將本文算法與FBP算法、基于TV先驗的懲罰加權最小二乘法(PWLS-TV)以及基于高斯MRF先驗的懲罰加權最小二乘法(PWLS-MRF)的重建圖像進行了比較。為了進一步評價重建圖像的質量，本文采用相關系數、信噪比、歸一化均方誤差對重建圖像進行定量描述分析，定義如下：

(1)相關系數(correlation coefficient，CORR)

(13)

(2)信噪比(signal to noise ratio，SNR)

(14)

信噪比越高表明重建圖像質量越好。

(3)歸一化均方誤差(root mean squared error，RMSE)

(15)

RMSE越小表明重建結果越接近原始圖像。

3.1 模擬數據實驗

圖2 Shepp-Logan模型在不同投影角度下 不同重建算法結果

圖2中第1列為FBP算法結果，第2列為PWLS-MRF算法結果，第3列為PWLS-TV算法結果，第4列為本文算法結果。為了更加清楚看清本文算法重建的圖像比其它重建算法重建的結果更好，選取圖2的兩個感興趣區(region of interest，ROI)，即ROI1和ROI2進行放大，結果如圖3所示。

由圖2、圖3可知，FBP重建結果中出現條形偽影，且角度越少越嚴重；采用PWLS-MRF算法進行重建，雖然重建圖像中的條形偽影得到了有效的抑制，但沒有很好保護圖像的邊緣信息，其邊緣模糊(如圖3箭頭所指邊緣)；PWLS-TV算法雖然保持了圖像的邊緣信息，但同時在平坦區域帶來了階梯效應(如圖3框中區域)，而本文算法在去除階梯偽影的同時保護了圖像的邊緣信息。

表1為各算法在不同角度下重建出的Shepp-Logan圖像相較于原始圖像的質量評價參數對比。從表1可知，本文算法和其它兩種重建算法相比，CORR、SNR值更大，RMSE更小。表明本文算法重建出的圖像質量更好，因此，

表1 Shepp-Logan模型在不同角度下 各算法的客觀評價參數

本文算法的改進是可行的。

3.2 臨床數據實驗

由于無法獲得實際的投影數據，所以本文實驗二使用了從The Cancer Image Archive(TCIA)中獲得的一張高劑量的實際肺部圖像切片(如圖4所示)進行實驗的仿真模擬來驗證本文所提算法的有效性。

圖4 高劑量肺部切片

實驗二采用仿真平行束掃描高劑量肺部切片來獲取投影數據，稀疏角度的個數分別為80、120、150，每個角度上有888個探測器，實驗二獲取原始投影數據和給投影數據中加噪聲來模擬低劑量CT稀疏投影數據的實驗參數設置同實驗一，然后對仿真的低劑量CT稀疏投影數據進行實驗，實驗結果如圖5和圖6所示。

圖5 肺部切片在不同投影角度下不同重建算法結果

圖6 肺部切片在不同角度下各重建圖像局部放大圖

圖6是將圖5中的ROI1和ROI2進行放大。從圖5、圖6可以看出，稀疏角度數越多重建結果越好，由于PWLS-MRF算法中權重是固定的，使得圖像邊緣較為模糊(如圖6箭頭所指邊緣)。PWLS-TV算法重建效果明顯改善，但是在欠采樣的情況下重建圖像產生了階梯偽影(如圖6框中區域)。而本文算法在平坦區域平滑的同時保持了圖像的邊緣信息，投影角度數越多越接近原始圖像。

表2為肺部切片采用PWLS-MRF、PWLS-TV和本文算法重建的圖像相對于模板圖像的RMSE、CORR、SNR等質量評價參數的對比。由表2可知，對比其它兩種重建算法，本文算法的CORR、SNR值更大，RMSE更小。表明本文算法重建出的圖像質量更好，而且去噪能力更強。

表2 肺部切片在不同角度下各算法的客觀評價參數

本文實驗一和實驗二使用向模板和高劑量CT圖像的投影數據加噪聲的方式來模擬低劑量CT投影數據，本文實驗三直接使用從The Mayo Clinic(USA)獲取的一張低劑量的實際腹部圖像切片(如圖7所示)，采用仿真平行束對其掃描來獲取投影數據，稀疏角度的個數分別為80、120、150，每個角度上有888個探測器，然后對仿真的低劑量CT稀疏投影數據進行實驗，實驗結果如圖8和圖9所示。

圖7 低劑量腹部切片

圖8 腹部切片在不同投影角度下不同重建算法結果

圖9 腹部切片在不同角度下各重建圖像局部放大圖

圖9是選取圖8的3個ROI，即ROI1、ROI2、ROI3進行放大。由圖8、圖9可以看出，第一列的FBP算法投影角度數越少，重建結果中噪聲越多，細節信息失去的越多。第二列的PWLS-MRF算法重建的圖像邊緣信息被過度平滑(如圖9箭頭所指邊緣)。第三列的PWLS-TV算法雖然較好地保持了圖像的邊緣，但是在噪聲多的地方產生階梯偽影(如圖9框中區域)。第四列的本文算法消除了PWLS-TV引入的階梯偽影，同時保持了圖像的邊緣細節信息。從視覺效果上來看，本文算法明顯優于其它3種重建算法。

表3為不同算法對投影角度數為80、120、150的噪聲投影數據進行重建結果的客觀評價參數，在條件相同的情況下，與其它兩種算法相比，本文算法重建圖像的SNR值和CORR值更大，RMSE值更小，表明本文算法的去噪能力更強且與原始圖像更相似。

4 結束語

本文提出了一種基于分區域處理的低劑量CT統計迭代重建算法，該算法能夠充分地利用高斯MRF和TV各自的優點對平坦區域和邊緣區域進行處理，使得重建出的圖像在平滑平坦區域的同時保留圖像的邊緣信息。將本文算法與其它重建算法比較可以發現，本文算法不僅去噪能力強，而且還能有效地保護重建圖像的邊緣細節信息。但是，本文所提的算法中有一些參數是憑經驗確定的，可能會影響重建圖像的質量。因此，可以將自適應參數選擇作為未來的研究方向。

表3 腹部切片在不同角度下各算法的客觀評價參數