基于局部空間變換網(wǎng)絡(luò)的醫(yī)學圖像配準

張 糾 劉曉芳* 楊 兵

1(中國計量大學計算機應用與技術(shù)研究所 浙江 杭州 310018) 2(中國計量大學電子信息與通信研究所 浙江 杭州 310018) 3(浙江省電磁波信息技術(shù)與計量檢測重點實驗室 浙江 杭州 310018)

0 引 言

基于計算機輔助技術(shù)的醫(yī)學影像分析為病灶識別與定位[1]、疾病診斷[2-3]等提供了可靠的信息基礎(chǔ)。其中，醫(yī)學圖像配準是指對浮動圖像(Moving image)進行空間變換，確定浮動圖像與固定圖像(Fixed image)的空間對應關(guān)系[4]，尋求最佳空間映射。

醫(yī)學圖像配準問題的核心在于如何確定圖像間的相似性以及如何對圖像進行空間映射[4]，在臨床應用中，通常將同一個人在不同時間段所獲得的圖像進行配準或?qū)⒍鄠€人配準到標準圖譜上。傳統(tǒng)基于灰度值的配準方法[5]將配準問題抽象為最小化相似度測量以尋找最優(yōu)參數(shù)的優(yōu)化問題，然而，此種方法存在效果較差，無法取得最優(yōu)解等局限性問題。近年來，由于神經(jīng)網(wǎng)絡(luò)具有能夠從圖像特征中學習復雜映射、自動優(yōu)化配準參數(shù)、建立優(yōu)化模型等優(yōu)點，提供了用于處理配準問題的新思路，已被廣泛用于醫(yī)學圖像配準中[6-8]。

基于編碼-解碼結(jié)構(gòu)的卷積神經(jīng)網(wǎng)絡(luò)如U-net[9]能夠自動對輸入圖像進行底層特征提取，其獨特的網(wǎng)絡(luò)結(jié)構(gòu)在處理醫(yī)學圖像解剖結(jié)構(gòu)相對固定、數(shù)據(jù)量較少的情況下具有較好的擬合性能。空間變換網(wǎng)絡(luò)(Spatial Transform Networks, STN)[10]是一種結(jié)合卷積運算、仿射變換、空間采樣的顯式特征空間變換網(wǎng)絡(luò)結(jié)構(gòu)，其常用于特征提取網(wǎng)絡(luò)之后，對配準場(Registration Field)進行空間變換得到配準圖像(Registered Image)，以確定圖像間的空間映射關(guān)系?；诳臻g變換網(wǎng)絡(luò)的編碼-解碼結(jié)構(gòu)是處理醫(yī)學圖像配準的常用網(wǎng)絡(luò)框架，然而，此類網(wǎng)絡(luò)具有沒有充分考慮到各通道特征重要性不同、基于單路的特征空間變換對局部圖像特征變換能力不足，以及特征與特征之間信息融合程度不高等問題，針對此問題，本文提出一種基于特征塊局部空間變換的雙路編碼-解碼網(wǎng)絡(luò)用于配準問題，本文方法在圖像特征塊的基礎(chǔ)上，進行了特征重要性加權(quán)和基于雙通道的局部特征空間變換。

1 相關(guān)工作

1.1 醫(yī)學圖像配準

醫(yī)學圖像配準是后續(xù)圖像融合，疾病診斷的基礎(chǔ)，其目標是根據(jù)浮動圖像與固定圖像之間的相似性關(guān)系，確定兩幅圖像間的最佳映射。根據(jù)空間變換的類型，醫(yī)學圖像配準方法可分為剛體配準與非剛體配準[11]。剛體配準方法具有變換參數(shù)較少、運算效率高等優(yōu)點，然而，大多數(shù)醫(yī)學圖像如心臟CT(Computed Tomography, CT)圖像都存在非剛體變換，因此對于基于非剛體變換的配準方法研究較多。

傳統(tǒng)基于灰度值的配準方法采用互相關(guān)系數(shù)(Cross Correlation)[12]等相似性測量指標，通過優(yōu)化參數(shù)模型確定最終映射關(guān)系，此外，互信息(Mutual Information)[13]及其改進的條件互信息(conditional Mutual Information)[14]、最大梯度距離(Maximum Distance-Gradient, MDG)[15]、熵圖像(Entropy Images)、流行學習(Manifold Learning)[16]等也常用于醫(yī)學圖像配準問題中。

神經(jīng)網(wǎng)絡(luò)本質(zhì)上也是一種參數(shù)優(yōu)化模型，其優(yōu)點是能夠通過最小化損失函數(shù)自動學習優(yōu)化參數(shù)，建立圖像間的配準模型。近年來，神經(jīng)網(wǎng)絡(luò)在處理配準問題上表現(xiàn)出較好的學習性能，如遞歸級聯(lián)網(wǎng)絡(luò)[17]通過網(wǎng)絡(luò)級聯(lián)的方式逐層對圖像特征進行空間變換，在肝臟配準問題上表現(xiàn)出較好的泛化性能，遞歸配準網(wǎng)絡(luò)[18]將遞歸神經(jīng)網(wǎng)絡(luò)引入圖像配準問題中，將圖像看作一系列序列模型。

醫(yī)學圖像配準常用的網(wǎng)絡(luò)框架如圖1所示，采用基礎(chǔ)特征提取網(wǎng)絡(luò)得到配準場，通過最小化損失函數(shù)得到浮動圖像與固定圖像間的空間映射參數(shù)。

圖1 醫(yī)學圖像配準框架

1.2 基于空間變換網(wǎng)絡(luò)的編碼-解碼配準模型

基于編碼-解碼結(jié)構(gòu)的網(wǎng)絡(luò)模型如U-net[9]、SegNet[19]等，利用跳躍連接(Skip Connection)結(jié)構(gòu)融合同級之間不同分辨率的特征，利用高分辨率特征具有輪廓、梯度、顏色等信息，低分辨率特征具有抽象化特征信息，將兩者的語義信息融合，使融合后的特征同時具有高層與底層特征的語義信息，能夠?qū)υ驾斎雸D像做出精細化特征表達。

圖2 編碼-解碼網(wǎng)絡(luò)結(jié)構(gòu)

圖2為編碼-解碼網(wǎng)絡(luò)基本結(jié)構(gòu)，其整體可分為三個組成部分：特征編碼，特征解碼，特征融合。一般來說，輸入圖像與輸出圖像具有相同或相近大小分辨率，輸入圖像經(jīng)過卷積、激活函數(shù)、下采樣等運算后得到編碼后的圖像特征，下采樣的次數(shù)一般為3到4次。過多的下采樣次數(shù)導致參數(shù)數(shù)量增加，計算速度減慢；過少的下采樣次數(shù)導致提取到的特征不足以表征圖像的深層細節(jié)。最大池化(max-pooling)是常見的下采樣運算，可看作正則化的一種方式。卷積運算是實現(xiàn)圖像特征提取的關(guān)鍵，其數(shù)學表達式如下：

(1)

式中:x為輸入圖像；k為卷積核；z為輸出圖像?；诰幋a-解碼結(jié)構(gòu)的卷積網(wǎng)絡(luò)通常采用跳躍連接的方式進行特征融合，其數(shù)學表達如下：

(2)

式中:X表示編碼層級的特征；Z表示解碼層級的特征；O表示通道拼接后的融合特征。

圖3 空間變換網(wǎng)絡(luò)基本結(jié)構(gòu)

圖3為空間變換網(wǎng)絡(luò)結(jié)構(gòu)，空間變換網(wǎng)絡(luò)是一種針對圖像特征變換的通用網(wǎng)絡(luò)模塊，能對輸入圖像進行仿射變換，由于其采用卷積、全連接等運算，使得空間變換網(wǎng)絡(luò)可以擬合復雜函數(shù)映射，輸出變換參數(shù)，得到空間變換后的圖像特征?？臻g變換網(wǎng)絡(luò)整體由定位器、網(wǎng)格生成器、采樣器三部分組成。

1) 定位器通過卷積等運算，輸出輸入特征與輸出特征之間的變換參數(shù)：

(3)

2) 網(wǎng)格生成器根據(jù)仿射變換矩陣以及輸出特征空間坐標，計算與之對應的輸入特征空間坐標。

3) 采樣器根據(jù)網(wǎng)格生成器得到的空間坐標點，利用雙線性插值對輸出特征進行填充：

(4)

2 基于特征塊空間變換的配準模型

2.1 網(wǎng)絡(luò)結(jié)構(gòu)

基于空間變換網(wǎng)絡(luò)的編碼-解碼網(wǎng)絡(luò)沒有考慮到不同通道的特征重要性不同，此外，空間變換網(wǎng)絡(luò)針對整個輸入特征，局部特征空間變換能力不足。針對此問題，本文提出了一種基于特征塊的雙通道空間變換網(wǎng)絡(luò)(Feature Patch Based Dual-channel Spatial Transform Network, FPDC-STN)。結(jié)合編碼-解碼網(wǎng)絡(luò)，本文的配準模型結(jié)構(gòu)如圖4所示。

圖4 基于特征塊的雙通道空間變換網(wǎng)絡(luò)配準模型

本文模型采用兩個編碼-解碼結(jié)構(gòu)的卷積神經(jīng)網(wǎng)絡(luò)分別對浮動圖像和固定圖像進行特征提取，同時將固定圖像的特征編碼部分與浮動圖像的特征解碼部分對齊，組成對稱網(wǎng)絡(luò)。為提高特征利用率，借鑒特征金字塔網(wǎng)絡(luò)[20]思想，從固定圖像和浮動圖像的編碼-解碼網(wǎng)絡(luò)中分別抽取了第2、第4、第6層特征，兩者特征融合之前用二維卷積進行了特征空間轉(zhuǎn)換。圖4中，F(xiàn)W-Block為特征加權(quán)(Feature Weighting)模塊，F(xiàn)PDC-STN為本文提出的特征塊雙通道變換模塊，兩者的網(wǎng)絡(luò)結(jié)構(gòu)如圖5、圖6所示。

特征加權(quán)模塊對原始輸入特征的通道加權(quán)可以分為4步：1) 對于輸入特征，采用兩個不同尺寸的分組卷積核得到兩組不同的輸出特征，并設(shè)置分組數(shù)K對兩組輸出特征進行分組，分別得到K個特征塊；2) 利用式(5)對特征塊進行全局平均池化(Global Average Pooling, GAP)，得到特征加權(quán)參數(shù)；3) 使用全連接和激活函數(shù)等對特征塊進行通道加權(quán)，如式(6)所示；4) 分別對K個特征塊進行合并，采用通道拼接的方式得到原始輸入特征通道加權(quán)后的輸出特征。

(5)

式中：H表示輸入特征塊的高；W表示輸入特征塊的寬；UC表示某個通道的特征塊；ZC表示某個通道的輸出特征塊。

O=σ(g(Z,w))=σ(w2δ(w1Z))

(6)

式中：Z表示輸出特征塊；w1、w2表示全連接參數(shù)；δ表示非線性激活函數(shù)如ReLU等；σ表示Sigmoid函數(shù)。

借鑒空間變換網(wǎng)絡(luò)結(jié)構(gòu)的思想，本文采用了基于特征塊空間變換的方式對原始輸入特征進行空間變換。此外，針對原始空間變換網(wǎng)絡(luò)感受野大小固定的局限性問題，本文還提出了一種基于雙通道的空間變換網(wǎng)絡(luò)，基于權(quán)值共享方式將特征塊空間變換與雙通道策略結(jié)合形成最終的基于特征塊的雙通道空間變換網(wǎng)絡(luò)(FPDC-STN)。FPDC-STN對輸入特征的空間變換可分為以下3步：1) 采用兩個不同尺寸大小的卷積核(3×3，5×5)對原始輸入特征進行分組卷積，然后對特征進行分塊；2) 對每個特征塊，分別進行空間定位-網(wǎng)格生成-采樣操作，得到空間變換后的輸出特征塊；3) 對輸出特征塊進行合并，基于通道拼接的方式得到最終空間變換后的輸出特征。

2.2 損失函數(shù)

本文所用的損失函數(shù)包括兩部分，分別為衡量浮動圖像與固定圖像之間的相似性指標Lsin(φ)，配準場的正則化項Lsmooth(φ)。兩者分別定義如下:

Lsim(F,R,φ)=-CC(F,R(φ))=

(7)

(8)

式中:CC(F,M(φ))表示局部相關(guān)系數(shù)(Correlation Coefficient)，衡量配準圖像R(φ)與固定圖像F之間的局部相關(guān)性；φ表示配準場(Registration Field);pi表示索引為i的像素?？倱p失函數(shù)定義為相似性指標Lsim(φ)與正則化項Lsmooth(φ)之和，如式(9)所示，λ為正則化參數(shù)。

L(F,M,φ)=Lsin(F,M,φ)+λLsmooth(φ)

(9)

2.3 整體算法流程

本文提出的基于局部空間變換網(wǎng)絡(luò)的醫(yī)學圖像配準算法整體流程如圖7所示。

圖7 本文整體算法流程

Step1輸入數(shù)據(jù)集劃分參數(shù)train_rate、validate_rate、test_rate，將數(shù)據(jù)集劃分為訓練集、驗證集、測試集。

Step2對訓練數(shù)據(jù)進行數(shù)據(jù)增強，包括圖像旋轉(zhuǎn)、圖像扭曲等。

Step3對待配準圖像進行預配準、去噪等預處理操作。

Step4初始化網(wǎng)絡(luò)參數(shù)、迭代次數(shù)Itr、正則化系數(shù)lamda等。

Step5更新網(wǎng)絡(luò)參數(shù)，以及分別建立特征塊雙通道變換模塊、特征加權(quán)模塊、特征編碼-解碼網(wǎng)絡(luò)。特征塊雙通道變換模塊(其中split_num表示劃分的特征塊數(shù)量)主要偽代碼表示如下(使用Pytorch框架編程)：

Apply 3×3 convolution and 5×5 convolution to input features;

Split input features into split_num;

# 空間定位

Location_module=nn.Squential(nn.Conv2d(1,8,3),

nn.MaxPool2d(2,2),

nn.Relu(True),

nn.Conv2d(8,10,5),

nn.MaxPool2d(2,2),

nn.ReLU(True))

#3*2 仿射變換回歸器

fc_loc=nn.Sequential(

nn.Linear(10 * 3 * 3, 32),

nn.ReLU(True),

nn.Linear(32, 3 * 2)

)

# 組合特征塊

Combine input features

特征加權(quán)模塊表示如下；

# 劃分特征塊

Split input features into split_num;

# 平均池化

avg_pool = nn.AdaptiveAvgPool2d(1)

# 特征加權(quán)以及特征賦值

fc=nn.Sequential(

nn.Linear(channel, channel // reduc tion, bias=False),

nn.ReLU(inplace=True),

nn.Linear(channel // reduction, channel, bias=False),

nn.Sigmoid()

)

# 組合特征塊

Combine input features

Step6網(wǎng)絡(luò)前向傳播獲得預測值。

Step7用式(9)計算損失函數(shù)，并計算配準評價指標。

Step8使用梯度下降算法，如式(10)所示，網(wǎng)絡(luò)反向傳播。

(10)

Step9迭代次數(shù)或損失函數(shù)是否滿足要求，若滿足要求，則執(zhí)行Step 10，若不滿足要求，則執(zhí)行Step 5。

Step10配準算法結(jié)束。

3 實 驗

實驗采用的硬件設(shè)備是顯存大小為12 GB的Nvidia GeForce GT1080ti GPU以及Intel i7 8700 CPU，軟件環(huán)境采用Python 3.6，機器學習處理包sklearn以及醫(yī)學圖像處理工具SimpleITK，配準實驗基于深度學習框架Pytorch[27]實現(xiàn)。初始學習率設(shè)置為0.000 01，并使用初始動量為0.9的Adam優(yōu)化器。迭代次數(shù)設(shè)置為2萬次，并采用隨機梯度下降(Stochastic Gradient Descent，SGD)更新網(wǎng)絡(luò)參數(shù)，正則化系數(shù)λ設(shè)置為1.0。

3.1 數(shù)據(jù)集與實驗參數(shù)

本文分別在腦部數(shù)據(jù)以及肝臟數(shù)據(jù)上進行了配準實驗。肝臟數(shù)據(jù)集來自MSD[21]以及SLIVER[22]，腦部數(shù)據(jù)集來自ADNI[23]、ABIDE[24]、LBPA[25]。所有的數(shù)據(jù)重采樣為132×132×132分辨率大小，并用醫(yī)學圖像處理工具FreeSurfer[26]進行了包括初步配準、腦部區(qū)域提取等標準預處理步驟。提取肝臟掃描數(shù)據(jù)、腦部掃描數(shù)據(jù)總計1 000個，其中800個用于訓練配準模型，100個用于驗證模型，100個用于測試模型。此外，為了增加實驗數(shù)據(jù)，所有數(shù)據(jù)都進行一系列數(shù)據(jù)增強操作，如二維旋轉(zhuǎn)、中心裁剪、數(shù)據(jù)扭曲等。實驗樣例圖像如圖8所示。

圖8 實驗樣例圖像

3.2 評價指標

本文采用Dice系數(shù)評價配準與固定圖像之間的重合程度：

(11)

式中:Dice系數(shù)范圍為[0,1]。當Dice系數(shù)為1時表示配準圖像與固定圖像完全匹配，當Dice系數(shù)為0時表示配準圖像與固定圖像完全不匹配。

3.3 腦圖像配準實驗

本文評估了不同實驗參數(shù)對配準效果的影響，評價標準為驗證集Dice系數(shù)、一段時間后評價指標是否變化、網(wǎng)絡(luò)收斂時間。評估參數(shù)包括正則化系數(shù)λ、迭代次數(shù)、Adam優(yōu)化器初始動量、學習率。結(jié)果如表1-表4所示。

表1 迭代次數(shù)對配準精度的影響

表2 正則化系數(shù)對配準精度的影響

表3 Adam初始動量對配準精度的影響

表4 學習率對配準精度的影響

根據(jù)表1可知，迭代次數(shù)為20 000次時，一段時間后驗證集評價不再變化，故本文設(shè)置的迭代次數(shù)為20 000次。根據(jù)表2可知，正則化系數(shù)λ=1時，驗證指標Dice系數(shù)最高，故本文設(shè)置的正則化系數(shù)為1。根據(jù)表3可知，Adam初始動量為0.90時，驗證集Dice系數(shù)指標最高，故本文設(shè)置的Adam初始動量為0.90。根據(jù)表4可知，學習率為0.000 01時，配準網(wǎng)絡(luò)已經(jīng)收斂，并且驗證集Dice系數(shù)指標最高，當學習率為0.000 1時，雖然配準網(wǎng)絡(luò)已經(jīng)收斂，但是其驗證集Dice系數(shù)低于學習率為0.000 01時Dice系數(shù)，故本文設(shè)置的學習率為0.000 01。

圖9為本文配準方法對于腦部圖像的配準結(jié)果，此外，為了說明本文配準模型的有效性與準確性，與UtilzReg方法[28]以及FAIM方法[29]進行了對比。

圖9 腦部圖像配準結(jié)果

本文選取了兩個切片的配準結(jié)果如圖9第一行以及第二行所示，圖9中第三行是第二行配準結(jié)果的局部放大圖。從固定圖像與浮動圖像中可以看出，固定圖像的外形輪廓與浮動圖像有較大差別，兩者內(nèi)部組織結(jié)構(gòu)也不相同，固定圖像內(nèi)部組織如腦白質(zhì)、腦灰質(zhì)等紋理較為簡單，而浮動圖像內(nèi)部組織的紋理比較復雜，這對配準算法的空間變換性能提出了較高要求。

圖9(c)為本文提出的配準方法的配準結(jié)果，本文方法可以針對浮動圖像與固定圖像之間的差異，對浮動圖像的外部輪廓以及內(nèi)部組織紋理進行準確的空間變換，這表明本文提出的基于特征塊的雙通道空間變換網(wǎng)絡(luò)FPDC-STN以及特征通道加權(quán)能夠捕捉圖像細節(jié)特征，對圖像局部做出準確配準。圖9(d)為FAIM配準方法的配準結(jié)果，可以看出，F(xiàn)AIM方法能夠根據(jù)固定圖像與浮動圖像之間的整體輪廓差異對浮動圖像做出較為準確的空間變換。然而，相比本文配準方法，F(xiàn)AIM方法針對圖像局部細節(jié)處理仍顯不足，配準結(jié)果中局部細節(jié)與浮動圖像差異不大，這說明FAIM方法不能準確地根據(jù)局部圖像細節(jié)對浮動圖像進行空間形變。圖9(e)為UtilzReg方法的配準結(jié)果，整體來說，此方法與FAIM方法差距不大，在處理圖像配準細節(jié)上較差。

表5為本文方法與FAIM方法以及UtilzReg方法的評價結(jié)果，本文進行了10次重復實驗，對配準結(jié)果計算Dice指標的平均值以及標準差?？梢钥闯觯疚姆椒ㄋ〉玫钠骄鵇ice指標相比FAIM方法提高2.7%，相比UtilzReg方法提高1.4%，在幾種方法中取得最高的Dice系數(shù)。從配準穩(wěn)定程度上看，本文的Dice指標的標準差最小，也是三種配準方法中最穩(wěn)定的算法。

表5 各方法腦部配準評價

3.4 肝臟圖像配準實驗

圖10為本文在肝臟圖像上的配準結(jié)果，同時也與FAIM算法以及UtilzReg算法比較了肝臟圖像配準的有效性與準確性，結(jié)果如表6所示。

圖10 肝臟圖像配準結(jié)果

表6 各方法肝臟配準評價

圖10第一行與第二行為選取的部分切片圖像配準結(jié)果，第三行為第二行配準結(jié)果的局部放大圖?？梢钥闯觯疚呐錅史椒軌蜉^準確地處理圖像細節(jié)，如圖10(c)所示，針對固定圖像中的局部亮點，本文配準方法能夠?qū)Ω訄D像進行局部空間變換，配準結(jié)果與固定圖像差異不大。對于FAIM以及UtilzReg方法，兩者能較好地對肝臟整體輪廓作出形變，但是在處理局部細節(jié)時仍然存在較多孤立像素點以及局部空間變換不足等問題。表6表明，與FAIM配準方法相比，本文方法平均Dice系數(shù)提高了8.5%，與UtilzReg方法相比，本文方法提高了3.6%。本文方法在保持最高平均Dice評教指標的同時，具有較好的穩(wěn)定性。

4 結(jié) 語

本文在現(xiàn)有基于空間變換網(wǎng)絡(luò)的編碼-解碼配準模型基礎(chǔ)上，提出了針對特征融合的特征加權(quán)模塊以及基于特征塊的雙通道空間變換網(wǎng)絡(luò)FPDC-STN；針對空間變換網(wǎng)絡(luò)無法有效解決圖像局部空間變換、感受野有限，以及空間變換參數(shù)無法共享等問題，本文提出的基于特征塊的雙通道空間變換網(wǎng)絡(luò)能夠針對拆分的特征塊進行準確局部變換；同時，利用雙通道策略以及兩種尺寸大小不同的分組卷積策略，共享空間變換網(wǎng)絡(luò)，有利于減少網(wǎng)絡(luò)參數(shù)量和提高網(wǎng)絡(luò)的空間變換能力。此外，本文提出的特征加權(quán)模塊以及特征塊變換模塊通過兩個大小不同的卷積核，獲取不同大小的感受野，本文提出的配準算法可以對自然圖像和醫(yī)學圖像進行特征變換和特征加權(quán)，作為具體應用實例，本文在醫(yī)學圖像上驗證了其配準有效性和準確性，進一步提升了圖像的空間變換性能，提高了圖像的配準精度。