劉 薇,陳雷霆

1.電子科技大學(xué) 計(jì)算機(jī)科學(xué)與工程學(xué)院,成都611731; 2.成都工業(yè)職業(yè)技術(shù)學(xué)院,成都 610218)(*通信作者電子郵箱WeiLiu0405@gmail.com)

基于自適應(yīng)切空間的MRI圖像配準(zhǔn)

劉 薇1,2*,陳雷霆1

1.電子科技大學(xué) 計(jì)算機(jī)科學(xué)與工程學(xué)院,成都611731; 2.成都工業(yè)職業(yè)技術(shù)學(xué)院,成都 610218)(*通信作者電子郵箱WeiLiu0405@gmail.com)

微分同胚是一種光滑可逆的變換,在MRI圖像配準(zhǔn)中可以保證圖像形變后的拓?fù)浣Y(jié)構(gòu)保持不變，同時(shí)避免出現(xiàn)不合理的物理現(xiàn)象。為了在空間變換中獲得更合理的同胚映射,高維空間中數(shù)據(jù)的非線性結(jié)構(gòu)被考慮,基于流形學(xué)習(xí)方法提出一種自適應(yīng)切空間的MRI圖像配準(zhǔn)算法。首先,把MRI數(shù)據(jù)構(gòu)造成對稱正定(SPD)的協(xié)方差矩陣,然后形成李群;接著,利用樣本點(diǎn)鄰域的局部切空間來表示李群的幾何結(jié)構(gòu)的非線性;接下來,在流形上用自適應(yīng)鄰域選擇的方法形成的線性子空間去逼近局部切空間,提高切空間的局部線性化程度,從而最大限度地保留流形的局部非線性結(jié)構(gòu),得到最優(yōu)的同胚映射。仿真數(shù)據(jù)和臨床數(shù)據(jù)的實(shí)驗(yàn)結(jié)果顯示,與傳統(tǒng)的非參數(shù)微分同胚配準(zhǔn)算法相比,該算法在高維稠密形變場上獲得更高的拓?fù)浔３侄?最終提高圖像配準(zhǔn)精度。

微分同胚;切空間;李群;鄰域選擇;MRI圖像配準(zhǔn)

0 引言

正常情況下,配準(zhǔn)后的解剖結(jié)構(gòu)圖像與參考圖像應(yīng)具有相同的拓?fù)浣Y(jié)構(gòu),不會(huì)產(chǎn)生新的組織,原有的結(jié)構(gòu)不會(huì)消失,形變域也不會(huì)出現(xiàn)撕裂、折疊、空洞等不合理的物理結(jié)構(gòu)。微分同胚是一種平滑、連續(xù)和可逆的變換,能夠使解剖結(jié)構(gòu)圖像變形后的拓?fù)浣Y(jié)構(gòu)保持不變,在醫(yī)學(xué)圖像配準(zhǔn)中是一個(gè)重要的應(yīng)用。

Marsland等[1]把微分同胚變換看作時(shí)變速度場,用測地插值樣條基構(gòu)造變換;Beg等[2]提出在速度場空間利用基于歐拉-拉格朗日方程的變分方法求解大形變的微分同胚變換,采用梯度下降法更新形變參數(shù);Ashburner等[3]對文獻(xiàn)[2]進(jìn)行改進(jìn),任意時(shí)刻的速度場都用初始速度場來表示,配準(zhǔn)的每次迭代都利用初始形變值來計(jì)算,減少了對內(nèi)存和外存的要求;Ashburner[4]提出非時(shí)變速度場的配準(zhǔn)方法,同胚變換構(gòu)成了復(fù)合運(yùn)算下的李群,將大形變同胚變換分解為一系列小形變來處理,使可逆變換相對容易,并降低計(jì)算代價(jià);Janssens等[5]在微分同胚方法中利用累加位移場的可逆性配準(zhǔn)對比增強(qiáng)度不同的圖像;Arsigny等[6]借助李群理論提出Log-Euclidean框架,在微分同胚空間里可以對向量場進(jìn)行分析統(tǒng)計(jì)并且能保持變換的可逆性;Vercauteren等[7]基于Demons算法,提出一種非參數(shù)微分同胚的Demons配準(zhǔn)算法——DD-NP算法,利用李群理論在連續(xù)域上進(jìn)行空間變換,使空間變換具有微分同胚性,從而形變場具有拓?fù)浔３中?閆德勤等[8]在保持圖像的局部和全局流形拓?fù)浣Y(jié)構(gòu)的基礎(chǔ)上,提出快速微分同胚變換速度場更新方法。

但是,上述方法都沒有考慮高維微分同胚變換空間中數(shù)據(jù)的非線性結(jié)構(gòu),高維數(shù)據(jù)包含更多的結(jié)構(gòu)信息,影響配準(zhǔn)結(jié)果。針對這個(gè)問題,結(jié)合流形學(xué)習(xí)方法[9],本文改進(jìn)DD-NP配準(zhǔn)算法,提出了一種基于自適應(yīng)切空間的MRI圖像配準(zhǔn)方法。首先,構(gòu)造正定對稱的協(xié)方差矩陣,協(xié)方差的空間結(jié)構(gòu)是一個(gè)高維黎曼流形,一定條件下形成李群;其次,利用樣本點(diǎn)鄰域的局部切空間來表示李群流形的幾何結(jié)構(gòu)的非線性;接下來,由于局部切空間的線性化程度越高, 高維數(shù)據(jù)流形的非線性結(jié)構(gòu)信息就越豐富, 在空間變換中則可以最大限度上保留流形的局部非線性結(jié)構(gòu)，于是用自適應(yīng)鄰域選擇的方法形成的線性子空間去逼近局部切空間,更好地在空間變換中保持圖像的拓?fù)浣Y(jié)構(gòu),從而提高圖像配準(zhǔn)的精度。

1 相關(guān)背景

1.1 李群和切空間

流形是一個(gè)拓?fù)淇臻g,微分流形是把m維流形M加上微分結(jié)構(gòu)。李群既是一種微分流形也是一種群。在李群的單位元Id處的切空間構(gòu)成一個(gè)向量空間, 即李代數(shù)。為了方便描述,先作如下規(guī)定：

1)GL(m)表示m×m矩陣構(gòu)成的李群；

2)gl(m)表示m×m矩陣構(gòu)成的李代數(shù)(線性向量空間)；

3)log()表示對數(shù)映射；

4)exp()表示指數(shù)映射。

指數(shù)和對數(shù)映射可由式(1)進(jìn)行轉(zhuǎn)換：

Si=exp(si)；si=log(Si)

(1)

其中:Si∈GL(m),si∈gl(m)。指數(shù)映射使m維流形局部同胚于m維向量空間,即李代數(shù)零元素的鄰域通過指數(shù)函數(shù)映射到李群的單位元鄰域,指數(shù)映射及其逆映射對數(shù)映射都是光滑可逆的同胚映射,如圖1所示。但是在高維流形空間,影響切空間逼近的鄰域大小目前并不清楚[10]。

圖1 李群和李代數(shù)之間的映射關(guān)系

1.2 正定對稱矩陣和李群

正定對稱的形式有很多,比如區(qū)域協(xié)方差特征描述子[11]和結(jié)構(gòu)張量[12]。正定對稱結(jié)構(gòu)是包含了豐富結(jié)構(gòu)信息的黎曼流形，但不具有群結(jié)構(gòu),記Sym+(m)是m×m的正定對稱矩陣構(gòu)成的流形空間。根據(jù)文獻(xiàn)[13],附加式(2)約束則可以形成李群(Sym+,⊙)。

Si⊙Sj=exp(log(Si)+log(Sj))

(2)

2 DD-NP模型

minE(u)=Sim(Ir,If°φ)+Reg(φ)

(3)

其中:Sim(·) 和Reg(·)分別表示相似性函數(shù)和正則化函數(shù)。式(3)的目的就是尋找一個(gè)最優(yōu)的空間變換φ(Si)=φ(Si,t)。其中：t∈[0,1]；φ(Si,t)是速度流；u是依賴于時(shí)間的速度場。φ通過求式(4)的常微分方程可得:

(4)

(5)

在條件exp(u(0))=Id約束下,李代數(shù)零點(diǎn)鄰域和李群單位元鄰域的同胚對應(yīng)關(guān)系是唯一的,以保證空間變換的一一映射特性,配準(zhǔn)前后的浮動(dòng)圖像仍保持鄰接關(guān)系,使配準(zhǔn)后的圖像結(jié)構(gòu)保持拓?fù)湫浴?/p>

3 基于自適應(yīng)切空間的配準(zhǔn)算法

3.1 自適應(yīng)鄰域選擇

流形上任一點(diǎn)附近的局部線性結(jié)構(gòu)都能被流形在該點(diǎn)處的切空間描述,切空間是該局部線性化子空間,線性化程度越高,包含的流形的非線性結(jié)構(gòu)越多。算法的目的是在李群流形下,用自適應(yīng)的鄰域選擇方法[9]去逼近單位元處的切空間。鄰域大小與樣本點(diǎn)的曲率密切相關(guān),流形上曲率是高度變化的,曲率大的樣本點(diǎn)的鄰域應(yīng)該相對比較小, 而流形上曲率小的樣本點(diǎn)處的鄰域比較大,這樣能更精確地逼近流形的局部幾何結(jié)構(gòu)。

假設(shè)李群的單位元Id處切空間是TG,其d維正交基矩陣為T=[τ1,τ2,…,τd],τi∈Rm,i=1,2,…,d,在局部坐標(biāo)系中,任意樣本點(diǎn)切向量空間(形變場)u可以表示為式(6)的線性組合：

(6)

其中：ui是給定坐標(biāo)下u的分量,是線性組合系數(shù)。

(7)

設(shè)U=[u1,u2,…,uk], 且U∈Rd×k,是對應(yīng)正交基矩陣的局部坐標(biāo)系:

(8)

U=diag(σ1,σ2,…,σd)VT

(9)

從而有

(10)

以及

(11)

由式(10)和(11)得到鄰域選擇標(biāo)準(zhǔn)式(12):

(12)

3.2 基于自適應(yīng)切空間配準(zhǔn)算法

本文算法步驟如下:

輸入:參考圖像Ir和浮動(dòng)圖像If。

輸出:配準(zhǔn)圖像。

步驟1 對Ir和If中的像素提取對稱正定的圖像特征。

步驟4 設(shè)定閾值α∈[0,1], 如果r<α, 則算法結(jié)束；否則到下一步。

步驟7 計(jì)算李群單位元Id處鄰域所對應(yīng)的切空間的正交基矩陣T=[τ1,τ2,…,τd],獲得形變場u。

步驟8 通過最小化式(3)中的代價(jià)函數(shù)得到最優(yōu)空間變換φ。

4 實(shí)驗(yàn)與分析

為了驗(yàn)證本文算法的有效性,在仿真數(shù)據(jù)和臨床數(shù)據(jù)上對大腦白質(zhì)和大腦灰質(zhì)的MRI圖像分別進(jìn)行實(shí)驗(yàn),并與DD-NP算法對比。實(shí)驗(yàn)平臺(tái)為:3.4GHzCPU,8GBRAM內(nèi)存,Windows8操作系統(tǒng),Matlab2012b, 采用C++和Matlab混合編程。評價(jià)標(biāo)準(zhǔn)采用可視評估和均方根誤差(RootMeanSquaredError,RMSE)。

4.1 構(gòu)造正定對稱的圖像特征

對選取的每個(gè)樣本像素構(gòu)造一個(gè)136×136的正定對稱的協(xié)方差矩陣作為圖像特征。先提取136維局部特征fi,包括128維的SIFT(Scale-Invariant Feature Transform)描述子、坐標(biāo)位置(x,y)、灰度值Ii(xy)、灰度值的一階梯度、一階梯度的模和二階梯度,如式(13)所示:

(13)

然后用fi和fi的轉(zhuǎn)置作外積計(jì)算,得到式(14)所示的正定對稱的協(xié)方差矩陣:

(14)

為了減少噪聲影響和計(jì)算代價(jià),對協(xié)方差矩陣形成的高維流形Sym+(m)作PCA降維處理,并在實(shí)驗(yàn)中發(fā)現(xiàn)當(dāng)維數(shù)m=126,116,106時(shí),得到的效果最好。

4.2 仿真數(shù)據(jù)實(shí)驗(yàn)及分析

仿真數(shù)據(jù)來源于BrainWeb數(shù)據(jù)庫,灰度值在0～255, 圖像大小為181×217×181。先從數(shù)據(jù)庫中隨機(jī)抽取兩個(gè)主體(subject18和subject46),然后實(shí)驗(yàn)分為兩組:第一組從兩個(gè)主體的大腦白質(zhì)圖像中隨機(jī)選擇橫斷面,冠狀面和矢狀面各一張,一張作為參考圖像,另一張作為浮動(dòng)圖像;第二組從兩個(gè)主體的大腦灰質(zhì)圖像中隨機(jī)選擇橫斷面,冠狀面和矢狀面各一張,一張作為參考圖像,另一張作為浮動(dòng)圖像。兩組實(shí)驗(yàn)各進(jìn)行10次,最后計(jì)算每10次實(shí)驗(yàn)的平均RMSE。

圖2～3顯示了不同情況下本文算法最好的可視化結(jié)果和DD-NP算法結(jié)果的對比。從可視效果的角度,本文算法在這幾種情況下的可視效果明顯優(yōu)于DD-NP算法。

圖2 仿真數(shù)據(jù)在大腦白質(zhì)圖像上的配準(zhǔn)結(jié)果

圖3 仿真數(shù)據(jù)在大腦灰質(zhì)圖像上的配準(zhǔn)結(jié)果

表1～2對比了在大腦白質(zhì)和灰質(zhì)中不同維數(shù)和鄰域?qū)?yīng)的配準(zhǔn)誤差結(jié)果，可看出配準(zhǔn)的精度不僅和鄰域大小有關(guān)，也和流形的維數(shù)有關(guān)。從三個(gè)不同斷面的圖像配準(zhǔn)結(jié)果可以看到,盡管本文算法有的情況下實(shí)驗(yàn)誤差要比DD-NP算法的誤差大,但是在每一個(gè)斷面都能獲得誤差小于DD-NP算法的結(jié)果。例如在白質(zhì)圖像橫斷面的配準(zhǔn)中,DD-NP的誤差為2.26,本文算法當(dāng)維數(shù)m=126,k=22;m=116,k=18和m=116,k=20這三種情況下,誤差((表中加粗?jǐn)?shù)據(jù))均小于DD-NP，樣本點(diǎn)非線性結(jié)構(gòu)的投影方向和流形的曲率密切相關(guān)，流形上點(diǎn)切方向的曲率反映了曲率沿著切方向的彎曲程度,也就是局部線性化的程度。樣本點(diǎn)的曲率越大,線性化程度越低,所以為了達(dá)到高精度的切空間逼近,那么鄰域選取可以較小;在極小主曲率方向的線性化程度最高,鄰域選取可以較大。從圖2～3可看出,冠狀面樣本點(diǎn)的曲率要比橫斷面和矢狀面小,所以選擇較大的鄰域可以得到更好的實(shí)驗(yàn)效果，但是鄰域的選取也會(huì)受流形維數(shù)的影響。

表1 大腦白質(zhì)圖像上RMSE對比 mm

表2 大腦灰質(zhì)圖像上RMSE對比 mm

4.3 臨床數(shù)據(jù)實(shí)驗(yàn)及分析

數(shù)據(jù)采集來自西門子核磁共振儀,回波時(shí)間TE=2.98 ms, 回波鏈長度ETL=1,反轉(zhuǎn)時(shí)間TI=1 100 ms,重復(fù)時(shí)間TR=2 530 ms,視野FOV=87.5,圖像大小448×512×192,20個(gè)健康主體大腦圖像。隨機(jī)抽選4個(gè)主體的矢狀面圖像,然后分為兩組:一組用于分割大腦白質(zhì),另一組用于分割大腦灰質(zhì)。每組選一張作為參考圖像,另一張為浮動(dòng)圖像。

圖4(a)和(b)分別顯示了臨床數(shù)據(jù)的大腦白質(zhì)和灰質(zhì)配準(zhǔn)在兩種算法上的對比可視化結(jié)果,本文算法的灰質(zhì)圖像情況下的效果明顯優(yōu)于DD-NP。圖5顯示在m=116,k=16時(shí)，本文算法與DD-NP算法在大腦白質(zhì)和大腦灰質(zhì)圖像上的誤差與迭代關(guān)系。由于需要在高維空間和切空間之間作映射,本文算法需要更多的迭代次數(shù)才能達(dá)到收斂。DD-NP算法迭代30次后開始產(chǎn)生收斂,而本文算法需要60次,但是第60次迭代以后本文算法的誤差小于DD-NP算法。

圖4 臨床數(shù)據(jù)的配準(zhǔn)結(jié)果對比

圖5 臨床數(shù)據(jù)中RMSE的對比

5 結(jié)語

傳統(tǒng)的微分同胚配準(zhǔn)算法沒有考慮高維空間數(shù)據(jù)的非線性結(jié)構(gòu);然而, 高維數(shù)據(jù)流形的非線性結(jié)構(gòu)受局部線性化程度影響。本文利用流形學(xué)習(xí)方法提出一種基于切空間的配準(zhǔn)算法,通過自適應(yīng)的鄰域選擇更精確地逼近流形的線性切空間,保持流形的非線性結(jié)構(gòu),使圖像在微分同胚變換中更好地保持拓?fù)浣Y(jié)構(gòu)。實(shí)驗(yàn)顯示,比起傳統(tǒng)算法,本文算法配準(zhǔn)精度有明顯提高。

進(jìn)一步的工作將關(guān)注以下兩點(diǎn):1)如何減少由于高維空間中的迭代映射造成的計(jì)算代價(jià);2)目前在理論上并不清楚多大的鄰域和流形維數(shù)能獲得最優(yōu)結(jié)果,算法只能通過枚舉去搜索,如果能找到規(guī)律,將大大減少工作量。

References)

[1] MARSLAND S, TWINING C J. Constructing diffeomorphic representations for the groupwise analysis of nonrigid registrations of medical images [J]. IEEE Transactions on Medical Imaging, 2004, 23(8): 1006-1020.

[2] BEG M F, MILLER M I, TROUVA, et al. Computing large deformation metric mappings via geodesic flows of diffeomorphisms [J]. International Journal of Computer Vision, 2005, 61(2): 139-157.

[3] ASHBUMER J, FRISTON K J. Diffeomorphic registration using geodesic shooting and Gauss-Newton optimization [J]. Neuroimage, 2011, 55(3): 954-967.

[4] ASHBUMER J. A fast diffeomorphic image registration algorithm [J]. Neuroimage, 2007, 38(1): 95-113.

[5] JANSSENS G, JACQUES L, DE XIVRY J O, et al. Diffeomorphic registration of images with variable contrast enhancement [J]. Journal of Biomedical Imaging, 2011, 2011: 3.

[6] ARSIGNY V, COMMOWICK O, PENNEC X, et al. A Log-Euclidean framework for statistics on diffeomorphisms[C]// MICCAI 2006: Proceedings of the 9th International Conference on Medical Image Computing and Computer-Assisted Intervention, LNCS 4190. Berlin: Springer, 2006: 924-931.

[7] VERCAUTEREN T, PENNEC X, PERCHANT A, et al. Diffeomorphic demons: efficient non-parametric image registration [J]. Neuroimage, 2009, 45(1): S61-S72.

[8] 閆德勤, 劉彩鳳, 劉勝藍(lán), 等. 大形變微分同胚圖像配準(zhǔn)快速算法[J]. 自動(dòng)化學(xué)報(bào), 2015,4(8):1461-1470.(YAN D Q, LIU C F, LIU S L, et al. A fast image registration algorithm for diffeomorphic image with large deformation[J]. Acta Automatica Sinica, 2015,4(8):1461-1470.)

[9] ZHANG Z, WANG J, ZHA H. Adaptive manifold learning [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2012, 34(2): 253-265.

[10] ARSIGNY V. Processing data in Lie groups: an algebraic approach. application to non-linear registration and diffusion tensor MRI [EB/OL]. [2016- 03- 01]. https://tel.archives-ouvertes.fr/tel-00121162/.

[11] TUZEL O, PORUKI F, MEER P. Region covariance: a fast descriptor for detection and classification [C]// ECCV 2006: Proceedings of the 9th European Conference on Computer Vision, LNCS 3952. Berlin: Springer, 2006: 589-600.

[12] GOH A, LENGLET C, THOMPSON P M, et al. A nonparametric Riemannian framework for processing high angular resolution diffusion images (HARDI) [C]// CVPR 2009: Proceedings of the 2009 IEEE Conference on Computer Vision and Pattern Recognition. Washington, DC: IEEE Computer Society, 2009: 2496-2503.

[13] ARSIGNY V, FILLARD P, PENNEC X, et al. Geometric means in a novel vector space structure on symmetric positive-definite matrices [J]. SIAM Journal on Matrix Analysis and Applications, 2007, 29(1): 328-347.

[14] TENENBAUM J B, DE SILVA V, LANGFORD J C. A global geometric framework for nonlinear dimensionality reduction [J]. Science, 2000, 290(5500): 2319-2323.

This work is partially supported by the Sheng-Bu Industry-Academia-Research Joint Project of Education Department and Science & Technology Department of Guangdong Province (2012A090300001).

LIU Wei, born in 1971, Ph. D. candidate, associate professor. Her research interests include computer vision, medical image processing.

CHEN Leiting, born in 1966, Ph. D., professor. His research interests include image processing, computer graphics, virtual reality.

MRI image registration based on adaptive tangent space

LIU Wei1,2*, CHEN Leiting1

(1. School of Computer Science and Engineering, University of Electronic Science and Technology of China, Chengdu Sichuan 611731, China;2. Chengdu Vocational and Technical College of Industry, Chengdu Sichuan 610218, China)

The diffeomorphism is a differential transformation with smooth and invertible properties, which leading to topology preservation between anatomic individuals while avoiding physically implausible phenomena during MRI image registration. In order to yield a more plausible diffeomorphism for spatial transformation, nonlinear structure of high-dimensional data was considered, and an MRI image registration using manifold learning based on adaptive tangent space was put forward. Firstly, Symmetric Positive Definite (SPD) covariance matrices were constructed by voxels from an MRI image, then to form a Lie group manifold. Secondly, tangent space on the Lie group was used to locally approximate nonlinear structure of the Lie group manifold. Thirdly, the local linear approximation was adaptively optimized by selecting appropriate neighborhoods for each sample voxel, therefore the linearization degree of tangent space was improved, the local nonlinearization structure of manifold was highly preserved, and the best optimal diffeomorphism could be obtained. Numerical comparative experiments were conducted on both synthetic data and clinical data. Experimental results show that compared with the existing algorithm, the proposed algorithm obtains a higher degree of topology preservation on a dense high-dimensional deformation field, and finally improves the registration accuracy.

diffeomorphism; tangent space; Lie group; neighborhood selection; MRI image registration

2016- 08- 15;

2016- 12- 23。 基金項(xiàng)目:廣東省省部級(jí)產(chǎn)學(xué)研聯(lián)合項(xiàng)目(2012A090300001)。

劉薇(1971—),女,四川成都人,副教授,博士研究生,主要研究方向:計(jì)算機(jī)視覺、醫(yī)學(xué)圖像處理; 陳雷霆(1966—),男,四川成都人,教授,博士生導(dǎo)師,博士,主要研究方向:圖像處理、計(jì)算機(jī)圖形、虛擬現(xiàn)實(shí)。

1001- 9081(2017)04- 1193- 05

10.11772/j.issn.1001- 9081.2017.04.1193

TP391.4

A