PET主流散射校正方法的评估

【作 者】 孙智鹏，李明，马健，马瑾瑾，梁国栋

东软医疗系统股份有限公司，沈阳市，110167

0 引言

PET数据采集过程中通过时间窗和能量窗对探测到的γ光子对进行筛选，满足筛选条件的称为符合事件。符合事件分为随机符合、散射符合和真符合。散射符合和随机符合都是噪声数据，需要估计其数量和分布，从而对数据进行校正。

常见的散射校正方法分为四大类：能窗法[1]、卷积法[2]、物理模型法和AI模型法。前两种方法出现的时间比较早，且都存在比较多的近似。

能窗法通过扫描一个特定放射源，估算不同空间位置上不同能窗内散射符合和真符合的比例，从而估算出每一条LOR上真符合的计数。但是实际上不同位置上的散射与真符合的比例，与放射源和物质的空间分布有关，理论上只有在被扫描物体的放射源分布、物质空间分布、摆位等属性与前面提到的特定放射源的分布等完全一致时，能窗法才是准确的。临床条件下被扫描物体更加复杂，难以确保散射估计的准确性。

卷积法假设空间中每一个点的散射分布都呈高斯分布，所以真符合与高斯核函数的卷积被认为是散射分布。但实际上散射的空间分布比较复杂，和物体的密度分布、形状关系较大，所以卷积法也不够准确。

基于上面提到的问题，以单散射模拟[3]为代表的模型法被提出，通过建立物理模型和统计模型，单散射模拟可以较准确地计算出不同源分布、不同物质分布下的单散射分布。但实际上临床条件下约有15%的散射是二次散射或多散射，所以单散射模拟被推广到了双散射模拟。双散射模拟虽然可以估算出双散射，但是解决不了视野外散射、多散射问题。

此外，还有一些基于人工神经网络的方法。这类方法[4]一般基于卷积神经网络（或者扩展网络，如U-Net等）。使用核素和物质的空间分布信息作为神经网络的输入，前者可以是PET符合数据，也可以是PET未经散射校正的图像；后者指的是衰减分布图或者是该图的投影。网络的输出是物理模型法计算的结果。通过迭代调整神经网络的参数，使得网络的输出结果与预设的结果越来越接近，直到满足收敛条件。但是因为神经网络本身难以解释，且存在训练样本不能保证覆盖所有可能情形，导致某些情况下结果不可预知的问题，所以目前此类方法尚未用于临床，故暂不介绍。

下面按照技术发展的时间顺序，分别介绍单散射模拟、视野外散射处理方法以及完全蒙特卡罗计算法，并进行效果对比。

1 材料与方法

1.1 单散射模拟

单散射模拟（single scatter simulation,SSS）是一种结合了物理模型和统计模型的方法[5-6]。如图1所示，一个正电子湮灭事件在椭圆柱中发生，生成了一对γ射线。其中向下的一条在S点处发生了康普顿散射，损失一部分能量，行进方向也发生了改变，形成了一个散射事件。

图1 散射事件示意Fig.1 Schematic diagram of scattering events

晶体对AB的散射可以用式(1)～式(3)来表示：

式中：VS是某个散射点S对应的体素；是晶体A和晶体B的立体角，表示空间几何上的探测概率；μ是线性衰减系数；σC是康普顿截面，两者之商是单位体积内的原子数（其物理单位是atom/cm3），与发生康普顿作用的概率成正比；是康普顿差分截面，表示发生当前散射角的概率；是探测器效率，是发生散射，光子能量降低后的效率；s为光子传播路径。

SSS的总体思路是，先确定所有可能的散射点，计算每个散射点与任意两个晶体形成折线上的散射值，遍历所有散射点和晶体对的组合后，计算出每个晶体对上的散射值。因为SSS算法计算出的是与真实散射空间等比例的散射分布，所以通常会通过尾部拟合的方式对计算结果进行缩放，得到与实际散射匹配的结果。

在实际使用中，尾部拟合有时会不稳定，特别是当被扫描目标体积较大时，留给拟合用的数据量变少、噪声增加，容易出现欠估计或过估计的情况。图2中展示了数据中的过估计和欠估计情况。其中，实线是真符合和散射的分布图，点状波动较大的虚线是真实的散射分布。三角形、点划线和星号分别是散射过估计、准确估计和欠估计这3种情况。散射估计误差对图像的影响，如图3所示。其中散射过估计会使得一部分真符合被丢弃，造成图像中间部分偏低（见图3(c)）：较低虚线是背景部分的理论值，较高虚线是热区理论值。图像中心部分的背景部分偏低，图像对比度过高；而散射欠估计（见图3(a)）则会导致图像中心的背景部分偏高，图像对比度降低。

图2 散射过估计和欠估计Fig.2 Scatter over-fitting and under-fitting

图3 散射估计误差对图像的影响Fig.3 Effects in PET images caused by scatter estimation error

为了防止这样的情况发生，一种基于蒙特卡罗的散射分数计算方法被引入SSS拟合过程中[7]。该方法通过建立探测器和核素、物质分布的模型，通过光子传播模型，只需要比较小的计算量，就可以准确地计算出数据中散射事件的比例。这样拟合的结果就可以被这个比例系数约束，输出的散射分布估计结果更准确、稳定。

1.2 基于蒙特卡罗方法

蒙特卡罗方法通过对探测器、放射性核素和衰减系数图进行建模或像素化，结合采样技术和物理模型，模拟出正电子湮灭和γ光子传播过程。甚至可以适当地向轴向视野外扩展，将一部分来自视野外的散射符合事件也考虑进来。

蒙特卡罗方法是一类把概率现象作为研究对象的数值模拟方法，适用于对离线系统进行计算仿真。对于某个过程，只要我们理清了其背后的概率特性，数值仿真逻辑并不会非常复杂，但往往伴随着非常大的计算量。对于模拟正电子湮灭现象，主要分为以下几个步骤，其基本流程如图4所示。

可以看到图4中有一个环路，其每个循环代表一个湮灭事件的模拟，直到模拟的湮灭事件达到预设值。

图4 蒙特卡罗方法的基本流程Fig.4 Flow-charts of Monte Carlo method

早在2000年前后，就有学者意识到了基于蒙特卡罗方法是一种准确的散射估计方法[8]，但是受限于当时的硬件水平，其计算耗时无法被临床接受。但随着近20年来硬件算力的快速提升，特别是GPU计算性能的提升，使得基于蒙特卡罗方法越来越实用化[9]。以目前比较新的Nvidia RTX 3080为例，其计算能力可以达到两百万计数每秒的仿真速度或更高。假设PET系统灵敏度为10 cps/kBq，那么单个GPU完成蒙特卡罗计算需要约2200 s。如果使用专业的显卡，这个计算耗时可以缩短到几分钟甚至更短，基本上可以满足临床的要求。

1.3 视野外散射的处理

早在20世纪90年代，3D PET刚刚开始成为主流方向的时候，就有学者注意到了视野外部的放射性核素对于散射的影响[10]。

如1.1节、1.2节中所述，目前主流的散射估计方法都需要核素空间分布信息和衰减分布信息作为输入。受限于PET有限的轴向视野，往往无法获知视野外的信息。即便是多个床位的连续扫描，也不能得到下一个床位所在位置的核素空间分布（衰减分布可知，因为CT扫描可以先于PET）。所以，如何估计视野外部分的核素大致分布是估计视野外散射的关键环节。

ANDREYEV等[11]提出了一种以散射估计结果为目标的方法，其主要算法如图5所示。该方法使用扩展了轴向视野的核素分布图像和衰减分布图像作为单散射模拟的输入，其目标是通过不断更新视野外部分的核素分布图像，使基于该图像计算出的散射分布与实采数据接近。但是该方法存在一个比较明显的问题：估计结果的精度不可控。因为整个图像更新过程随机，无法保证结果收敛；同时受限于该方法的目标ü ü 只关心散射估计的结果，而不关心视野外图像的估计是否准确，所以其精度很难保证。

图5 视野外散射估计方法Fig.5 Estimation of out-of-FOV (field of view) scatter

因为估计散射的时刻，PET视野内、外的衰减信息已知，视野内的核素分布可以估计，所以在实验中可以通过拟合的方式（如多项式拟合、神经网络拟合等）计算出衰减系数与核素之间的大致映射关系，粗略估计出视野外一定距离内的核素分布，如式(4)所示。其中θ是拟合参数，μ和I分别是衰减分布图像和核素分布图像，下标in和out表示轴向视野内或外，F(·)是映射函数，L(·)是似然函数。

在计算出映射函数F的参数之后，我们可以用视野外的衰减分布图像μout作为函数输入，估算出视野外的核素分布图像Iout。这样就可以得到轴向视野内外的衰减分布图像和核素分布图像，然后使用单散射模拟或蒙特卡罗方法进行散射估计。这种方法可能不是最优，但稳定性有保障、计算量不高，比较适合在临床场景中使用。

1.4 结果对比

我们做了一组仿真实验，探测器模型是东软医疗的NeuWise Pro PET/CT产品，其轴向视野约为224 mm。放射源的轴向长度分别为200 mm和400 mm，半径都是100 mm，位于探测器视野中心。从GATE符合数据中取出散射符合，并进行单层重组，按层求和，得到如图6所示的视野外核素对视野内散射分布的影响曲线。

图6 视野外核素对视野内散射分布的影响Fig.6 Effect on scatter distribution from out-of-FOV nuclide

图6中曲线横坐标代表不同的轴向层位置，纵坐标代表每一层上的总散射符合计数。可以看出不同长度的放射源，虽然视野内部分的核素分布完全一致，视野外的核素对视野内的散射分布存在明显的影响。

所以，为了给出一个准确且全面的散射校正效果评估，我们对单散射模拟和蒙特卡罗方法进行了实验，并且这两种方法又分为考虑视野外散射和不考虑视野外散射两种情形。

图7(a)中的实线是真实的散射分布；长虚线和断虚线分别对应单散射模拟和蒙特卡罗方法；星号或三角标记的，表示考虑了视野外的散射。在散射的轴向分布上，如果不考虑视野外的核素衰变产生的散射，那么散射估计的精度会呈现从探测器轴向视野中心到两边递减的趋势。这一现象是由康普顿散射角的概率分布导致的：随着康普顿散射角度的增大，其发生康普顿散射的概率会降低。所以视野外部分的散射事件分布在探测器两端比较多，中心比较少。在图7(b)散射的径向分布对比中，可以看出无论是单散射模拟，还是蒙特卡罗方法，视野外部分的影响都非常小。在相同的视野条件下（即都考虑或都不考虑视野外散射），因为SSS模型是针对单次散射设计的简化模型，而蒙特卡罗方法则没有散射次数的限制，模型与真实情况更接近。故蒙特卡罗方法比SSS方法更准确一些。

图7 不同散射估计方法在径向/轴向上散射分布的估计结果Fig.7 Axial &radial scatter distribution estimation of different methods

在散射估计方法的实际使用中，通常会在散射拟合过程中使用散射分数对散射做归一化处理，所以这里使用KL距离来度量不同散射估计方法的结果SE与真实散射分布ST的差异。

表1是4种散射估计方法与真实散射分布的KL距离，越小越好。

表1 不同散射估计方法与真实散射分布之间的偏差Tab.1 Deviation between different scatter estimation methods and real scatter distribution

2 结果

从实验结果来看，视野外部分的核素产生的散射事件，如果没有被散射估计算法纳入模型，则会引起散射估计出现偏差。因为实验中的仿真模体在轴向上的形状、活性浓度变化不大，所以视野外的核素对散射的轴向分布影响较大，径向影响较小。但是在实际临床扫描条件下，如果视野外部分刚好是活性浓度较高的区域，如脑部、肝部、肾脏和膀胱等，则很可能会导致径向的散射也出现偏差。从理论上来讲，这种偏差会引起图像均匀性、图像定量精度等指标的下降。

此外，更准确的模型，如基于蒙特卡罗的散射估计算法，有助于提高散射分布的估计的精度。蒙特卡罗（含视野外）方法取得了最接近真实散射的结果。

3 结论

笔者简要介绍了PET散射校正的发展历程。通过对目前主流方法的原理进行剖析，引出了主流方法中存在的问题，如SSS尾部拟合稳定性问题、视野外核素散射的影响等。

实际上，目前主流的散射校正还是局限在单散射的处理上，双散射甚至多散射的处理技术尚未铺开。WATSON等[6]提出了将单散射模拟扩展到双散射模拟的方法，该方法可以比较准确地计算出散射次数不大于2的散射事件分布；相比单（双）散射模拟法，蒙特卡罗方法减少了空间抽样、数据近似，同时可以对探测器、物质、核素甚至电子线路进行建模，模型精度更高，并且可以支持大于2次的散射。所以，蒙特卡罗方法是目前业内公认的最接近真实情形的散射计算方法。随着算法的进步和硬件性价比的提升，相信这类方法会出现在商业化的PET产品中。

对于PET/CT或PET/MR来讲，只要下一个床位中存在核素分布，就需要考虑来自视野外的散射。在扫描完当前床位时，通常我们不清楚下一个床位位置的核素分布。如何准确地估计下一个床位位置的核素分布，是实现视野外散射估计的一个重要前提。随着人工智能的快速发展，相关技术已经被引入医疗领域，相信今后会有更优的方法来估算未知的核素分布。此外，超长轴向视野PET的逐步商业化也使得视野外散射的估计更加准确、可信。