CT约束的电阻抗医疗成像方法

摘要

一种CT约束的电阻抗医疗成像方法，包括：(1)在待测区域布设电极，所述电极包括注入电极和测量电极；(2)向注入电极注入微电流，通过测量电极获取测量电压，与此同时，对待测部位进行CT影像扫描，获取待测部位CT值；(3)按照经验，建立一个应用于成像的初始模型；(4)执行CT约束下的电阻抗联合反演过程；(5)获取电压残差曲线，判断电压残差曲线是否收敛，如果不收敛，将新的电阻抗模型定义为初始模型，计算电压残差，计算交叉梯度函数，再次进行迭代得到新的电阻抗模型；如果收敛，则输出当前反演所得电阻抗模型，作为最终反演结果。本发明所述方法能显著提高传统EIT成像的空间分辨率。

说明书

技术领域

本发明涉及电阻抗成像技术以及有限元正演建模技术，具体涉及一种用CT约束的电阻抗医疗成像方法。

背景技术

随着科学技术的发展，计算机断层成像技术(Computed Tomography， CT)为疾病的诊断提供了一种客观、准确、迅速的手段。但是传统的CT 是一种医学结构成像，在一些病症的诊断上存在局限性。近年来，又出现了一种基于生物组织电学特性的电阻抗断层成像技术(Electrical Impedance Tomography，EIT)。EIT是根据物体内部组织电特性参数(如电阻率、电容率)的不同，通过对其表面施加安全激励电流或电压，同时测量物体表面的电压或电流信号来获知物体内部电特性参数的分布，进而重建出反映物体内部结构的图像。这对生物体内部电特性的研究具有重要意义，因为不同组织和器官的电特性不同。这种图像不仅包含了丰富的解剖学信息，而且可以获得某些组织和器官的电特性随其病理、生理功能状态而改变的信息。

某些人体组织的生理功能变化能引起组织阻抗的变化，如组织充血和放电等，脂肪过多，某些组织病理改变也能引起组织阻抗的变化(如早期白血病显示骨髓导电率先发生变化等)。当加上正弦激励后，这些信息必然会改变点位参数的分布，这些信息将会在EIT图像中体现出来，所以 EIT具有功能成像的性质。利用EIT技术，可以显示人体内组织的阻抗分布图像、人体组织随频率变化图像、人体器官进行生理活动(如呼吸、心脏搏动)时的阻抗变化图像等，在临床上可用于检测和监护。该技术具有无损伤，低成本，操作简单和信息丰富等特点，近十几年受到国际学术界的广泛关注，并呈现出很好的应用前景。相对于CT，EIT是一种功能成像，对人体器官病变的诊断具有一定的优越性。但是，传统的EIT存在着空间分辨率低的问题。

发明内容

为了解决现有技术中存在的上述问题，本发明提出了一种高分辨率电阻抗成像方法，应用交叉梯度算子，将同一部位的CT成像的结构信息作为约束加入到EIT成像中，从而提高EIT成像的空间分辨率。

为了达到以上目的，本发明提出了一种CT约束的电阻抗医疗成像方法，包括：

(1)在待测区域布设电极，所述电极包括注入电极和测量电极；

(2)向注入电极注入微电流，通过测量电极获取测量电压，与此同时，对待测部位进行CT影像扫描，获取待测部位CT值；

(3)按照经验，建立一个应用于成像的初始模型；

(4)执行CT约束下的电阻抗联合反演过程，所述电阻抗联合反演过程：根据初始模型，通过有限元正演计算理论电压，与测量电压相减获得电压残差；根据初始模型以及CT模型，计算交叉梯度函数；通过高斯牛顿迭代法得到电阻抗残差，将电阻抗残差改正到初始模型得到新的电阻抗模型；

(5)获取电压残差曲线，判断电压残差曲线是否收敛，如果不收敛，将新的电阻抗模型定义为初始模型，计算电压残差，计算交叉梯度函数，再次进行迭代得到新的电阻抗模型；如果收敛，则输出当前反演所得电阻抗模型，作为最终反演结果。

优选地，如果不建立初始模型，初始模型设为默认均匀模型。

优选地，所述电阻抗联合反演过程利用最小化一个目标函数来求解最优的电阻抗分布，所述目标函数为：

Φ(m)＝||d_E-S(m_E)||²+||αLm_E||²+||βt(m_E，m_C)||²

其中d_E表示测量电压，S(m_E)表示正演计算得到的理论电压，m_E表示阻抗分布，m_C表示CT值分布，L表示正则化算子，α表示正则化权重因子，t(m_E，m_C)表示交叉梯度函数，β表示交叉梯度权重因子。

优选地，交叉梯度函数的定义为一个给定三角形单元的重心T(x₀，y₀)>1，y₁)、B(x₂，y₂)、C(x₃，y₃)分别为三个三角形的重心，交叉梯度函数t(m_E，m_C)的函数表现形式为：

其中，

其中m_E1、m_E2、m_E3分别为所述三个相邻三角形中的阻抗分布，m_C1、>C2、m_C3分别为所述三个相邻三角形中的CT值分布。

优选地，测量电压的步骤包括：在同一部位设计多点分别输入微量直流电流，在每一次输入电流的同时，测量该部位的多点电压值，改变电流输入点，直至测量充分覆盖该目标区域。

优选地，步骤(1)中，注入电极以及测量电极的布设遵从对向驱动模式。

优选地，步骤(2)中，设计多点分别输入微量直流电流，对于每一次输入电流的同时，测量该部位的多点电压值，改变电流输入点，直至测量充分覆盖该目标区域。

优选地，步骤(2)中，对于同一测量部位，在不同时间段重复获取相同点的测量电压。

优选地，步骤(3)中，所述初始模型通过参考CT影像扫描信息来建立。

优选地，步骤(4)中，所述交叉梯度函数采用在有限元下的表现形式。

与现有技术相比，本发明具有以下有益效果：

本发明利用交叉梯度函数开发出一种CT约束EIT成像方法。由于该方法结合了CT空间结构信息，所以能显著提高传统EIT成像的空间分辨率。此外，本发明将交叉梯度法推广到了有限元法中，这将有可能运用到基于有限元理论的各种成像方法。

附图说明

图1为反演方法流程图。

图2为对向驱动示意图。

图3(a)阻抗-CT值交叉梯度函数在矩形网格中的定义。对于一个二维有限差分网格。m_E表示电阻率，m_C表示CT值，交叉梯度函数t定义为一个给定单元的中心(用X表示)与其右侧单元(标定为r)以及下侧单元(标定为b)之间的关系。(b)阻抗-CT值交叉梯度函数在二维三角网格中的定义。对于一个向l₁，l₂以及l₃方向扩展的二维三角网格。m_E表示电阻率，m_C表示CT值，交叉梯度函数t定义为一个给定三角形单元的重心(用T(x0，y0)表示)与其第一条边相邻三角形(标定为1)以及第二条边相邻三角形(标定为2)以及第三条边相邻三角形(标定为3)的对应关系。A(x1，y1)、B(x2，y2)、C(x3，y3)分别为三个三角形的重心。

图4为一个简单的块状分布模型：(a)EIT模型；(b)CT模型。

图5为不同参数α选择下标准EIT成像方法的重建图像：(a)α＝1； (b)α＝0.1；(c)α＝0.05；(d)α＝0.01。

图6为不同参数α和β选择下CT约束EIT成像方法的重建图像：(a) α＝0.05，β＝1；(b)α＝0.05，β＝10；(c)α＝0.05，β＝100；(d)α＝1，β＝50； (e)α＝0.1，β＝50；(f)α＝0.01，β＝50。

图7为误差收敛曲线：(a)标准EIT成像和CT约束EIT成像的电压误差随迭代次数变化曲线；(b)标准EIT成像和CT约束EIT成像的电阻率模型误差随迭代次数变化曲线。

图8为电阻率剖面：(a)水平电阻率剖面y＝0；(b)垂直电阻率剖面x＝0。

图9为CT约束初始模型对应的标准EIT成像结果。(a)初始模型；(b) 标准EIT重建图像；(c)电阻率误差收敛曲线。

图10为CT约束初始模型对应的标准EIT成像结果。(a)初始模型； (b)标准EIT重建图像；(c)电阻率误差收敛曲线。

图11为一个二维胸腔模型：(a)原始CT图像；(b)离散化的CT 图像，X和Y方向的坐标单位为cm；(c)根据离散CT图像建立的离散 EIT图像，X和Y方向的坐标单位为cm。

图12为成像结果：(a)原始CT图像；(b)离散化的CT图像，X和 Y方向的坐标单位为cm；(c)根据离散CT图像建立的离散EIT图像，X 和Y方向的坐标单位为cm。

图13为误差收敛曲线：(a)标准EIT成像和CT约束EIT成像的电压误差随迭代次数变化曲线；(b)标准EIT成像和CT约束EIT成像的电阻率模型误差随迭代次数变化曲线。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明作进一步的详细说明。

交叉梯度法最初被提出用来解决地球物理问题，交叉梯度法利用两个或多个不同物理量梯度的叉乘来实现对多种不同地球物理参数的约束。在本发明中，将交叉梯度法应用于EIT成像技术，将CT成像信息与EIT成像信息通过交叉梯度函数联系起来，利用联合约束提高传统EIT成像的空间分辨率。该本发明提出的一种CT约束的电阻抗医疗成像方法以图1的方法流程图进行展示。所述方法包括：

在待测区域布设电极，包括：注入电极和测量电极；

向注入电极注入微电流，通过测量电极获取测量电压。与此同时，对待测部位进行CT影像扫描，获取待测部位CT值；

按照经验，建立一个应用于成像的初始模型；

执行CT约束下的电阻抗联合反演过程，所述电阻抗联合反演过程：根据初始模型，通过有限元正演计算理论电压，与测量电压相减获得电压残差，根据初始模型以及CT模型，计算交叉梯度函数。通过高斯牛顿迭代法得到电阻抗残差，将电阻抗残差改正到初始模型得到新的电阻抗模型；

获取电压残差曲线，判断电压残差曲线是否收敛，如果不收敛，将新的电阻抗模型定义为初始模型，进行迭代；如果收敛，则输出当前反演所得电阻抗模型，作为最终反演结果。

下面，本发明将具体阐述该方法。

A.数据获取

CT数据的获取：通过对待测区域进行CT扫描获得CT数据。

电阻抗数据的获取：对向驱动模式

对向驱动以相向的一对电极作为驱动电极，是一种比较成熟的驱动模式。对向驱动能对中心区域提供较大的电流密度、提高信噪比和降低测量信号的动态范围。通过测量电极获得待测区域测量电压数据。

B.正演过程

正演过程是电阻抗成像的基础问题，正演过程是在已知待测区域电阻抗分布的情况下，通过注入一定量的电流，进行数值模拟得到测量电极的测量电压。通常应用于理论电压的求解。

场域电位分布函数φ与电阻抗分布γ满足拉普拉斯方程

其边界条件为

φ(x，y)＝f(x，y)，(x，y)∈Ω(B-2)

其中，表示场域Ω的边界；f表示已知边界电位；j表示流入场域Ω的电流密度；n表示场域Ω的外法向向量。

通过有限元法，(B-1)～(B-3)可以转化为

KΦ＝I(B-4)

K定义为系数矩阵，和阻抗分布密切相关，Φ表示节点电势矢量，I表示注入电流矢量。

方程(B-4)揭示了阻抗分布，节点电势与注入电流之间的定量关系，也就是说，在已知阻抗分布以及注入电流的情况下，通过(B-4)即可求解节点电势，进而获得测量电压。

C.反演过程

电阻抗成像的反演过程是指通过测量数据计算得到阻抗分布的过程，这个过程可以等效为最优化一个目标函数。现在本发明提出一个函数，利用最小化目标函数来求解最优的电阻抗分布。

Φ(m)＝||d_E-S(m_E)||²+||αLm_E||²+||βt(m_E，m_C)||²(C-1)

其中d_E表示测量电压，S(m_E)表示正演计算得到的理论电压，m_E表示阻抗分布，m_C表示CT值分布，L表示正则化算子，α表示正则化权重因子，t(m_E，m_C)表示交叉梯度函数，β表示交叉梯度权重因子。

方程(C-1)指示了本发明设计的目标函数，通过高斯牛顿法，最小化目标函数(C-1)可以转化为求解：

其中，

通过方程(C-2)可以计算出Δm，通过方程(C-4)可以求出新的电阻抗分布：

m_E＝m_E+Δm(C-4)

为了确保点阻抗分布和CT值分布之间的结构一致性，在目标函数中引入了交叉梯度函数t(m_E，m_C)来进行约束反演，交叉梯度函数定义如下：

本发明使用交叉梯度函数的理念是，即使在不同的医学成像技术使用的不同模型参数物理性质之间，没有定量的解析关系的情况下，仍然预期它们所提供的模型反演结果具有一定程度的结构相似性。但本发明不强调它们的数值之间有任何联系。使用这种交叉梯度函数，如果t(A(x，y)，B(x，y)) 在任何地方为零，则代表物性A和物性B的模型参数在结构上是相同的，即具有相同的结构变化边界。因为它意味着物性A和物性B变化梯度的方向完全相同或相反。

D.交叉梯度函数在有限元法下的表现形式

交叉梯度函数在二维下的函数表现形式为：

将它展开：

在解决地球物理问题时，方程(D-2)可以通过有限差分方法来近似表示为：

这里Δx，Δy，m_Ec，m_Eb等变量的定义如图3（a）所示。

本发明将交叉梯度函数推广至有限元方法中，方程(D-2)在有限元模型中的定义如图3（b）所示。

根据图3（b），m_E的梯度可以被定义为如下的形式：

l₁，l₂，l₃分别是三个方向的单位向量(在图3(b)中展示)，其中：

i，j是x，y方向的单位向量。

同样的，

m_E和m_C的梯度同样可以被写成

对比方程(D-11)和方程(D-13)，

对比方程(D-12)和方程(D-14)，

将方程(D-15)～方程(D-18)带入方程(D-2)，可以得到交叉梯度函数在有限元下的离散表达形式。

实施例1

在该实施例中，建立一个简单块状分布来检测本发明的方法。模拟的电阻率和CT值分布如图4所示。背景电阻率和CT值为ρ＝5(任意单位) 和s＝20(任意单位)；模拟“肺部”为ρ＝10(任意单位)和s＝-600(任意单位)；模拟“心脏”为ρ＝1(任意单位)和s＝50(任意单位)。整个模型被划分为1385个三角网格，节点数为N＝737，电极个数为L＝16。采用一个阻抗为ρ＝5的均匀模型作为初始模型执行标准的EIT成像以及 CT约束下的EIT成像。

合适的正则化参数对获取一个准确的成像结果十分重要。正则化参数α控制着数据误差和正则化参数的权重。在此次实验中，测试了一系列的α。图5描述了标准EIT的重建图像，其中α＝1，0.1，0.05，0.01。标准的 EIT成像可以重建一个可靠的阻抗图像，但是分辨率较低。

为了更好的理解参数α和β在CT约束下的EIT成像中的作用，采用固定一个参数改变另外一个参数的方法。固定参数α改变参数β，在图6（a）， 6（b），6（c）中α＝0.05，β分别为1，10，100。可以发现，在α固定的情况下，β越大，阻抗图像中各个“器官”的边界越明显。当然，需要指出的是，β如果选取太大，会导致成像过程极度不稳定。固定参数β改变参数α，在5d， 5e，5f中β＝50，α分别为1，0.1，0.01。可以发现，在β固定的情况下，α越大，图像的分辨率就越低。当然，需要指出的是，α如果选取太小，同样会造成成像过程极度不稳定。

为了进一步验证和对比实验结果，计算了每次迭代的数据(测量电压) 误差和模型(重建电阻率)误差。图7（a）和图7（b）指示了标准EIT成像和 CT约束EIT成像方法的数据误差和模型误差收敛曲线。在两种方法的数据误差收敛到同一个水平的前提下，相对于标准EIT成像，CT约束EIT 成像方法的模型误差收敛到了更小的值，这说明CT约束EIT成像方法能提供一个更精确的阻抗分布模型。

为了进一步比较图5和图6的成像质量，绘制了标准EIT成像和CT 约束EIT成像的电阻率剖面图。图8（a）显示了y＝0对应的水平电阻率剖面图。图8（b）显示了x＝0对应的垂直电阻率剖面图。可以发现，标准EIT成像的电阻率剖面是光滑的，相对的，CT约束EIT成像提高了电阻率分布的准确度。

在迭代反演过程中，初始模型的选取可能直接影响着最后的成像结果，本发明尝试将CT信息作为先验信息加入到初始模型来测试标准EIT成像方法是否能得到高分辨率的成像结果。图9（a）、10（a）描述了用作标准EIT成像方法的初始模型；图9（b）、10（b）描述了标准EIT成像方法的重建图像；图 9（c）、10（c）描述了电阻率误差收敛曲线。从上述结果中可以发现，在CT约束初始模型的前提下，标准EIT成像并不能提供一个更高分辨率的成像结果。

实施例2

在该实施例中，建立一个二维胸腔模型来检测本发明的方法。模拟的电阻率和CT值分布如图11所示。一个轴向的CT影像图片(图11（a）)被分为三种类型的组织：肺，心脏，软组织(图11（b）)。其中对应的CT值为：肺部为-700HU，心脏为50HU，软组织为40HU。根据CT影像建立对应的电阻率模型：肺部为6.6Ω.m，心脏为1.43Ω·m，软组织为3Ω·m(图 10（c）)。整个模型被划分为1770个三角网格，节点数为N＝935，电极个数为L＝32。采用一个阻抗为ρ＝1.43Ω·m的均匀模型作为初始模型执行标准的EIT成像以及CT约束下的EIT成像。为了更加接近实际情况，向测量电压加入了信噪比为80dB的高斯白噪声。

图12描述了不同成像方法的重建图像。图12（a）展示了标准EIT成像结果。标准EIT成像方法可以重建一个可靠的电阻抗图像，但是分辨率不足。图12（b）展示了CT约束EIT成像结果，可以发现，在CT信息的约束下，肺部的轮廓更加清晰。

为了进一步验证和对比实验结果，计算了每次迭代的数据(测量电压) 误差和模型(重建电阻率)误差。图13（a）和图13（b）指示了标准EIT成像和 CT约束EIT成像方法的数据误差和模型误差收敛曲线。在两种方法的数据误差收敛到同一个水平的前提下，相对于标准EIT成像，CT约束EIT 成像方法的模型误差收敛到了更小的值，这说明CT约束EIT成像方法能提供一个更精确的阻抗分布模型。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。