基于動態優化的X射線CT低劑量重建

王 康，趙 琦，李 銘

(1.南京理工大學 電子工程與光電技術學院，江蘇 南京 210094；2.中國科學院 蘇州生物醫學工程技術研究所，江蘇 蘇州 215163)

1 引 言

計算機斷層成像(Computed Tomography,CT)是通過無損方式獲取人體內部解剖信息的一種醫學影像技術。CT掃描的普及，給臨床診斷帶來了極大方便，但也引發了受檢者的輻射損傷。由于X射線輻射劑量的累加效應，接受過量的X光照射可能會顯著增加病人誘發癌癥和白血病等疾病的概率。因此，國際輻射防護委員會建議，在滿足臨床診斷需求的情況下，CT輻射劑量越低越好[1-2]。

低劑量CT成像技術是減少輻射損傷的有效方法。臨床應用中的CT低劑量掃描成像可大致分為兩類：降低管電流或管電壓和降低采樣數[3-4]。前者多適用于投影數據采集完備的螺旋CT設備，后者則更適用于投影數據欠采樣的CBCT(Cone-beam CT,CBCT)設備。目前，國內主流的CT設備仍采用濾波反投影法(Filtered Back-projection,FBP)進行重建，采用CT低劑量掃描協議將導致重建圖像中偽影和噪聲的增加以及對比度的下降，直接影響醫生對病變組織的診斷分析。因此，如何保證低劑量掃描下的圖像質量是CT成像領域的熱點問題之一，而研究迭代重建算法[5-14]則是抑制噪聲和偽影，提高CT重建圖像質量的關鍵。Idris等人[5]首先采用復合泊松分布模型對投影數據的噪聲特性進行建模，并提出基于這一模型的統計迭代重建算法(Statistical Iterative Reconstruction，SIR)。Lasio等[6]分析了基于泊松分布模型的統計迭代重建算法在實際CT系統中的性能，并驗證了泊松分布模型對于實際CT設備采集投影數據能取得較好的重建效果。在此基礎上，為了進一步提升重建效果，需改進目標函數中的先驗約束項。Zhang等[7]和Wu等[8]將L1范數字典學習(DL)正則化項引入統計迭代重建模型(SIR)，減輕了L2最小化帶來的過平滑效果并保留更多圖像細節。Gou等[9]測試了一種基于Lp字典學習(DL)稀疏約束的新穎重建算法，結論是較低的p值會呈現出更好的重建性能。基于字典學習的CT低劑量重建對于圖像細節有所改善，但是其重建時間較長，且由于低劑量投影數據會導致迭代前期存在偽影，這些偽影會被字典所表示從而影響最終的重建結果。Xu等[10]則結合全變分最小化(Total Variation，TV)，提出了SIR-TV算法。該算法能較好地抑制噪聲和偽影，但對于稀疏性表征差的重建目標，重建圖像中易出現塊狀偽影。此后，Sidky等[11]又提出分數階全變分算法(TpV)。該算法通過調整p值可有效改進待重建目標的稀疏性表征，但是TpV算法易導致重建圖像中出現散點噪聲。Lei等[12]對TV結構進行改進，通過引入自適應的加權權重，提出了自適應變權全變分算法(Adaptively Reweighted TV,ARWTV)。ARWTV算法對于噪聲和偽影有很好的抑制效果，但重建結果中弱對比度組織邊緣模糊。近年來，隨著深度學習的發展，許多學者將基于深度學習的方法應用于醫學領域[15-23]。與傳統方法相比，深度學習的圖像質量有所提高，但這需要大量的訓練數據，且重建質量取決于數據集的質量。目前大多數深度學習中使用的數據是通過噪聲插入來模擬的，這可能無法反映低劑量CT掃描的實際噪聲分布。相反，傳統的基于壓縮感知(Compressed Sensing ,CS)的重建方法不需要大量的訓練數據，并且可以很容易地集成到傳統的迭代重建框架中。

針對CT低劑量掃描下出現的噪聲和偽影問題，文中提出一種基于分數階全變分動態優化(TpV dynamic optimization,TpV-DO)的CT低劑量成像算法。首先通過一族雙曲正切函數集構造TpV的動態復合函數模型，以增強待重建目標的稀疏性表征；其次，通過引入動態優化設計以實現對噪聲和偽影的抑制以及對弱對比度組織邊緣的保護；最后，結合統計迭代重建框架，使用交替優化算法實現對目標函數的優化求解。采用數值模型和動物掃描獲取的投影數據對新算法進行了實驗，并與FBP算法、TpV算法、ARWTV算法進行了比較。

2 算法及理論分析

2.1 統計迭代成像模型

基于泊松分布模型來模擬投影數據采集過程，將統計迭代模型與正則項結合的目標函數可以表示為[10]：

(1)

2.2 基于TpV的動態優化模型描述

在2D圖像空間中，基于TpV模型的稀疏性變換可以用雙索引的方式表示如下[11]：

(2)

其中：

(3)

其中：j=(m-1)×W+n，m=1,…,H，n=1,…,W，H和W分別代表重建圖像x的寬度和高度；p為彈性的范數模型控制參數，通過p值調整可以實現范數模型在L0范數與L2范數之間選擇。

壓縮感知理論指出[24]：重建信號的稀疏性表示越好，越有利于精確地恢復原始信號。為了應用TpV模型增強重建圖像的稀疏表征，并建立動態優化模型。現考慮如下一族雙曲正切函數集，其具體表示形式如下：

(4)

其對應的范數求和形式如下：

(5)

其中：σ為動態調整參數，用于控制圖像表征和相鄰迭代過程的近似性。

如圖1所示，當σ參數取值較大時，對重建圖像的稀疏性表征較差，重建結果容易導致過平滑，圖像紋理細節模糊；當σ參數取值較小時，對重建圖像的稀疏性表征較好，重建結果容易陷入局部極值，圖像噪聲水平顯著。基于上述雙曲正切函數集構造TpV動態復合函數模型，其具體形式如下：

圖1 雙曲正切函數曲線動態演化圖Fig.1 Dynamic evolution of hyperbolic tangent function curve

(6)

其中：p和σ均為動態優化參數，對式(6)分析，當p趨近于2，σ取較大值時，R(x)接近于L2范數，對于偽影和噪聲抑制效果較好；當σ趨近于0，p取較小值時，R(x)接近于理想的L0范數，對圖像紋理細節的保護效果較好。因此通過對p和σ的合理選擇從而在保持圖像邊緣和細節的同時能有效克服TV正則項存在的階梯效應問題。

為了使上述的正則化估計同時表現出較好的紋理保護和噪聲抑制性能，本文提出動態最優化過程，即初始迭代優化時，采用較大的p和σ值以更好地抑制偽影和噪聲，接下來的優化過程中依次降低p和σ值以更好地保存圖像紋理信息。

實驗中，為了保證算法初始時具有較好的偽影和噪聲抑制效果，初始設定：σ0=0.8，p0=1.2；為了保證相鄰兩次迭代的近似性，相鄰迭代的遞減因子Δσ和Δp分別設定為0.016和0.01；同時為了避免算法收斂于局部極值，實驗中設定：σmin≥0.1，pmin≥0.8。

2.2 交替最小化DO-TpV重建算法

如式(1)為典型的統計重建目標函數，該目標函數由對數似然數據保真項和正則化項兩部分構成，本文采用交替最小化的方法對其進行優化求解，具體算法流程如下：

步驟一、最小化對數似然數據保真項：采用典型的可分離拋物面算法[10]進行優化，其迭代更新公式為：

(7)

式中：t表示迭代次數，x0為初始迭代圖像。

步驟二、式(6)的動態最優化求解：采用梯度下降法對其進行優化。其中式(6)的梯度公式可以表示為：

(8)

其中：TpV模型的梯度公式可表示為：

(9)

式中：ζ為擾動約束因子，通常取ζ=1×10-8。

基于上述的梯度公式，采用梯度下降法進行迭代更新，其計算公式為：

(10)

上述算法具體實施中，首先進行參數x、σ、p的初始化，接下來重復執行步驟一與步驟二直至滿足預先設定的迭代終止條件。下文實驗中將迭代終止條件設定為t不大于預先設定的迭代總次數Tmax或相鄰兩次迭代誤差小于設定閾值err，如果t>Tmax，則迭代終止；否則，t=t+1 ，更新σ和p值，并啟動下一輪迭代。本文提出的基于分數階全變分動態優化算法(DO-TpV)的CT低劑量成像重建流程如下所示：

輸入:A,l^,β,σ,p,err,初始化x0,σ0,p0,Δσ,Δp,t=1;

輸出:重建圖像x.While(t‰Tmaxor‖xt-xt-1‖/‖xt‖≥err){步驟一:優化SIR過程xtj=xt-1j-∑Ndi=1ai,jyi[Axt-1]i-^li ∑Ndi=1ai,jyi∑Jk=1ai,k ,j=1,...,J;步驟二:分數階全變分動態模型優化(TpV-DO)d=‖xt-xt-1‖2;dx=?R(x)?xj;xtj=xtj-β·d·dx;步驟三:參數更新t=t+1;ifσ≥0.1,p≥0.8σ=σ-Δσ,p=p-Δp.}

3 實驗與結果

為了驗證本文算法在噪聲抑制和軟組織細節保護方面的有效性，分別采用數值仿真實驗和動物掃描實驗對本文提出的CT低劑量成像算法(DO-TpV)進行實驗驗證，并對其與解析重建算法(FBP)、TpV算法、字典學習(DL)、ARWTV算法進行了比較。其中，算法的優化參數設定：本文所提算法采用2.2節所述方式；TpV算法正則化參數β=0.2，動態優化參數p=0.8～1.2；ARWTV算法正則化參數β=0.2，變全次數設定為5；字典學習(DL)算法正則化參數β=0.2，圖像塊尺寸設定為7×7，過完備字典的稀疏度設定為4。

3.1 數值模擬實驗

本部分實驗采用數值模型生成圓軌跡扇形束掃描投影數據，投影參數設定如下：X射線源到旋轉中心的距離為54.1 cm，探測器到旋轉中心的距離為40.8 cm，CT圖像的像素數為512×512，尺寸大小為20 cm×20 cm，探測器數為642，相鄰探測器間距為0.672 mm。掃描范圍為0°～360°，步長為1°和2°，以生成180°和360°兩種采樣角度數下的投影數據。實驗中采用Patrick等[10]提出的Poisson分布模擬光子探測的隨機過程，其中：實驗總光子數設定為8.0×105。通過對探測器記錄信號進行對數化處理后，獲得含噪聲的投影數據。分別使用FBP算法、TpV算法、字典學習(DL)算法、ARWTV算法和DO-TpV算法進行CT圖像重建。圖2和圖3分別給出了數值模型在180°和360°兩種采樣數下的實驗結果，窗口顯示范圍[0.185,0.205]。

圖3 數值模型實驗結果(360采樣數)Fig.3 Numerical model results (360 views)

如圖2和圖3所示，在低劑量掃描條件下，FBP算法重建圖像中出現嚴重的偽影和噪聲，且軟組織細節完全被噪聲和偽影遮蓋，總體圖像質量極差。對比解析算法(FBP)結果，4種迭代算法重建圖像的質量有了明顯提升。對比4種迭代算法結果，ARWTV算法取得了很好的偽影和噪聲抑制效果，但重建圖像軟組織細節模糊嚴重，如圖2(e)紅色箭頭所示的放大區域。相比ARWTV算法，字典學習算法(DL)和TpV算法在軟組織細節保護方面有了顯著的提升，但字典學習(DL)算法的部分軟組織產生變形模糊，如圖2(d)紅色箭頭所示，且字典學習(DL)從整體看尤其在180°較低采樣數下會存在條紋狀偽影，而對于TpV算法部分軟組織細節結構不清晰，如圖2(c)紅色箭頭所示的放大區域。由圖2和圖3可以看出，DO-TpV算法在軟組織細節保護方面表現最好，且取得較好的偽影和噪聲抑制效果，總體圖像質量的提升顯著。

為了評價不同算法的重建結果對高對比度邊緣細節的保護能力，本文提取圖2紅色實線所示位置，并采用剖面密度曲線顯示TpV算法、字典學習(DL)算法、ARWTV算法和DO-TpV算法的重建結果與金標準圖像的接近程度。圖4實驗結果表明，與其他算法的重建結果相比，DO-TpV算法的重建結果更好地擬合了金標準圖像的密度曲線。

圖4 重建結果的剖面灰度值比較圖(圖2)Fig.4 Line profiles of reconstructed results of Fig.2

本文采用歸一化平均絕對偏差(Normalized mean absolute deviation,NMAD)和信噪比(Signal noise ratio,SNR)2個標準定量評價重建結果，NMAD和SNR的定義如下：

(11)

(12)

由表 2 可以看出，DO-TpV算法的NMAD小于FBP算法、TpV算法、字典學習(DL)算法、ARWTV算法的NMAD，表明DO-TpV算法重建的圖像與理想圖像更接近。在180°視角和360°視角的數值模型實驗中，DO-TpV算法重建圖像的信噪比比FBP算法、TpV算法、字典學習(DL)算法、ARWTV算法要高，表明DO-TpV算法能更好地保護圖像紋理信息和抑制圖像噪聲，重建圖像的質量更高。

表2 數值模型定量分析Tab.2 Quantitative evaluation of the numerical results

3.2 動物數據實驗

為了進一步驗證DO-TpV算法的性能，采用真實CT系統掃描小鼠來生成投影數據。系統掃描參數設定如下：X光管電壓50 kVp，X光管電流1 mA，曝光時間0.466 9 s，X射線源到掃描中心距離為22.188 cm，X射線源到探測器中心的距離為65.85 cm，平板探測器陣列數1 536×880，探測單元尺寸0.15 mm×0.15 mm，掃描間隔為2°，共采集180個等間距投影圖像。分別使用FBP算法、TpV算法、字典學習(DL)算法、ARWTV算法和DO-TpV算法對中心片層投影數據進行重建實驗。重建像素大小為512×512，重建像素單元尺寸為0.1 mm×0.1 mm。圖5給出小鼠掃描數據重建結果，窗口顯示范圍為[-0.09，1.19]。

圖5 小鼠掃描數據重建結果Fig.5 Reconstructed results of scanned mouse datasets

從圖5可以看出，低劑量掃描協議下，解析算法(FBP)重建圖像中出現嚴重的偽影和噪聲，且軟組織邊緣細節極其模糊，如圖5(b)放大區域所示。對比解析算法重建結果，4種迭代算法重建結果在抑制偽影和噪聲，保護軟組織結構方面取得了明顯改善。對比4種迭代算法結果，ARWTV算法在抑制偽影和噪聲方面效果顯著，但重建圖像部分細節模糊不清晰，如圖5(c)放大區域紅色箭頭所示位置。相比ARWTV算法，TpV算法能夠重建出圖像的邊緣細節，但重建結果中存在一些散點噪聲，如圖5(c)放大區域所示。而對于字典學習(DL)方法，偽影和噪聲得到了一定程度的抑制，但結果仍然存在由于過完備字典而導致的不規則偽影和噪聲，如圖5(d)放大區域所示。對比ARWTV算法和TpV算法的重建圖像，DO-TpV算法重建的圖像細節更為清晰，且圖像的噪聲水平也顯著降低，總體圖像質量最佳。

4 結 論

本文提出一種基于動態優化的CT低劑量成像算法。該算法通過雙曲正切函數集構造分數階變分的動態復合函數模型，使基于分數階變分的重建逐步接近理想的L0范數重建。通過動態最優化過程逐次增強當前重建圖像的稀疏表征，改善了DO-TpV算法在圖像紋理細節保護和噪聲抑制方面的性能，提升了重建圖像的質量。數值模型實驗結果表明，在180個采樣角度下，文中算法重建圖像的信噪比分別比FBP算法、TpV算法、字典學習(DL)方法、ARWTV算法高出29.51，8.03，9.15，6.81 dB。動物數據實驗結果表明，DO-TpV算法重建結果有效抑制了噪聲和偽影，清晰重建出小鼠軟組織細節，極大地提高了低劑量采集數據重建圖像的質量。