基于共軛梯度算法的EIT仿真

侯雪，王超 ，劉準(zhǔn)

（1.天津大學(xué) 自動化學(xué)院，天津 300072；2.天津輕工職業(yè)技術(shù)學(xué)院，天津 300350；3.中儲糧油脂有限公司，天津 300350）

電阻抗斷層成像（EIT）作為一種獨特的醫(yī)療成像技術(shù)，根據(jù)不同生物組織特性、相關(guān)物理參數(shù)及其變化的分布，反映目標(biāo)的生理、病理狀態(tài)［1］。通過配置于目標(biāo)斷層外圍的電極陣列向斷面施加安全電流，對目標(biāo)域進行無創(chuàng)激勵，在體外測量電極響應(yīng)電信號以重建人體內(nèi)部的電阻抗特性分布變化的圖像［2］。是一種理想并具有廣泛應(yīng)用前景的醫(yī)學(xué)成像技術(shù)。

由生物體模型內(nèi)部的阻抗分布及邊界激勵信號，求物體內(nèi)部或表面的電壓和電流分布，是EIT正問題［3］。正問題求解是利用給定的邊界條件和阻抗分布模型求解拉普拉斯方程：

σ為目標(biāo)內(nèi)部的電導(dǎo)率分布；φ0為給定邊界上的點位，JP是給定邊界上外加激勵的電流密度，無注入電流時為零。

1 有限元剖分

由于有限元法不受分析區(qū)域邊界形狀的限制，因此EIT中的正問題通常采用有限元方法進行求解，使用最為普遍。Comsol軟件與傳統(tǒng)編程方法相比計算耗時短，可以與Matlab對接，完成最終的圖像重構(gòu)，計算結(jié)論還可以通過繪圖表現(xiàn)出來，隨著剖分單元數(shù)目增加，近似解收斂于精確解。但是，模型不能無限細(xì)分，剖分越細(xì)，計算周期越長。仿真實驗采用相鄰電極激勵、相鄰電極測量的方式，根據(jù)成像對象的位置和大小不同分為三種情況，得到三種情況下圓形場域剖分疏密程度不同時，測量電極電壓的變化曲線，為逆問題求解和實際測量提供依據(jù)。三種情況下測量電極電壓的變化曲線如圖1所示。

圖1 三種情況下測量電極電壓的變化曲線Fig.1 Electrode voltage curve of three Situation

從仿真結(jié)果可以看出：在三種情況下，一次細(xì)分和二次細(xì)分的結(jié)果相差并不明顯，所以在實際正問題計算中，可以采用一次細(xì)分來剖分網(wǎng)絡(luò)以節(jié)省求解時間。

2 圓形場共軛梯度重建算法研究

由表面測量電壓、電流分布及邊界激勵信號，求生物體模型內(nèi)部的阻抗分布，稱之為EIT逆問題［4］。逆問題的求解也稱為阻抗圖像的重構(gòu)。阻抗圖像的重構(gòu)是EIT技術(shù)的核心也是關(guān)鍵技術(shù)所在。本論文采用共軛梯度算法進行圖像重建。

共軛梯度法最早用于解正定系數(shù)矩陣的線性方程組，克服了最速下降法收斂慢的缺點，又避免了牛頓法需要存儲和計算Hesse矩陣并求逆的缺點，共軛梯度法是一個典型的共軛方向法，它的每一個搜索方向是互相共軛的，而這些搜索方向d僅僅是負(fù)梯度方向與上一次迭代的搜索方向的組合［5］。

共軛梯度算法適用于求解系數(shù)矩陣為對稱正定矩陣的方程，對于非對稱正定矩陣，可以對原方程進行轉(zhuǎn)換，即轉(zhuǎn)換成ATAx=ATb，通過正則化處理，上式轉(zhuǎn)化為(ATA+λI)x=ATb，此時，系數(shù)矩陣為對稱正定矩陣，可以利用共軛梯度算法進行求解。

本實驗利用共軛梯度算法成像，采用100次迭代，j代表迭代次數(shù)，第一次迭代結(jié)果為A1［N］，第二次迭代結(jié)果為 A2［N］，第 j次迭代結(jié)果為Aj［N］，N為剖分單元數(shù)。隨著迭代次數(shù)的增加，目標(biāo)函數(shù)將收斂于極值點。這將體現(xiàn)在相鄰的迭代結(jié)果差范圍變小，為定量的表示這一變化趨勢，定義

得到三種情況下100步迭代過程中的ε變化趨勢曲線如圖2所示。

ε然后利用作為停止迭代條件進行圖像重建，圖像重建結(jié)果如圖3所示。

圖2 ε變化趨勢曲線Fig.2 The trend curve of ε

圖3 圓形場的重建圖像Fig.3 Reconstruction image of the round field

從仿真結(jié)果可知：ε的變化趨勢在開始迭代10步左右的過程中，ε變化較大，20步之后，變化基本平穩(wěn)，根據(jù)ε變化趨勢迭代20步左右的ε大概在10e-10數(shù)量級，ε變化整體趨勢是遞減的，與理論推導(dǎo)一致。結(jié)果表明ε在10e-011數(shù)量級時（迭代次數(shù)在50～70次之間）重建圖像效果較好。

采用ε=10e-011作為迭代停止條件，進行圓形場重建實驗，針對單個成像對象進行分辨能力的評估。r為成像對象半徑相對于場域半徑的比值。單個成像對象圖像重建結(jié)果如圖4所示。

當(dāng)成像對象在中心位置，r＜0.4時，圖像重建結(jié)果并不理想，這說明場域中心的靈敏度較低。

采用ε=10e-011作為迭代停止條件，進行了圓形場重建實驗，針對兩個成像對象進行分辨能力的評估。d為兩個成像對象的中心距相對于場域半徑的比值。兩個成像對象圖像重建結(jié)果如圖5所示。

當(dāng)兩個成像對象d＜0.4時，圖像重建結(jié)果并不能對兩個成像對象進行區(qū)分，兩個成像對象的間距越大，成像效果相對越好。

3 基于胸部模型的仿真研究

圓形場域難以較好地模擬人體胸部的形狀特性。為了使研究結(jié)果更接近真實情況，有必要建立逼近真實研究對象的模型。根據(jù)人體胸部的CT掃描圖，構(gòu)建胸部先驗?zāi)Ｐ?，心臟電導(dǎo)率設(shè)為0.67S/m，脊椎電導(dǎo)率設(shè)為0.006S/m，皮下組織電導(dǎo)率設(shè)為0.037S/m，為了體現(xiàn)肺部的呼吸過程，根據(jù)生理數(shù)據(jù)，肺部電導(dǎo)率在0.042～0.138 S/m范圍內(nèi)變化［6］。

圖4 單個成像對象重建圖像Fig.4 Reconstruction image of single image object

圖5 兩個成像對象的重建圖像Fig.5 Reconstruction image of two imaging objects

圖6 模型1Fig.6 Model 1

圖7 模型2Fig.7 Model 2

圖像重建時，利用共軛梯度算法，分別采用圓形場模型、胸部廓形模型（模型1）如圖6所示和胸部先驗信息模型（模型2）如圖7所示。進行圖像重建。圖像重建結(jié)果如圖8所示。

圖8 三種模型的圖像重建Fig.8 Reconstruction image based on the three models

5 結(jié)論

仿真及實驗結(jié)果表明，采用基于肺部先驗知識的共軛梯度算法所重建的圖像質(zhì)量明顯優(yōu)于基于圓形場域的重建圖像，所重構(gòu)的圖像分辨率高，且胸部組織結(jié)構(gòu)輪廓更為清晰.應(yīng)用胸部廓形模型能夠?qū)崿F(xiàn)在呼吸過程中對人體胸腔進行監(jiān)測和實時成像的目的。有可能為醫(yī)學(xué)的臨床應(yīng)用和研究提供有價值的診斷信息。

［1］唐化勇，牟永閣，楊常運.生物電阻抗成像技術(shù)在醫(yī)學(xué)中的應(yīng)用［J］.第四軍醫(yī)大學(xué)吉林軍醫(yī)學(xué)院學(xué)報，2000，22（3）：174-176.

［2］王超.醫(yī)學(xué)電阻抗斷層成像的研究［D］.天津：天津大學(xué)，2001.

［3］任超世，王妍，鄧娟，電阻抗斷層成像應(yīng)用技術(shù)研究，中國醫(yī)療器械雜志，2007，31（4）：235-238.

［4］侯文生，楊東，彭承琳.阻抗斷層成像中的圖像重建技術(shù)［J］.生物醫(yī)學(xué)工程學(xué)雜志，2000，17（2）：214-217.

［5］董秀珍，秦明新，湯孟興.電阻抗斷層成像系統(tǒng)及重構(gòu)算法.第四軍醫(yī)大學(xué)學(xué)報，1999，20（3）：21.

［6］侯衛(wèi)東，莫玉龍.電阻抗成象中求解反問題的新方法［J］.計算機工程，2001，27（5）：71-7.