论文阅读：迭代CT重建中线性插值模型的比较

A comparison of linear interpolation models for iterative CT reconstruction

目标：
最近的研究表明，基于模型的迭代重建方法可以改善计算机断层扫描(CT)中的图像质量。现在所面临的一个困难是，如何在各类迭代重建方法中根据实际应用进行选择，包括前向投影模型和惩罚项目的选择。
目前，关于这一选择步骤的指导文献相当缺乏，而这些选项将会影响图像质量。在这里，作者研究了三种依赖于线性插值的正向投影模型的优点：距离驱动方法（the distance-driven method）、约瑟夫方法（Joseph’s method）和双线性方法（the bilinear method）。
本文的选择是由三个因素推动的：(1)在CT中，线性插值通常被视为权衡离散化误差和计算成本之间的最合适方法，(2)前两种方法受到制造商的欢迎，(3)第三种方法可以评估其他方法中一个关键假设的重要性。

方法：
评价正投影模型的一种方法是考察它们对离散图像的影响，以及它们的转置对数据集的影响，但这类研究的意义尚不明确，因为在迭代重建中总是联合使用矩阵及其转置。另一种方法是在模型使用的环境中研究它们，即，连同统计权重和惩罚条款。不幸的是，这种方法需要选择一个首选的目标函数，并且不能提供关于模型固有特性的明确信息。作者采用了以下两阶段的方法。首先，作者分析了在没有统计权重和惩罚项的重建图像中逐步包含模型奇异值分解成分的图像。接下来，作者研究了权重和惩罚对观察到的差异的影响。

结果：
我们研究了16种不同的扇形束成像场景的图像质量指标，以探测所有模型的各个方面。该指标包括计算成本的替代，以及偏差、噪声和估计任务，所有这些都具有匹配的分辨率。分析揭示了偏差和噪声方面的根本差异。似乎需要基于任务的评估来欣赏噪音的差异;作者选择的估计任务表明，这些差异抵消了产生相似的性能。一些场景强调了距离驱动方法在偏差方面的优点，但增加了计算成本。统计权重和惩罚项的三种组合表明，观察到的差异保持不变，但强保边惩罚可以显著降低这些差异的量级。

讨论：
在许多情况下，Joseph的方法似乎在成本和计算工作量之间提供了一种有趣的折中。
距离驱动方法提供了减少偏差的可能性，但增加了计算成本。
双线性方法表明，其他两种方法中的一个关键假设具有很强的鲁棒性。

最后，强边缘保持惩罚可以作为前向投影模型的不足的补偿器，使所有模型在最具挑战性的成像场景中达到相似的水平。此外，作者发现，他们的评估方法有助于欣赏模型，统计权重和处罚期限是如何进行相互作用的。

1、介绍

非线性迭代重建方法已成为x射线计算机断层扫描(CT)研究的热点。这一努力的动机是越来越多的证据表明，这种方法可以在同等剂量下改善图像质量或在同等图像质量下节省辐射剂量如文献1-7中的临床研究。与经典的滤波反投影不同，迭代重建方法允许在重建过程中包括额外的信息，如测量的光子统计(参考文献8，第2章)，x射线束的物理性质，9-11和图像约束这些优点都有助于观察到的图像质量改善。
利用半渐近迭代的希望是，在发散发生之前，用相对较少的迭代次数，获得与目标图像很好的近似。所有算法的开发都需要快速实现，根据算法的具体情况使用CPU或GPU。
另一个重要的研究方向是对重建过程中涉及的众多因素对图像质量的依赖进行定量分析，我们的贡献就属于这一类。目前关于这一方向的研究文献相对匮乏，我们的贡献填补了一个重要的知识空白。
考虑到非线性迭代重建方法具有难以置信的大量自由度，分析是特别具有挑战性的。因素列表从选择定义最终图像的目标函数开始，然后迅速扩展到包含该函数中出现的各种参数，以及正向投影模型。这个方向的研究由于需要进行多次重构来解释不同的噪声实现和几何变化而变得更加复杂，这对于有意义的观测是必不可少的。缓解这种困难的一种方法是专注于2D而非3D重建，希望获得的知识至少部分扩展到3D。
在本文中，我们研究了前向投影模型对图像质量的影响。尽管我们只调查了这一个因素，但由可用的的前向预测模型很多，我们不得不做出选择；我们决定重点研究以下三种二维扇束CT几何模型：

Joseph方法，距离驱动方法以及双线性方法

这个决定的主要原因是：
（1）这三种方法都基于线性插值，这通常被视为离散误差和诊断CT图像重建计算工作量之间的合适权衡。依赖于最近邻插值的算法，如流行的Siddon方法，在诊断CT中往往缺乏鲁棒性，除非采取对策，（2）这三种方法都可以很容易地扩展到3DCT几何形状。
（3）这些方法很受制造商的欢迎，尤其是前两种，它们似乎被用于临床产品。
（4）对这些重要方法的结果认识目前相当有限的仅在CT几何诊断。
研究前向投影模型对图像质量的影响并不简单。一个关键的问题是，确定模型是应单独研究，还是应在迭代CT重建的背景下研究，该背景具有特定的目标函数，包括统计权重和惩罚项等典型成分。我们选择了前一个选项，因为我们认为首要目标应该是首先了解模型本身可以提供什么。统计权重和处罚条款的影响在调查的第二阶段可能更容易理解。研究模型本身的价值的一种方法是使用奇异值分解。这样的方法可以提供可恢复特性和相关稳定性方面的完整图景，但这幅图景可能太丰富了，无法完全理解。我们确定了一种中间方法，该方法逐渐将奇异值分解的分量注入到重建图像中，并依赖于经典的图像质量度量。

本文的组织结构如下：
1>第2节：描述了用于比较的扫描仪几何形状，对三种正向投影模型的简要描述，强调它们之间基本差异的简短讨论;
2>第3节：解释了实验比较条件，包括图像的形成过程、几何设置，以及用于测试的幻影；
3>第4节：详细描述了用于定量比较的基本和基于任务的度量标准；
4>第5节：介绍并讨论了使用这些度量标准获得的结果；
5>第6节：研究了图像重建的统计权重和惩罚条款的例子；
6>第7节：总结和结论。

2、背景

2.A 扇形波束几何

数据模拟采用带有曲面探测器的第三代CT扫描仪的几何结构，如图1所示。在这种几何结构中，x射线源与探测器围绕一个固定的中心同步旋转，从而产生扇束数据。

图1 数据采集几何模型

设λ为相对于x轴测量的极角，设eu(λ) =[−sin(λ)，cos(λ)]T和ew(λ) = [cos(λ)，sin(λ)]T为两个正交的单位向量，其中上标T表示转置算子。然后，将圆扫描路径定义为：

其中，R0为扫描半径，即从光源到旋转中心的距离，λs为起始角，λe为结束角，δ(λ)解释了数据采集过程中可能出现的飞行焦点(FFS)特征。
当δ(λ) = 0时，从源a(λ)到探测器上任意点的距离为d，探测器上的每个点使用曲线坐标u确定，它在eu(λ)方向的探测器上被正测量。原点u = 0与源到旋转中心的连线重合;这条线在−ew(λ)方向。设L(λ，u)是从a(λ)开始的，指向探测器上u位置的射线。沿L(λ，u)的测量值为g(λ，u)， L(λ，u)和连接a(λ)到旋转中心的直线之间的角度称为γ。角度γ在逆时针旋转时被测得正，如图1所示，因此u = Dγ。

CT测量也使用平行束坐标进行识别，如图2所示。

图2 平行束模型 ，线L(θ，s)表示与图1中 L(λ，u)相同的测量值

即在这些坐标中，每条射线都用角度θ和带符号的距离s来描述。

：表示射线从源到探测器的方向；
s：表示射线在θ = (cos θ， sinθ)方向上离原点的距离；
(θ，s)与(λ，γ)之间的联系为：

被扫描物体在源轨迹所定区域内以原点为中心的圆盘内。这个磁盘称为视场(FOV)，其半径为rFOV < R0视场内的点由笛卡尔坐标x和y确定。重构是在有限的点集合上进行的，称为(xk, yl)，其中k和l是指标。这些点在x和y中均匀分布，在全文中，x和y的采样距离相同，这个距离记为∆，使xk = k∆，yl = l∆。x中的样本数为Nx, y中的样本数为Ny。每个点(xk, yl)可以关联宽度为∆的正方形像素;这个像素以(xk, yl)为中心。在保持l不变的情况下，改变k得到的样本称为一列样本。同样，保持k不变，改变l得到的样本组成一列样本。
最后，x射线的线性衰减记为f (x, y)。章节2.B - 2.D展示了我们比较的三种线性插值模型中f (xk, yl)和CT测量值之间的联系。

2.B Joseph’s Method

约瑟夫方法将每一个测量作为一个直线积分模型，并沿每一个直线的近似积分可以是x方向的和，也可以是y方向的和，取决于直线的方向。（看射线与探测器点位置的连线是更平行于X轴还是Y轴，可以看斜率）
这个概念允许所有的计算都可以用一个简单的一维插值来完成。当和在x(或者y)上时，插补在y上(分别为x)，采用线性插补，并根据式(2)中的角度θ选择求和方向。当| cos θ| > 1/√2时，直线比x轴更平行于y轴，求和在y上；反之则在x中。

图3 Joseph方法的几何描述，沿着这条线的测量是通过在x上的线性插值和在y上的求和创建的。

图3说明了在y中求和的情况。对于每个L(小)，L(λ，u)和坐标yl之间有一个交点，并且有一个指标k(l)，使这个交点x在xk(l)和xk(l)+1之间。wl是交点和xk(l)在x中的距离。那么，Joseph方法的基本公式可以表示为：

当xk(l)超出可用样本范围时，f (xk(l)，yl) = 0, f (xk(l)+1, yl)也适用此公式。前面的比例因子是雅可比矩阵，它使得和是在y上而不是沿着直线进行的。
当求和值在x而不是y时，公式与式(3)非常相似。即求和值在k（X轴方向）而不是l（Y轴方向），l变成了依赖于k的指标，wl转化为y中的距离hk，求和值前的比例因子为∆/|sinθ|。
在我们对约瑟夫方法的描述中，一个重要的方面是，在任何给定的投影中，求和的方向都可能改变。另一种方法是对视图中的所有光线使用相同的求和方向。由于Joseph的方法实际上只描述了平行波束数据，所以对于扇形束数据，Joseph的方法是什么还没有达成共识。直观地说，如果为每条射线单独选择求和的方向，特别是对于具有高放大系数的几何图形，可以期望获得更好的精度。

2.C 距离驱动方法

距离驱动方法与Joseph的方法非常相似，除了一个关键的方面是，距离驱动法是判断光线路径占检测器像素的大小为多少。在约瑟夫的方法中，每一次测量都是使用x或y的总和创建的，但总和的方向对同一视图内的所有射线是相同的。如果|sinλ| > 1/√2，求和在y中;否则，它在x中。一旦选择了一个求和方向，比如在y中，使用对检测器像素在每一行样本上的投影长度的积分来计算检测器像素的有限大小。如果和的方向为x，则在样本的柱上可以识别出投影的长度。

图4 距离驱动方法|sin(λ)| > |cos(λ)|。(左)测量是通过y中的求和来创建的，每个行位置的贡献是通过对探测器像素投影长度的积分来定义的，表示为[q1, q2]。(右)参与积分的行内样本数量随∆、q1、q2变化。上面的图显示了涉及两个样本的情况的放大视图，而下面的图涉及三个样本。

图4说明了求和为y的情况下的概念。要创建与给定检测器像素相关的测量值，我们首先画两条线，这两条线连接源位置和检测器像素的边缘。对于每一个l，这两条线在x的两个位置与坐标yl相交，称为q1(l)和q2(l)， q1(l) < q2(l)。这两个位置定义了要对其进行集成的探测器像素的投影长度。这个积分在x中使用最近邻插值模型来表示f (x, yl)。结果是：

式中，dk,l为区间[q1(l)，q2(l)]与区间[xk−∆/2,xk +∆/2]重叠的长度。如果没有重叠，那么dk,l = 0。l求和前的比例因子和Joseph的方法是一样的，意义也是一样的。根据定义，在式(4)中对k求和只涉及少数项，但项的确切数目取决于∆，也取决于l、λ、u和探测器像素大小。图4中右侧显示了特定例子，物体的数量如何随着x中的采样而变化。
从实现的观点来看，距离驱动方法的首选方法可能是一次对每个λ和每个l执行所有u值的计算，而不是在对u值进行循环之前对l求和，这种方法在使用缓存内存方面具有优势，因为内存中彼此接近的值可以同时使用。同样的概念也可以用于Joseph的方法。

2.D 双线性方法

就像约瑟夫的方法一样，双线性方法将每个测量模型作为直线积分。然而，它不需要选择首选的笛卡尔方向，这可能对某些光线更好。其概念为：
(1)用双线性插值法估计测量线上每个位置f的值；
(2)整合估计值。双线性方法由Kak和Slaney(参考文献66，第7.4节)提出，他们提出了一种近似实现方法:双线性插值，在沿直线的有限个数点处得到f，然后求和。据我们所知，双线性方法从未在没有近似值的情况下被检验过。
双线性方法中的参数仍然是位置(xk, yl)处f的值。为方便起见，这些值在这里记作ckl，令：

为1阶b样条，设b(x, y) = b1(x) b1(y)，由b(x, y)我们可以把双线性插值得到的f (x, y)的估计f a(x, y)写成

根据公式(2)的平行束坐标，利用f a的这个方程，可以估计出每个线积分为：

采用式(7)，利用r(θ，s)的表达式实现无近似的双线性方法。如图5(左)所示。

图5 (左)双线性方法示意图: 直线被分割成几段，每段位于四个样本之间;利用这些样本进行双线性插值，得到线段上每个位置的函数;然后对各部分的结果进行整合。(右)Joseph方法(灰色点)和双线性方法(灰色和黑点)使用样本的比较。

概念上，r(θ，s)是b(x, y)的拉东变换（Radon transform），这种变换可以表示为：

其中单侧幂函数[x]3+的定义为：


注意，当θ是π/2的倍数时，g(θ，s)降为b1(s)。
对于每一个λ，双线性方法可能是最好的实现使用环在k和l之前。一个重要的计算负担在于r(θ，s)的计算，但使用查找表插值可以大大减少这个负担。然而，我们没有使用这样的表，因为我们的主要目标是用双线性方法评估可能的最佳性能。还要注意的是，使用(参考66，章节7.4)中建议的近似实现可能会快速使用GPU特性。

2.E 比较讨论

这三种方法均可归为线性插值方法。然而，线性插值的应用方式与其他方法有着本质的不同，这会影响迭代重建的图像质量和计算成本。
第一个基本区别是像素值的利用。图6(上)说明了像素大小相对于光线之间的距离较小的情况。这种情况在Joseph的方法中产生了不平衡，即某些像素值与其他像素值相比没有得到充分利用。在图中，可以看到在这个特殊的视图中，像素b和c之间的插值是不需要的。如果∆非常小，或者如果光线更分开，甚至可能有些像素在此视图中的测量建模中不起作用。这个观察结果可能被认为是约瑟夫方法的一个不良特征。考虑到它的定义，双线性方法表现出相同的特征，尽管对给定的像素大小不具有相同的重要性，因为在该方法中，x和y中的像素大小都很重要。另一方面，距离驱动的方法不会受到类似的影响。
第二个基本区别是插值过程中涉及的像素数量。不像其他的有两种方法，在距离驱动的方法中，这个数字会有很大的变化。图6(底部)说明了λ接近π/4的高倍几何。由此可以看出，插补长度q2−q1可以变得相当大。这种不平衡可能会对计算工作以及图像质量产生负面影响。另外，Joseph方法与双线性方法所涉及的像素数也有显著差异，如图5(右)所示。

图6 影响正投影模型的几何特征。
(上图)Joseph的方法可能没有充分利用像素;例如，在这个视图中，像素b和c之间不需要插值。
(下图)插值长度可以变得相当大的距离驱动方法;例如，在本视图中，u为负值时的q2−q1比u为正值时的q2−q1大得多。

第三个基本差异是 |tanλ|接近1（180度）的视图插值的潜在精度在这些视图中，在x或y上求和而不是沿着射线求和的概念显得不够理想。因为它们在x或y上求和，约瑟夫的方法和距离驱动方法可能不如双线性方法精确，双线性方法沿着射线求和。但是请注意，通过允许在视图中从一条射线到另一条射线的求和方向改变，Joseph的方法可能比距离驱动方法更少暴露在这个问题上。

3、实验比较的设置

4、定量比较:指标

5、定量比较:结果

6、统计权重和处罚期限的影响

以上几部分此时对我意义不大，没看。

7、摘要和结论

我们提出了三个前向投影模型之间的比较，依赖于线性插值的迭代重建CT。其中两个模型，即Joseph方法和距离驱动方法，是常用的，第一个是因为它的简单，第二个是因为它能够建模探测器像素大小。第三种方法是首选的，因为它不需要选择首选的笛卡尔方向进行测量，不像其他两种方法。我们的比较包括大量不同的扫描几何形状，以及每个几何形状中相当数量的重复扫描，以评估对光子统计的依赖，这些几何图形是当代CT扫描的代表。结合使用两种不同尺寸的图像像素和具有挑战性的幻影，他们允许对模型的行为进行彻底的调查。

首先，当小像素与测量次数较少的几何图形一起使用时，距离驱动方法比其他模型更稳健。然而，请记住，在这种情况下，距离驱动方法提供的改进可能会带来大量的计算工作量的增加;在我们的研究中，观察到系统基质大小的变化高达50%。另外，我们认为Joseph的方法在RMSE和计算成本之间提供了一个有吸引力的平衡。何时使用距离驱动方法而不是Joseph的方法取决于扫描的解剖、重建和扫描参数以及归因于RMSE和计算成本的重要性。

其次，三种模型的像素噪声和像素之间的相关性存在较大差异，因此需要基于任务的评估对模型进行排序。对于我们选择的估计任务，我们观察到效果的显著平衡，因此所有模型在有足够测量的几何形状中表现相似。

距离驱动的方法可能会比其他模型表现得更好，因为它考虑到了有限的探测器尺寸。在这方面，需要注意的是，本工作中的数据模拟不仅包括检测器模糊;它还包括模糊，由于有限的焦点大小和由于连续的x射线发射。此外，所有模糊效果的建模使用一个现实的非线性线积分组合。因此，从物理建模的观点来看，距离驱动方法不一定是一个更好的模型。

此外，我们还比较了三种前向投影模型在使用三种不同的惩罚加权最小二乘目标函数形成图像时的性能。额外的研究还表明，惩罚项可以显著地补偿一个弱RMSE，所以Joseph的方法在使用边缘保持惩罚时，在所有研究的几何图形的成本和图像质量之间的权衡方面显得具有很强的竞争力。

我们报道的研究为诊断性CT的迭代重建带来了大量关于重要的前向投影模型行为的信息。然而，由于研究是在二维扇形束几何上进行的，因此有一定的局限性。我们预计，大多数观测结果在三维多层CT中都是成立的，就像解析重建方法中的观测结果一般从二维到三维一样。