基于總變差正則化算法的三維肺部呼吸過程電阻抗成像方法

王 琦，陸紀(jì)璇，李秀艷，段曉杰

（1.天津工業(yè)大學(xué)生命科學(xué)學(xué)院，天津 300387；2.天津工業(yè)大學(xué)天津市光電檢測(cè)技術(shù)與系統(tǒng)重點(diǎn)實(shí)驗(yàn)室，天津300387；3.天津工業(yè)大學(xué)電子與信息工程學(xué)院，天津 300387）

電阻抗斷層成像技術(shù)（electrical impedance tomography，EIT）是一種通過電極給予人體微弱的電流激勵(lì)，以此測(cè)量人體表面電壓來計(jì)算人體電阻抗分布的技術(shù)，具有實(shí)時(shí)、安全、便攜的優(yōu)點(diǎn)，并對(duì)電導(dǎo)率變化的捕捉十分敏感。由于肺部充滿空氣，肺部內(nèi)區(qū)域與周圍組織之間的電導(dǎo)率差異較大，且肺部容量不斷變化，導(dǎo)致其電導(dǎo)率分布處于持續(xù)變化的狀態(tài)，因此，EIT 成像技術(shù)非常適用于肺部成像[1-2]，可用于肺功能的連續(xù)監(jiān)測(cè)。現(xiàn)已證明EIT 在肺部功能測(cè)試中有顯著效果，在保護(hù)肺部通氣過程中能夠進(jìn)行可視化實(shí)時(shí)反饋[3-4]，在治療慢性阻塞性肺病（chronic obstructive pulmonary disease，COPD）[5]、哮喘[6]、囊性纖維化[7]、氣胸[8]和肺水腫[9]等方面使用EIT 成像已經(jīng)取得很大進(jìn)展[10-12]，目前已有臨床EIT 產(chǎn)品問世。使用EIT 技術(shù)能夠?qū)Ψ尾考膊∵M(jìn)行安全有效的指導(dǎo)[13]，不僅對(duì)成人患者有效，對(duì)新生兒和幼兒患者均有益處。由于人體的三維特性，傳統(tǒng)的二維EIT 成像方法觀察范圍有限，不能很好地滿足診斷、監(jiān)測(cè)需要，因此，實(shí)現(xiàn)三維成像是EIT 技術(shù)發(fā)展的趨勢(shì)。

總變差（total variation，簡稱TV）正則化算法采用總變差函數(shù)作為約束項(xiàng)，可以有效克服Tikhonov 正則化算法導(dǎo)致的過度平滑問題，并具有運(yùn)行速度快、保邊性好等優(yōu)點(diǎn)。因此，TV 算法對(duì)于肺部三維圖像重建是一個(gè)很好的選擇。

本文利用EIT 重建技術(shù)對(duì)肺部呼吸運(yùn)動(dòng)進(jìn)行監(jiān)測(cè)。首先，基于CT 切片構(gòu)建人體三維胸腔模型，即呼氣末模型和吸氣末模型，進(jìn)行仿真研究；其次，采用TV 算法與Tikhonov 正則化算法對(duì)肺部仿真數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行三維圖像重建，并對(duì)2 種算法的圖像重建結(jié)果進(jìn)行對(duì)比與評(píng)價(jià)。

1 EIT 重建基礎(chǔ)

EIT 通常采用差分成像方法，將感應(yīng)域電導(dǎo)率相對(duì)于參考電導(dǎo)率圖像的變化與測(cè)量的邊界電壓差V的關(guān)系表示為：

式中：F 為非線性正演算子；e 為測(cè)量噪聲。在實(shí)際情況下，式（1）可線性化，表示為[14]：

式中：V 為歸一化電壓測(cè)量；J 為雅可比矩陣或靈敏度矩陣；σ 為物場(chǎng)與空?qǐng)鲋g的電導(dǎo)率變化向量。

EIT 技術(shù)包含正問題與逆問題。正問題求解過程為通過“原因”求“結(jié)果”，此過程一般是常態(tài)的、適定的，具有唯一且穩(wěn)定的解。逆問題求解過程即圖像重建過程，由邊界電極測(cè)得的邊界電壓聯(lián)合正問題中求得的靈敏度矩陣，反向推導(dǎo)被測(cè)區(qū)域內(nèi)的電導(dǎo)率分布。逆問題求解過程即為由“結(jié)果”推“原因”，因此，該過程是病態(tài)的、不適定的，它的解不是唯一存在且穩(wěn)定的，微小的變化都會(huì)引起解的極大變化。

由于電導(dǎo)率分布的變化對(duì)最終邊界電壓測(cè)量值的變化影響較小，導(dǎo)致逆問題的求解計(jì)算過程變得十分不穩(wěn)定。對(duì)于邊界電壓測(cè)量值，微小的改變便會(huì)對(duì)最終解產(chǎn)生巨大影響，導(dǎo)致計(jì)算求得的電導(dǎo)率分布與實(shí)際電導(dǎo)率分布可能會(huì)相差甚遠(yuǎn)[15]，這是EIT 逆問題求解過程中面臨的最大難點(diǎn)。

2 基于EIT 正問題的人體三維胸腔模型構(gòu)建

本文基于CT 掃描圖片，提取人體真實(shí)胸部輪廓與肺部輪廓，進(jìn)行組合之后形成模型。三維胸腔模型的構(gòu)建過程如圖1 所示，基于最大類間方差算法[16]，首先提取出胸腔和肺部輪廓并進(jìn)行疊加，得到模擬胸腔的初步模型；然后利用SolidWorks 軟件對(duì)三維模型進(jìn)行平滑處理，使正演問題得到穩(wěn)定求解；最后，在COMSOL Multiphysics 軟件中將肺部、胸廓模型進(jìn)行組合并安裝電極，得到最終實(shí)驗(yàn)?zāi)Ｐ汀?/p>

圖1 建立人體三維胸腔模型的步驟Fig.1 Steps of establishing 3D thoracic model of human body

為了研究肺部EIT 的呼吸運(yùn)動(dòng)，本文在吸氣末和呼氣末2 種特定狀態(tài)下建立了人體三維胸腔模型，對(duì)人體EIT 數(shù)據(jù)進(jìn)行重構(gòu)，如圖2 所示。根據(jù)人體電導(dǎo)率特性[17]，將肺組織的電導(dǎo)率設(shè)置為σexp=120 mS/m來模擬呼氣末端狀態(tài)，胸腔電導(dǎo)率設(shè)置為σbkg=480 mS/m，吸氣時(shí)會(huì)導(dǎo)致正常肺部區(qū)域的電導(dǎo)率變小，通過改變肺部電導(dǎo)率至σinsp=60 mS/m 來模擬吸氣末端狀態(tài)[17-18]。

圖2 2 種呼吸狀態(tài)下的3D 胸腔模型Fig.2 3D thoracic models under two breathing states

通過對(duì)模型進(jìn)行仿真實(shí)驗(yàn)，得到5 層電極測(cè)量值，對(duì)測(cè)量值進(jìn)行圖像重建，得到三維重建結(jié)果。由于內(nèi)存限制和速度要求，將三維重建圖像設(shè)為111 個(gè)重建層，即重建層間距為0.25 cm。建立三維逆問題網(wǎng)格，如圖3 所示，圖3 中所示的感知域縱向劃分為111層，每層分辨率為32×32。

圖3 用于逆問題計(jì)算的網(wǎng)格剖分Fig.3 Grid division for inverse problem calculation

根據(jù)電場(chǎng)的三維效應(yīng)，通過獲得單層電極平面的實(shí)測(cè)數(shù)據(jù)，可以得到重建平面圖像。首先用單層電極環(huán)模型測(cè)試重建圖像，如圖4 所示，該模型電極層由16 個(gè)電極組成，每個(gè)電極直徑為1 cm。采用相鄰激勵(lì)相鄰測(cè)量模式，在3 kHz 的激勵(lì)頻率下注入10 mA 電流進(jìn)行仿真，用COMSOL 軟件求解三維模型的正問題。在測(cè)量中加入隨機(jī)噪聲，以此模擬測(cè)量當(dāng)中可能出現(xiàn)的噪聲。采用TV 算法和Tikhonov 正則化算法對(duì)每個(gè)平面上的圖像進(jìn)行單獨(dú)重建。

圖4 一層EIT 電極的胸腔模型Fig.4 Thoracic models with a layer of EIT electrodes

本文選擇相鄰電極層間的20 個(gè)重構(gòu)層，即每5 cm一層進(jìn)行測(cè)量。因此，5 層電極環(huán)（每層電極環(huán)有16 個(gè)電極）可以覆蓋整個(gè)胸腔區(qū)域的測(cè)量，如圖5 所示，此種方式可以提高z 軸方向的軸向分辨率。

圖5 5 層電極的3D 胸腔模型Fig.5 3D thoracic models with 5 layers of electrodes

3 Tikhonov 正則化算法和總變差正則化算法

3.1 Tikhonov 正則化算法

Tikhonov 正則化方法是常用的正則化方法，它可以求解很多在測(cè)量中出現(xiàn)的不適定問題，如擬合推估問題、病態(tài)問題等，通過選擇合適的正則化參數(shù)，可以得到較好的結(jié)果。

Tikhonov 正則化方法的基本思想是將目標(biāo)函數(shù)最小化：

經(jīng)證明，STS+μI 的逆是存在的。因此，基于標(biāo)準(zhǔn)Tikhonov 正則化方法的解可以表示為：

正則化質(zhì)量主要取決于正則化因子μ 的大小。μ值較小，能夠?qū)υ紗栴}給出好的近似，但誤差會(huì)影響解，使其不穩(wěn)定；μ 值較大，會(huì)削弱解對(duì)誤差的敏感性，但估計(jì)解一般會(huì)偏離真實(shí)值。因此，在實(shí)際應(yīng)用中μ 值可以取中位數(shù)。理論上存在μ 的最優(yōu)值，但由于其計(jì)算量過大，對(duì)于實(shí)時(shí)性要求較高的場(chǎng)合，選取μ 值主要還是依靠實(shí)際經(jīng)驗(yàn)。

由于Tikhonov 正則化方法能夠得到逆問題真實(shí)解的近似穩(wěn)定解，并且重建速度較快，目前已被廣泛應(yīng)用。但由于此類方法采用2 范數(shù)形式的正則化函數(shù)（即（‖L（g-g）2‖））進(jìn)行求解，即通過對(duì)重建圖像進(jìn)行一定程度的平滑后使重建過程趨于穩(wěn)定，喪失了邊緣信息，因此，圖像的分辨能力不高[19]。

3.2 總變差正則化算法

由于Tikhonov 正則化算法得到的圖像分辨率較差，Rudin 等[20]于1992 年首次提出一種基于總變差（total variation）正則化的圖像重建算法，簡稱TV 算法。TV算法是一種應(yīng)用于圖像處理中的圖像除噪技術(shù)，由于該方法在圖像重建時(shí)具有較好的保邊緣性[21]，因此，已經(jīng)在圖像處理領(lǐng)域得到廣泛認(rèn)可與應(yīng)用。該算法能夠提高逆問題解的穩(wěn)定性，并能夠提高非連續(xù)分布介質(zhì)的區(qū)域成像分辨能力，具有良好的保邊緣性，并且成像速度較快，可以滿足EIT 系統(tǒng)對(duì)成像實(shí)時(shí)性的要求。

TV 算法在應(yīng)用于EIT 圖像重建時(shí)，采用歸一化的介電常數(shù)σ 的總變差作為正則化函數(shù)，即：

式中：（xσ）、（yσ）和（zσ）分別為三維模型在x、y、z方向上的一階局部差；指數(shù)i、j、k 為有限元網(wǎng)格中的單元。在傳統(tǒng)的TV 算法求解中，通常引入平滑參數(shù)β[22]來近似原始TV 算法，使其在原點(diǎn)處可微。

4 圖像重建

利用Tikhonov 正則化、TV 正則化方法對(duì)肺部EIT進(jìn)行二維與三維圖像重建工作，并根據(jù)重建結(jié)果對(duì)肺部呼吸運(yùn)動(dòng)進(jìn)行分析。

4.1 二維重建

將仿真數(shù)據(jù)導(dǎo)入對(duì)應(yīng)的Matlab 程序當(dāng)中運(yùn)行，通過不同算法對(duì)肺部進(jìn)行二維分層圖像重建，結(jié)果如圖6 所示，左圖部分為仿真模型和成像區(qū)域，右圖部分為二維重建圖像結(jié)果。

由于呼氣末與吸氣末的肺部狀態(tài)不同，吸氣末肺部狀態(tài)充盈，呼氣末肺部狀態(tài)癟小，因此，造成兩種狀態(tài)的有效肺部圖像層數(shù)不一致，呼氣末肺部狀態(tài)有效圖層為60 張，吸氣末肺部狀態(tài)有效圖層為72 張。由圖6 能夠看出，Tikhonov 正則化算法和TV 算法均能在二維空間上反映肺部區(qū)域的變化。

圖6 2 種算法重建二維肺部圖像結(jié)果Fig.6 Reconstruction of two-dimensional images of lung with two kinds of algorithms

4.2 三維重建

4.2.1 評(píng)價(jià)參數(shù)

從理論上講，重建圖像可以反映不同呼吸狀態(tài)下肺容積的變化。為有效評(píng)估2 種算法對(duì)重建圖像的影響，評(píng)估不同方法構(gòu)成的胸部模型與人體胸部真實(shí)結(jié)構(gòu)的差異，本文采用肺部通氣區(qū)域占胸部截面除肺部以外區(qū)域的比例即肺區(qū)域比（lung regional ratio，簡寫為LRR，用rLR表示）[23]及其誤差對(duì)成像結(jié)果進(jìn)行評(píng)價(jià)，其定義如式（10）所示：

式中：t（x，y）為像素點(diǎn)（x，y）的像素值；T 為肺部區(qū)域阻抗閾值，表示肺部區(qū)域最低阻抗值。那么SRV和SSV分別可由式（13）和式（14）表示：

本文采用最大類間方差算法[24]確定區(qū)分肺部區(qū)域和非肺部區(qū)域的閾值。設(shè)f（x，y）為重建圖像位置（x，y）處的灰度值，具有L 個(gè)灰度級(jí)（本文選取L 為20，即將歸一化的電導(dǎo)率范圍平均分為20 份），用變量b 表示灰度級(jí)。設(shè)灰度級(jí)bi的所有像素個(gè)數(shù)為fi，其中i=0，1，…，L-1，且b0＜b1＜…＜bL-1。重建圖像像素?cái)?shù)為N，則第bi級(jí)灰度出現(xiàn)的概率為：

式中：rLR1為基于仿真數(shù)據(jù)或?qū)嶒?yàn)數(shù)據(jù)的肺區(qū)域比；rLR0為基于真實(shí)數(shù)據(jù)的肺區(qū)域比。

4.2.2 基于仿真數(shù)據(jù)的三維圖像重建

分別采用Tikhonov 正則化算法和TV 算法對(duì)肺部圖像進(jìn)行三維重建，如圖7 所示。

圖7 基于模擬數(shù)據(jù)的肺部三維重建結(jié)果Fig.7 Three-dimensional lung reconstruction results based on simulated data

由圖7 可以看出，2 種算法生成的三維重建結(jié)果在不同位置均具有一致性。對(duì)于不同的肺部情況，三維重建圖像在肺葉位置和體積變化方面表現(xiàn)出的趨勢(shì)一致。

Tikhonov 正則化算法能夠生成一個(gè)平滑解，并且從圖像上看可能會(huì)由于過于平滑造成圖像邊界模糊，與TV 方法相比，Tikhonov 正則化方法提高了反演結(jié)果的清晰度與表面光滑度，但由于過于平滑，并不適用于對(duì)被測(cè)物進(jìn)行圖像重建。

根據(jù)式（10）計(jì)算2 種呼吸狀態(tài)下基于三維重建圖像仿真數(shù)據(jù)的LRR 值，如表1 所示。在CT 掃描圖像的基礎(chǔ)上建立仿真模型，得到真實(shí)的LRR。

表1 基于仿真數(shù)據(jù)的LRR 比較Tab.1 Comparison of LRR based on simulation data

由表1 數(shù)據(jù)可以得出，基于仿真數(shù)據(jù)的圖像重建結(jié)果中，與傳統(tǒng)Tikhonov 正則化算法相比，TV 算法的eLRR降低了58.86%，LRR 值準(zhǔn)確率提高了14.58%，更加接近真實(shí)值，不存在顯著性差異。由此說明，TV 算法在降低成像誤差方面有顯著優(yōu)勢(shì)，能夠使圖像重建結(jié)果更精確。

4.2.3 基于實(shí)驗(yàn)數(shù)據(jù)的三維圖像重建

通過提取人體CT 的胸腔以及肺部圖像，利用3D打印技術(shù)，構(gòu)建不同呼吸狀態(tài)下胸腔、肺部的實(shí)驗(yàn)?zāi)Ｐ停捎霉饷魳渲牧系玫酵该鞯男厍荒Ｐ秃桶咨姆尾磕Ｐ停鶕?jù)測(cè)算之后安裝5 層電極。兩種呼吸狀態(tài)下的胸腔與肺部模型如圖8 所示。

圖8 3D 打印胸腔與肺部模型Fig.8 Models of chest and lung by 3D printing

肺部模型內(nèi)部處于中空狀態(tài)，通過導(dǎo)柱向其中填充質(zhì)量濃度為9 g/L 的NaCl 溶液（即生理鹽水）進(jìn)行導(dǎo)電，以此模擬真實(shí)人體肺部。將激勵(lì)頻率設(shè)為50 kHz，激勵(lì)電流設(shè)置為5 mA，采用相鄰激勵(lì)相鄰測(cè)量的激勵(lì)測(cè)量方法，通過TJPU-EIT 測(cè)量系統(tǒng)內(nèi)的模塊進(jìn)行數(shù)據(jù)采集，得到相鄰電極電壓數(shù)據(jù)。測(cè)量過程如圖9 所示.

圖9 電壓數(shù)據(jù)的測(cè)量過程Fig.9 Measuring process of voltage data

TJPU-EIT 系統(tǒng)為每層電極環(huán)提供一整套測(cè)量值Vj，k（j，k=1，…，16）。由于電流驅(qū)動(dòng)電極上的電壓容易產(chǎn)生誤差，并受到其他未知接觸阻抗的影響，因此，只能使用測(cè)量數(shù)據(jù)進(jìn)行圖像重建工作。與仿真數(shù)據(jù)實(shí)驗(yàn)相似，5 層電極覆蓋了111 層胸腔切片。利用實(shí)驗(yàn)數(shù)據(jù)建立三維成像結(jié)果，如圖10 所示。

圖10 基于實(shí)驗(yàn)數(shù)據(jù)的三維肺部成像結(jié)果Fig.10 Three-dimensional imaging results of lung based on experimental data

根據(jù)三維圖像重建結(jié)果計(jì)算LRR 及其誤差，結(jié)果如表2 所示。

表2 實(shí)驗(yàn)數(shù)據(jù)的LRR 比較Tab.2 LRR comparison of experimental data

由表2 數(shù)據(jù)可以得出，基于實(shí)驗(yàn)數(shù)據(jù)的圖像重建結(jié)果中，與傳統(tǒng)Tikhonov 正則化算法相比，TV 算法的LRR 值準(zhǔn)確率提高了16.33%，誤差eLRR降低了53.81%更接近真實(shí)值，不存在顯著性差異，說明TV 算法成像更精確，能夠有效降低成像誤差，提高成像效果。

5 結(jié) 論

本文利用三維肺部EIT 對(duì)呼吸運(yùn)動(dòng)進(jìn)行研究，首先建立呼氣末與吸氣末狀態(tài)的三維胸部模型，并根據(jù)非電極所在層的圖像重建結(jié)果，討論電極層間距的設(shè)計(jì)；其次，介紹了用于EIT 圖像重建的Tikhonov 正則化算法和TV 算法，同時(shí)基于兩種算法，對(duì)模型進(jìn)行仿真成像與實(shí)驗(yàn)成像，并比較兩種算法在圖像重建結(jié)果上的優(yōu)劣，結(jié)果表明：

（1）仿真數(shù)據(jù)中，與Tikhonov 正則化算法相比，TV算法使三維肺部EIT 成像誤差降低了58.86%，肺區(qū)域比（LRR）準(zhǔn)確率提高了14.58%；

（2）實(shí)驗(yàn)數(shù)據(jù)中，與Tikhonov 正則化算法相比，TV算法使三維肺部EIT 成像誤差降低了53.81%，LRR準(zhǔn)確率提高了16.33%；

（3）三維EIT 在肺部呼吸運(yùn)動(dòng)成像方面具有巨大潛力，TV 算法可以有效提高三維圖像質(zhì)量，從而提高肺部呼吸運(yùn)動(dòng)參數(shù)監(jiān)測(cè)的準(zhǔn)確性。