姜?杉，張紅運，楊志永，張國彬

基于無監(jiān)督學習的三維肺部CT圖像配準方法研究

姜?杉，張紅運，楊志永，張國彬

(天津大學機械工程學院，天津 300350)

三維肺部電子計算機斷層掃描(computed tomography，CT)圖像非剛性配準是醫(yī)學圖像配準領(lǐng)域中最重要的任務(wù)之一．但是，肺部組織受呼吸運動影響而產(chǎn)生的非線性形變與大尺度位移給三維肺部CT圖像的非剛性配準帶來巨大挑戰(zhàn)．針對這一難題，設(shè)計開發(fā)了一種基于無監(jiān)督學習端到端的配準方法．通過改進現(xiàn)有U-Net神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)，在跳接之間引入Inception模塊，充分融合多尺度深層特征生成高精度的稠密位移向量場．為保證位移向量場光滑，在損失函數(shù)中加入雅可比正則化項，以達到訓練中顯式懲罰位移向量場中奇點的目的．另外，為緩解現(xiàn)有公開數(shù)據(jù)資源有限導致的過擬合問題，提出了一種基于三維薄板樣條(3D-thin plate spline，3D-TPS)變換的數(shù)據(jù)增強方法實現(xiàn)對訓練數(shù)據(jù)的擴充，將具有60套三維肺部CT圖像的訓練數(shù)據(jù)集EMPIRE10擴充為6060套以滿足卷積神經(jīng)網(wǎng)絡(luò)訓練的需要．設(shè)計驗證實驗，通過與基于學習的Voxelmorph方法和兩個包含傳統(tǒng)方法配準工具包ANTs和Elastix進行比較．實驗結(jié)果表明：在公開可用的DIR-Lab 4DCT數(shù)據(jù)集上，所提出的方法在目標配準誤差(target registration error，TRE)上達到次優(yōu)的2.09mm，平均Dice得分達到最優(yōu)的0.987，同時所生成的扭曲圖像中幾乎不存在折疊體素．

非剛性配準；卷積神經(jīng)網(wǎng)絡(luò)；無監(jiān)督學習；數(shù)據(jù)增強；雅可比正則化

三維肺部CT圖像已廣泛應(yīng)用于圖像引導下的放射治療[1]．肺作為典型的運動器官，在面向肺癌的放射治療之前，為估計所需的精確劑量必須追蹤肺呼吸運動[2]．術(shù)中需要將不同呼吸狀態(tài)的三維肺部CT圖像與參考狀態(tài)的三維肺部CT圖像進行對齊，跟蹤每個體素的劑量．三維肺部CT圖像非剛性配準是準確跟蹤肺部運動的一種可行方法．

非剛性配準是指在一對待配準圖像間建立非線性稠密映射關(guān)系．現(xiàn)有的傳統(tǒng)算法大多通過幾何方法解決每個圖像對的非剛性配準優(yōu)化問題，如SyN (standard symmetric normalization)[3]和Diffeomorphic Demons[4]．傳統(tǒng)方法計算開銷很大，每當配準一對新圖像時，需重新迭代運算．近年來，越來越多的學者傾向于利用基于學習方法完成配準任務(wù)．訓練后的卷積神經(jīng)網(wǎng)絡(luò)模型可以在幾秒鐘內(nèi)完成一對三維醫(yī)學圖像的配準，且其精度足以媲美較成熟的傳統(tǒng)方法．深度學習下的圖像配準大致可分為兩類：基于監(jiān)督學習的配準方法和基于無監(jiān)督學習的配準方法．

基于監(jiān)督學習的配準方法在訓練網(wǎng)絡(luò)時，需要提供與訓練樣本相對應(yīng)的標簽(ground truth，GT)．Fan等[5]使用傳統(tǒng)方法對圖像進行配準，將獲取到位移向量場作為GT，然而此方法獲得的GT與真實GT有所偏差．Hu等[6]提出用分割的解剖結(jié)構(gòu)掩膜作為GT來訓練卷積神經(jīng)網(wǎng)絡(luò)．在該方法中，卷積神經(jīng)網(wǎng)絡(luò)以固定和運動圖像對作為輸入，通過對解剖結(jié)構(gòu)掩膜進行對齊達到配準圖像的目的．Miao等[7]將原始圖像作為浮動圖像，將被模擬位移向量場扭曲的原始圖像作為固定圖像，將模擬位移向量場作為GT．盡管監(jiān)督學習在圖像配準領(lǐng)域有很大的潛力，但是通過傳統(tǒng)的配準工具獲取GT過程較為繁瑣．

基于無監(jiān)督學習的配準無需額外手工生成的GT，在訓練過程中利用固定圖像和扭曲圖像之前的差異來指導網(wǎng)絡(luò)參數(shù)的優(yōu)化．VoxelMorph為一種典型的基于無監(jiān)督學習的三維醫(yī)學圖像端到端配準方法，將配準中發(fā)生的扭曲變換定義為參數(shù)函數(shù)，并對其參數(shù)進行了優(yōu)化[8]．Zhao等[9]采用級聯(lián)配準子網(wǎng)絡(luò)預(yù)測位移向量場，取得了較好的效果．上述方法在腦、肝圖像上取得了良好的效果，而神經(jīng)網(wǎng)絡(luò)在肺部CT圖像配準中應(yīng)用較少．相比較于腦、肝圖像，不同時刻采集到的肺部圖像往往會由于人體的呼吸運動產(chǎn)生大范圍復(fù)雜非線性形變和位移，這將導致配準肺部圖像較為困難．

本文提出一種基于無監(jiān)督學習的卷積神經(jīng)網(wǎng)絡(luò)方法用于三維肺部CT圖像非剛性配準．在現(xiàn)有U-Net神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)上進行改進，將Inception模塊引入至跳接之間以提升網(wǎng)絡(luò)的特征檢測能力，充分融合多尺度深層特征生成高精度的稠密位移向量場．為了抑制不可逆變形，對稠密位移向量場中的具有負雅可比行列式值的變換進行懲罰．實驗表明，在變形后的圖像中幾乎沒有折疊體素．訓練后的卷積神經(jīng)網(wǎng)絡(luò)模型，可實現(xiàn)端到端的三維肺部CT圖像非剛性配準，同時具有較強的泛化能力．

1?方?法

圖1 基于無監(jiān)督學習的三維肺部CT圖像配準流程

Fig.1 Flow chart of 3D lung CT image registration based on unsupervised learning

1.1?卷積網(wǎng)絡(luò)結(jié)構(gòu)

如圖2所示，使用的卷積網(wǎng)絡(luò)是基于U-Net的改進．通常情況下，為生成高精度的稠密位移向量場，可在U-Net上堆疊卷積層加深網(wǎng)絡(luò)來實現(xiàn)．但是隨著網(wǎng)絡(luò)深度的增加，可能會出現(xiàn)梯度消失現(xiàn)象，進而導致模型過擬合．另外，不同個體的肺部CT圖像中信息位置存在較大差異，因此選取合適尺寸的卷積核變得十分困難．針對上述問題，本文不再通過盲目加深U-Net來優(yōu)化網(wǎng)絡(luò)性能，而是通過引入Inception模塊加寬網(wǎng)絡(luò)，在保證U-Net網(wǎng)絡(luò)一定深度的基礎(chǔ)上，通過Inception模塊所特有的多尺度卷積操作，促進網(wǎng)絡(luò)捕捉并學習更加豐富的多尺度特征．配合U-Net中跳躍連接所融合的上、下文信息，以達到生成高精度稠密位移向量場的目的．

圖2?3D卷積網(wǎng)絡(luò)模型結(jié)構(gòu)

1.2?損失函數(shù)

損失函數(shù)整體可以表示為

2?實?驗

2.1?數(shù)據(jù)集

本文使用3個包含同一病人不同時刻的三維肺部圖像數(shù)據(jù)集：EMPIRE10[12]、DIR-Lab 4DCT[13]和POPI[14]．DIR-Lab 4DCT數(shù)據(jù)集包含10個四維肺部CT圖像，每個四維肺部CT圖像中包含一個完整呼吸周期中的10個呼吸相位．選擇具有最大相對變形的兩個階段數(shù)據(jù)：呼氣末時刻和吸氣末時刻，兩階段分別標注有300個地標點．EMPIRE10數(shù)據(jù)集由30對三維肺部CT組成．POPI數(shù)據(jù)集由6個四維肺部CT圖像組成，使用最大相對變形的兩個階段：呼氣末時刻和吸氣末時刻，兩階段分別標注有100個地標點.

2.2?數(shù)據(jù)增強

3D-TPS插值函數(shù)的參數(shù)能夠通過求解下面線性方程組得出．

(9)

圖3?3D-TPS生成圖像示例

2.3?評價指標

由于CT圖像中肺實質(zhì)區(qū)域往往會占據(jù)大量的體素，因此使用地標點間的TRE來評估肺內(nèi)部的配準質(zhì)量．

Dice雖然不能很好地評價肺實質(zhì)區(qū)域的配準質(zhì)量，但可用于評估肺邊界的對齊程度．Dice分數(shù)表達式為

此外，圖像折疊在解剖學上是不合理的．第1.2節(jié)中提到的折疊體素的數(shù)量也作為算法的評價指標之一．

2.4?對比實驗

首先，將所提出的算法與ANTs[15]配準工具包中提供的SyN算法進行了比較．實驗表明，使用SyN默認參數(shù)對肺部CT圖像進行配準，達不到最佳配準性能．通過在多個數(shù)據(jù)上開展大量實驗，可得到更優(yōu)參數(shù)設(shè)置：SyN步長為0.25，高斯參數(shù)為(3，0)，在4種尺度下，每一種最多進行219次迭代．同時，將本文算法與基于Elastix[16]工具包中提供的B樣條(B-Spline)配準方法進行了比較．使用歸一化交叉相關(guān)作為相似性度量，進行了5個分辨率的1000次迭代．另外，還與基于學習的VoxelMorph算法進行了比較．VoxelMorph提出了兩種以U-Net為網(wǎng)絡(luò)框架的Vm1和Vm2網(wǎng)絡(luò)結(jié)構(gòu)，本文分別在Vm1和Vm2上進行了對比實驗．最初的VoxelMorph算法被應(yīng)用于腦組織配準，為保證公平對比，在EMPIRE10數(shù)據(jù)集上重新對VoxelMorph進行了訓練．最后，為驗證Inception模塊的有效性，筆者將本文方法與未添加Inception模塊的U-Net進行了對比實驗．

2.5?實現(xiàn)細節(jié)

3?實驗結(jié)果

通過使用測試數(shù)據(jù)集中包含的3000對地標點計算TRE．TRE結(jié)果如表1所示，括號中為標準偏差.本文提出的方法取得的平均配準誤差為2.09mm，標準偏差為1.55mm，取得了次優(yōu)結(jié)果．

表1?不同算法在DIR-Lab4DCT數(shù)據(jù)集上的TRE

Tab.1?TRE of different algorithms on the DIR-Lab 4DCT dataset

注：括號內(nèi)為標準偏差，無量綱．

如圖5所示，通過計算測試數(shù)據(jù)集上固定圖像掩膜與扭曲圖像掩膜間的平均Dice得分，相比于VoxelMorph、SyN和B-Spline，本文所提出的方法可實現(xiàn)最優(yōu)的配準性能(平均Dice＝0.987)．

能夠進行主觀性評價的實驗結(jié)果如圖6所示.稠密位移向量場以RGB圖像表示，3個通道對應(yīng)位移場的3個維度，越亮位置表示該位置發(fā)生的變形越大．圖中也給出被稠密位移向量場扭曲的形變網(wǎng)格圖像．值得注意的是，本文提出的方法與VoxelMorph使用的采樣方法相同，與ANTs和Elastix的采樣方式不同．

在稠密位移向量場的雅可比行列式圖像中，比0大的位置映射為藍色，而小于等于0的位置映射為紅色，即紅色標記了出現(xiàn)奇點的位置．在雅可比圖像行列式圖像中，本文方法并未顯現(xiàn)出折疊體素．在RGB圖像中和網(wǎng)格圖像中，能夠看出本文方法獲得了更加豐富的變形．

表2總結(jié)了所有的配準實驗結(jié)果，并給出了所有方法在GPU和CPU上的運行時間對比．到目前為止，還沒有針對ANTs和Elastix的GPU實現(xiàn)．結(jié)果表明，本文所提方法的配準速度相比傳統(tǒng)方法(SyN和BSpline)有著指數(shù)級提升，與VoxelMorph相近．值得注意的是，本文方法在獲得高配準精度的同時幾乎不產(chǎn)生折疊體素．另外，表2中的差異和是絕對差異圖像中所有體素值的總和，配準效果越好則差異和越?。畬嶒灲Y(jié)果表明本文方法在差異圖像的定量評價上取得了最優(yōu)．

表2?各算法在DIR-Lab4DCT數(shù)據(jù)集上的實驗結(jié)果

Tab.2?Experimental results of different algorithms on the DIR-Lab 4DCT dataset

注：括號內(nèi)為各變量的標準偏差，無量綱．

4?討?論

本文提出了一種基于無監(jiān)督的三維肺部CT圖像配準方法．利用浮動圖像和固定圖像對之間的NCC來優(yōu)化卷積神經(jīng)網(wǎng)絡(luò)參數(shù)，而不需要額外GT.在損失函數(shù)中引入雅可比正則項，有效地減少了扭曲圖像的折疊體素數(shù)量．采用基于3D-TPS變換的數(shù)據(jù)增強方法人工生成大規(guī)模訓練數(shù)據(jù)，滿足了卷積神經(jīng)網(wǎng)絡(luò)模型的訓練需求．實驗結(jié)果表明，該方法取得TRE值為2.09mm，平均Dice得分為0.987，且扭曲圖像中僅包含極少量的折疊體素．測試集上的實驗結(jié)果也表明該方法在三維肺部CT圖像配準中具有良好的魯棒性．

5?結(jié)?語

本文提出了一種基于無監(jiān)督學習的卷積神經(jīng)網(wǎng)絡(luò)方法實現(xiàn)三維肺部CT圖像的非剛性配準．基于3D-TPS變換的數(shù)據(jù)增強能夠解決數(shù)據(jù)資源有限導致的過擬合問題．在跳接之間引入Inception模塊，充分融合多尺度深層特征生成高精度的位移向量場，以獲得更佳配準效果．實驗結(jié)果表明，相較于傳統(tǒng)方法，本文算法可實現(xiàn)可靠的配準性能，且配準效率得到了指數(shù)級提升．

［1］Kaus M R，Brock K K，Pekar V，et al. Assessment of a model-based deformable image registration approach for radiation therapy planning[J]. International Journal of Radiation Oncology Biology Physics，2007，68(2)：572-580.

［2］Yu Z H，Lin S H，Balter P，et al. A comparison of tumor motion characteristics between early stage and locally advanced stage lung cancers[J]. Radiotherapy & Oncology，2012，104(1)：33-38.

［3］Avants B B，Epstein C L，Grossman M，et al. Symmetric diffeomorphic image registration with cross-correlation：Evaluating automated labeling of elderly and neurodegenerative brain[J]. Medical Image Analysis，2008，12(1)：26-41.

［4］Thirion J P. Image matching as diffusion process：An analogy with Maxwell’s demons[J]. Medical Image Analysis，1998，2(3)：243-260.

［5］Fan J，Cao X，Yap P T，et al. BIRNet：Brain image registration using dual-supervised fully convolutional networks[J]. Medical Image Analysis，2018，54(1)：193-206.

［6］Hu Yipeng ，Marc M，Eli G，et al. Label-driven weakly-supervised learning for multimodal deformable image registration[EB/OL]. https://arxiv.org/abs/1711. 01666v2，2017-11-05.

［7］Miao S，Wang Z J，Liao R. A CNN regression approach for real-time 2D/3D registration[J]. IEEE Transactions on Medical Imaging，2016，35(5)：1352-1363.

［8］Balakrishnan G，Zhao A，Sabuncu M R，et al. VoxelMorph：A learning framework for deformable medical image registration[J]. IEEE Transactions on Medical Imaging，2019，38(8)：1788-1800.

［9］Zhao S，Dong Y，Chang E，et al. Recursive cascaded networks for unsupervised medical image registration[C]// IEEE/CVF International Conference on Computer Vision. Seoul，Korea，2019：10599-10609.

［10］ Max J，Karen S，Andrew Z，et al. Spatial transformer networks[EB/OL]. https：//arxiv.org/abs/1506.02025，2016-02-04.

［11］ Balakrishnan G，Zhao A，Sabuncu M R，et al. An unsupervised learning model for deformable medical image registration[C]// IEEE/CVF Conference on Computer Vision and Pattern Recognition. Salt Lake City，USA，2018：9252-9260.

［12］ Murphy K，van Ginneken B，Reinhardt J M，et al. Evaluation of registration methods on thoracic CT：The EMPIRE10 challenge[J]. IEEE Transactions on Medical Imaging，2011，30(11)：1901-1920.

［13］ Castillo E，Castillo R，Martinez J，et al. Four-dimensional deformable image registration using trajectory modeling[J]. Physics in Medicine and Biology，2010，55(1)：305-327.

［14］ Vandemeulebroucke J，Rit S，Kybic J，et al. Spatiotemporal motion estimation for respiratory-correlated imaging of the lungs：Spatiotemporal motion estimation for 4D CT[J]. Medical Physics，2010，38(1)：166-178.

［15］ Avants B B，Tustison N J，Song G，et al. A reproducible evaluation of ANTs similarity metric performance in brain image registration[J]. NeuroImage，2011，54(3)：2033-2044.

［16］ Klein S，Staring M，Murphy K，et al. Elastix：A toolbox for intensity-based medical image registration[J]. IEEE Transactions on Medical Imaging，2010，29(1)：196-205.

Research on a 3D Lung Computed Tomography Image Registration Method Based on Unsupervised Learning

Jiang Shan，Zhang Hongyun，Yang Zhiyong，Zhang Guobin

(School of Mechanical Engineering，Tianjin University，Tianjin 300350，China)

Deformable registration of 3D lung CT images is crucial in medical image registration. However，nonlinear deformation and large-scale displacement of lung tissues caused by respiratory motion pose great challenges in the deformable registration of 3D lung CT images. Thus，we present a fast end-to-end registration method based on unsupervised learning. We optimized the classic U-Net model and added Inception modules between skip connections. The Inception module aims to capture and merge information at different spatial scales for generating a high-precision dense displacement vector field. To ensure a smooth displacement vector field，we introduced the Jacobian regularization term into the loss function to directly penalize the singularity of the displacement field during training. The existing publicly available datasets cannot implement model training. To address over-fitting caused by limited data resources and to expand the training data，we proposed a data augmentation method based on a 3D thin plate spline transform. Moreover，6060 CT scans will be generated based on the EMPIRE10 dataset，which contains 60 original CT scans to meet the requirement of convolution neural network training. Regarding the DIR-Lab 4DCT dataset，we achieved a target registration error of 2.09mm，an optimal Dice score of 0.987，and almost no folding voxels in comparison with the experimental results obtained using the deep learning method Voxelmorph and registration packages，such as advanced normalization tools (ANTs) and Elastix.

deformable registration；convolution neural network；unsupervised learning；data augmentation；Jacobian regularization

TP242

A

0493-2137(2022)03-0247-08

10.11784/tdxbz202010040

2020-10-21；

2020-12-10.

姜?杉（1973—??），女，博士，教授.

姜?杉，shanjmri@tju.edu.cn.

國家自然科學基金資助項目(51775368，81871457，51811530310)；天津市科技資助項目(18YFZCSY01300)；天津市津南區(qū)科技計劃資助項目(20200110)．

Supported by the National Natural Science Foundation of China(No.51775368，No.81871457，No.51811530310)，Tianjin Science and Tech-nology Project(No.18YFZCSY01300)，Tianjin Jinnan District Science and Technology Planning Project(No. 20200110).

(責任編輯：王曉燕)