基于張量分解的電阻抗三維成像方法

李秀艷，尉 蕾，王 琦，3，段曉杰，馬會(huì)梅，陸紀(jì)璇

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

電阻抗層析成像（EIT）是近年來發(fā)展起來的一種新興功能成像技術(shù)。EIT技術(shù)是通過對(duì)被測(cè)物場(chǎng)施加激勵(lì)電流，檢測(cè)物場(chǎng)邊界電壓信息變化，進(jìn)而獲得被測(cè)物場(chǎng)內(nèi)部阻抗信息分布的非侵入方法。EIT成像系統(tǒng)[1]由于其無輻射、高速、非侵入性的傳感能力，在生物醫(yī)學(xué)成像領(lǐng)域得到了迅速發(fā)展。根據(jù)不同介質(zhì)具有不同阻抗的理論，可以通過同時(shí)測(cè)量邊界上的電壓或電流來恢復(fù)內(nèi)部阻抗。EIT是一種極具潛力的醫(yī)學(xué)成像技術(shù)，但是在圖像重建過程中對(duì)圖像領(lǐng)域信息的充分利用仍然是一個(gè)挑戰(zhàn)。

根據(jù)人體肺部三維結(jié)構(gòu)特點(diǎn)，EIT成像問題本身就是一個(gè)3D問題。而目前大多數(shù)EIT重建算法集中在二維成像，其橫截面切片只能顯示真實(shí)的三維物體的部分信息，所以限制了EIT的圖像重建效果。此外，由于電流的路徑遍布整個(gè)三維域，平面外電導(dǎo)率的變化通常會(huì)影響整個(gè)電極平面，并在重建的2D圖像中產(chǎn)生相當(dāng)大的畸形，這些在2DEIT重建產(chǎn)生的問題都促進(jìn)了3DEIT圖像重建的發(fā)展[2]。為此，研究者設(shè)計(jì)了多電極實(shí)驗(yàn)?zāi)Ｐ?，增加了縱向的測(cè)量，以達(dá)到三維成像的效果。王海剛等[3]提出通過獲取圖像的灰度值，構(gòu)造小尺寸的等值面，使用圖像處理方式將等值面進(jìn)行拼接，最后以三維形式顯示三維圖像的構(gòu)造方法。趙進(jìn)創(chuàng)等[4]提出將二維斷層的圖像進(jìn)行重建，提取圖像輪廓特征，在此基礎(chǔ)上重建三維圖像。目前的三維重建方法一般首先對(duì)每個(gè)二維切片進(jìn)行圖像處理[5]，之后羅列為三維圖像，沒有考慮每層圖像間的相關(guān)性，從而實(shí)驗(yàn)數(shù)據(jù)所包含的空間信息沒有被充分利用。本文提出一種基于張量的電阻抗三維成像方法，該方法充分利用了圖像數(shù)據(jù)的多路性質(zhì)，保留了圖像的3D結(jié)構(gòu)。通過圖像的張量表示方法，能夠保留圖像的高階特征，可以避免圖像信息的損失?；趶埩康膱D像視為三維張量，其在時(shí)間、空間尺度上的展開矩陣分別表示圖像幀在時(shí)間序列上的強(qiáng)關(guān)聯(lián)，以及單幀圖像在空間尺度上像素間的弱關(guān)聯(lián)，可以有效地提高圖像的重建質(zhì)量。同時(shí)TV算法[6]可以保持清晰的界面，有保護(hù)圖像邊緣的效果[7]，本文將其與基于張量的電阻抗三維成像方法進(jìn)行對(duì)比實(shí)驗(yàn)。

1 背景和相關(guān)工作

1.1 三維胸腔輪廓EIT正問題模型構(gòu)建

為了獲得三維重建圖像的結(jié)果，設(shè)計(jì)了一種呼氣末的三維EIT模型的電極，在人體三維胸部模型周圍均勻布置5層電極陣列，每層電極陣列上有16個(gè)電極用于測(cè)量整個(gè)胸腔[8]，如圖1所示。

圖1 設(shè)置5個(gè)電極環(huán)的3D胸部模型Fig.1 3D thorax modelswith designed 5 electric ring

考慮到內(nèi)存限制和速度要求，選擇相鄰環(huán)間的20個(gè)重構(gòu)層，即5cm進(jìn)行測(cè)量；有111層重建三維重建圖像，即重建層分離為0.25 cm。建立反問題三維網(wǎng)格如圖2所示。圖中感知域縱向劃分為111層，每層的分辨率為32×32。

圖2 用于反問題計(jì)算的體素網(wǎng)格Fig.2 Voxel mesh for inverseproblem calculation

1.2 EIT圖像重建基本原理

上述模型的構(gòu)建屬于EIT正問題，那么EIT的逆問題[9]是EIT圖像的重建過程。倘若被測(cè)量的物場(chǎng)和空?qǐng)龅南鄬?duì)電導(dǎo)率對(duì)的變化較小，那么可以通過線性方程來求解逆問題：

式中：δU為邊界電壓的變化值；δσ是物場(chǎng)和空?qǐng)鲋g的電導(dǎo)率變化向量；J為雅可比矩陣或靈敏度矩陣[10]。

為了從邊界電壓測(cè)量中估計(jì)電導(dǎo)率分布，常用方法是將其重新表述為Tikhonov正則化[11]問題：

為了幫助清晰的界面的結(jié)構(gòu)，TV正則化[12]被納入EIT模型，產(chǎn)生：

式中：λ為超參數(shù)。

2 基于張量分解的EIT圖像重建

2.1 張量的基礎(chǔ)知識(shí)

張量的mode-n矩陣化是對(duì)張量進(jìn)行降階，稱為張量的展開或者張量的平整化。張量I1×I2×…×IN的mode-n矩陣表示為矩陣

為了捕捉信號(hào)中的潛在模式和重要特征，并且移除原始數(shù)據(jù)的冗余成分，經(jīng)常對(duì)張量進(jìn)行Tucker分解[13]，Tucker分解可以看做是矩陣的奇異值分解的推廣。Tucker分解[14]是將一個(gè)N階張量分解為一個(gè)核心張量與因子矩陣A1，A2，…，AN的n-模乘積形式，即

式中：Ak∈RIk×Jk（k=1，…，N）是正交矩陣。

2.2 EIT張量模型構(gòu)建

由于人體的胸腔本身就是三維結(jié)構(gòu)，根據(jù)1.1節(jié)建立的三維胸腔輪廓模型，本文根據(jù)圖像的分辨率以及層數(shù)構(gòu)建張量模型其中p×q為圖像的分辨率，m為圖像的總層數(shù)。由于像素?cái)?shù)會(huì)隨著圖像輪廓形狀的變化而發(fā)生變化，所以肺部三維EIT圖像是不規(guī)則的，從而導(dǎo)致二維圖像也是不規(guī)則的。為此本文使用零填充的方法來獲取有效像素[15]?？紤]到每張圖片的有效像素?cái)?shù)為L(zhǎng)，每張圖片的重建圖像的總像素?cái)?shù)為n，總共m層圖片，定義2個(gè)分解算子分別為ΡEi∈RL×n和ΡLi∈Rn×m，ΡEi用來提取有效像素，ΡLi用來選擇對(duì)應(yīng)的圖片層數(shù)，則圖片中的Δσi∈RL×m，表示為Δσi=ΡEiΡLiΔσ，圖3所示為一層EIT像素?cái)U(kuò)展原理圖。

圖3 EIT像素?cái)U(kuò)展原理Fig.3 Schematic of EIT pixel expansion

通過1.1節(jié)構(gòu)建的胸部仿真模型獲得的仿真電壓數(shù)據(jù)，通過式（6）可以進(jìn)一步估計(jì)電導(dǎo)率的變化Δσ∈RN×1。通過零填充方法，提取每層圖片的有效像素值，最終將電導(dǎo)率變化值Δσ擴(kuò)展為32×32的矩陣，最后將111個(gè)擴(kuò)展之后的數(shù)據(jù)，在z軸方向羅列，最后整合為三階張量32×32×111。

2.3 基于張量分解的EIT重建算法

在實(shí)際應(yīng)用中，由于人體胸腔的三維結(jié)構(gòu)以及電場(chǎng)的三維分布特性，EIT三維圖像在各圖層間具有一定空間相關(guān)性，因此2.2節(jié)中得到的EIT張量模型可進(jìn)行線性-低秩近似，提取圖像三維空間相關(guān)特征，去除冗余信息。因此一個(gè)EIT信號(hào)可以寫為[16]：

本文利用張量的多路結(jié)構(gòu)，通過將N階張量的每個(gè)模式乘以感知矩陣，來獲取壓縮測(cè)量集。根據(jù)公式（4）中對(duì)mode-n乘積的定義，將數(shù)據(jù)張量乘以在mode-n的感知矩陣。

其中Φn∈RRn×In是對(duì)應(yīng)mode-n的感知矩陣[17]，論文中設(shè)置為單位矩陣。注意公式（8）是原始張量乘除moden以外的所有模式的感知矩陣。

結(jié)合公式（8）和（9）得到重建公式為：

式中：“+”代表矩陣的MP偽逆；Zn是壓縮多路測(cè)量集的mode-n展開（n=1，2，3）。

圖4 三階張量的展開示意圖Fig.4 Unfolding diagram of the 3rd-order tensor W

當(dāng)τ→0時(shí)，W*τ→W+，“+”代表矩陣的MP偽逆，τ為截?cái)嚅撝怠?/p>

對(duì)角矩陣S*τ∈RI2×I1定義如下：

因此重建公式修改為如下公式：

圖5為多線性一秩一（R1，R2，R3）三階張量的恢復(fù)。圖5說明了原始張量可以從多路測(cè)量集Z(n)和W中恢復(fù)。

圖5 多線性-秩-（R1，R2，R3）三階張量的恢復(fù)Fig.5 Restoration of 3rd-order tensor of multilinear-r and（R1，R2，R3）

本文中誤差?以及閾值τ由經(jīng)驗(yàn)確定誤差，其中?=8.5×10-20，τ=4.15×10-20。

基于張量的EIT圖像重建算法總結(jié)如下：

2.4 仿真實(shí)驗(yàn)

為了評(píng)估新算法的性能，對(duì)模擬數(shù)據(jù)進(jìn)行了實(shí)驗(yàn)。運(yùn)用第2.2節(jié)中Comsol建立基于人體結(jié)構(gòu)先驗(yàn)信息的胸部仿真模型，進(jìn)行EIT正問題計(jì)算。仿真模型設(shè)置5個(gè)電極環(huán)，每個(gè)電極環(huán)均勻分布16個(gè)電極，用于測(cè)量整個(gè)胸腔區(qū)域[19]，其中肺部的電導(dǎo)率是0.1 S/m，胸腔電導(dǎo)率可看作脂肪電導(dǎo)率0.024 S/m。為了模擬實(shí)際測(cè)量系統(tǒng)中的典型噪聲水平，將高斯零均值隨機(jī)噪聲添加到模擬電壓中。噪聲幅度是模擬電壓幅度的0.5%。采用相鄰激勵(lì)相鄰測(cè)量方法[20]，每次激勵(lì)可測(cè)得13個(gè)電壓數(shù)據(jù)，共需進(jìn)行16次電流激勵(lì)，因此共得到208個(gè)獨(dú)立電壓值用于重建一幅EIT圖像。所有的實(shí)驗(yàn)都是在英特爾i5處理器的筆記本上使用Matlab軟件進(jìn)行的。為了評(píng)價(jià)重建質(zhì)量，使用峰值信噪比，定義為：

為了反映重建圖像和真實(shí)圖像在空間結(jié)構(gòu)上的相似程度，采用結(jié)構(gòu)相似度（SSIM）參數(shù)作為評(píng)價(jià)指標(biāo)，定義為：

式中：μI、μK分別表示原圖像與重建圖像的均值分別表示原圖像與重建圖像的方差；σIK表示原圖像和重建圖像之間的協(xié)方差；C1、C2均為常數(shù)，通常取值非常小，目的是避免分母為零。

圖6 3種EIT三維重建結(jié)果Fig.6 EIT 3D reconstruction resultsusing three methods.

選取上述三維EIT重建結(jié)果對(duì)應(yīng)的第49、53、61、77這4層進(jìn)行二維圖像顯示，可以很明顯地發(fā)現(xiàn)，進(jìn)行閾值截?cái)?，即?4.15×10-20后的重建圖像能夠較為清楚地顯示邊緣信息，與真值圖的結(jié)構(gòu)相似度最高，成像質(zhì)量高。圖7所示為基于3種方法的3PEIT重建中的2D切片成像圖。

圖7 基于3種方法的3D EIT重建中典型的2D切片成像圖Fig.7 Typical 2D slicesfrom 3D EIT reconstruction resultsbased on threemethods

為了將本文所提出的三維重建方法和直接二維成像效果進(jìn)行對(duì)比，進(jìn)行直接二維成像實(shí)驗(yàn)，將測(cè)量得到的電壓數(shù)據(jù)，通過TV的方法計(jì)算電導(dǎo)率值，通過零填充的方式提取有效像素，并將其擴(kuò)展為大小為32×32的矩陣X，由于矩陣是張量的一種特殊形式，所以類似于三維張量重建方法，將公式（14）轉(zhuǎn)換為二維公式進(jìn)行EIT圖像重建，其中感知矩陣是隨機(jī)生成的，最后選取相應(yīng)位置的像素值進(jìn)行圖像顯示。由1.1節(jié)得知人體胸腔模型在感知域縱向劃分為111層，所以需要對(duì)111個(gè)橫截面進(jìn)行單獨(dú)處理。本文選擇了較為明顯的三層圖像進(jìn)行顯示，如圖8所示。圖8中依次顯示了原圖、TV法、閾值τ=0和τ=τ0的EIT直接二維重建圖像，其中τ0依舊根據(jù)實(shí)驗(yàn)經(jīng)驗(yàn)選取，即τ0=4.15×10-20。為了便于與上述三維張量重建結(jié)果做對(duì)比，本文將111張重建圖像的PSNR和SSIM進(jìn)行了平均值計(jì)算，3種重建方法的PSNR平均值分別是29.6、9.8和173 dB，SSIM分別是0.656、0.432和0.783。通過實(shí)驗(yàn)結(jié)果很明顯地發(fā)現(xiàn)，在二維直接成像中τ0=4.15×10-20的張量分解重建方法較好，但是由于沒有利用每層圖像間的關(guān)聯(lián)性，所以結(jié)構(gòu)相似度顯然不如三維張量重建算法，重建圖像質(zhì)量也相對(duì)較差。

圖8 基于3種方法的直接二維EIT重建圖像Fig.8 2D EIT reconstruction results by three methods

3 實(shí)驗(yàn)結(jié)果

為了驗(yàn)證基于張量的電阻抗圖像重建方法的可行性，采用TJPU-EIT系統(tǒng)進(jìn)行胸形模體實(shí)驗(yàn)，如圖9所示。

圖9 基于肺部模型的EIT實(shí)驗(yàn)系統(tǒng)Fig.9 Experimental equipment for thoracic model imaging

根據(jù)第2.4節(jié)討論的有5個(gè)電極環(huán)仿真模型幾何尺寸，采用樹脂材料3D打印出清晰的胸腔模體以及白色肺部模體。呼氣末下胸腔邊界和肺的三維模型如圖10所示。模體填充導(dǎo)電性為0.32 S/m的鹽水溶液。在相鄰頻率50 kHz處注入5 mA電流，測(cè)量相鄰電壓數(shù)據(jù)。

圖10 胸腔以及肺部模型Fig.10 3D printed phantom of thorax

基于Tivkhonov正則化、TV法以及張量法重建圖像的實(shí)驗(yàn)成像結(jié)果如圖11所示。

圖11 基于3DEIT實(shí)驗(yàn)結(jié)果的圖像重建對(duì)比Fig.11 Comparison of reconstruction results based on 3D EIT

圖11中，計(jì)算了重建圖像的PSNR和SSIM參數(shù)，PSNR分別是12.39、19.39和189.1 dB，SSIM分別是0.607、0.721和0.835。TV法和Tivkhonov正則化重建法都是基于三維模型的橫截面進(jìn)行單獨(dú)圖像重建，然后羅列為三維達(dá)到三維重建圖像的目的?；跈M截面的重建方式，忽略了人體胸腔三維結(jié)構(gòu)的特點(diǎn)，平面化致使縱向信息相關(guān)聯(lián)的特性無法利用，最終導(dǎo)致成像效果較差?；趶埩糠纸獾娜S重建方式是將EIT數(shù)據(jù)作為整體進(jìn)行圖像重建，充分利用了圖像的縱向關(guān)聯(lián)性，所以基于張量分解的EIT圖像重建效果最佳。

為了證明基于張量分解的重建算法充分考慮了每層EIT數(shù)據(jù)層間的空間結(jié)構(gòu)關(guān)系，成像質(zhì)量明顯改善，選取111張EIT圖像中的第64、90張進(jìn)行單獨(dú)顯示，如圖12所示。從實(shí)驗(yàn)結(jié)果中可以發(fā)現(xiàn)，基于張量分解重建算法中每一層肺部成像的細(xì)節(jié)較為清晰，并且SSIM參數(shù)明顯提高，提高了空間信息的利用率，成像質(zhì)量好。

圖12 基于3D EIT實(shí)驗(yàn)重建結(jié)果的2D切片F(xiàn)ig.12 Typical 2D slices from 3D EIT experimental reconstruction results

4 結(jié)論

本文提出了一種基于張量的電阻抗三維成像方法，利用實(shí)驗(yàn)中測(cè)得的數(shù)據(jù)建立張量模型，進(jìn)行基于張量的三維電阻抗成像。

（1）圖像的張量表示可以保留圖像的固有結(jié)構(gòu)，更好地利用層與層之間的結(jié)構(gòu)信息，避免圖像信息的損失。

（2）解決了通常由于平面外電導(dǎo)率的變化會(huì)影響電極平面，并在產(chǎn)生的2D重建圖像中產(chǎn)生較大的偽影的問題。

（3）通過計(jì)算基于張量分解的三維重建算法與TV、Tikhonov以及二維張量重建算法的信噪比可知，基于張量分解的圖像重建質(zhì)量最好，相比于基于二維張量的重建算法信噪比提高了13.8%。