梯度與曲率變形力相結(jié)合的圖像配準(zhǔn)方法
——大形變非剛性圖像配準(zhǔn)

賈義亭，祝 軒，藺恩標(biāo)，雷文娟

西北大學(xué) 信息科學(xué)與技術(shù)學(xué)院，西安710127

1 引言

圖像配準(zhǔn)就是對(duì)在不同時(shí)間、不同傳感器等不同條件（天候、照度、攝像位置和角度等）下獲得的兩幅或多幅圖像根據(jù)某種配準(zhǔn)測(cè)度函數(shù)，尋找一種空間變換關(guān)系，使經(jīng)過(guò)空間變換后兩幅圖像間的相似性最大或者差異性最小。圖像配準(zhǔn)（Image registration）被廣泛地應(yīng)用于遙感數(shù)據(jù)分析、計(jì)算機(jī)視覺、醫(yī)學(xué)圖像分析和軍事目標(biāo)自動(dòng)識(shí)別等領(lǐng)域，是圖像處理領(lǐng)域的一個(gè)重要研究方向。圖像配準(zhǔn)方法通常分為兩類：一類是剛性配準(zhǔn)，另一類是非剛性配準(zhǔn)。剛性配準(zhǔn)適合于剛性物體形變的配準(zhǔn)，其只能對(duì)圖像進(jìn)行全局線性變換，使它在圖像配準(zhǔn)方面存在很大的局限性，目前剛性配準(zhǔn)技術(shù)已經(jīng)發(fā)展得比較成熟。自然界中多數(shù)物體的形變都屬于非剛性形變，非剛性配準(zhǔn)[1]能夠更精確地實(shí)現(xiàn)局部變形，其配準(zhǔn)技術(shù)有待進(jìn)一步提高，基于此圖像的非剛性配準(zhǔn)算法研究受到相關(guān)人員的高度關(guān)注。

非剛性配準(zhǔn)方法主要分為兩大類：基于像素的方法和基于特征的方法。其中基于像素的方法不需要對(duì)圖像進(jìn)行分割，可避免因分割帶來(lái)的邊緣效應(yīng)和誤差，是目前公認(rèn)的較好的方法之一。Hellier 等人[2]對(duì)配準(zhǔn)方法的比較研究表明Demons 配準(zhǔn)算法對(duì)非剛性形變具有較高的配準(zhǔn)精度。文獻(xiàn)[3]則指出Demons 算法只適用于小形變配準(zhǔn)，對(duì)大形變問(wèn)題Demons 算法將不能保持圖像的拓?fù)浣Y(jié)構(gòu)，并產(chǎn)生物理上的不合理形變。此外，Demons算法還存在僅靠梯度信息來(lái)確定浮動(dòng)圖像的位移量，當(dāng)梯度非常小時(shí)，圖像形變的方向不能被確定，往往會(huì)導(dǎo)致錯(cuò)誤的配準(zhǔn)變換[3]。針對(duì)Demons 算法的弊病，Wang 等人[4]提出了改進(jìn)的Active Demons 算法，將牛頓作用力與反作用力原理引入Demons算法，假設(shè)擴(kuò)散是雙向的，認(rèn)為不但可以用參考圖像的梯度信息作為驅(qū)動(dòng)力使浮動(dòng)圖像向參考圖像擴(kuò)散，同時(shí)施加反向作用力使參考圖像向浮動(dòng)圖像擴(kuò)散。Active Demons算法能夠部分克服Demons 算法的缺陷，特別是在提高大形變配準(zhǔn)的正確性和一致性方面具有一定優(yōu)勢(shì)。

Active Demons 配準(zhǔn)算法[5]是一種將形變配準(zhǔn)視作擴(kuò)散問(wèn)題的基于灰度的非剛性配準(zhǔn)方法，其利用參考圖像的灰度梯度信息作為驅(qū)動(dòng)力實(shí)現(xiàn)浮動(dòng)圖像與參考圖像的配準(zhǔn)。

隨著對(duì)圖像局部特征研究的深入，人們認(rèn)識(shí)到僅用一階微分量——梯度來(lái)表征圖像局部特征是不夠的，二階微分量中含有更豐富的信息。Sapiro 和Caselles 指出圖像中物體的形狀信息取決于等照度線（isophoto）的形態(tài)學(xué)特征[6-7]。在文獻(xiàn)[3]中將圖像等照度線的曲率作為一個(gè)校正形變的驅(qū)動(dòng)力因素引入Demons 擴(kuò)散方程，建立了一個(gè)基于梯度和曲率變形力相結(jié)合的擴(kuò)散模型（簡(jiǎn)稱G&C model），實(shí)現(xiàn)了小形變圖像的精確配準(zhǔn)。更進(jìn)一步將等照度線的曲率引入Active Demons 算法，提出改進(jìn)的Active G&C 模型，并在Active G&C 模型應(yīng)用于大形變圖像配準(zhǔn)的算法實(shí)現(xiàn)過(guò)程中加入多分辨率策略，使梯度與曲率變形力相結(jié)合的Active G&C 配準(zhǔn)模型不僅能夠克服Active Demons算法僅依靠梯度信息易導(dǎo)致錯(cuò)誤配準(zhǔn)的弊病，而且在提高大形變圖像的配準(zhǔn)精度方面明顯優(yōu)于Active Demons 算法。本文通過(guò)實(shí)驗(yàn)對(duì)上述內(nèi)容進(jìn)行驗(yàn)證，結(jié)果表明提出的方法行之有效。

2 Active Demons算法

Demons 算 法 是Thirion 等 人[5]根 據(jù) 經(jīng) 典 的Maxwell’s demons 熱力學(xué)實(shí)驗(yàn)，提出的一種基于光流場(chǎng)的非剛性圖像配準(zhǔn)模型，其配準(zhǔn)思想是將配準(zhǔn)過(guò)程看作是浮動(dòng)圖像的像素向參考圖像逐步擴(kuò)散的過(guò)程，擴(kuò)散強(qiáng)度由參考圖像的梯度信息確定，擴(kuò)散方程如下：

其中m、s分別表示浮動(dòng)圖像和參考圖像，?s表示參考圖像的梯度。

文獻(xiàn)[4]提出改進(jìn)的Active Demons 算法，該算法根據(jù)牛頓作用力與反作用力的原理，假設(shè)擴(kuò)散是雙向的，不僅利用參考圖像的灰度梯度信息?s使浮動(dòng)圖像向參考圖像擴(kuò)散，而且利用浮動(dòng)圖像的灰度梯度信息?m使參考圖像向浮動(dòng)圖像擴(kuò)散，定義如下基于浮動(dòng)圖像的灰度梯度信息驅(qū)動(dòng)力um：

其中，“-”表示反作用力。

結(jié)合式（1）和（2），得Active Demons算法驅(qū)動(dòng)力：

其中，前項(xiàng)促使浮動(dòng)圖像向參考圖像移動(dòng)，后項(xiàng)促使參考圖像向浮動(dòng)圖像移動(dòng)。

具體在算法實(shí)現(xiàn)過(guò)程中引入歸一化參數(shù)α[4]，用來(lái)自適應(yīng)地調(diào)節(jié)每次迭代中的驅(qū)動(dòng)力，以加快收斂速度，得式（4）：

3 梯度與曲率變形力相結(jié)合的大形變非剛性圖像配準(zhǔn)方法

記灰度圖像u(x,y)，其梯度是一階量，可表示為?u=(ux,uy)，其中ux=?u/?x，uy=?u/?y，梯度信息反映圖像中像素灰度值的變化強(qiáng)弱。圖像的等照度線定義為Xλu={[x,y],u(x,y)=λ}，等照度線是指圖像中光順的曲線，反映圖像中物體的形狀。文獻(xiàn)[6-7]指出圖像中物體的形狀信息取決于等照度線的形態(tài)學(xué)特征。鑒于真實(shí)圖像中的等照度線足夠光順，相應(yīng)地描述等照度線形態(tài)學(xué)特征的二階微分量——曲率的變化應(yīng)該較小，然而因?yàn)楦?dòng)圖像存在形變，造成其中的等照度線發(fā)生畸變，使曲率發(fā)生顯著變化。基于上述分析，認(rèn)為Demons 和Active Demons 算法僅靠梯度信息決定配準(zhǔn)驅(qū)動(dòng)力是不夠的，為了使畸變的等照度線也變得比較光順，將等照度線曲率作為另一個(gè)驅(qū)動(dòng)力因素引入Active Demons 方程，使浮動(dòng)圖像在梯度和曲率雙重驅(qū)動(dòng)作用力下形變，定義如下兩個(gè)曲率驅(qū)動(dòng)力uk(m)和uk(s)：

其中，f(k(s))和f(k(m))分別是以參考圖像的曲率k(s)和浮動(dòng)圖像的曲率k(m)為自變量的曲率驅(qū)動(dòng)函數(shù)，其作用是保持圖像的結(jié)構(gòu)信息，使等照度線變的光滑，“-”表示反作用力。

結(jié)合式（5）和（6），并引入歸一化參數(shù)α，得：

將該模型稱為梯度和曲率變形力相結(jié)合的配準(zhǔn)模型（簡(jiǎn)稱Active G&C model）。

本文選取如下曲率驅(qū)動(dòng)函數(shù)：

以參考圖像的曲率ks計(jì)算為例，如式（9）：

浮動(dòng)圖像的曲率km計(jì)算只需將式（9）中的s 換成m。

為了更進(jìn)一步提高大形變圖像的配準(zhǔn)精度，在迭代求解Active G&C model 的過(guò)程中加入金字塔多分辨率策略，對(duì)浮動(dòng)圖像和參考圖像進(jìn)行多級(jí)金字塔分解，將低級(jí)金字塔圖像計(jì)算得到的位移偏移量超采樣后作為上一級(jí)金字塔圖像的初始變換矩陣（如圖1 所示），按照由粗到精逐級(jí)迭代的方式對(duì)配準(zhǔn)結(jié)果優(yōu)化。多分辨率策略不僅有利于提高配準(zhǔn)算法的執(zhí)行速度和精度，更重要的是可以避免配準(zhǔn)過(guò)程中出現(xiàn)局部極值[8]。

圖1 多分辨率策略示意圖

Active G&C 算法實(shí)現(xiàn)過(guò)程如下：

步驟1設(shè)定初始偏移矩陣T0，設(shè)置金字塔級(jí)數(shù)K ，對(duì)浮動(dòng)圖像m 和參考圖像s 進(jìn)行K 級(jí)金字塔分解。

步驟2當(dāng)前金字塔級(jí)數(shù)k ＜K 時(shí)，根據(jù)式（7）進(jìn)行n 次迭代，得浮動(dòng)圖像m 相對(duì)參考圖像s 的位移偏移矩陣更新偏移矩陣Tk=Tk-1+ukn。

步驟3用Tk使浮動(dòng)圖像形變，用均方誤差MSE（式（10））作為形變后的浮動(dòng)圖像與參考圖像的相似性測(cè)度，當(dāng)相似性測(cè)度小于某設(shè)定閾值時(shí)，認(rèn)為第k 級(jí)分辨率圖像已收斂，進(jìn)入第k+1 級(jí)，執(zhí)行步驟4，否則繼續(xù)迭代，計(jì)算ukn+1。

步驟4對(duì)前一級(jí)圖像的偏移矩陣Tk-1超采樣，作為當(dāng)前圖像配準(zhǔn)的初始位移偏移量 執(zhí)行步驟2。

步驟5K 級(jí)金字塔迭代完成后，將最終的位移偏移矩陣作用于浮動(dòng)圖像，得到最終的配準(zhǔn)結(jié)果。

4 實(shí)驗(yàn)結(jié)果及評(píng)價(jià)

4.1 實(shí)驗(yàn)結(jié)果

在Matlab7.6 環(huán)境下，分別用Demons、Active Demons和本文方法分別對(duì)大形變Woman 圖像進(jìn)行配準(zhǔn)。實(shí)驗(yàn)結(jié)果如圖2 所示。實(shí)驗(yàn)迭代次數(shù)為200 次，進(jìn)行三級(jí)金字塔分解，Active Demons和本文Active G&C 算法中α 取2.5。

一般情況下金字塔級(jí)數(shù)越高，配準(zhǔn)效果越好，但是配準(zhǔn)時(shí)間越長(zhǎng)。可根據(jù)所處理對(duì)象的不同，權(quán)衡配準(zhǔn)精度和時(shí)間效率之間的利弊加以選擇。

觀察圖2可以看出本文方法在配準(zhǔn)大形變圖像時(shí)具有更好的視覺效果（如眼睛和嘴的部位），明顯優(yōu)于其他兩種方法。

4.2 結(jié)果評(píng)價(jià)

為了進(jìn)一步客觀評(píng)價(jià)上述三種配準(zhǔn)方法的配準(zhǔn)質(zhì)量，以大形變Elaine 圖為例，采用不同的配準(zhǔn)方法得圖3 所示配準(zhǔn)結(jié)果。通過(guò)在Elaine 圖中的相同位置劃線，繪制相應(yīng)的灰度偏移曲線，比較不同配準(zhǔn)方法配準(zhǔn)大形變Elaine 圖像后與參考圖像之間的灰度偏移量，得圖4（a）和（b），其中特別給出了本文方法和Active Demons 方法的灰度偏移比較曲線，如圖4（b）。

觀察對(duì)比圖4 中灰度偏移曲線可知，采用多分辨策略的Active G&C 方法在處理大形變圖像配準(zhǔn)問(wèn)題時(shí)優(yōu)越性明顯，配準(zhǔn)結(jié)果中的各個(gè)像素點(diǎn)的灰度值與參考圖像基本一致，可實(shí)現(xiàn)圖像局部的精確配準(zhǔn)，視覺效果更好。

分別用均方誤差、峰值信噪比、相關(guān)系數(shù)、歸一化互信息和結(jié)構(gòu)相似度[9]五種客觀評(píng)價(jià)方法比較Demons、Active Demons 和Active G&C 方法配準(zhǔn)Elaine 圖像的配準(zhǔn)質(zhì)量，得表1。

圖2 Woman 圖像配準(zhǔn)

圖3 Elaine圖像配準(zhǔn)

圖4 灰度偏移比較

表1 質(zhì)量評(píng)價(jià)比較

相關(guān)系數(shù)可以反映兩幅圖像之間的相似度，按式（11）計(jì)算：

Scc越大配準(zhǔn)效果越好。

結(jié)構(gòu)相似度（SSIM）按式（12）計(jì)算：

其中，s為參考圖像，T(m)為配準(zhǔn)結(jié)果，S(s,T(m))是參考圖像s與配準(zhǔn)結(jié)果T( )m之間的結(jié)構(gòu)相似度度量函數(shù)，α＞0,β＞0,γ＞0，用來(lái)調(diào)整圖像亮度、對(duì)比度和結(jié)構(gòu)信息所占的比重，l(·),c(·),s(·)分別是浮動(dòng)圖像和參考圖像之間亮度、對(duì)比度和結(jié)構(gòu)的比較函數(shù)。

從表1 中可以看出，Active G&C 方法相比其他兩種配準(zhǔn)方法優(yōu)勢(shì)更明顯。

5 結(jié)束語(yǔ)

本文將等照度線曲率引入Active Demons 算法，同時(shí)在算法迭代過(guò)程中采用金字塔多分辨率策略，建立了一個(gè)具有曲率和梯度雙重驅(qū)動(dòng)力相結(jié)合的配準(zhǔn)擴(kuò)散模型，該模型克服了僅靠梯度信息易導(dǎo)致錯(cuò)誤圖像配準(zhǔn)以及對(duì)大形變圖像配準(zhǔn)精度不高的弊病。實(shí)驗(yàn)表明，本文方法是一種行之有效的圖像配準(zhǔn)算法。

[1] Crum W R，Hartkens T.Non-rigid image registration：theory and practice[J].British Journal of Radiology，2004，77：140-153.

[2] Hellier P，Barillot C，Corouge I，et al.Retrospective evaluation of inter-subject brain registration[J].IEEE Transactions on Medical Imaging，2003，22（9）：1120-1130.

[3] 雷文娟.非剛性圖像配準(zhǔn)方法研究[D].西安：西北大學(xué)信息科學(xué)與技術(shù)學(xué)院，2008.

[4] Wang He，Dong Lei.Validation of an accelerated“demons”algorithm for deformable image registration in radiation therapy[J].Physics in Medicine Biology，2005，50（12）：2887-2905.

[5] Thirion J P.Image matching as a diffusion process：an analogy with Maxwell’s demons[J].Medical Image Analysis，1998，2（3）：243-260.

[6] Caselles V.Shape preserving local contrast enhancement[C]//Proceedings of the 1997 International Conference on Image Processing（ICIP’97）.Washington：IEEE Computer Society，1997：314-317.

[7] Sapiro G，Caselles V.Histogram modification via differential equations[J].Journal of Differential Equations，1997，135（2）：238-268.

[8] 張紅穎，張加萬(wàn)，孫濟(jì)洲.改進(jìn)Demons 算法的非剛性醫(yī)學(xué)圖像配準(zhǔn)[J].光學(xué)精密工程，2007，15（1）：145-150.

[9] Wang Zhou，Bovik A C.Modern image quality assessment（synthesis lectures on image，video & multimedia processing）[M].[S.l.]：Morgan & Claypool Publishers，2006.