1.   引言

电阻抗层析成像(Electrical Impedance Tomography, EIT)技术是继形态、结构成像之后出现的一种新型可视化检测技术^[1]，以电磁场理论为基础，其敏感场的数学模型是基于Maxwell方程组的电磁场模型^[2]，具有响应速度快、结构简单、成本低、非侵入性、无辐射性、实时无损、功能成像等优势^[3]，在管道工业^[4]、生物医学辅助诊断^[5]、生物技术^[6]等领域具有广阔的应用前景。

到目前为止，在医学应用领域，EIT重建图像分辨率尚无法与计算机成像^[7](Computed Tomography, CT)等主流影像技术相比，存在一定的局限性。在进行图像重建时，由于软场效应及逆问题的不适定性，现有的成像算法只能得到一个近似的最优解。重建算法不精确以及物场剖分的低分辨率，使得目标边缘轮廓存在较大的伪迹。这些伪迹降低了图像的质量，甚至会掩盖真正需要重建的目标信息，在医学成像领域也会极大降低重建图像的临床诊断价值。因此，如何有效提高EIT技术的空间分辨率、减少重建图像伪迹一直是该领域的研究重点和难点。

近年来，随着对EIT成像空间分辨率要求的不断提高，基于优化理论的迭代类算法被应用于EIT问题。2003年，李英^[8]提出了一种将Tikhonov正则化方法和代数重建法(Algebra Reconstruction Technique, ART)相结合的组合型图像重建算法，使得重建图像分辨率得到提高；2009年罗辞勇等人^[9]结合计算机断层成像(CT)中的均值反投影理论，对等位线反投影算法进行改进，减少了重构图像的星状伪迹；2016年，李东晔^[10]针对传统Tikhonov正则化图像分辨力低及边界模糊的问题，提出一种基于全变差的正则化优化模型，提高了重构图像质量，减少了伪迹；2017年，岳士弘等人^[11]针对电学层析成像技术空间分辨率较低且稳定性差的问题，提出了一种基于灵敏度系数更新及伪迹辨识的电学层析成像算法；2019年，陈瑞娟等人^[12]提出一种基于结构先验信息的磁探测电阻抗成像算法，在重建图像结构和电导率重建数值精度方面得到显著提高。上述算法的研究有效提高了重建图像的空间分辨率，在一定程度上减少了重建图像的伪迹，但是伪迹的不利影响依然存在，仍然具有进一步改善的空间。定量地评价重建图像质量高低是进行伪迹修正的前提。经典的图像质量评价指标(如相关系数、相对误差等)都是基于有监督的评价方式，即需要原始参考图像(真实电导率分布)。但在实际应用中，被测物场的真实电导率分布是未知的，没有一个预定义的参考模型。这就需要设计一种客观的无监督的图像质量评价方法，即不需要参考图像而只依靠重建图像结果进行图像质量评价，反映图像质量。同时，一种有效的无监督图像质量评价方法也可用于图像处理领域，如衡量图像的失真程度以及预测失真的图像质量等。

基于以上研究分析，本文将快速模糊聚类算法引入EIT重建图像质量评价中，并结合邻域信息提出了一种无监督的图像质量评价指标，再利用该指标进行图像重建。通过COMSOL和MATLAB联合仿真建立位置、数量和大小不同的病变模型，结合电导率的改变，将该算法与Tikhonov正则化算法的重建图像进行对比来验证基于邻域信息和快速模糊C均值聚类的重建算法的伪迹修正能力，最后通过仿体实验进一步验证该算法在重建图像精度上的有效性和鲁棒性，为提高重建图像分辨率提供一种新的思路。

2.   重构图像伪迹优化算法

2.1   Tikhonov正则化算法

针对EIT逆问题求解中的不适定性，采用正则化技术是一种有效的手段。Tikhonov正则化算法通过引入先验信息和约束项、压缩逆问题的解空间，进而降低EIT逆问题求解的不适定性，已广泛应用于EIT图像重建^[13]。

Tikhonov正则化基本思想即最小化目标函数：

其中，$\lambda$为正则化因子；R(g)称为正则化函数。R(g)通常采用L2范数形式，即

其中，$\overline {\boldsymbol{g}}$是根据先验信息进行空场重置后得到的解的估计值；L是正则化矩阵。通常情况下，不考虑先验信息($\overline {\boldsymbol{g}}$为零)，即为标准Tikhonov正则化方法，它等价于求解：

可以证明${{\boldsymbol{S}}^{\rm{T}}}{\boldsymbol{S}} + \lambda {\boldsymbol{I}}$的逆总是存在的。于是标准Tikhonov的解可以表示为

在后续介绍中，将Tikhonov正则化算法简记为TR。

2.2   无监督图像质量评价指标

快速模糊C均值聚类算法(fast Fuzzy C-Means clustering, f-FCM)也称为1维像素聚类算法，是FCM算法在图像处理领域的扩展应用。该算法用于处理1维灰度值，若所有的像素构成集合X，像素集X被分为c个类且每一个像素分别为L个灰度级中的一个。假设H(l)是属于第l个灰度级的像素数量，即频数，模糊隶属度函数u_il表示第l个灰度水平对第i个类的隶属度，则f-FCM算法的目标函数为

其中，$m$为模糊指数，v_i为类中心。 f-FCM算法的迭代过程需要按不断更新式(6)中的模糊隶属度u_il及式(7)类中心v_i：

当式(8)中的误差E小于一定的数值，f-FCM算法停止迭代，输出模糊隶属度矩阵和聚类中心。根据隶属度的大小判断像素所属的类及对应的类中心，类中心就是每个类的代表灰度值。

本文利用邻域信息和f-FCM定义一种无监督的图像质量评价指标。图1(a)中，x1-x8为任意像素V的8个邻域像素，模型中含有两个圆形目标，其电导率为3 S/m，背景电导率为1 S/m。EIT重建图像中的所有像素都可以根据f-FCM划分为3类，分别分布在两个圆的内部(目标内)、外部(背景)和边界上(伪迹)。图1(b)是利用TR算法对图1(a)中的模型重建的图像，图1(c)是利用f-FCM对图1(b)进行的二分类聚类图像。

图 1  空间邻域V的图解

下载: 全尺寸图片 幻灯片

如果EIT算法能正确无失真地重建模型中的所有对象，那么V₁和V₂中的任意像素及其8个相邻像素的灰度值肯定是相同的。V₃中8个邻域像素可以分成两个部分，每个部分由相同灰度值的连通像素组成。然而，重建的图像通常与真实的图像不一致。一般情况下，EIT图像中位于V₁和V₂的像素分别具有较高和较低的灰度级。为了评估EIT图像分辨率高低，我们利用f-FCM将所有像素灰度聚成两类，并将灰度二值化为0和1。具有高空间分辨率的EIT图像应具备以下特征：(1)目标内部像素和背景像素。任何内部像素及其8个相邻像素都有相同的灰度值1，而任何背景像素都有相同的灰度值0。(2)目标边界像素。任意边界像素的8个邻域像素由两个连通部分组成：其中一部分的灰度值都是1，另一部分的灰度值都是0。

根据以上特征，当EIT图像中任意像素的8个邻域像素具有以下特征时，就认为该像素的灰度值是最不可靠的：

(1)一个像素的8个连通邻域交替分布灰度值0和1。事实上，任何EIT算法难以识别像素大小的目标，因此这些交替分布的灰度值可能是随机分布的噪声、伪迹等。

(2)一个像素的8个连通邻域具有相同的灰度值，但是与该像素的灰度值正好相反。

根据上述特征，利用f-FCM对EIT重建图像进行二值化后，将评价EIT图像的指标描述为

其中，p_i是V_i的8个邻域中灰度为1的最大连通的像素个数，q_i是V_i的8个邻域中灰度值为0的最大连通的像素个数。如果V_i的灰度值为1，则${\rm{num}}({p_i}T{q_i},{V_i})$等于p_i；如果V_i的灰度值为0，则${\rm{num}}({p_i}T{q_i},{V_i})$等于q_i。式中两项的分母为6和8的目的是将NC值归一化到[0,1]区间。

图2显示了12种不同0,1分布的邻域类型及其NC值，领域类型分别记为V₁～V₁₂。其中，V₁～V₆的灰度分布对应了低分辨率的情况，V₇～V₁₂的灰度分布对应了高分辨率的情况。

图 2  12个不同灰度分布的邻域类型及NC值

下载: 全尺寸图片 幻灯片

NC值越大，表明像素的灰度值越可靠，即图像的分辨率越高。由图2结果看出，V₁-V₆的NC值都较小，这意味着它们可能具有较低的空间分辨率，这与我们视觉上直观的判断是一致的；V₆-V₁₂的NC值都较大，这意味着它们可能具有较高的空间分辨率，也与我们前面分析的结果一致。因此，NC值可以反映EIT图像的空间分辨率，NC值越大，图像的空间分辨率越高。该方法根据重建图像灰度进行聚类分割，利用邻域信息判断像素可靠性，不需要参考真实电导率分布，因此NC是一种是无监督的图像质量评价指标。

2.3   伪迹修正算法

医学EIT重建图像中的伪迹会极大降低其临床诊断价值。因此，需要研究相应的伪迹修正算法，减少伪迹对诊断结果的影响。1.2节提出的无监督图像质量评价指标(NC)可以反映目标像素的可靠性，进而也能有效判断哪些像素是不可靠的，即能够准确地反映出伪迹的位置和覆盖区域。Tikhonov正则化方法虽然能够得到近似稳定的图像重建效果，但算法的解过于光滑，导致对目标边缘重建效果并不理想^[14]。本文利用NC指标对Tikhonov正则化算法进行修正，该修正算法记作TR-NC。如图3所示，利用NC指标重建的图像能很好地反映伪迹的位置和大小，同时TR-NC算法能够有效地减少重构图像伪迹，使目标边界更加清晰，提高成像质量。

图 3  伪迹修正流程

下载: 全尺寸图片 幻灯片

TR-NC算法的基本思想：先用TR算法进行一次成像，然后利用f-FCM将1次成像灰度聚成两类，将两类灰度二值化为0和1，计算其NC值。根据NC生成的对角阵D替换TR求解公式中的单位矩阵，即得

根据式(10)求解结果进行成像，实现的具体步骤如表1所示。

表 1  TR-NC的计算流程

输入：边界测量值U，灵敏度矩阵S
输出：灰度向量g
步骤：1．利用TR得1次成像灰度g0；
2．利用f-FCM将灰度g0聚为2类并将灰度二值化为0和1；
3．计算NC值，并据此生成对角矩阵D；
4．将D代入式(12)进行求解，得到灰度值g。

下载: 导出CSV 

| 显示表格

3.   实验结果与分析

3.1   仿真实验分析

为了验证TR-NC算法在解决EIT图像伪迹修正中的有效性和可行性，我们开展了两种类型的仿真实验。一种是验证病变的数量及位置不同时，TR-NC算法能否有效地修正重建图像的伪迹；另一种则是验证在不同电导率的情况下，相比于传统的TR算法，TR-NC算法是否能够更早、更准确地检测出病变的大小、位置及数量。

仿真建模：本文利用COMSOL 3.5进行建模，圆形敏感场域半径设置为16cm，采用16个电极、相邻电流激励-相邻测量的模式，激励电流设为1mA，所测电压数为208。在逆问题求解中，被测物场被剖分成812个网格来重建电导率分布。

3.1.1   不同位置及数量下的病变模型仿真

为了验证病变数量及位置不同时，TR-NC算法能否有效地修正重建图像的伪迹，我们构建了5个仿真模型，分别对应病变位于不同位置(左肺上、下叶，右肺上、中、下叶)，病变数量依次递增。同时，作为性能对比，基于相同的测量数据选取TR算法进行对比参考。借助相关系数(Correlation Coefficient, CC)^[15]和相对误差(Relative Error, RE)^[16]来定量评估f-FCM算法的性能。其中，背景(蓝色)电导率设为1 S/m，肺部(绿色)电导率设为2 S/m，病变目标(红色)电导率设为3 S/m。

重建结果如图4所示。第1行是原始模型，第2行是TR算法得到的重建图像，第3行是TR-NC算法得到的重建图像。从图像质量上可以看出，TR算法虽然可以重建目标的数量，但是存在较大的伪迹，不能精确呈现出目标的真实尺寸及位置。TR-NC算法成像效果较好，在目标数量较少且分布分散时，可以准确地呈现出目标的位置、边界和数量，特别是模型1和模型2；随着目标数量的增加且位置相距较近时，病变部分重建图像的形状与原始模型稍有差异，但与TR算法相比，重建图像伪迹较少，目标区域空间分辨率较高，特别是模型4和模型5。

图 4  5个仿真模型的重建图像

下载: 全尺寸图片 幻灯片

图4的对比结果可以用CC和RE进一步量化解释，如图5所示。与TR算法相比，TR-NC算法所得图像CC值平均提高了17.8%，RE值平均降低了16.8%，图像质量得以改善，说明该算法能够有效地提高重建图像的空间分辨率。

图 5  重建图像评价指标

下载: 全尺寸图片 幻灯片

3.1.2   不同电导率的病变模型仿真

根据分子病理学的发展机制，肺组织由正常向癌变发展的病程是一个逐渐变化的过程。课题组前期研究也表明，肺癌组织的病变过程与组织电导率成正相关关系，即当肺组织发生癌变后，其电导率通常会增加。

为了验证在不同电导率的情况下，相比于传统的TR算法，TR-NC算法能否更早、更准确地检测到病变的大小、位置及数量，即当组织电导率变化多大时，病变能够被准确检测出来，从而有利于尽早进行针对性治疗。对此我们通过改变肺部内某些像素块(红色标记)的电导率来模拟较早期的病变组织进行仿真试验。

为了测试病变组织尺寸大小、位置等因素对重建结果的影响，我们构建了病变位置及尺寸有差异的仿真模型，分别如图6所示。其中，病变组织对应图6的红色部分，其电导率变化增量为Δ σ，肺部组织对应于图中绿色部分，整个人体胸腔对应于图中蓝色部分。成像算法采用TR算法与TR-NC算法进行对比，第1行对应于不同的算法，第2～5行分别对应于Δσ取0.1%，0.3%，1% 和3% 时的重建图像。

图 6  原始图像

下载: 全尺寸图片 幻灯片

从图7的仿真结果可以看出，两个病变位于两个肺部组织的中心位置，且病变尺寸不同。基于相同的电导率及其他测量数据，在各种Δσ值下，两个算法都能准确地识别出病变组织的位置。但由于两个病变的尺寸不同，相应的成像效果也存在较大差异，左侧病变尺寸较大，成像质量较好；而右侧病变尺寸较小，当Δσ= 0.3%时，TR算法不能够较为准确地呈现病变组织的位置及尺寸，相比较而言，TR-NC算法可以较为准确地判断病变组织的位置和大小，且边界清晰，只是尺寸有些许偏差。随着电导率增量的不断增加，当Δσ≥1%时，TR算法所呈现的图像位置及尺寸较好，但仍然存在较大的伪迹，而TR-NC算法已经能够准确检测出病变的位置及尺寸，成像效果大为改观。从算法角度分析，TR-NC算法比TR算法有更高的空间分辨率，不论是对病变组织的位置还是形态大小的判断都要优于后者。

图 7  重建图像

下载: 全尺寸图片 幻灯片

图7的对比结果通过CC和RE进一步量化，如图8所示。从两个定量评价指标的直方图上可以更清晰地看出不同算法之间的性能差异。当Δ σ = 3%时，TR-NC算法的相关系数CC值由TR算法的0.31增大到0.38，相对误差RE值由0.29降低到0.20，验证了TR-NC算法在重建图像的精确度方面相比于TR算法得到了显著提高。

图 8  重建图像评价指标

下载: 全尺寸图片 幻灯片

3.2   仿真体验实验分析

为了验证本文所提算法在医学应用中的可行性，选用16电极EIT系统进行实际实验，进一步评估该算法的性能，采用相邻激励和相邻测量的模式，为了在胸腔模型表面放置电极且保持与人体测量的一致性，电极采用心电电极，该电极具有贴附紧密、导电均匀的优势。

实验容器和测试目标采用3D技术打印的导电实体模型，测量容器的内径为160 mm，每个电极直径为5 mm。为了模拟背景脂肪的电导率，采用向自来水中添加盐的方法调制盐水，利用电导率仪不断测量盐水电导率，在逐渐加盐的过程中充分搅拌并读取示数，当达到脂肪电导率0.037 S/m时停止该操作。肺部实体模型的电导率约为0.1 S/m，左肺中有一个直径10 mm的空心圆柱，通过注水孔向空心圆柱内注入不同电导率的盐水来模拟病变组织。空场时，空心圆柱内注入盐水的电导率为0.1 S/m，即和肺部正常组织电导率一致；满场时，不同电导率的盐水(σ S/m)被注入到空心圆柱内。肺部病变组织相比正常组织的电导率变化率为：α=(σ–0.1)/0.1 。

由于统计误差和系统误差的存在，EIT成像技术在实际应用中的测量数据可能与仿真数据存在较大差异。我们测试了5组不同α值的静态实验。成像过程中基于相同的测量数据，采用TR算法作为参考对比。

实验结果图9所示，第1行是不同的α值，第2行是TR算法对应不同α值下的重建图像，第3行是TR-NC算法对应不同α值下的重建图像。从重建结果可以看出，当α较小时(3% , 6%)，病变组织能被TR算法检测到，但是重建目标不是很清晰，TR-NC算法的修正效果较差。当α较大(大于10%)时，TR算法已经能够清晰地检测出病变组织，但伪迹较大，边界较为模糊；TR-NC算法能够清晰地检测出病变组织的位置及尺寸，目标边界清晰，伪迹较小。随着α的增大，TR-NC算法的重建目标的位置和大小更接近于真实情况。

图 9  不同α值对应的TR-NC重建图像

下载: 全尺寸图片 幻灯片

由此可见，静态实验和仿真实验的结果一致，基于快速模糊C均值和邻域信息的TR-NC算法具有较好的伪迹修正能力和较高的空间分辨率，在肺癌等疾病的早期诊断检测中有很大潜力。

4.   结论

为了解决电阻抗成像技术的图像重建过程中存在的空间分辨率较低的问题，提出了基于快速模糊C均值聚类和空间邻域信息的无监督图像质量评价指标，并利用该指标设计的修正策略对TR算法进行修正。通过建立不同位置、数量和大小的病变模型，我们开展了两类仿真实验，重建结果都表明TR-NC算法能有效减少重建图像的伪迹，使得病变目标的位置更加精确，边界更加清晰，相比于TR算法，该算法能够更早、更准确的检测到病变。仿体实验结果也表明TR-NC算法具有良好的图像重构质量及伪迹修正能力，在医学领域具有较高的实用性及优越性。