改進(jìn)的差分搜索算法的醫(yī)學(xué)圖像配準(zhǔn)*

桂 鵬，邵黨國(guó)，祝曉紅，相 艷，王 碩，馬 磊

昆明理工大學(xué) 信息工程與自動(dòng)化學(xué)院，昆明 650504

1 引言

圖像配準(zhǔn)是通過(guò)一種或一系列的空間變換，使得兩幅或多幅圖像的特征點(diǎn)或?qū)?yīng)點(diǎn)在空間上達(dá)到一致，這種一致指的是匹配圖像的位置一致、角度一致、大小一致[1]。在醫(yī)學(xué)圖像配準(zhǔn)上可以理解為，人體上的同一解剖點(diǎn)在兩張匹配圖像上有相同的空間位置。醫(yī)學(xué)圖像配準(zhǔn)作為圖像處理的一項(xiàng)關(guān)鍵技術(shù)，在臨床診斷和手術(shù)治療等實(shí)際應(yīng)用中具有重要意義。Viergever等人比較了近20年來(lái)醫(yī)學(xué)圖像配準(zhǔn)方法的差異，并得出近20年醫(yī)學(xué)圖像配準(zhǔn)的主流算法并沒(méi)有特別大的改變，大多數(shù)經(jīng)典的方法現(xiàn)今仍舊適用[2]。當(dāng)前，醫(yī)學(xué)圖像配準(zhǔn)的兩類主要方法有基于灰度的方法和基于特征的方法[3]。其中，基于灰度的配準(zhǔn)算法因?yàn)椴恍枰M(jìn)行人工交互和圖像分割處理，可以直接利用圖像的灰度信息進(jìn)行配準(zhǔn)，使得該方法成為了醫(yī)學(xué)圖像配準(zhǔn)的常用方法之一。相似性測(cè)度和優(yōu)化算法作為基于灰度的圖像配準(zhǔn)的兩個(gè)關(guān)鍵問(wèn)題，一直是圖像配準(zhǔn)領(lǐng)域的研究熱點(diǎn)，近年來(lái)互信息理論更是被廣泛應(yīng)用于圖像配準(zhǔn)方法中[4]。交叉累計(jì)剩余熵（cross cumulative residual entropy，CCRE）就是眾多互信息中的一種，Wang等人通過(guò)研究發(fā)現(xiàn)，CCRE比傳統(tǒng)互信息（mutual information，MI）和歸一化互信息（normalized mutual information，NMI）曲線更光滑、魯棒性更強(qiáng)[5]。醫(yī)學(xué)圖像配準(zhǔn)的實(shí)質(zhì)是相似性測(cè)度的尋優(yōu)問(wèn)題，即尋找互信息最大時(shí)的幾個(gè)空間變換參數(shù)值?；诨バ畔⒌尼t(yī)學(xué)圖像配準(zhǔn)常用的優(yōu)化算法有Powell算法[6]、單純形法（simplex method，SM）[7]、遺傳算法（genetic algorithm，GA）[8]和粒子群算法（particle swarm optimization，PSO）[9]等。這些優(yōu)化算法各有優(yōu)點(diǎn)和不足之處：Powell法極易陷入局部極值，而單純形法收斂速度過(guò)慢。遺傳算法能夠跳出局部最優(yōu)，但計(jì)算時(shí)間較長(zhǎng)，且不時(shí)會(huì)搜索到錯(cuò)誤的方向而不能得到最優(yōu)解；PSO能夠有效地按照概率進(jìn)行全局搜索，從而可以克服局部極值并獲得全局最優(yōu)解，但迭代過(guò)程中容易出現(xiàn)早熟，對(duì)于離散的優(yōu)化問(wèn)題處理不佳，從而陷入局部最優(yōu)。

Civicioglu在觀察自然界生物群體從一個(gè)食物匱乏的棲息地持續(xù)不斷遷徙到其他給養(yǎng)充足的棲息地這一自然規(guī)律后，提出了一種新的元啟發(fā)式算法，稱為差分搜索算法（differential search algorithm，DSA）。DSA已經(jīng)證明了該算法在無(wú)約束優(yōu)化問(wèn)題上相對(duì)于粒子群算法、遺傳算法、人工蜂群算法（artificial bee colony algorithm，ABC）等的優(yōu)越性[10]，DSA也已在有約束全局優(yōu)化、無(wú)功規(guī)劃、電化學(xué)顯微加工等問(wèn)題上得到應(yīng)用[11-13]。但是，Civicioglu提出的DSA模型存在生物群搜尋的方向過(guò)于發(fā)散，從而導(dǎo)致算法不易沿著最佳路徑進(jìn)行優(yōu)化，且DSA在迭代后期不易找到峰值等局限性，將其直接應(yīng)用于醫(yī)學(xué)圖像配準(zhǔn)會(huì)導(dǎo)致配準(zhǔn)精度不理想等情況。針對(duì)DSA的缺點(diǎn)，文獻(xiàn)[14]將DSA模型應(yīng)用到了功率通量問(wèn)題上，調(diào)整了DSA模型使之適應(yīng)了文章目標(biāo)參數(shù)的求解，但并沒(méi)有對(duì)DSA進(jìn)行優(yōu)化。文獻(xiàn)[15]提出了一種CDSA（composite differential search algorithm）的解決方法，但該方法僅通過(guò)增加算法的復(fù)雜度來(lái)求取最優(yōu)候選解，并未對(duì)算法的尋優(yōu)策略進(jìn)行優(yōu)化。文獻(xiàn)[16]使用了一種BDSA（binary differential search algorithm）的方法用于解決多維度參數(shù)的背包問(wèn)題，但該方法由于其尋優(yōu)的特性并不適用于醫(yī)學(xué)圖像配準(zhǔn)中低維度的優(yōu)化問(wèn)題。對(duì)此，本文基于原始的DSA模型，調(diào)整了原模型中布朗運(yùn)動(dòng)的搜索模式和活動(dòng)范圍，并優(yōu)化了虛擬超級(jí)生物體尋優(yōu)的方向，更新了算法迭代的條件，基于原始的DSA模型提出了改進(jìn)的差分搜索算法（modified differential search algorithm，MDSA），且成功使用MDSA作為搜索策略進(jìn)行相似性測(cè)度尋優(yōu)。實(shí)驗(yàn)部分以CCRE為相似性測(cè)度，分別使用DSA和MDSA作為優(yōu)化算法，對(duì)兩組不同的醫(yī)學(xué)圖像進(jìn)行了優(yōu)化曲線對(duì)比、單模和多譜配準(zhǔn)實(shí)驗(yàn)。實(shí)驗(yàn)結(jié)果證明，MDSA相比DSA，不易陷入局部極值，迭代速度快，魯棒性強(qiáng)，能高效地完成醫(yī)學(xué)圖像配準(zhǔn)中的尋優(yōu)步驟。

2 基于交叉累計(jì)剩余熵的圖像配準(zhǔn)理論

在圖像配準(zhǔn)中，雖然兩幅圖像來(lái)源于不同的設(shè)備或不同的時(shí)間，但它們都是基于共同個(gè)體的同一位置的解剖信息，當(dāng)兩幅圖像的空間位置完全一致時(shí)，重疊部分所對(duì)應(yīng)像素的灰度互信息達(dá)到最大值。但是由于互信息對(duì)重疊區(qū)域的變化比較敏感，灰度互信息曲線往往不平滑，尋優(yōu)算法難以找到最大值甚至陷入局部極值，因此引入一種高效的相似性測(cè)度是非常有必要的。

Wang和Rao等人基于隨機(jī)變量的累計(jì)分布原理研究出一種與香農(nóng)熵（Shannon entropy）相類似的全新的相似性測(cè)度，稱之為累計(jì)剩余熵（cumulative residual entropy，CRE）。CRE 中的隨機(jī)變量 X，是由累計(jì)剩余分布函數(shù)P(x>u)中x所定義的，可以表示為：

累積剩余分布函數(shù)P(x>u)描述了X的值大于u的概率。采用Wang等人對(duì)CRE的定義[5]，計(jì)算CCRE在浮動(dòng)圖像F和參考圖像R之間灰度值的方法如下：

式中，ε(F|R)表示先驗(yàn)設(shè)定的CRE，E[·]為期望算子。根據(jù)Wang等人提供的函數(shù)定義，CCRE表示為累積和邊緣分布函數(shù)的公式。

式中，P(t>u,r=u)表示在圖像F的像素點(diǎn)中，比灰度值u更高的聯(lián)合概率；將這些像素點(diǎn)的坐標(biāo)置于圖像R中，則可得到各像素點(diǎn)的灰度值v。圖像F的累計(jì)剩余分布函數(shù)用PF(t>u)表示。

Wang等人[5]已經(jīng)證明了相對(duì)于互信息（MI），CCRE在醫(yī)學(xué)圖像配準(zhǔn)中展現(xiàn)出了更加優(yōu)異的魯棒性，CCRE作為相似性測(cè)度能夠得到更好的配準(zhǔn)質(zhì)量。CCRE之所以能夠優(yōu)于MI（基于香農(nóng)熵），其最主要原因是CCRE在處理參數(shù)化轉(zhuǎn)換領(lǐng)域時(shí)具有更好的噪聲穩(wěn)健性和更廣的收斂區(qū)間。

3 本文算法

3.1 差分搜索算法原理

差分搜索算法的創(chuàng)作靈感來(lái)源于自然界的各類生物體的遷徙過(guò)程。自然界中許多群體生物需要不斷地遷徙到更富饒的土地上以獲取更豐富的食物資源，當(dāng)大量個(gè)體同時(shí)進(jìn)行遷徙時(shí)，這些個(gè)體就形成一個(gè)超級(jí)生物體。Civicioglu在DSA中使用了布朗隨機(jī)運(yùn)動(dòng)的搜索方式來(lái)模擬超級(jí)生物體在大自然中遷徙的過(guò)程。在DSA中，算法產(chǎn)生的虛擬超級(jí)生物體在經(jīng)過(guò)隨機(jī)初始點(diǎn)的首次尋優(yōu)迭代后產(chǎn)生的位置被稱作棲息點(diǎn)，借此來(lái)模擬超級(jí)生物體的出發(fā)點(diǎn)。虛擬的超級(jí)生物體會(huì)不斷地迭代尋優(yōu)，每次尋優(yōu)找到的位置被稱為停留點(diǎn)，此時(shí)算法會(huì)判斷該位置的停留點(diǎn)是否優(yōu)于上一個(gè)位置，如果答案是肯定的，則把此時(shí)的停留點(diǎn)替換為最優(yōu)點(diǎn)。算法會(huì)不斷更新停留點(diǎn)的值，直到達(dá)到算法的迭代終止條件為止。算法流程如下：

（1）初始化DSA并設(shè)定初始參數(shù)，種群的規(guī)模popm,m={1,2,…,S},其中S表示個(gè)體總數(shù)，總維度popn,n={1,2,…,D},其中D表示所優(yōu)化問(wèn)題的維數(shù)。

（2）根據(jù)初始化參數(shù)求出棲息點(diǎn)的值。

此時(shí)的popmn即為棲息點(diǎn)的值，由式（4）得出的隨機(jī)值決定。

（3）遷徙過(guò)程中根據(jù)布朗隨機(jī)運(yùn)動(dòng)尋找停留點(diǎn)，能有效避免陷入局部最優(yōu)。

其中，Stopover即為停留點(diǎn)，表示虛擬超級(jí)生物體當(dāng)前位置；R代表虛擬個(gè)體搜索的范圍，用于模擬布朗隨機(jī)運(yùn)動(dòng)；map是一個(gè)由問(wèn)題維數(shù)構(gòu)成的由0和1組成選擇器，決定問(wèn)題維度是否更新，0代表該問(wèn)題方向在此輪迭代中不更新，1反之；Pos表示超級(jí)生物體目前棲息的位置；Dir-Pos表示超級(jí)生物體遷徙的方向。

（4）限制生物個(gè)體的邊界，Sitemn為虛擬種群的位置，當(dāng)Sitemnupn時(shí)，對(duì)Sitemn按照下式進(jìn)行隨機(jī)位置分配，保證個(gè)體的位置不超出搜索范圍。

（5）當(dāng)個(gè)體得到各自的Sitem后，要與當(dāng)前個(gè)體位置區(qū)域進(jìn)行比較，判斷所尋找的位置是否比當(dāng)前位置更優(yōu)。如果Sitem的位置更加優(yōu)越，則替換為該個(gè)體當(dāng)前位置。根據(jù)此時(shí)Sitem的值即可計(jì)算出CCRE的值。

（6）虛擬超級(jí)生物體不斷向全局最優(yōu)位置遷徙并持續(xù)更新CCRE的值，當(dāng)優(yōu)化算法滿足迭代終止條件或滿足迭代終止目標(biāo)時(shí)，結(jié)束優(yōu)化并輸出CCRE值最大時(shí)的幾個(gè)空間變換參數(shù)值；否則返回步驟（2）繼續(xù)探索。

DSA與人工蜂群算法[17]和蟻群算法[18]等仿生智能優(yōu)化算法不同，每個(gè)虛擬超級(jí)生物體中都包含了若干個(gè)虛擬生物體，算法會(huì)考慮每一個(gè)虛擬生物體的值，不完全依賴于當(dāng)前全局最優(yōu)解，很好地結(jié)合了搜索的目的性與隨機(jī)性。

3.2 改進(jìn)的差分搜索算法

如何平衡MDSA的全局搜索能力與局部搜索能力[19]是改進(jìn)算法的關(guān)鍵。DSA使用了較大的搜尋范圍結(jié)合單一迭代方向的迭代模型，這一點(diǎn)充分利用了算法的全局搜索性能，使得DSA難以陷入局部極值，但此迭代模型會(huì)在后期的局部搜索中變得乏力，難以找到峰值。而由于搜尋范圍過(guò)大，DSA會(huì)頻繁地使用到式（6）進(jìn)行邊界限制，而經(jīng)式（6）隨機(jī)產(chǎn)生的種群位置很難優(yōu)于現(xiàn)有的已更新過(guò)的停留點(diǎn)，這種尋優(yōu)模式極大地降低了算法的更新效率，與算法優(yōu)化的目標(biāo)相悖[20]。MDSA調(diào)整了式（5）中R×map的搜尋模型，改進(jìn)了這兩個(gè)特征的尋優(yōu)模式。相對(duì)于DSA對(duì)單一方向進(jìn)行大規(guī)模搜索的迭代模型，MDSA調(diào)整成為了對(duì)多個(gè)方向進(jìn)行小規(guī)模搜索的迭代模型。此模型的改進(jìn)使得MDSA整體偏向于積極的迭代更新，配合小規(guī)模的布朗運(yùn)動(dòng)搜索，規(guī)避掉了一部分算法因?yàn)殡S機(jī)運(yùn)動(dòng)帶來(lái)的不穩(wěn)定性，使得算法在迭代前期能夠穩(wěn)步更新，迭代后期更加容易找到峰值，更好地平衡了算法前期的廣度探索能力和后期的深度探索能力。DSA與MDSA在三維度優(yōu)化問(wèn)題上，同一種群中的9個(gè)粒子在迭代中的示意圖如圖1與圖2所示。

Fig.1 Particle trajectory by DSA圖1 DSA種群粒子運(yùn)動(dòng)軌跡

圖1和圖2中所有種群粒子均按照式（4）隨機(jī)分布初始位置，其迭代搜索方式按照式（5）中R·map·(Dir-Pos)進(jìn)行。圖1中，DSA的種群粒子尋優(yōu)范圍較長(zhǎng)且迭代方向單一，圖2中MDSA的種群粒子的尋優(yōu)范圍較短且約有1/3的種群粒子會(huì)進(jìn)行兩個(gè)方向的迭代（用虛線表示），多方向迭代的種群粒子會(huì)疊加這兩個(gè)方向的尋優(yōu)結(jié)果（用實(shí)線表示）。

在DSA的原始代碼中，式（5）中的R給出了5種預(yù)定的設(shè)置可供選擇，分別是R1=4×randn；R2=4×randg;R3=lognrnd(rand,5×rand);R4=1/gamrnd(1,0.5);R5=1/normrnd(0,5)，上述 R1～R5 的取值都存在波動(dòng)范圍較大且不穩(wěn)定的特點(diǎn)。MDSA迭代模型對(duì)R值的要求是數(shù)值更小且變化較小，考慮到布朗運(yùn)動(dòng)的特性，即隨機(jī)性與不穩(wěn)定性，故需要調(diào)整MDSA中的R參數(shù)，使其搜索范圍縮小且保持穩(wěn)定即可。在本文中，R被賦值為R=2×rand。

Fig.2 Particle trajectory by MDSA圖2 MDSA種群粒子運(yùn)動(dòng)軌跡

算法改進(jìn)的差分搜索算法（MDSA）

要求：N，種群數(shù)量，此處i={1,2,…,N}；D，問(wèn)題的維度；G，最大迭代次數(shù)。

如上述偽代碼所示，MDSA在激活虛擬個(gè)體的模式上與原始DSA完全不同[10]。首先，MDSA移除了DSA中人為設(shè)定的用于調(diào)節(jié)map激活狀態(tài)的p1、p2兩個(gè)控制參數(shù)，加強(qiáng)了算法的自動(dòng)化程度且提高了算法的更新效率；其次，MDSA中虛擬個(gè)體的激活方式也被重新調(diào)整了（偽代碼7～20行）。在DSA中，所有虛擬個(gè)體的活性，即式（5）中map是被隨機(jī)激活的，未被激活的虛擬個(gè)體在此次迭代中不進(jìn)行搜索。在MDSA中改進(jìn)了虛擬個(gè)體的激活選擇模式，圖3表示了在1 000次模擬下，3個(gè)維度的種群粒子激活比例，即在初始化時(shí) popm,m={1,2,…,S},S=1 000,popn,n={1,2,…,D},D=3的條件下，式（5）中map在3個(gè)方向上的虛擬個(gè)體是否激活。

如圖3所示，在DSA中激活單個(gè)虛擬個(gè)體（001,010,100）的比例非常大（90%左右），即DSA在大多數(shù)迭代過(guò)程中，僅會(huì)執(zhí)行某一個(gè)問(wèn)題維度的求解，導(dǎo)致其尋優(yōu)效率偏低。此外DSA還包括了000這種不激活的模式，不激活即代表在3個(gè)問(wèn)題維度上此輪迭代都不會(huì)進(jìn)行更新，即在此輪迭代中，算法并未進(jìn)行更新。MDSA則加大了多個(gè)體（110,101，011）的激活比例（33%左右）。因?yàn)镸DSA算法對(duì)偽代碼中7～20行代碼的改進(jìn)，使得MDSA中的多問(wèn)題維度能夠同時(shí)更新，使得MDSA能夠更快、更高效地更新停留點(diǎn)，且舍棄掉了000這種不更新的選擇，使得MDSA能夠有效地規(guī)避掉某些無(wú)意義的隨機(jī)搜索。

Fig.3 Comparison of DSAand MDSAof activation 1 000 operating map mode圖3 DSA與MDSA激活1 000次種群粒子（map）比較

Fig.4 Flow chart of medical image registration based on MDSA圖4 基于MDSA的醫(yī)學(xué)圖像配準(zhǔn)流程圖

3.3 基于改進(jìn)的差分搜索算法的醫(yī)學(xué)圖像配準(zhǔn)

以二維醫(yī)學(xué)圖像配準(zhǔn)為例，可以將圖像的X軸和Y軸以及旋轉(zhuǎn)的角度θ等變換參數(shù)看作解空間中的一個(gè)三維向量，并利用基于CCRE的不同互信息測(cè)度，將根據(jù)該變換參數(shù)計(jì)算出的互信息值作為目標(biāo)函數(shù)值，這樣便可以使用優(yōu)化算法在解空間中進(jìn)行搜索。以MDSA為優(yōu)化算法，CCRE為相似性測(cè)度的醫(yī)學(xué)圖像配準(zhǔn)流程如圖4所示。

上述流程圖中具體步驟如下所示：

（1）初始化MDSA并設(shè)定初始參數(shù)，種群的規(guī)模popm,m={1,2,…,S}，其中S表示個(gè)體總數(shù)，總維度popn,n={1,2,…,D}，其中D表示所優(yōu)化問(wèn)題的維數(shù)。在本文實(shí)驗(yàn)中，個(gè)體總數(shù)S=10,優(yōu)化問(wèn)題的維數(shù)D=3;個(gè)體位置的上限與下限值分別設(shè)定為upn=10,lown=-10，最大迭代次數(shù)為200次。

（2）根據(jù)初始化參數(shù)求出棲息點(diǎn)的值，由此棲息點(diǎn)可得到一個(gè)CCRE值，并將此CCRE值設(shè)為全局最大值。

（3）更新停留點(diǎn)的值，由此停留點(diǎn)可得到一個(gè)CCRE值，比較此時(shí)的CCRE值與全局最大值，若CCRE的值優(yōu)于全局最大值，則將此時(shí)的CCRE值替換為全局最大值。

（4）判斷此時(shí)算法是否滿足迭代終止條件或滿足迭代終止目標(biāo)，若為否則返回步驟（3）。

（5）輸出此時(shí)的全局最大值，即CCRE的最大值，此外還有此時(shí)的3個(gè)空間變換參數(shù)值。

4 實(shí)驗(yàn)與討論

4.1 單模圖像配準(zhǔn)

本文使用的單模配準(zhǔn)圖像數(shù)據(jù)集是美國(guó)Vanderbilt大學(xué)的“回顧性圖像配準(zhǔn)評(píng)估（Evaluation of Retrospective Image Registration）”項(xiàng)目里patient_001和patient_002兩位病人的腦部CT圖像（512×512×29 pixels）和T1、T2圖像（256×256×26 pixels）。所有的圖像都選擇第15片斷層作為實(shí)驗(yàn)圖像，如圖5所示。實(shí)驗(yàn)將經(jīng)過(guò)空間變換后的CT、T1、T2圖像當(dāng)作浮動(dòng)圖像，實(shí)驗(yàn)使用CCRE作為相似性測(cè)度，DSA和MDSA作為優(yōu)化算法進(jìn)行尋優(yōu)。

Fig.5 Single mode image registration圖5 單模圖像配準(zhǔn)

實(shí)驗(yàn)在X軸和Y軸以及旋轉(zhuǎn)的角度θ上，用在某一區(qū)域（6～7）隨機(jī)產(chǎn)生的空間變換參數(shù)來(lái)模擬浮動(dòng)圖像與參考圖像之間的空間變換。由于元啟發(fā)式算法的不穩(wěn)定性質(zhì)，對(duì)每一幅圖像DSA和MDSA都會(huì)進(jìn)行50次實(shí)驗(yàn)。計(jì)算出來(lái)的誤差值，即真實(shí)參數(shù)和配準(zhǔn)參數(shù)在X軸和Y軸以及旋轉(zhuǎn)角度θ上的誤 差值用ΔX、ΔY和Δθ來(lái)表示，此外本實(shí)驗(yàn)加入了一個(gè)綜合性能指標(biāo)實(shí)驗(yàn)將patient_001和patient_002的CT、T1、T2圖像配準(zhǔn)結(jié)果取平均值來(lái)當(dāng)作最后的實(shí)驗(yàn)結(jié)果。6幅實(shí)驗(yàn)圖像配準(zhǔn)誤差參數(shù)值如表1所示。

Table 1 Statistical analysis of DSAand MDSAin 6 single mode medical image registrations表1 6組單模醫(yī)學(xué)圖像配準(zhǔn)中DSA和MDSA的性能統(tǒng)計(jì)

表1中6組數(shù)據(jù)分別代表了圖5中（a）～（f）的圖像配準(zhǔn)結(jié)果。在上述所有實(shí)驗(yàn)中，MDSA的最差配準(zhǔn)值都是最小的，說(shuō)明MDSA相對(duì)于其他兩種算法更加穩(wěn)定，魯棒性更強(qiáng)。就平均誤差而言，MDSA在3個(gè)性能指標(biāo)上的表現(xiàn)（ΔX，ΔY，Δθ）大部分優(yōu)于DSA，就整體誤差Delta而言，在上述所有實(shí)驗(yàn)中，MDSA均為最小，此兩點(diǎn)說(shuō)明MDSA的配準(zhǔn)精度更高。就配準(zhǔn)完成所需時(shí)間而言，MDSA的配準(zhǔn)速度也更快。實(shí)驗(yàn)證明MDSA相對(duì)于DSA而言，魯棒性也更強(qiáng)，配準(zhǔn)精度更高，配準(zhǔn)耗時(shí)更少。

為了進(jìn)一步展現(xiàn)MDSA和DSA在優(yōu)化過(guò)程中的差異，下面引入了兩種算法在patient_002 T1圖像中，空間變換參數(shù)ΔX、ΔY和Δθ取值均為7時(shí)的迭代優(yōu)化曲線，如圖6所示。

圖6中，MDSA比DSA更早完成迭代，MDSA在17次迭代時(shí)就完成了尋優(yōu)，DSA在40次迭代時(shí)才完成尋優(yōu)；MDSA的空間變換參數(shù)ΔX、ΔY、Δθ最終配準(zhǔn)結(jié)果分別為7.94、6.23、7.08，相對(duì)于DSA的最終配準(zhǔn)結(jié)果，即8.18、5.94、6.82，MDSA的3項(xiàng)空間變換參數(shù)都更接近預(yù)設(shè)值7，配準(zhǔn)完成時(shí)MDSA相似性測(cè)度CCRE的值（67.55）也優(yōu)于DSA的相似性測(cè)度CCRE的值（65.99）。

Fig.6 Iterative optimization process of DSAand MDSA圖6 DSA與MDSA的配準(zhǔn)優(yōu)化曲線

4.2 多譜圖像配準(zhǔn)

本文使用的多譜配準(zhǔn)圖像數(shù)據(jù)集是Brain Web數(shù)據(jù)集里的T1、T2和PD圖像，如圖7所示。圖7中（a）～（c）的圖像的噪聲等級(jí)為0%，強(qiáng)度非均勻等級(jí)分別為0%，切片厚度均為1 mm，切面層數(shù)均為第90斷層。

實(shí)驗(yàn)將T1和T2圖像作為參考圖像，將PD圖像作為浮動(dòng)圖像，在X軸和Y軸以及旋轉(zhuǎn)角度θ上，用3組不同的空間變換參數(shù)（0～1，3～4，6～7）來(lái)模擬浮動(dòng)圖像與參考圖像之間的空間變換。DSA和MDSA的各種參數(shù)與單模實(shí)驗(yàn)配置相同，每一次配準(zhǔn)仍舊進(jìn)行50次實(shí)驗(yàn)后取均值，表2所示數(shù)據(jù)為6組配準(zhǔn)結(jié)果的平均值。實(shí)驗(yàn)設(shè)定了6項(xiàng)評(píng)價(jià)指標(biāo)，分別是交叉累計(jì)剩余熵（CCRE）、傳統(tǒng)互信息（MI）、歸一化互信息（NMI）、歸一化積相關(guān)（normalized cross correlation，NCC）、最小均方差（least mean-square error，LMSE）以及時(shí)間（Time），文中的時(shí)間單位統(tǒng)一以秒為單位。

Fig.7 Multi-spectral image registration圖7 多譜圖像配準(zhǔn)

表2列出了多譜醫(yī)學(xué)圖像配準(zhǔn)中DSA和MDSA的性能統(tǒng)計(jì)，大部分的實(shí)現(xiàn)結(jié)果表明MDSA是優(yōu)于DSA的。表3為表2所示所有數(shù)據(jù)取平均值的結(jié)果，從整體上看，MDSA在多譜圖像配準(zhǔn)實(shí)驗(yàn)中，所有指標(biāo)都優(yōu)于DSA。相對(duì)于DSA而言，MDSA在CCRE這項(xiàng)指標(biāo)上提高了近1個(gè)單位；在MI上提高了0.5個(gè)單位；在NMI這項(xiàng)指標(biāo)上，DSA劣于MDSA 0.06個(gè)單位；在NCC這項(xiàng)指標(biāo)上，DSA也略低于MDSA；LMSE數(shù)值越低性能越優(yōu)秀，相對(duì)于DSA，MDSA在此項(xiàng)的表現(xiàn)也更為優(yōu)異。從時(shí)間上看，MDSA整體迭代的時(shí)間也更快，從表3能夠看出MDSA比DSA平均配準(zhǔn)速度提高了近1 s，證明了雖然MDSA迭代更加頻繁，但是由于MDSA較DSA會(huì)更少地用到邊界限制，MDSA尋優(yōu)的速度相對(duì)于DSA反而變快了。上述結(jié)果說(shuō)明了MDSA在單模和多譜醫(yī)學(xué)圖像配準(zhǔn)中均展現(xiàn)出了更高的配準(zhǔn)精度、更快的配準(zhǔn)速度以及更好的魯棒性。

Table 2 Statistical analysis of DSAand MDSAin 6 multi-spectral medical image registrations表2 6組多譜醫(yī)學(xué)圖像配準(zhǔn)中DSA和MDSA的性能統(tǒng)計(jì)

Table 3 Average statistical performances between DSAand MDSAfor multi-spectral medical image registration表3 多譜圖像配準(zhǔn)整體DSA和MDSA的平均值統(tǒng)計(jì)

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

MDSA重新定義了停留點(diǎn)的更新方法，改進(jìn)了虛擬個(gè)體的激活方式，配合小規(guī)模的迭代更新，減少了邊界限制的使用頻率，從而提高了算法的尋優(yōu)效率，相對(duì)于DSA更好地平衡了算法的全局搜索和局部搜索性能。此外，相對(duì)于原始DSA，MDSA并沒(méi)有增加算法的復(fù)雜度或增添更多的特征值，且簡(jiǎn)化了控制參數(shù)。在實(shí)驗(yàn)部分，本文使用了MDSA和DSA對(duì)Vanderbilt大學(xué)以及BrainWeb數(shù)據(jù)集進(jìn)行了實(shí)驗(yàn)比較。結(jié)果表明，本文提出的算法在進(jìn)行二維單模剛體配準(zhǔn)時(shí)，配準(zhǔn)精度高，配準(zhǔn)速度快，不易陷入局部極值，魯棒性強(qiáng)；進(jìn)行多譜剛體配準(zhǔn)時(shí)，各項(xiàng)指標(biāo)均表現(xiàn)優(yōu)異，且配準(zhǔn)耗時(shí)較短。本次實(shí)驗(yàn)在一定程度上提高了醫(yī)學(xué)圖像配準(zhǔn)算法在臨床醫(yī)學(xué)中的實(shí)際應(yīng)用價(jià)值。