基于改進的差分進化算法的電阻抗圖像重建

劉海，閆丹丹，荊會敏

中國計量大學(xué)信息工程學(xué)院，浙江杭州 310018；

生物電阻抗斷層成像（electrical impedance tomography，EIT）是繼形態(tài)、結(jié)構(gòu)成像后的新一代無創(chuàng)功能性成像技術(shù)。它通過向體表注入較小的安全交變電流，由體表電極陣列測量得到體外邊界電壓（電勢），并采用重建算法重建組織內(nèi)部電阻率分布[1]。與正常的生物阻抗相比，生物組織在發(fā)生病變時，其阻抗值變化十分顯著[2]。EIT技術(shù)的實質(zhì)是利用反映組織生理及功能狀態(tài)的電阻抗信息進行阻抗重建。因此，EIT技術(shù)可用于腫瘤、癲癇等多種疾病的早期發(fā)現(xiàn)與診斷；并且EIT技術(shù)不產(chǎn)生電離輻射、對人體無害、相對成本較低、可多次重復(fù)使用、便攜，能對患者進行長期、實時動態(tài)監(jiān)護?；谝陨咸攸c，盡管目前應(yīng)用EIT技術(shù)與主流醫(yī)學(xué)成像技術(shù)在圖像分辨率上還有一定差距，但其在醫(yī)學(xué)工程及相關(guān)領(lǐng)域仍引起極大的關(guān)注[3-4]。

目前主要的重建算法有修正的牛頓-拉夫遜法、牛頓單步誤差重構(gòu)法、層析法、雙限定法等，上述方法收斂速度較快，但計算量大且過程復(fù)雜。近年來基于粒子群算法、神經(jīng)網(wǎng)絡(luò)算法和遺傳算法等相關(guān)的智能尋優(yōu)算法逐漸應(yīng)用于EIT圖像重建，這些智能算法與傳統(tǒng)算法相比，計算量較小、較易修正局部最優(yōu)值，具有較快的收斂速度和全局搜索能力等優(yōu)勢[5]。

本文采用改進的差分進化（differential evolution，DE）算法對頭部電阻抗分布進行重建，在三維頭球模型上進行仿真實驗，并與其他進化類算法進行對比。

1 EIT的基本原理

在EIT研究中，電阻抗分布與電流之間的關(guān)系由麥克斯韋方程給出。EIT測量中采用較低（10～100 kHz）的激勵源頻率，可忽略介電常數(shù)的影響，此時電流場可當(dāng)作穩(wěn)態(tài)電流場。為簡化計算，忽略電極與皮膚間的接觸阻抗，電壓與電流間的非線性關(guān)系可由滿足相應(yīng)邊界條件的拉普拉斯方程確定，其求解場域的數(shù)學(xué)模型[6-7]見公式（1），應(yīng)滿足的第一類或第二類邊界條件（、）見公式（2）。

其中，表示散度，表示梯度，ρ為目標(biāo)內(nèi)部電阻率分布，φ為場域中電位分布函數(shù)，φ0為邊界電位，n為方向向量，為邊界電位分布的法向?qū)?shù)，Jn為電流密度，無電流注入時為0。

EIT正問題是指由人體內(nèi)部阻抗分布及邊界條件，求解體表或內(nèi)部電壓分布，這是求解逆問題的基礎(chǔ)。EIT逆問題是指由體表電壓分布及邊界條件，求解人體內(nèi)部阻抗分布。經(jīng)過有限元法離散后，EIT逆問題可轉(zhuǎn)化為公式（3）中求電位分布殘差的相對均方誤差最小值[RMS(ρ)]問題。

其中，上標(biāo)的T表示矩陣的轉(zhuǎn)置，ρ為目標(biāo)內(nèi)部電阻率分布，φ(ρ)是邊界電位的計算值，φ0是邊界電位“測量值”。

因此，EIT逆問題是找到合適的電阻抗分布ρ使目標(biāo)函數(shù)RMS(ρ)達到最小，即搜尋一個m維參數(shù)矢量使RMS(ρ)最小，其中m為有限元剖分單元數(shù)，即電阻率分布自由度。從而，EIT逆問題即轉(zhuǎn)化為多參數(shù)函數(shù)優(yōu)化問題。

2 基于改進DE的EIT重建算法

2.1 基本 DE的算法 DE算法是一種基于種群進化的全局啟發(fā)式算法[8]，它在許多優(yōu)化問題上的效果均優(yōu)于遺傳算法[9]、粒子群算法[10]、模擬退火算法[11]等。該算法一般包括變異、交叉、選擇3個基本操作。首先在整個搜索空間隨機生成一組初始種群，然后對種群中每個個體進行變異、交叉操作得到臨時種群，最后用臨時種群中適應(yīng)性更好的個體替代原始種群中相應(yīng)的個體。整個算法的步驟見圖1。

圖1 DE基本算法流程

在整個算法迭代過程中，變異是其中的核心過程，它采用隨機偏差擾動生成變異向量的方式能夠增加種群多樣性，增強算法全局搜索能力。根據(jù)生成變異向量方法的不同，能夠組合成不同的變異策略，其中基本變異策略DE/rand/1/bin方程見公式（4）。

此外，其他典型的變異策略 DE/best/1/bin和DE/rand-to-best/1/bin的方程分別見公式（5）、（6）。

其中，νi,G+1表示生成的變量向量，i表示種群個體序列號，G表示種群進化代數(shù)，r1≠r2≠r3≠i，且i、r1、r2、r3∈[1,NP]，NP表示種群規(guī)模，xbest,G為第G代種群中的最優(yōu)個體，F(xiàn)、λ是縮放因子，為常數(shù)。

2.2 改進DE算法 本文在基本DE算法的基礎(chǔ)上，引入隨機最佳變異和局部增強算子，提出一種能在保證種群多樣性前提下加速收斂、提高局部搜索能力的改進DE算法。具體實現(xiàn)如下：

2.2.1 隨機最佳變異 在基本 DE算法中，為克服DE/rand/*/*變異策略收斂速度較慢的缺點，普遍采用DE/best/*/*變異策略，在變異過程中，種群中個體均采用相同的最佳個體，對加快算法收斂性有良好的引導(dǎo)作用，但最佳個體周圍集中大量的尋優(yōu)個體，種群多樣性遭到破壞，容易導(dǎo)致算法局部收斂。

為使變異向量具有較好的多樣性和引導(dǎo)性，本文采用一種隨機最佳變異策略，見公式（7）。

其中，r1≠r2≠r3≠i，且i、r1、r2、b∈[1,NP]，RMS(xb,G)≤min{RMS(xr1,G), RMS(xr2,G)}。

隨機最佳變異策略將3個隨機選擇的個體中的局部最優(yōu)個體當(dāng)作變異算子的基向量，而另外2個向量作為差向量對，這樣既保證種群多樣性，又為種群優(yōu)化提供了一定的引導(dǎo)作用，避免種群大量聚集于最佳個體周圍，從而避免了算法的“早熟”。

2.2.2 局部增強算子 引入隨機最佳變異策略，保證了算法收斂速度和種群多樣性之間的相對平衡，但算法收斂速度相對較慢，為進一步加快算法收斂，引入局部增強算子，對除最佳個體外的其他種群個體以比例因數(shù)CP（0＜CP＜1）的概率進行重新賦值，使部分個體向最佳個體靠攏，以增強這些個體的貪婪性，加快算法收斂速度，局部增強算子見公式（8）。

其中，r1≠r2≠i，且i、r1、r2∈[1,NP]，?i,G+1是增強后的新個體，用來代替?zhèn)€體xi,G+1，xbest,G+1是新種群中的最佳個體，g為算法迭代次數(shù)，改進算法引入新的比例因數(shù)CP，當(dāng)CP=1時，除最佳個體外，其余所有個體向量都需重新賦值；當(dāng)CP=0時，未選中個體進行局部增強。

局部增強算子的目的是在增加部分個體貪婪性的同時保證種群多樣性，以使算法達到既快又好地逼近到全局最優(yōu)解。算法在增加局部增強算子后，每次迭代僅按一定比例使部分個體變得更加貪婪，從而在一定程度上限制了算法的貪婪性，避免影響算法的收斂性。需注意的是，隨著迭代過程的進行，種群中這部分個體隨機擾動范圍也在動態(tài)縮小，這種方式對算法局部搜索能力有較大的提高，尤其是能減少收斂到全局最優(yōu)解所需的迭代次數(shù)，但算法也有陷入局部最優(yōu)的可能。

3 仿真實驗與結(jié)果

3.1 創(chuàng)建仿真模型 在ANSYS 10.0軟件上采用有限元法對頭部進行建模，并對EIT正問題進行求解。三層有限元頭球模型如圖 2A，由 14 332個單元和20 193個節(jié)點組成。頭部模型分為頭皮、顱骨、大腦及病變組織4個部分[12]。為簡化計算，假定頭部模型為各向同性均質(zhì)導(dǎo)體，且電阻率值已知。根據(jù)先驗知識，將模型 4個部分的電阻率按比例設(shè)定為：頭皮∶顱骨∶大腦∶病變組織=1∶15∶1∶2[13]，其目標(biāo)電阻率分布的剖面圖如圖2B，其中淺色區(qū)域為病變組織。

圖 2 采用有限元法對頭部進行建模。A.頭球模型；B.目標(biāo)電阻率分布

3.2 仿真參數(shù)設(shè)定 本文對三維頭球模型進行仿真實驗，并采用相對電極注入方式向模型左右注入±5 mA的電流。在逆問題重建實驗中，將頭球模型的頭皮、顱骨、大腦及病變組織4個部分的上限設(shè)為[1.20，17.40，1.20，2.25]，下限設(shè)為[0.81，13.60，0.81，1.75]，初始種群需從此約束邊界范圍內(nèi)隨機生成。由于大腦和頭皮電阻率的設(shè)定值相等，為簡化計算，假定變量個數(shù)D=3，種群規(guī)模NP=60，目標(biāo)函數(shù)為公式（3）。算法終止條件為迭代次數(shù)達到100次。根據(jù)文獻將比例因數(shù)F、交叉因數(shù)CR分別設(shè)定為0.8、0.9[13]。局部增強算子的比例因數(shù)CP等到下文分析此參數(shù)對重建結(jié)果的影響后再另行確定。

3.3 仿真實驗與結(jié)果 為了對算法圖像重建質(zhì)量進行定量分析，分別引入相對誤差RE和誤差總和TE準則，見公式（9）、（10）。

其中，ρ*是指模型各個部分的重建電阻率值，ρ0是指相應(yīng)的模型各個部分的目標(biāo)電阻率值。

其中，m為電阻率分布自由度，ρitar和ρirec分別代表單元i的目標(biāo)電阻率值和重建電阻率值。

圖3為改進DE算法的重建圖像誤差總和TE隨局部增強因子的比例因數(shù)CP變化的曲線，當(dāng)CP=5/60時，重建圖像的整體誤差最小，算法效果最佳；當(dāng)CP為1/60～5/60時，在相同的迭代次數(shù)下，隨著CP取值的增加，重建圖像的誤差逐漸減小，說明當(dāng)聚集在最優(yōu)個體周圍的個體在一定范圍內(nèi)時，會在一定程度上加快算法收斂的速度，提高重建圖像的精度；當(dāng)CP為5/60～15/60時，隨著CP取值的增加，重建圖像的誤差反而逐漸增加，說明當(dāng)聚集在最優(yōu)個體周圍的個體超出一定范圍后，算法的貪婪性過大，反而使算法陷入“早熟”。因此建議采用此改進方法進行三維腦部阻抗圖像重建時，將CP的取值范圍確定在 3/60～7/60，即可能獲得更好的圖像重建效果。在本文的仿真對比實驗中選擇局部增強因子CP=4/60。

圖3 誤差總和隨局部增強因子的比例因數(shù)變化的關(guān)系

采用上述實驗參數(shù)設(shè)置，在相同的初始條件下，將本文提出的改進算法與其他算法進行仿真對比實驗。不同算法的電阻抗重建圖像結(jié)果見圖4。

圖4 不同算法電阻率重建分布

其中算法DE-rand、DE-best分別為采用文獻[13]中的DE/rand/1/bin變異策略和DE/best/1/bin變異策略，算法PSO采用文獻[14]中的粒子群優(yōu)化算法。比較圖4及圖2B發(fā)現(xiàn)，這些算法均能在電阻抗重建過程中得到較為清晰的圖像，且能對病變組織進行準確定位，但不同算法的電阻率重建圖像與目標(biāo)電阻率的分布圖像僅在顏色深淺上有細微差別，用肉眼無法看出明顯差異，因此需要對阻抗重建結(jié)果進行定量分析。

表 1是經(jīng)過 100輪迭代后不同算法的圖像重建結(jié)果，其中t表示各個算法單次循環(huán)所需平均時間，其單位為分鐘。由表1可見，改進算法相對于算法DE-rand、PSO、DE-best而言，其模型各個部分整體的重建誤差總和TE分別降低了66.996%、32.723%、49.834%；其單次循環(huán)所需平均時間分別降低了 8.331%、2.017%、4.793%，表明在相同條件下改進算法圖像重建精度更高、迭代效率更高，算法整體性能更好。

表1 不同算法的圖像重建實驗結(jié)果

由表2可見，與其他部分相比，各個算法對病變組織部分的重建誤差更大，表明在求解EIT正問題時，大腦顱骨會阻礙一部分電流的流入，導(dǎo)致病變部分的重建誤差變大，大體上符合理論情況。從總體上看，改進算法相比算法 DE-rand、PSO、DE-best，其模型的大腦、顱骨、頭皮、病變部分電阻率重建誤差分別降低86.794%、68.690%、86.720%、86.874%，-16.084%、33.431%、-18.571%、81.164%，80.378%、53.955%、80.238%、70.146%，表明改進算法有效降低了圖像重建的誤差，提高了圖像重建的精度，算法性能得到一定的提升。

表2 頭球模型不同部分的圖像重建相對誤差（%）

由圖 5可見，所有算法在迭代初期收斂速度較快，隨著迭代次數(shù)的增加，收斂速度慢慢趨于平緩。當(dāng)?shù)螖?shù)＜18時，改進DE算法較其他算法無明顯優(yōu)勢，說明改進算法在種群多樣性和算法收斂速度上保持了較好的水平；當(dāng)?shù)螖?shù)為 18～27時，改進算法的局部增強算子逐漸發(fā)揮作用，使改進算法的重建效果逐漸優(yōu)于算法1和算法2；當(dāng)算法迭代次數(shù)為18～100時，改進算法明顯優(yōu)于其他算法的效果和精度。

圖5 種群平均目標(biāo)函數(shù)值隨迭代次數(shù)變化關(guān)系

4 結(jié)束語

電阻抗圖像重建算法是 EIT研究的關(guān)鍵技術(shù)之一。本文提出了一種新的用于電阻抗重建的改進 DE算法。該算法通過引入隨機最佳變異策略和局部增強算子，在算法收斂速度和種群多樣性之間取得良好的平衡，并在一定程度上加快了算法收斂速度，增強了算法局部搜索能力。結(jié)果表明，改進算法能夠獲得較為清晰的電阻抗重建圖像，且能對病變組織進行準確定位。同時，與其他算法相比，改進算法的圖像重建誤差減小30%以上，有效提高了算法的整體性能，為EIT臨床研究提供了一定的參考。

由于DE算法是一種相對較新的全局類智能尋優(yōu)算法，仍有較多的相關(guān)理論問題亟需完善。改進算法由于每次均需對所有種群個體進行目標(biāo)函數(shù)評價，仍存在測量數(shù)據(jù)、計算時間較長的缺點。因此，如何在一定的實驗條件下增加測量數(shù)據(jù)并減少仿真計算時間是下一步研究的重點。

[1] Silvera-Tawil D, Rye D, Soleimani MA. Electrical impedance tomography for artificial sensitive robotic skin: a review. IEEE Sens J, 2015, 15(4): 2001-2016.

[2] 陳會, 閆丹丹. 基于磁共振技術(shù)的生物電阻抗成像研究進展. 中國醫(yī)學(xué)影像學(xué)雜志, 2014, 22(9): 714-717.

[3] Wagenaar J, Adler A. Electrical impedance tomography in 3D using two electrode planes: characterization and evaluation. Physiol Meas, 2016, 37(6): 922-937.

[4] Schullcke B, Gong B, Moeller K. Steps towards 3D electrical impedance tomography. Conf Proc IEEE Eng Med Biol Soc, 2015: 5323-5326.

[5] 萇飛霸, 張和華, 顏樂先, 等. 電阻抗斷層成像技術(shù)研究.中國醫(yī)療器械雜志, 2016, 40(1): 52-54.

[6] Trepte CJ, Phillips CR, Solà J, et al. Electrical impedance tomography (EIT) for quantification of pulmonary edema in acute lung injury. Crit Care, 2016, 20(1): 18.

[7] 張輝, 李穎, 王西明, 等. 電阻抗斷層成像的MPSOMNR算法研究. 計算機工程與應(yīng)用, 2013, 49(9): 29-32.

[8] Ribeiro RR, Feitosa AS, Souza RD, et al. A modified differential evolution algorithm for the Reconstruction of electrical impedance tomography images//Biosignals and Biorobotics Conference. Salvador: IEEE, 2014: 1-6.

[9] Rashid A, Kim S, Liu D, et al. A dynamic oppositional biogeography-based optimization approach for time-varying electrical impedance tomography. Physiol Meas, 2016,37(6): 820-842.

[10] Martin S, Choi CT. Nonlinear electrical impedance tomography reconstruction using artificial neural networks and particle swarm optimization. IEEE Trans Magn, 2016,52(3): 1-4.

[11] Zhang XL, Wang W, Sze G, et al. An image reconstruction algorithm for 3-D electrical impedance mammography.IEEE Trans Med Imaging, 2014, 33(12): 2223-2241.

[12] Hamilton SJ, Lassas M, Siltanen S. A direct reconstruction method for anisotropic electrical impedance tomography.Inverse Probl, 2014, 30(7): 1375-1378.

[13] 閆丹丹, 陳會. 基于電阻抗成像技術(shù)的腦病變組織檢測仿真研究. 中國醫(yī)學(xué)影像技術(shù), 2014, 30(7): 1113-1116.

[14] 陳民鈾, 楊艷利, 何為, 等. 基于粒子群優(yōu)化算法的電阻抗圖像重建. 重慶大學(xué)學(xué)報, 2011, 34(1): 82-87.