基于4D CT的肺部運(yùn)動呼吸模型

梁 睿，張建勛，馮昌利，趙汝琛

(1. 南開大學(xué)機(jī)器人與信息自動化研究所，天津 300071；2. 廣州雅敦微創(chuàng)科技有限公司，廣州 510110)

基于4D CT的肺部運(yùn)動呼吸模型

梁 睿1，張建勛1，馮昌利1，趙汝琛2

(1. 南開大學(xué)機(jī)器人與信息自動化研究所，天津 300071；2. 廣州雅敦微創(chuàng)科技有限公司，廣州 510110)

肺部的呼吸運(yùn)動會在圖像導(dǎo)航穿刺手術(shù)中引入誤差，并且該誤差是穿刺導(dǎo)航手術(shù)的主要誤差．為解決該問題，提出了一個接近亞體素級精度、方便易用、運(yùn)算簡便的呼吸模型．首先根據(jù)肺內(nèi)部各標(biāo)記點(diǎn)的真實運(yùn)動信息推算運(yùn)動狀態(tài)參數(shù)；然后根據(jù)運(yùn)動狀態(tài)參數(shù)的運(yùn)動趨勢，提取各點(diǎn)運(yùn)動的總體規(guī)律；最后結(jié)合已獲得的樣本數(shù)據(jù)，得到一個一般化的呼吸模型．交叉實驗結(jié)果表明，運(yùn)用該模型可以將由呼吸引起的原始誤差降低 75%左右，在呼吸幅度較小的情況下，模型可以達(dá)到亞體素級的精度，并且具有運(yùn)算簡便、運(yùn)用簡單的特點(diǎn)．

呼吸模型；曲線擬合；運(yùn)動狀態(tài)參數(shù)；誤差

隨著手術(shù)導(dǎo)航技術(shù)的發(fā)展成熟，手術(shù)導(dǎo)航定位越來越精準(zhǔn)．但是，與此同時，人們對手術(shù)導(dǎo)航的精度要求也越來越高．

手術(shù)導(dǎo)航系統(tǒng)的誤差主要來自于兩方面：一方面來源于圖像與病體組織的配準(zhǔn)誤差，主要原因是病體標(biāo)記點(diǎn)圖像掃描產(chǎn)生的偽影、位置傳感器的定位精度不高，目前已經(jīng)有很多非常成熟的偽影消除方法[1]、補(bǔ)償方法[2]等解決上述問題；另一方面來自于病體自然的生理運(yùn)動，尤其是人體的呼吸．呼吸是人體內(nèi)部器官最大的生理運(yùn)動，根據(jù)文獻(xiàn)記載，其運(yùn)動幅度一般有0.5～4.0,cm[3]，對于4.0,cm以下的病灶，呼吸的研究顯得尤為重要．呼吸一刻不能停止，它引起的肺部移動量也最大，連帶著其他臟器做較大運(yùn)動，所以研究呼吸時肺部產(chǎn)生的三維變化和位移，建立相應(yīng)的運(yùn)動模型預(yù)測其運(yùn)動趨勢，對手術(shù)導(dǎo)航系統(tǒng)精確度的改進(jìn)具有重大意義．

此外，CT掃描技術(shù)也日益發(fā)展，逐漸由 2D CT發(fā)展為3D CT，進(jìn)而出現(xiàn)了4D CT．4D CT能同時監(jiān)控呼吸信號，為人體呼吸模型的研究開辟了一個新途徑[4]．

以往呼吸模型的研究根據(jù)指導(dǎo)思想的不同分為4個類別：基于生物學(xué)的模型、基于信號的模型、基于圖像配準(zhǔn)的模型[5]和基于統(tǒng)計學(xué)的模型[6]．

基于生物學(xué)的模型是根據(jù)生理學(xué)和解剖學(xué)的理論建立力學(xué)模型，大都利用有限元分析的方法．但是該方法的假設(shè)過于理想化，無法真正詮釋肺的運(yùn)動，實驗結(jié)果的準(zhǔn)確度不高．

基于信號的模型，其外貼標(biāo)記點(diǎn)無法測知內(nèi)部運(yùn)動，而內(nèi)嵌標(biāo)記點(diǎn)由于個數(shù)有限，只能感知局部個別點(diǎn)．該方法雖然比生物學(xué)模型準(zhǔn)確，但是研究對象有限的特點(diǎn)限制了該模型的研究進(jìn)展．

基于圖像配準(zhǔn)的模型是近年來的研究熱點(diǎn)．但是，圖像配準(zhǔn)尤其是 3D圖像配準(zhǔn)耗時長，通常達(dá)到20～30,h，并且其配準(zhǔn)誤差不可知．利用帶有誤差的數(shù)據(jù)進(jìn)行研究，容易產(chǎn)生累計誤差．

基于統(tǒng)計學(xué)的模型是提取若干個樣本的呼吸運(yùn)動信息，再歸納出一個普遍的規(guī)律，因此，該方法普遍應(yīng)用性強(qiáng)，但沒有真正表達(dá)其個性特征．根據(jù)個性特征調(diào)整呼吸模型以適應(yīng)每一個人的呼吸運(yùn)動，是目前研究發(fā)展的方向．

本文提出了一種基于4D CT的肺部運(yùn)動呼吸模型，結(jié)合基于圖像配準(zhǔn)模型和基于統(tǒng)計學(xué)模型的思想，同時避免了基于圖像配準(zhǔn)模型由圖像配準(zhǔn)所帶來的缺點(diǎn)．首先對4個樣本的4D CT的圖像數(shù)據(jù)進(jìn)行處理，為了避免配準(zhǔn)帶來的誤差，利用醫(yī)學(xué)專家針對整個肺部特征做的標(biāo)記點(diǎn)進(jìn)行運(yùn)動信息分析；然后獲得每一個呼吸相位各標(biāo)記點(diǎn)的位置信息，歸納每一個呼吸相位的運(yùn)動特征；最后，綜合 4個樣本的肺部運(yùn)動特征，擬合出一個普遍的統(tǒng)計學(xué)模型．

1 呼吸模型

1.1 研究材料

本文所用第 3方數(shù)據(jù)以及信息來自文獻(xiàn)[7]．該論文及其網(wǎng)站上提供了研究圖像配準(zhǔn)算法時所用的實驗材料、實驗過程以及實驗結(jié)果．筆者利用其中 4個樣本的4D CT圖像以及醫(yī)學(xué)專家針對這4個樣本所做的特征標(biāo)記點(diǎn)，其中3個樣本有100個標(biāo)記點(diǎn)，1個樣本有 41個標(biāo)記點(diǎn)(標(biāo)記點(diǎn)位置分布如圖 1所示)，每個樣本在整個呼吸周期的每個相位均做了標(biāo)記點(diǎn)[8]．4個樣本均患有非小細(xì)胞肺癌，4D CT掃描采用螺旋模式，分辨率為 512×512×141，所有切片間距均為 2,mm，每張切片兩個方向的像素間距為0.976,6,mm．呼吸信號偵測利用 Lafayette儀器公司生產(chǎn)的氣動胸部壓力帶進(jìn)行．標(biāo)記過程采用了一個半自動化的軟件，該軟件具有良好的交互界面，能夠自動計算點(diǎn)與點(diǎn)的相似度，協(xié)助觀察者找到對應(yīng)點(diǎn)的位置．

圖1 41個標(biāo)記點(diǎn)空間分布投影Fig.1 Projection of 41 fiducials’ space distribution

1.2 研究方法

由于呼吸運(yùn)動的速度和幅度具有個體內(nèi)和個體間的差異性，病患的每個呼吸周期都沒有重復(fù)性，不同病患的呼吸運(yùn)動都不同，因此無法研究時間與呼吸運(yùn)動的關(guān)系，但是多項研究表明，呼吸周期的相位與呼吸運(yùn)動之間有一定的規(guī)律[9]．目前確定呼吸相位的方法有很多[10]，其中應(yīng)用比較普遍的有兩種：一種是根據(jù)呼吸通氣設(shè)備或者腹部壓力傳感器得到信號的振幅來對應(yīng)呼吸相位；另一種是根據(jù)這些設(shè)備得到信號的相位來對應(yīng)呼吸相位．針對不同呼吸相位對 CT圖像的肺部呼吸進(jìn)行建模．臨床應(yīng)用時，可以一方面利用相關(guān)設(shè)備監(jiān)測呼吸相位，另一方面運(yùn)用模型對圖像進(jìn)行調(diào)整，使其實時對應(yīng)該病患體內(nèi)的肺部運(yùn)動狀態(tài)．

該模型的建立基于第1.1節(jié)的CT圖像．每一個樣本的4D CT圖像包含10組，分別是呼吸相位0、10%、20%、30%、40%、50%、60%、70%、80%和90%．其中，按照4D CT的慣例標(biāo)準(zhǔn)，0是吸氣末端，50%或者 60%是呼氣末端．為了更好地找出肺部每個點(diǎn)的普遍運(yùn)動規(guī)律，根據(jù)肺部運(yùn)動規(guī)律，假設(shè)從 0到50%(或 60%)的變化是肺部運(yùn)動的最大幅度，并且肺部的每一個點(diǎn)都是從0到50%的位置上做往復(fù)運(yùn)動，那么該點(diǎn)某個時刻的位置可以描述為

式中：t為相位；k為標(biāo)記點(diǎn)；i為運(yùn)動狀態(tài)參數(shù)，表示呼吸運(yùn)動幅度的百分比，取值為 0～100；為 t相位下、i運(yùn)動幅度時，第 k個標(biāo)記點(diǎn)的理論估計坐標(biāo)值；為0相位下第k個標(biāo)記點(diǎn)的坐標(biāo)；為50%相位下第k個標(biāo)記點(diǎn)的坐標(biāo)；ti為t相位下各點(diǎn)的運(yùn)動幅度．

根據(jù)式(1)需要確定 t與 i的對應(yīng)關(guān)系．本文采用誤差最小化法找到這個對應(yīng)關(guān)系，即 t相位下，找到一個i值令誤差最小，且有

式中,tkP 為 t相位下第 k個標(biāo)記點(diǎn)的坐標(biāo)值．該誤差的標(biāo)準(zhǔn)差為

呼吸模型研究過程流程如圖2所示．

圖2 呼吸模型研究過程流程Fig.2 Flow chart of breathing model research process

1.3 呼吸模型的研究

根據(jù)上述研究過程，針對4個樣本的材料分別進(jìn)行研究．每一個相位都會得到一個誤差最小時的運(yùn)動狀態(tài)參數(shù)，稱為最優(yōu)運(yùn)動參數(shù)．得出 4個樣本的整個呼吸周期中的最優(yōu)運(yùn)動參數(shù)的散點(diǎn)圖，如圖 3所示.

圖3 4個樣本的最優(yōu)運(yùn)動參數(shù)散點(diǎn)圖Fig.3 Plot of optimal motion parameters for the four samples

2 多樣本一般化的呼吸模型

根據(jù)以上4個樣本的各相位最優(yōu)運(yùn)動參數(shù)，分別采用不同方法進(jìn)行曲線擬合，對擬合結(jié)果進(jìn)行比較，得出最優(yōu)的一般化呼吸運(yùn)動趨勢．選取幾種效果較好的結(jié)果進(jìn)行對比，如圖4所示．其擬合方法有傅里葉曲線擬合、高斯曲線擬合、多項式曲線擬合和正弦曲線擬合．

傅里葉曲線擬合式為

擬合結(jié)果為

高斯曲線擬合式為

擬合結(jié)果為

多項式曲線擬合式為

擬合結(jié)果為

正弦曲線擬合式為

擬合結(jié)果為

何為優(yōu)質(zhì)？韓光曙分析稱，醫(yī)院是救死扶傷的場所，必須尊重生命，對生命有敬畏感?！捌蜇づc王子平等”的濟(jì)世情懷在醫(yī)院傳承百余年?！暗魞H強(qiáng)調(diào)平等，機(jī)械地完成工作也可以。鼓樓醫(yī)院希望的是，在提供治病救人服務(wù)的同時，讓患者感到溫馨?！?/p>

圖4 4種方法的擬合結(jié)果Fig.4 Fitting results with four kinds of methods

運(yùn)用不同函數(shù)進(jìn)行曲線擬合，得出結(jié)果如表1所示，誤差單位為該點(diǎn)運(yùn)動幅度的百分比．其中SSE為和方差，RMSE為均方根，則有

SSE、RMSE越小，說明模型選擇和擬合越好．

R-square為確定系數(shù)，其計算式為

R-square越接近1，方程變量對y的解釋能力越強(qiáng)．

R-square-a為調(diào)整后的確定系數(shù)，它是嵌套模型的最優(yōu)參考指標(biāo)，其計算式為

表1 4種方法結(jié)果比較Tab.1 Comparison of fitting results with four kinds of methods %

經(jīng)過比較發(fā)現(xiàn)高斯曲線擬合和正弦曲線擬合產(chǎn)生了局部畸變，這種運(yùn)動狀態(tài)在肺的自然生理運(yùn)動中是不存在的，所以只在傅里葉曲線擬合和多項式曲線擬合的方法中進(jìn)行取舍．從表1中，可以看到傅里葉曲線擬合算法的精度較高，故選擇傅里葉擬合結(jié)果來進(jìn)行驗證．

3 實驗及模型誤差

由于4D CT輻射很大，國內(nèi)該類型儀器很少，而且實驗成本高昂，樣本取材很受限制，所以采取樣本間的交叉驗證方法對僅有的 4組樣本數(shù)據(jù)來進(jìn)行如下4組數(shù)值實驗，分別計算了其在模型作用下的各相位各點(diǎn)的平均誤差、標(biāo)準(zhǔn)差以及分解到x、y、z方向上的平均誤差．

該交叉驗證實驗結(jié)果如表2和表3所示．

表2 交叉驗證實驗結(jié)果1Tab.2 Result 1 of the cross validation test

表3 交叉驗證實驗1的xyz方向上的誤差Tab.3 Comparison of the errors in xyz directions in result 1

計算的整個呼吸運(yùn)動中各點(diǎn)平均運(yùn)動位移為7.670,8,mm，標(biāo)準(zhǔn)差為 5.054,7，所以，根據(jù)模型計算的結(jié)果誤差減小為無模型計算結(jié)果的 16.87%．從表3可以看出，該實驗結(jié)果誤差達(dá)到了亞體素級．

3.2 用第1、2、4個樣本擬合，用第3個樣本驗證

該交叉驗證的實驗結(jié)果如表4和表5所示．

表4 交叉驗證實驗結(jié)果2Tab.4 Result 2 of the cross validation test

表5 交叉驗證實驗2的xyz方向上的誤差Tab.5 Comparison of the errors in xyz directions in result 2

計算的整個呼吸運(yùn)動各點(diǎn)平均運(yùn)動位移為14.044,5,mm，標(biāo)準(zhǔn)差為 7.204,3，所以，根據(jù)模型計算的結(jié)果誤差減小為無模型計算結(jié)果的 26.47%．該組數(shù)據(jù)的呼吸運(yùn)動幅度較大，所以誤差較大．

3.3 用第1、3、4個樣本擬合，用第2個樣本驗證

該交叉驗證的實驗結(jié)果如表6和表7所示．

表6 交叉驗證結(jié)果3Tab.6 Result 3 of the cross validation test

表7 交叉驗證實驗3的xyz方向上的誤差Tab.7 Comparison of the errors in xyz directions in result 3

計算的整個呼吸運(yùn)動各點(diǎn)平均運(yùn)動位移為5.904,8,mm，標(biāo)準(zhǔn)差為 2.732,7，所以，根據(jù)模型計算的結(jié)果誤差減小為無模型計算結(jié)果的 22.87%．從表7可以看出，該組實驗的誤差也幾乎達(dá)到了亞體素級別．

3.4 用第2、3、4個樣本擬合，用第1個樣本驗證

該交叉驗證的實驗結(jié)果如表8和表9所示．

表8 交叉驗證結(jié)果4Tab.8 Result 4 of the cross validation test

表9 交叉驗證實驗4的xyz方向上的誤差Tab.9 Comparison of the errors in xyz directions in result 4

計算的整個呼吸運(yùn)動各點(diǎn)平均運(yùn)動位移為5.747,4,mm，標(biāo)準(zhǔn)差為 2.773,9，所以，根據(jù)模型計算的結(jié)果誤差減小為無模型計算結(jié)果的 25.81%．從表9中可以看出，該組實驗的誤差幾乎達(dá)到了亞體素級別．

綜合 4個實驗可以看出，運(yùn)用該呼吸模型進(jìn)行估計可以把呼吸運(yùn)動引起的誤差降低為 25%左右，而且在呼吸運(yùn)動較小的情況下可以達(dá)到亞體素別．

3.5 討 論

根據(jù)表10的無模型下各相位的平均誤差可以看出，本文模型雖然在整體的誤差消除作用比較顯著，達(dá)到了 25%左右，但在首相位(10%)和末相位(90%)時，誤差消除作用比較微弱，甚至偶有增長．這是因為，呼吸運(yùn)動在運(yùn)動開始和末端時候的位移相對于整個呼吸過程來講比較小，而且運(yùn)動增長速率也比較小，這時的呼吸運(yùn)動存在遲滯效應(yīng)，運(yùn)動模式也比較混沌[11]．該階段下肺臟各點(diǎn)運(yùn)動具有非線性不規(guī)則的特點(diǎn)，每個點(diǎn)運(yùn)動差異較大，這種現(xiàn)象的分析和建模非常復(fù)雜．而本文為了增強(qiáng)臨床的實用性和運(yùn)算的簡便性，對該現(xiàn)象的建模進(jìn)行了統(tǒng)一簡化，使這兩個相位下模型的運(yùn)動參數(shù)幾乎為 0，也就是弱化了此時模型的作用，從而使模型對于這兩個相位下的誤差消除作用比較?。怯捎谠撓辔幌赂鼽c(diǎn)的運(yùn)動幅度相對較小，在正常的呼吸幅度較小的情況下，呼吸誤差已經(jīng)達(dá)到了亞體素級別，在圖像的呈現(xiàn)上看不出具體差別．即使該相位下對呼吸誤差消除作用較弱，也保證了較小的誤差，使其在呼吸幅度較小的情況下也達(dá)到了亞體素級別．所以此相位下的誤差消除作用較弱并不影響整體的誤差消除效果．這樣既保證了文中模型的有效性，也保證了其實用性和簡便性．

表10 無模型下每個呼吸相位各點(diǎn)的平均誤差Tab.10 Average movement of each point in each phase without the model

4 結(jié) 語

為了處理穿刺手術(shù)中呼吸運(yùn)動引入的誤差，本文根據(jù)已有的肺部呼吸數(shù)據(jù)，找到普遍呼吸運(yùn)動規(guī)律，簡化為一個運(yùn)算簡捷、方便易用的呼吸模型．首先根據(jù)肺部標(biāo)記點(diǎn)運(yùn)動數(shù)據(jù)進(jìn)行參數(shù)尋優(yōu)，然后針對參數(shù)進(jìn)行多項式擬合，最后根據(jù)多個樣本得出一般化的呼吸模型．

該模型很容易運(yùn)用到穿刺手術(shù)場景，而且運(yùn)算量很小，不會影響實時性．在呼吸運(yùn)動較小的情況下，該模型的誤差可以達(dá)到亞體素級別．

該模型也有其局限性．該模型假設(shè)同一個呼吸相位時每個點(diǎn)的運(yùn)動參數(shù)相同，這只是比較一般的規(guī)律，如果能夠根據(jù)位置的不同具有不同的運(yùn)動參數(shù)變化規(guī)律，就會使精度更好，這需要進(jìn)一步進(jìn)行研究，同時需要增加樣本參數(shù)進(jìn)行驗證．

［1］ Bal M，Spies L. Metal artifact reduction in CT using tissue-class modeling and adaptive prefiltering[J]. Medical Physics，2006，33(8)：2852-2859.

［2］ Yu Hengyong，Zeng Kai，Bharkhada D K，et al. A segmentation-based method for metal artifact reduction[J]. Academic Radiology，2007，14(4)：495-504.

［3］ Li Xing，Thorndyke B，Schreibmann E，et al. Overview of image-guided radiation therapy[J]. Medical Dosimetry，2006，31(2)：91-112.

［4］ Keall P，Yamamoto T，Suh Y. Introduction to 4D Motion Modeling and 4D Radiotherapy[M]. Berlin Heidelberg：Springer Berlin Heidelberg，2013.

［5］ McClelland J R，Blackall J M，Tarte S，et al. A continuous 4D motion model from multiple respiratory cycles for use in lung radiotherapy[J]. Medical Physics，2006，33(9)：3348-3358.

［6］ George R，Vedam S S，Chung T D，et al. The application of the sinusoidal model to lung cancer patient respiratory motion[J]. Medical Physics，2005，32(9)：2850-2861.

［7］ Vandemeulebroucke J. Lung motion modelling and estimation for image guided radiation therapy[EB/OL]. http：//www. creatis. insa-lyon. fr/rio/popi-model，2012-12-20.

［8］ Murphy K，van Ginneken B，Klein S，et al. Semiautomatic construction of reference standards for evaluation of image registration[J]. Medical Image Analysis，2011，15(1)：71-84.

［9］ Yamamoto T，Kabus S，Lorenz C，et al. 4D CT lung ventilation images are affected by the 4D CT sorting method[J]. Medical Physics，2013，40(10)：101907.

［10］ Heinrich H P，Jenkinson M，Brady S，et al. MRF-based deformable registration and ventilation estimation of lung CT[J]. IEEE Transactions on Medical Imaging， 2013，32(7)：1239-1248.

［11］ Rahni A A A，Lewis E，Wells K. Characterisation of respiratory motion extracted from 4D MRI[C]//SPIE Medical Imaging International Society for Optics and Photonics. Florida，USA，2013：86692Z.

(責(zé)任編輯：金順愛)

Model of Lung Respiratory Movement Based on 4D CT

Liang Rui1，Zhang Jianxun1，F(xiàn)eng Changli1，Zhao Roger2
(1. Institute of Robotics & Automatic Information System，Nankai University，Tianjin 300071，China；2. Acton Medical Device Corporation of Guangzhou，Guangzhou 510110，China)

Lung respiratory movement can cause errors in the operation of image navigation surgery and they are the main errors in the navigation system. To solve this problem，a convenient，simple and easy-to-use breathing model whose precision was close to the sub-voxel was proposed. Motion parameters were first calculated according to the actual lung movement information of each point. The general movement law of each point was then obtained from the movement trend of the motion parameters. Finally a generic lung breathing model was constructed based on all the sample data. Results of cross experiments show that the proposed model can reduce the original errors caused by lung breath by about 75%. Besides，in circumstances of lung breath with small amplitude，the proposed model is simple and easy-to-use and can reach the sub-voxel precision.

breathing model；curve fitting；motion parameter；error

TP391.4

A

0493-2137(2015)04-0298-07

10.11784/tdxbz201311065

2013-11-22；

2014-03-10.

國家科技支撐計劃資助項目(2012BAI14B00).

梁 睿（1987— ），女，博士研究生，lrcalr@163.com.

張建勛，zhjx@robot.nankai.edu.cn.

時間：2014-03-13.

http：//www.cnki.net/kcms/doi/10.11784/z201311065.html.