基于變分光流估計的肺部4D-CT圖像超分辨率重建

耿鳳歡 劉 慧 郭 強(qiáng) 尹義龍

1(山東財經(jīng)大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院 濟(jì)南 250014) 2(山東省數(shù)字媒體技術(shù)重點實驗室 濟(jì)南 250014) 3 (山東大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院 濟(jì)南 250013) (gfh0928@outlook.com)

基于變分光流估計的肺部4D-CT圖像超分辨率重建

耿鳳歡1,2劉 慧1,2郭 強(qiáng)1,2尹義龍1,3

1(山東財經(jīng)大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院 濟(jì)南 250014)2(山東省數(shù)字媒體技術(shù)重點實驗室 濟(jì)南 250014)3(山東大學(xué)計算機(jī)科學(xué)與技術(shù)學(xué)院 濟(jì)南 250013) (gfh0928@outlook.com)

由于受到掃描時間和照射劑量的限制，肺部4D-CT數(shù)據(jù)中縱向采樣率遠(yuǎn)小于面內(nèi)采樣率.為了得到更高質(zhì)量的肺部圖像，從醫(yī)學(xué)圖像固有的自相似性出發(fā)，提出了一種基于局部和全局相結(jié)合的變分光流估計的圖像序列超分辨率重建技術(shù)，用于提高4D-CT圖像重建質(zhì)量.首先，構(gòu)建了一個用于求解肺部4D-CT不同相位圖像之間的光流場的變分光流模型；然后，利用快速交替方向乘子法求解該模型，得到不同相位圖像之間的光流場；最后，基于光流場，并利用非局部迭代反投影超分辨率重建算法，實現(xiàn)了高分辨率肺部圖像的重建.實驗結(jié)果表明：與已有算法相比，本方法在增強(qiáng)圖像紋理結(jié)構(gòu)的同時更好地保留了圖像的輪廓.

4D-CT圖像；超分辨率重建；光流估計；交替方向乘子法；迭代反投影

放射治療是醫(yī)學(xué)上治療肺癌最常用的手段.肺部4D-CT圖像因其能夠提供必要的呼吸運動信息來引導(dǎo)醫(yī)生進(jìn)行精確地放射治療，從而得到越來越多的關(guān)注.肺部4D-CT圖像是由呼吸運動周期中的多個相位各自對應(yīng)的3D-CT圖像按照獲取時間的先后順序排序得到的，能夠包含呼吸周期內(nèi)肺部器官以及腫瘤區(qū)域的空間結(jié)構(gòu)(可通過圖像分割算法[1]得到)以及運動信息.這些信息可以輔助醫(yī)生精確定位靶區(qū)，有助于病人的放射治療.但是，眾所周知CT圖像的獲取伴隨著高劑量照射[2-3]，為了降低對病人的輻射量，往往只能降低沿縱向的采樣，從而導(dǎo)致肺部4D-CT圖像的層間分辨率遠(yuǎn)低于層內(nèi)分辨率.在觀察3D肺部圖像的冠狀圖時，為了還原真實肺部形狀，通常會將3D圖像沿Z軸方向進(jìn)行插值放大.傳統(tǒng)插值算法會使圖像變模糊.因此，如何提高4D-CT圖像超分辨率重建圖像質(zhì)量是亟待解決的問題.

目前，研究者們常利用3類算法實現(xiàn)4D-CT圖像的超分辨率重建：傳統(tǒng)插值算法、投影重建以及基于字典學(xué)習(xí)的方法.

1) 傳統(tǒng)插值算法.如雙三次插值(Bicubic)算法[4]，是根據(jù)插值公式，利用鄰域內(nèi)已知像素值得到待插值點像素值，其實現(xiàn)簡單且運行速度快，但由于插值過程中未考慮圖像的各向異性特性，導(dǎo)致放大后圖像的輪廓和紋理比較模糊，降低了圖像質(zhì)量.

2) 投影重建算法.大部分是利用4D-CT不同相位間的互補信息來提高當(dāng)前相位圖像的分辨率，如基于全搜索(full search)的運動估計超分辨率重建方法[5]利用迭代反投影思想及全搜索得到的不同相位圖像間的仿射變換信息重建出高分辨率圖像，但是該方法的缺點是對不同相位圖像間運動估計的精度不高；此外，還有基于配準(zhǔn)的凸集投影算法[6]，這類算法的精度很大程度取決于光流估計誤差的大小.

3) 基于字典學(xué)習(xí)的超分辨重建方法[7-9].它們需要大量的訓(xùn)練數(shù)據(jù)集和較長運算時間，很難滿足實時性的需求.

綜合比較上述3類方法，其中投影重建類算法能夠在較短的計算時間內(nèi)得到可靠的超分辨率圖像，較好地迎合了4D-CT超分辨率重建的需求.常用的投影包括迭代反投影算法[10]、凸集投影算法[11]以及最大后驗估計算法[12].其中迭代反投影算法因收斂速度快得到了廣泛的應(yīng)用.

因此，為了提高重建圖像的質(zhì)量和速度，我們采用迭代反投影方法，充分利用4D-CT不同相位間的互補信息實現(xiàn)快速的CT圖像超分辨率重建，并在2個方面進(jìn)行了創(chuàng)新：1)鑒于影響超分辨率重建算法效果的主要因素是圖像配準(zhǔn)的精度，本文采用了局部和全局相結(jié)合(combined local and global, CLG)的變分光流估計模型獲得更加精準(zhǔn)稠密的光流場，并采用快速交替方向乘子法[13-14]實現(xiàn)光流模型的快速求解；2)采用基于邊緣的非局部迭代反投影算法，有效增強(qiáng)超分辨率圖像中的邊緣和紋理細(xì)節(jié).

本文算法主要包含2個步驟：1)將待重建相位圖像作為參考相位，采用CLG變分光流估計算法獲取其他相位與參考相位之間的光流場；2)利用其他相位圖像的互補信息，采用改進(jìn)的非局部迭代反投影算法迭代求解重建圖像.

1 光流模型及求解

CLG變分光流估計模型是用于計算視頻中相鄰幀間光流場的模型，4D-CT數(shù)據(jù)中不同相位間光流場亦可用此模型求解.本節(jié)主要介紹針對肺部4D-CT圖像的CLG變分光流估計模型，以及利用快速交替方向乘子算法求解該光流模型的過程.

1.1CLG變分光流模型

光流估計[15]是利用圖像序列中參考幀與當(dāng)前幀之間像素值數(shù)據(jù)的時域相關(guān)性來確定當(dāng)前幀中各像素位置相對于參考幀的“運動”.光流計算模型都是基于光流基本等式，通過引入光流約束條件來實現(xiàn)光流求解的.CLG變分光流估計[16]是一種對噪聲具有魯棒性且更加準(zhǔn)確的光流計算模型.首先，該模型中引入了對光流的l1范數(shù)約束，從而保證了模型對噪聲的魯棒性；其次，模型引入了對數(shù)據(jù)項的雙邊濾波操作，以減小相鄰像素之間的相互影響；最后，該模型基于圖像局部特性實現(xiàn)了光流的各向異性傳播，有效地抑制光流過傳播帶來的負(fù)影響.下面構(gòu)建適用于肺部4D-CT數(shù)據(jù)不同相位圖像間光流估計的計算模型.

首先，基于光流場基本理論，將肺部4D-CT圖像看作肺部圖像I隨時間發(fā)生形變而生成的圖像序列，也稱為圖像一致性假設(shè)，即呼吸運動過程中肺部圖像各個像素位置亮度沿著運動軌跡保持不變，即：

I(x(t),y(t),t)=C.

(1)

其次，引入對光流的l1范數(shù)約束和各項異性擴(kuò)散約束，以及對數(shù)據(jù)項的雙邊濾波操作，得到肺部4D-CT圖像的CLG變分光流估計模型：

(2)

(3)

其中，X=(x,y)為像素點坐標(biāo)，U0是光流場的初始估計.數(shù)據(jù)殘差項的局部約束使得鄰域像素的光流會影響中心像素點光流的計算，為了確保在得到稠密光流場的同時弱化鄰域像素的影響，最有效的方法是引入對數(shù)據(jù)殘差項的濾波操作.濾波算法有雙邊濾波[18-19]以及基于奇異值分解的濾波算法[20-21]等，其中雙邊濾波算法實現(xiàn)簡單且能滿足需求，故本文采用雙邊濾波對數(shù)據(jù)殘差項進(jìn)行約束.此外，在CLG變分光流估計模型的求解過程中，光流場的各向同性傳播會造成光流的過傳播，從而降低光流場的精度.本文中CLG變分光流模型通過加入光流的各向異性濾波操作，有效地抑制了光流的過傳播.

1.2求解光流模型

為了實現(xiàn)3D-CT序列圖像間光流場的快速求解，我們采用快速交替方向乘子法(alternating direc-tion method of multipliers, ADMM)算法求解本文中的變分光流模型.ADMM是一種求解優(yōu)化問題的計算框架，通過將全局問題分解成多個相互獨立、較小且易于求解的局部子問題，之后協(xié)調(diào)子問題的解得到全局問題的解.ADMM算法結(jié)合了增廣拉格朗日方法以及求解凸優(yōu)化問題的對偶分解策略的優(yōu)勢[13]，實現(xiàn)了對凸優(yōu)化問題的快速求解.定義Ψi(X)為張量X在第i維度上的融合算子，光流模型可轉(zhuǎn)化成等價問題：

s.t.U=Vi, 1≤i≤2，

(4)

其中，Vi(1≤i≤2)是松弛變量.為了求解式(4)中的條件極值問題，本文將該約束問題轉(zhuǎn)換成無約束的增廣拉格朗日函數(shù)：

(5)

為了求解式(5)，本文采用了ADMM算法，通過依次迭代更新U,Vi,Wi(i=1,2)來最小化式(5)中的目標(biāo)函數(shù)，得到光流U的局部最優(yōu)值.U,Vi,Wi(i=1,2)的更新規(guī)則：

具體迭代步驟為

Step1. 固定Vi和Wi，更新Uk+1：

(6)

最優(yōu)解為

(7)

Step2. 固定U和Wi，求解Vi,i=1,2

(8)

式(8)是可分解的，我們可以獨立地計算Vi,i=1,2，且可等價寫成:

(9)

(10)

其中，ν是張量Vi所包含的一維張量，t是Γi中與ν相對應(yīng)的一維張量，k是ν的維數(shù)，di是D中與ν相對應(yīng)一維張量的第i個元素.式(10)是一維全變分問題有很多算法可以快速求解.并行計算的引入使得Vi的計算可以在毫秒時間內(nèi)完成.

Step3. 計算Wi,i=1,2

(11)

該迭代求解過程在原始?xì)埐詈蛯ε細(xì)埐頪13]滿足特定條件時結(jié)束.

文獻(xiàn)[14]將多維變分問題轉(zhuǎn)換成標(biāo)準(zhǔn)的l1范數(shù)問題，并證明了轉(zhuǎn)換后的問題滿足收斂的約束條件.而本文中光流求解算法僅是該算法中的特例(二維變分問題)，因此，可證明本文光流求解過程是收斂的.

2 改進(jìn)的非局部迭代反投影超分辨率重建

2.1迭代反投影算法

(12)

(13)

其中，Tk(·)是指圖像的仿射變換；h是點擴(kuò)散函數(shù)；↓s為下采樣算子.反投影重建過程可表示為

(14)

其中，T-1k是Tk的逆；p是反投影算子，p的取值將影響算法迭代時間；↑s表示上采樣算子.

2.2改進(jìn)的非局部迭代反投影算法

為了更好地利用不同幀圖像間的互補信息，我們結(jié)合迭代反投影算法和非局部迭代反投影算法提出了改進(jìn)的非局部迭代反投影，實現(xiàn)4D-CT圖像超分辨率重建.本文算法在2個方面具有創(chuàng)新性：1)對IBP算法中的投影重建誤差進(jìn)行非局部均值濾波，來保留圖像中的高頻細(xì)節(jié)；2)對反投影重建過程中的圖像融合策略進(jìn)行了改進(jìn).

IBP算法在反投影重建過程中采用了灰度平均的像素級圖像數(shù)據(jù)融合策略，該策略簡單且便于實現(xiàn)，但在輸入低分辨率圖像間光流場估計誤差較大的情況下會嚴(yán)重影響超分辨率重建結(jié)果的質(zhì)量.基于此，本文利用數(shù)據(jù)加權(quán)平均法進(jìn)行融合，根據(jù)低分辨率圖像光流估計誤差確定該圖像的重建誤差在反投影過程中的權(quán)重.這樣便可以很好地抑制光流估計誤差對重建結(jié)果的影響，從而得到高質(zhì)量的高分辨率圖像.修改后的反投影過程可表示為

(15)

(16)

其中,t是控制核衰變速度的參數(shù).

我們將4D-CT圖像所包含的所有相位圖像作為已知的低分辨率圖像集合，從中選出某一相位圖像作為浮動圖像.首先，用CLG變分光流估計算法得到浮動圖像與其余相位圖像間的光流場；其次，結(jié)合該光流場以及圖像放大倍數(shù)，得到圖像間的仿射變換；最后，采用改進(jìn)的非局部迭代反投影算法重建出高分辨率圖像.算法1給出了本文方法的實現(xiàn)過程.

算法1. 基于變分光流估計的4D-CT圖像超分辨率重建算法.

輸出：高分辨率圖像Ih.

①U←利用式(2)得到Ik與參考幀Il間的光流場；

②Tk←由光流場U和倍數(shù)l得到仿射變換；

③ do

④e(n)←0； /*e(n)為第n次迭代中的誤差*/

⑤ fork=1 toK

⑨ end for

3 光流模型中參數(shù)問題

當(dāng)待配準(zhǔn)圖像與基準(zhǔn)間存在較大位移時，式(1)圖像灰度一致性假設(shè)會造成估計結(jié)果誤差較大.為此，我們在求解過程中，利用由粗到細(xì)的圖像金字塔方法提高光流估計精度.式(2)中各項異性擴(kuò)散因子D的參數(shù)設(shè)置為典型值α=5,β=1/2，權(quán)重系數(shù)λ=1 000.增廣拉格朗日函數(shù)中的懲罰項參數(shù)設(shè)置為ρ=10.由光流求解過程可知，雙邊濾波器的鄰域大小的選取對光流估計結(jié)果影響最大.因此，我們用6種不同大小濾波器并分別對2幅圖像間存在大位移和小位移的情況進(jìn)行了對比實驗，如圖1所示.

Fig. 1 The influence of filter size on the accuracy of optical flow圖1 濾波器大小對光流精度的影響

圖1中，橫軸表示濾波器鄰域大小，縱軸表示圖像插值誤差[24].插值誤差是一種光流估計效果評價標(biāo)準(zhǔn)，其值為真實圖像與基于光流補償?shù)玫焦烙媹D像間的均方差.由圖1可知，當(dāng)濾波器大小為3×3，5×5時，對于大位移圖像間的光流估計效果較差；當(dāng)濾波器大小為9×9,11×11,13×13時，雖然對圖像間存在大位移的情況下光流估計結(jié)果較精確，但計算時間較長且對圖像間存在小位移情況估計偏差變大了，這是由于濾波器的鄰域變大，使得影響當(dāng)前中心像素的鄰域像素變多，從而影響估計結(jié)果；當(dāng)濾波器大小為7×7時，能夠很好地處理大位移光流估計問題，且處理時間也能達(dá)到要求.故實驗中將濾波器的大小設(shè)置為7×7.

圖2顯示了使用基于CLG變分光流估計算法對不同相位中的3幅冠狀面圖像進(jìn)行光流估計的2個實例.這2幅冠狀面圖像是從同一幅4D-CT圖像中選取不同相位上的冠狀面圖像圖2(a)、圖2(b)和圖2(f)，其中圖2(a)作為參考幀圖像，圖2(b)和圖2(f)分別作為圖像，使用CLG變分光流估計算法分別得到圖2(a)和圖2(b)之間的光流場圖2(c)以及圖2(a)和圖2(f)之間的光流場圖2(g).利用得到的光流場圖2(c)對圖2(a)進(jìn)行運動補償后得到圖2(d)；利用得到的光流場圖2(g)對圖2(a)進(jìn)行運動補償后得到圖2(h).圖2(e)是圖2(b)和圖2(d)之間的差圖像，圖2(i)是圖2(f)和圖2(h)之間的差圖像.觀察可知，圖2(b)和圖2(f)中肺部相對于圖2(a)整體向下運動且圖2(a)和圖2(b)之間的位移較小、圖2(a)和圖2(f)之間的位移較大，這與得到光流場一致.結(jié)果可見，圖2(b)和圖2(d)以及圖2(f)和圖2(h)的差異很小，說明CLG變分光流估計能準(zhǔn)確估計圖像.間的位移場，有助于更好地進(jìn)行超分辨率重建.

Fig. 2 The solved optical flow field and the difference between the original image and the compensated according to the optical flow field圖2 光流場、補償圖像以及誤差圖像

4 實驗結(jié)果對比與分析

實驗所用數(shù)據(jù)采用的是德克薩斯大學(xué)安德森癌癥中心DIR-lab實驗室*https://www.dir-lab.com/Referenceta.html提供的公開的標(biāo)準(zhǔn)肺部4D-CT數(shù)據(jù)集[25]，文獻(xiàn)[6-7]均采用了該數(shù)據(jù)集.我們從中選取了5組4D-CT數(shù)據(jù)(數(shù)據(jù)1～數(shù)據(jù)5).每組數(shù)據(jù)含有呼吸運動過程中的10個不同相位，包括呼氣末端相位和吸氣末端相位.數(shù)據(jù)的層內(nèi)像素尺寸范圍從(0.97×0.97)mm2到(1.16×1.16)mm2，層厚為2.5mm.下面將分別展示CLG變分光流估計和最終的超分辨率重建結(jié)果.

4.1實驗對比

實驗將本文算法與雙三次插值(Bicubic)算法、基于全搜索運動估計(Full Search)的超分辨算法以及非局部迭代反投影(NLIBP)重建算法進(jìn)行了對比，以評估本文算法在邊界保持方面的效果，并論證光流估計對超分辨重建效果的影響.

Bicubic算法利用待插值像素點鄰域內(nèi)的16個像素值得到插值點的像素值，實現(xiàn)了圖像的放大.但是，該方法未考慮圖像的各向異性特性，導(dǎo)致放大后的圖像輪廓和紋理比較模糊.Full Search算法通過基于圖像塊的全搜索算法對低分辨率序列進(jìn)行運動估計，再利用迭代反投影算法重建出高分辨率圖像，但是，由于全搜索運動估計精度較低，不能精確引導(dǎo)投影重建過程，造成重建結(jié)果質(zhì)量較差.NLIBP算法利用低分辨率圖像自身的非局部冗余信息引導(dǎo)重建誤差的反向投影，很好地抑制了高對比度邊緣的“振鈴”現(xiàn)象，但是，該方法僅利用了單幅低分辨率圖像的自身信息，導(dǎo)致重建圖像缺失部分細(xì)節(jié)信息.相比而言，本文算法在提高光流估計精度的基礎(chǔ)上，改進(jìn)了非局部反投影算法，更好地利用了不同相位圖像間的互補信息，從而得到具有更多細(xì)節(jié)信息且結(jié)構(gòu)清晰的高分辨率圖像.我們利用圖3和圖4中的對比結(jié)果對上述方法進(jìn)行驗證.

Fig. 3 Visual comparison of reconstruction results of coronal images圖3 不同算法得到的冠狀圖超分辨率圖像視覺效果對比

Fig. 4 Visual comparison of reconstruction results of sagittal images圖4 不同算法得到的矢狀面超分辨率圖像視覺效果對比.

圖3展示了采用不同算法對數(shù)據(jù)2進(jìn)行冠狀面圖像重建的結(jié)果圖像.由于數(shù)據(jù)2中層內(nèi)像素尺寸為(1.16×1.16)mm2且層厚2.5mm，為了使軸向采樣率和層內(nèi)采樣率相同，本文將數(shù)據(jù)2冠狀圖上采樣倍數(shù)設(shè)置為2.15.圖3中列圖像分別是利用Bicubic算法、Full Search算法、NLIBP算法以及本文算法得到的冠狀面圖像超分辨率結(jié)果.前3行圖像是分別對應(yīng)數(shù)據(jù)2中不同相位(相位5、相位8和相位1)的超分辨率重建結(jié)果.為了更好地對比不同算法重建效果，第4行給出了第3行圖像中局部細(xì)節(jié)的放大圖.對比前3行圖像可以發(fā)現(xiàn)，當(dāng)以不同相位圖像作為基準(zhǔn)圖像時得到的超分辨率圖像間存在差異，但在圖像質(zhì)量方面相差無幾.對比放大后的肺實質(zhì)周圍的血管與組織邊緣處的細(xì)節(jié)可以看出，本文算法重建出的圖像具有更清晰的結(jié)構(gòu)，同時邊緣和細(xì)節(jié)信息也得到加強(qiáng).

圖4則給出采用不同算法對數(shù)據(jù)2進(jìn)行矢狀面圖像重建得到的結(jié)果圖像.圖4的圖片排列方式與圖3一致.對比放大后的血管與組織邊緣處的細(xì)節(jié)可以看出，Bicubic算法重建效果最差；較Full Search方法，NLIBP算法與本文算法重建出的圖像具有更清晰的結(jié)構(gòu)；較NLIBP算法，本文方法重建出的圖像中細(xì)節(jié)信息更加明顯.

為了更好地評價重建結(jié)果圖像的清晰度，我們采用平均梯度作為量化標(biāo)準(zhǔn).圖像邊緣附近的灰度值變化率較大，這種變化率大小可以用來表示圖像的清晰度.圖像平均梯度就是用圖像多維方向上灰度變化率之和來表征圖像的相對清晰度，值越大圖像越清晰：

(17)

表1列出了分別利用Bicubic算法、Full Search算法、NLIBP算法以及本文算法對5幅4D-CT圖像選擇不同相位圖像為浮動圖像，對不同平面進(jìn)行超分辨率重建，然后計算得到所有重建圖像平均梯度值的平均值.對比表1中數(shù)據(jù)可知，較前3種方法，本文方法得到的平均梯度值更高，圖像清晰度明顯增強(qiáng).

Table 1 Average Gradient Using Four Different Algorithms

Fig. 5 Comparison of average gradient of super reconstruction results圖5 不同算法得到的高分辨率結(jié)果對比

此外，利用山東省千佛山醫(yī)院提供的10組4D-CT數(shù)據(jù)對相關(guān)算法進(jìn)行了對比分析.圖5給出了不同算法得到的高分辨率圖像平均梯度，可以看出，對于10組實驗數(shù)據(jù)，本文算法得到的超分辨率重建結(jié)果均好于前3種算法.

4.2結(jié)果分析

從圖3和圖4中的細(xì)節(jié)放大圖中可以看出，全搜索算法由于得到的光流場不夠精確，導(dǎo)致重建圖像中出現(xiàn)一些噪聲.本文算法利用基于CLG變分光流估計算法計算得到比全搜索估計更加精確和稠密的光流場，所以超分辨結(jié)果要優(yōu)于全搜索算法.由于非局部迭代反投影只利用低分辨率圖像自身的冗余信息，導(dǎo)致超分辨率結(jié)果圖像在具有較少信息冗余的器官組織邊緣出現(xiàn)了模糊.而本文算法則是很好地利用不同幀圖像間的互補信息，在保留細(xì)節(jié)信息的同時提高了重建圖像的清晰度.對比表1中5份4D-CT數(shù)據(jù)重建結(jié)果圖像的平均梯度值可知，本文算法得到的平均梯度值都是最大的，這表明本文算法更好地保留了圖像的細(xì)節(jié)信息.

5.3算法時間效率分析

本文利用ADMM算法求解光流模型的過程中，主要工作是迭代更新U,Vi,Wi(i=1,2).設(shè)實驗圖像包含的像素數(shù)為n，迭代次數(shù)為k1，則在光流求解過程中，更新U的時間復(fù)雜度為O(n)，更新Vi(i=1,2)的時間復(fù)雜度為O(m)(其中m為實驗圖像像素的列數(shù))，更新Wi(i=1,2)的時間復(fù)雜度為O(n).由于m

實驗用C++語言實現(xiàn)了本文中光流求解算法，并在配備英特爾i5處理器(主頻為3.2 GHz)的主機(jī)上進(jìn)行了測試.對于分辨率為256×99的浮動圖像和基準(zhǔn)圖像，平均光流計算時間為0.323 s，極大地縮短了光流計算時間.因此，對于10幅圖像組成的實驗圖像序列，使用本文超分辨率算法可以在30 s內(nèi)得到超分辨結(jié)果.

5 總 結(jié)

本文提出了一種基于CLG變分光流估計的超分辨率重建算法，用于提高4D-CT冠狀面圖像分辨率.文中，我們將不同相位的3D-CT圖像看作不同“幀”圖像，基于CLG變分光流模型，并利用快速ADMM算法迅速收斂的性質(zhì)，實現(xiàn)了浮動“幀”圖像與基準(zhǔn)“幀”圖像間光流場的快速計算；然后，基于得到的光流場，改進(jìn)了非局部迭代反投影算法，并利用不同相位間的互補信息增強(qiáng)重建了圖像的結(jié)構(gòu)細(xì)節(jié)信息.實驗結(jié)果表明：本文算法不僅有效增強(qiáng)了圖像的紋理結(jié)構(gòu)，而且能夠更好地保留圖像中的輪廓信息.

在未來的研究工作中，我們將研究實現(xiàn)光流模型及其求解過程中參數(shù)的自適應(yīng)設(shè)定，以進(jìn)一步提高光流估計的精度和超分辨率重建質(zhì)量.此外，我們還將應(yīng)用GPU并行計算技術(shù)來縮短本文算法的運行時間，用于滿足實時性需求.

[1] Liu Hui, Zhang Caiming, Deng Kai, et al. Research on fast FCM pulmonary nodule segmentation algorithm using improved self-adaption[J]. Journal of Computer-Aided Design and Computer Graphics, 2014, 26(10): 1727-1736 (in Chinese)(劉慧, 張彩明, 鄧凱, 等. 改進(jìn)局部自適應(yīng)的快速FCM肺結(jié)節(jié)分割方法[J]. 計算機(jī)輔助設(shè)計與圖形學(xué)學(xué)報, 2014, 26(10): 1727-1736)

[2] Khan F, Bell G, Antony J, et al. The use of 4D-CT to reduce lung dose: A dosimetric analysis[J]. Medical Dosimetry, 2010, 34(4): 273-278

[3] Li Tianfang, Schreibmann E, Thorndyke B, et al. Radiation dose reduction in four-dimensional computed tomography[J]. Medical Physics, 2005, 32(12): 3650-3660

[4] Carey W, Chuang D, Hemami S. Regularity-preserving image interpolation[J]. IEEE Trans on Image Processing, 1997, 8(9): 1293-1297

[5] Su Heng, Zhou Jie, Zhang Zhihao. Survey of super-resolution image reconstruction methods[J]. Acta Auto-matica Sinica, 2013, 39(8): 1202-1213 (in Chinese)(蘇衡, 周杰, 張志浩. 超分辨率圖像重建方法綜述[J]. 自動化學(xué)報, 2013, 39(8): 1202-1213)

[6] Zhang Yu, Wu Xiuxiu, Yang Wei, et al. Super-resolution reconstruction for 4D computed tomography of the lung via the projections onto convex sets approach[J]. Medical Physics, 2014, 41(11): No.111917

[7] Zhang Yu, Wu Guorong, Yap P, et al. Hierarchical patch-based sparse representation—A new approach for resolution enhancement of 4D-CT lung data[J]. IEEE Trans on Medical Imaging, 2012, 31(11): 1993-2005

[8] Li Min, Chen Jian, Le Xiang, et al. Super-resolution based on sparse dictionary coding[J]. Journal of Software, 2012(5): 1315-1324 (in Chinese)(李民, 程建, 樂翔, 等. 稀疏字典編碼的超分辨率重建[J]. 軟件學(xué)報, 2012(5): 1315-1324)

[9] Yang Shuaifeng, Zhao Ruizhen. Image super-resolution reconstruction based on low-rank matrix and dictionary learning[J]. Journal of Computer Research and Development, 2016, 53(4): 884-891 (in Chinese)(楊帥鋒, 趙瑞珍. 基于低秩矩陣和字典學(xué)習(xí)的圖像超分辨率重建 DLOF[J]. 計算機(jī)研究與發(fā)展, 2016, 53(4): 884-891)

[10] Irani M, Peleg S. Motion analysis for image enhancement: resolution, occlusion, and transparency[J]. Journal of Visual Communication and Image Representation, 1993, 4(4): 324-347

[11] Stark H, Oskoui P. High-resolution image recovery from image-plane arrays using convex projections[J]. Journal of the Optical Society of America Optics and Image Science, 1989, 6(11): 1715-1726

[12] Levitan E, Herman G. A maximum a posteriori probability expectation maximization algorithm for image reconstruction in emission tomography[J]. IEEE Trans on Medical Imaging, 1987, 6(3): 185-192

[13] Boyd S, Parikh N, Chu E, et al. Distributed optimization and statistical learning via the alternating direction method of multipliers[J]. Foundations and Trends in Machine Learning, 2011, 3(1): 89-97

[14] Yang Sen, Wang Jie, Fan Wei, et al. An efficient ADMM algorithm for multidimensional anisotropic total variation regularization problem[C] //Proc of the 19th ACM SIGKDD Int Conf on Knowledge Discovery and Data Mining. New York: ACM, 2013: 641-649

[15] Horn B, Schunck B. Determining optical flow[J]. Artificial Intelligence, 1981, 17(1/2/3): 185-203

[16] Marius D, Sergiu N. Total variation regularization of local-global optical flow[C] //Proc of the 14th Int IEEE Conf on Intelligent Transportation Systems. Piscataway, NJ: IEEE, 2011: 318-323

[17] Perona P, Malik J. Scale space and edge detection using anisotropic diffusion[J]. IEEE Trans on Pattern Analysis and Machine Intelligence, 1990, 7(12): 629-639

[18] Tomasi C, Manduchi R. Bilateral filtering for gray and color images[C] //Proc of the 6th Int Conf on Computer Vision(ICCV). Los Alamitos, CA: IEEE Computer Society, 1998: 839-846

[19] Durand F, Dorsey J. Fast bilateral filtering for the display of high-dynamic-range images[J]. ACM Trans on Graphics, 2002, 21(3): 257-266

[20] Guo Qiang, Zhang Caiming, Zhang Yunfeng. An Efficient SVD-Based Method for Image Denoising[J]. IEEE Trans on Circuits and Systems for Video Technology, 2016, 26(5): 868-880

[21] Guo Qiang, Zhang Caiming, Zhang Yunfeng, et al. Low-rank image denoising based on minimum variance estimator[J]. Journal of Computer-Aided Design and Computer Graphics, 2015, 27(12): 2237-2246 (in Chinese)(郭強(qiáng), 張彩明, 張云峰, 等. 基于最小方差估計的圖像低秩去噪[J]. 計算機(jī)輔助設(shè)計與圖形學(xué)學(xué)報, 2015, 27(12): 2237-2246)

[22] Liu Jun, Yuan Lei, Ye Jieping. An efficient algorithm for a class of fused lasso problem[C] //Proc of the 19th ACM SIGKDD Int Conf on Knowledge Discovery and Data Mining. New York: ACM, 2010: 323-332

[23] Dong Weisheng, Zhang Lei, Shi Guangming, et al. Nonlocal back-projection for adaptive image enlargement[C] //Proc of the 16th IEEE Int Conf on Image. Piscataway, NJ: IEEE, 2009: 348-352

[24] Baker S, Scharstein D, Lewis J, et al. A dataset and evaluation methodology for optical flow[J]. International Journal of Computer Vision, 2011, 92(1): 1-31

[25] Castillo R, Castillo E, Guerra R, et al. A framework for evaluation of deformable image registration spatial accuracy using large landmark point sets[J]. Physical in Medicine and Biology, 2009, 54(7): 1849-1870

VariationalOpticalFlowEstimationBasedSuper-ResolutionReconstructionforLung4D-CTImage

Geng Fenghuan1,2, Liu Hui1,2, Guo Qiang1,2, and Yin Yilong1,3

1(CollegeofComputerScienceandTechnology,ShandongUniversityofFinanceandEconomics,Jinan250014)2(ShandongProvincialKeyLaboratoryofDigitalMediaTechnology,Jinan250014)3(CollegeofComputerScienceandTechnology,ShandongUniversity,Jinan250013)

The acquisition of lung 4D computed tomography (4D-CT) data is limited by the scanning time and radiation dose, which leads to the sampling rate in the axial direction is much less than that in the in-plane direction. In order to get better quality of 4D-CT images, based on the inherent self-similarity of medical images, a new method of image sequence super-resolution reconstruction is proposed in this paper. This method uses the local and global variational optical flow estimation to improve the quality of enlarged 4D-CT image. Firstly, we present a combined local and global variational optical flow model, in order to estimate the motion fields (i.e., the optical flow fields) between different phases in the corresponding positions. Then, the optical flow field is obtained by solving the model with the fast alternating direction method of multiplier. Finally, according to the calculated motion fields, we employ the improved non-local iterative back projection (NLIBP) algorithm to reconstruct high resolution lung images. The experimental results have shown that, in both quantification standard and visual perception, this method outperforms non-local iterative back projection algorithm and full search block matching based iterative back projection technique. Furthermore, our method can generate clear edges while enhancing the texture of images.

4D-CT image; super-resolution reconstruction; optical flow estimation; alternating direction method of multipliers; iterative back projection

Geng Fenghuan, born in 1992. Master candidate. His main research interests include image processing and machine learning.

Liu Hui, born in 1978. PhD, professor and master supervisor. Member of CCF. Her main research interests include image processing, computer vision, and machine learning.

Guo Qiang, born in 1979. PhD, associate professor and MSc supervisor. Member of CCF. His main research interests include image processing and computer vision.

Yin Yilong, born in 1972. PhD, professor and PhD supervisor. Senior member of CCF. His main research interests include machine learning, data mining, and computational medicine.

2017-05-23；

：2017-06-20

國家自然科學(xué)基金項目(61572286，61472220)；山東省重點研發(fā)項目(2014GGX101037)；濟(jì)南市高校自主創(chuàng)新項目(201401216)；山東省高校優(yōu)勢學(xué)科人才團(tuán)隊培育項目 This work was supported by the National Natural Science Foundation of China (61572286,61472220), the Project of Shandong Provincial Key Research and Development (2014GGX101037), the Project of Independent Innovation of Universities in Jinan (201401216), and the Fostering Project of Dominant Displine and Talent Team of Shandong Province Higher Education Institutions.

劉慧(liuh_lh@sdufe.edu.cn)

TP391