一种PET成像中卡尔曼滤波图像重建方法

摘要

本发明的PET成像中卡尔曼滤波图像重建方法，通过PET正电子发射断层扫描仪得到原始投影线的正弦图，然后建立状态空间体系，通过基于状态空间的卡尔曼滤波法得出放射性活度分布，重建图像。由于本发明利用基于状态空间体系的卡尔曼滤波法重建PET图像，有效地提高了重建图像的质量；通过与现有重建方法的实验比较，重建结果的偏差和方差都有了很大地提高。

说明书

技术领域

本发明涉及一种医学断层影像重建方法，尤其是涉及一种在PET成像中的卡尔曼滤波图像重建方法。

背景技术

正电子发射断层成像(positron emission tomography，PET)作为一种生物医学研究技术和临床诊断手段是核医学成像装置中最为重要的应用之一。PET的基本思想：向生物体内部注射正电子同位素标记的化合物，在体外通过图像重建技术测量它们的空间及时间分布。如今许多临床医学领域已经开始广泛使用PET图像进行肿瘤、心脏疾病(心肌存活的检测、心脏移植、介入治疗前后监测等)、神经和精神系统疾病(如癫痫病灶的确定、帕金森病、老年痴呆和精神病等)的诊断，部分发达国家已经把部分疾病的PET检查列入到医疗保险范围。

PET探测人体内发出的放射性信号，经过符合和采集系统处理，形成投影线，并以正弦图(Sinogram图)方式存放于计算机硬盘中。计算机以Sinogram图为输入，调用图像重建模块，计算得出人体横切断层图像，反映体内某组织器官的新陈代谢情况，其实这类重建逆问题都可以用方程Y(t)＝DX(t)+v(t)来表示，Y为PET的测量数据，X是要重构的图像，D是系统模型。基于含有噪声的数据Y估计X的方法，可以划分为两类：解析法和统计迭代法。相对于解析法而言，统计迭代算法包括基于泊松模型的最大似然法、最大后验概率法和基于高斯模型的最小二乘法，是以统计规律为基础，对不完全数据的适应性好，能够得到更精确的结果，也因此受到研究人员的广泛关注。

迭代法是从一幅假设的初始图像出发，采用逐步逼近的方法，将理论投影值同实际测量投影值进行比较，在某种最优化准则指导下寻找最优解。在PET中常用的迭代法包括最大似然法(Maximum LikelihoodExpectation-Maximization，MLEM)和OSEM(有序子集最大似然法，Ordered Subset Expectation Maximization)算法。迭代法优点之一是可以根据具体成像条件引入与空间几何有关的或与测量值大小有关的约束条件，如可进行空间分辨不均匀性的校正，物体几何形状约束，平滑性约束等控制迭代的操作，在某些场合下，比如在相对欠采样、低计数的核医学成像中可发挥其高分辨的优势。迭代法最大的缺点是计算量大，计算速度慢。

迭代法求解过程是：

a.假定一初始图像；

b.计算该图像投影；

c.同测量投影值对比；

d.计算校正系数并更新初始图像值；

e.满足停步规则时，迭代中止，否则以新的重建图像作为初始图像从b步开始。

OSEM是近年来发展完善的快速迭代重建算法，每次重建时只使用一个子集对投影数据进行校正，重建图像更新一次，这样所有的子集都对投影数据校正一次，称为一次迭代。和传统的迭代算法MLEM相比，大大加快了图像重建速度，缩短了重建时间。

发明内容

本发明提供一种正电子断层扫描中的基于状态空间理论的卡尔曼滤波图像重建方法，大幅度提高PET重建图像的质量。

本发明的卡尔曼滤波图像重建方法，包括以下步骤：

(1)输入原始投影线的正弦图；

Y(t)=DX(t)+v(t)

(2)建立状态空间体系

X(t+1)=AX(t)+w(t)，

其中，t表示时间；Y就是正弦图数据；D为系统矩阵，表示发射的光子被探测器接收到的概率；A是状态转移矩阵，在瞬时稳定状态下是一个单位矩阵；X为放射性浓度状态变量，为需要重建的对象；w是过程噪声，服从正态高斯分布ω(t)～N(0，Q(t))；v为测量噪声是各种噪声的集合，服从正态高斯分布v(t)～N(0，R(t))；

(3)利用卡尔曼滤波，根据下列方程得出放射性活度分布，重建图像；

X(t+1，t)＝AX(t)                                ①

P(t+1，t)＝AP(t)A^T+Q(t)                     ②

K(t+1)＝P(t+1，t)D^T[DP(t+1，t)D^T+R(t)]^-1  ③

$\tilde{Z}(t+1)=Z(t+1)-\mathrm{DX}(t+1,t)$④

$X(t+1)=X(t+1,t)+K(t+1)\tilde{Z}(t+1)$⑤

P(t+1)＝[I-K(t+1)D]P(t+1，t)                           ⑥

其中，Z(t)为测量值，X₀为浓度初始值，P₀为初始浓度误差协方差，迭代从初始X₀，P₀，Q₀，R₀出发，通过测量值Z(t)，不断修正估计出来的浓度值X，最终给出放射性活度分布，重建图像。

上述的正弦图(即Sinogram图)通过PET正电子发射断层扫描仪得到。

PET正电子发射断层扫描仪进行透射扫描和发射扫描，透射扫描得到图像的衰减校正系数，发射扫描在180度内选取48个采样角度和34个径向单元。

在上述步骤(3)中，具体地迭代重建包括以下步骤：

1)首先设定状态的初始值和初始协方差X₀，P₀；

2)利用时间更新方程方程①和②及时向前推算出当前状态和误差协方差的先验估计值；

3)利用离散采集时间点测得的观测值Z(t)，根据状态更新方程⑤、⑥来更新先验值，得到这一时刻的最优估计；

4)重复交替步骤2)和步骤3)直至获得最优重建结果。

本发明的优点是：

通过基于状态空间体系的卡尔曼滤波法重建PET图像，有效地提高了重建图像的质量；通过与现有重建方法的实验比较，重建结果的偏差和方差都有了很大地提高。

附图说明

附图1是扫描图像数据；

附图2a是MLEM法对附图1的重建结果示意图；

附图2b是本发明对附图1的重建结果示意图；

附图3是模型中间一列像素对应真值的放射性浓度分布值示意图；

附图4是体模示意图；

附图5是MLEM法对附图4的重建结果示意图；

附图6是本发明对附图4的重建结果示意图。

具体实施方式

正电子发射断层扫描仪探测人体内发出的放射性信号，经过符合和采集系统处理，形成投影线，并以Sinogram图(正弦图)方式存放于计算机硬盘中。计算机以Sinogram图为输入，调用图像重建模块，计算得出人体横切断层图像。

在应用PET正电子发射断层扫描仪时，进行透射扫描和发射扫描。透射扫描得到图像的衰减校正系数。发射扫描在180度内选取48个采样角度和34个径向单元，系统矩阵由Fessler教授研究小组提供的软件包生成。

对原始采集到的Sinogram图数据进行各类校正，从而建立一个离散的测量方程：

Y(t)＝DX(t)+v(t)

其中，Y是Sinogram图数据，D是系统矩阵，v表示观测噪声。

建立发展方程：人体内的同位素分布在成像过程中暂时稳定时有X(t)＝⁰，通过对发展方程的离散化，并且加入过程噪声，可得到更一般的表达：

X(t+1)＝AX(t)+w(t)

其中，A是状态转移矩阵，在同位素分布暂时稳定时是一个单位矩阵。

构建状态空间体系：把离散的测量方程和发展方程联合起来，建立一个状态空间体系：

Y(t)＝DX(t)+v(t)

X(t+1)＝AX(t)+w(t)

其中，v服从正态高斯分布v(t)～N(0，R(t))，w服从正态高斯分布ω(t)～N(0，Q(t))。

利用卡尔曼滤波，根据下列方程重建图像：

X(t+1，t)＝AX(t)①

P(t+1，t)＝AP(t)A^T+Q(t)②

K(t+1)＝P(t+1，t)D^T[DP(t+1，t)D^T+R(t)]^-1

$\tilde{Z}(t+1)=Z(t+1)-\mathrm{DX}(t+1,t)$④

$X(t+1)=X(t+1,t)+K(t+1)\tilde{Z}(t+1)$⑤

P(t+1)＝[I-K(t+1)D]P(t+1，t)⑥

其中，Z(t)为测量值，X₀为浓度初始值，P₀为初始浓度误差协方差，迭代从初始X₀，P₀，Q₀，R₀出发，通过测量值Z(t)，不断修正估计出来的浓度值X，最终得出放射性活度分布，重建图像。

采用Kalman滤波进行图像重建时，主要分成四个步骤：

1)首先设定状态的初始值和初始协方差X₀，P₀；

2)利用时间更新方程①和②及时向前推算出当前状态和误差协方差的先验估计值；

3)利用离散采集时间点测得的观测值Z(t)，根据状态更新方程⑤、⑥来更新先验值，得到这一时刻的最优估计；

4)重复交替步骤2)和步骤3)直至获得最优重建结果。

本发明技术的实验结果如下：

应用本发明技术进行计算机模拟实验，且和MLEM法的重建结果作比较。采用Zubal胸腔体模合成发射扫描数据，如图1所示。图像的原始分辨率为32×32像素，正弦图模拟生成过程中采用180度旋转角度下48个采样角度和34个径向采样的投影数据。为了模拟真实的采集条件，S inogram数据中加入了10％的随机噪声，以此作为PET探测到的数据。使用MLEM法和本发明的方法来进行图像重建。实验利用J.A.Fessler教授研究小组提供的MATLAB工具箱生成系统矩阵。

图2a和图2b分别是用MLEM法和本发明的方法得到的重建结果，模型中间一列像素对应的真值及不同方法重建的放射性浓度分布值如图3所示。进一步我们计算了两种算法的重建结果的偏差和方差。MLEM法的偏差达到0.1320，而本发明的偏差为0.1054；MLEM法的方差达到0.0468，而本发明的方差为0.0254。可以看出两种算法基本上都能很好地还原出原来的图像，各个区域的轮廓都很清晰，中间列像素的重建值与真实曲线基本符合。但是能从本发明重建的图像2b中仍能清楚地分辩出中间两个小黑点，而MLEM的结果太过平滑以至于不能分辨出这两个黑点。重建结果的偏差和误差也表明本发明优于传统的MLEM法。

为进一步说明本发明提出的在PET成像中卡尔曼滤波图像重建方法的有效性，我们给出应用不同的算法对真实体模在临床PET系统采集得到的数据进行重建研究的结果。图4是所用的体模示意图，在六个小球体内关注F18溶液，周围充满水，一共采集25帧，分别为5×10s，5×30s，5×60s，5×120s，和5×180s。图5a、5b、5c、5d是利用MLEM法对22，23，24，25帧数据的重建结果；图6a、6b、6c、6d是本发明对22，23，24，25帧数据的重建结果。可以看出，本发明在真实实验条件下的重建效果也较好。