基于分段反投影张量退化特征编码的牙科锥形束计算机断层扫描运动伪影校正

A segmented backprojection tensor degradation feature encoding model for motion artifacts correction in dental cone beam computed tomography

摘要

目的 为了去除患者在牙科锥形束计算机断层扫描（CBCT）扫描过程中发生躯体运动导致的伪影，提升重建图像质量，提出一种基于分段反投影张量退化特征编码的运动伪影校正模型（SBP-MAC）。 方法 该模型由一个生成器和一个退化编码器构成。将分段有限角度重建的子图像堆叠成张量并作为模型输入；用退化编码器提取张量中空间变化的运动信息，自适应调制生成器的各级跳跃连接特征，从而指导模型校正不同运动波形导致的伪影；最后设计伪影一致性损失来简化生成器的学习任务。 结果 该模型能有效地去除运动伪影，提升重建图像质量。在仿真数据上的峰值信噪比提升了8.28%，结构相似度提升了2.29%，均方根误差降低了23.84%；在真实数据上的专家评分最高4.4221（5分制），与所有对比方法之间的评分差异具有统计学意义（P<0.05）。 结论 本文提出的SBP-MAC模型能够有效提取张量中空间变化的运动信息，实现从张量域到图像域的自适应伪影校正，提升牙科CBCT图像质量。

Abstract

Objective We propose a segmented backprojection tensor degradation feature encoding (SBP-MAC) model for motion artifact correction in dental cone beam computed tomography (CBCT) to improve the quality of the reconstructed images. Methods The proposed motion artifact correction model consists of a generator and a degradation encoder. The segmented limited-angle reconstructed sub-images are stacked into the tensors and used as the model input. A degradation encoder is used to extract spatially varying motion information in the tensor, and the generator's skip connection features are adaptively modulated to guide the model for correcting artifacts caused by different motion waveforms. The artifact consistency loss function was designed to simplify the learning task of the generator. Results The proposed model could effectively remove motion artifacts and improve the quality of the reconstructed images. For simulated data, the proposed model increased the peak signal-to-noise ratio by 8.28%, increased the structural similarity index measurement by 2.29%, and decreased the root mean square error by 23.84%. For real clinical data, the proposed model achieved the highest expert score of 4.4221 (against a 5-point scale), which was significantly higher than those of all the other comparison methods. Conclusion The SBP-MAC model can effectively extract spatially varying motion information in the tensors and achieve adaptive artifact correction from the tensor domain to the image domain to improve the quality of reconstructed dental CBCT images.

Graphical abstract

关键词

运动伪影校正 / 牙科锥形束计算机断层扫描 / 分段反投影张量

Key words

motion artifact correction / dental cone beam computed tomography / segmented backprojection tensor

锥形束计算机断层扫描（CBCT）在牙科诊断和治疗评估中得到了广泛的应用^［1，2］。然而，由于牙科CBCT的扫描时间长达5.4~40 s，因此患者更容易自主或不自主地发生运动^［3］。相当一部分牙科检查存在患者运动，其中18%为儿科患者，24%为老年患者^［4］。这些运动来源包括患者害怕龙门旋转、患有帕金森病的老人，或者患者无法配合（如儿科患者）等^［5-7］。即使这些运动仅在毫米级别，成像系统与投影数据的几何对齐性也会被破坏，导致重建图像出现条纹伪影、双轮廓、解剖结构的模糊和失真^［8］；而且，二次召回重新扫描会增加辐射剂量，危害患者身体健康^［9］。所以，运动伪影校正在牙科CBCT成像发展中至关重要。

近年来，研究人员提出了一系列牙科CBCT运动伪影校正的方法，可以归类为直接抑制运动来源、运动估计与补偿、图像域自聚焦技术，以及基于深度学习的运动伪影校正。

直接抑制运动来源是从根本上避免严重的运动伪影。其中一类方法是使用特定装置来固定患者，例如头带、固定咬合块、下巴托等^［10］。这类方法足以限制患者的剧烈运动，但无法避免轻微运动。另一类方法是将患者摆位在相对稳定的位置。研究表明牙科CBCT坐位扫描对图像质量的退化影响最小^［11］，但这种摆位的抑制效果也有限。

运动估计与补偿通过前瞻地记录或者回顾地估计来获取运动波形，然后将它们整合到系统投影矩阵中，以恢复投影数据和牙科CBCT成像系统之间的几何对齐性。现有基于硬件的方法通常使用外部运动跟踪系统或基准标记来前瞻地记录运动波形，随后进行运动补偿^{［12， 13］}。然而，这些方法不仅增加了额外的成本，而且难以落地牙科检查的临床应用。相反，无硬件的运动估计与补偿才是主流。比如，一些研究人员借用了实时几何校正中的2D/3D配准思想来估计运动波形^［14-17］，但运动估计的准确性和迭代执行的时间需要得到权衡；其他研究人员在运动补偿期间用数据一致性条件来保持投影一致性^［18-20］，但这些手工设计的一致性条件带来了不少设计的难度，且运动估计的准确性有待提升。

图像域自聚焦技术的概念借鉴于光学成像技术^［21］，其使用高效的迭代算法来搜索运动参数或模糊核，以最大化或最小化运动失准伪影度量（MMAMs）。例如，总变分、图像熵或基于图像梯度的指标等^［22-24］。然而，这些MMAMs需要研究人员精心设计，且大多数优化算法通常会陷入局部最优解。

深度神经网络（DNN）能够凭借其优良的特征学习能力弥补上述传统伪影校正方法的不足。其中第1类改进策略是用可学习的2D/3D配准提升运动估计与补偿的精度^［25］。但对于牙科CBCT的运动配准，其先验3D图像含有运动伪影，研究人员往往要用精心设计的正则方法改良先验3D图像^［26］。第2类改进策略是用可学习的MMAMs构造图像域自聚焦方法的目标函数^［27］，在其他条件不变的情况下试图用更好的目标函数获得更优解。但即使网络学到了不错的MMAMs，也无法避免运动参数搜索中迭代算法的加速问题，迭代时间和图像质量需要得到权衡^［28］。第3类改进策略是用配对的数据直接让网络学习有伪影图像到无伪影图像的映射^{［29， 30］}。第3类方法虽然避免了以上两类改进策略的问题，但这种端到端的图像域DNN存在以下两个问题。首先，这类模型仅以CT图像作为输入，而在这之前的中间数据，比如分段有限角度重建子图像的信息丰富性和高时间分辨率没有得到分析和利用，传统的CT图像输入可能会限制模型的伪影校正性能。其次，这类模型仅学到有伪影图像到无伪影图像的简单映射，试图用固定的网络权重处理复杂多变的退化过程，缺少具体的CBCT运动退化建模。

鉴于此，本文分析了这些子图像的特性，借鉴了自然图像去模糊中的可学习退化编码^［31］，提出一种新的基于分段反投影张量退化特征编码的牙科CBCT运动伪影校正模型（SBP-MAC）。该模型由生成器和退化编码器构成，实现了从张量域到图像域的伪影校正，其主要创新在于从输入端提升传统图像域DNN的伪影校正性能，引入退化编码器指导模型自适应伪影校正，设计伪影一致性损失提升模型整体性能。

1 方法

1.1 分段反投影张量

图1A展示了分段反投影张量的由来：首先，对原始投影加权并滤波，得到预处理后的投影数据；然后，按照投影角度将其均分为 N_{seg} 份并分别反投影，得到 N_{seg} 张分段有限角度重建的子图像；最后，将它们沿通道方向堆叠，构成大小为 C\times H\times W 的张量，其通道大小C 恰好为N_{seg}。

张量中的子图像具有更丰富的信息。观察图1A，扇形面积表示CT图像重建中反投影的求和范围，求和操作越多则压缩的原始信息越多，显然分段有限角度重建的范围小于全角度重建，而所以这些子图像保留了更多原始信息。张量中的子图像具有更高的时间分辨率。观察图1A，全角度重建图像表示患者状态的时长t_{full}=\mathit{\tau }，对应的时间分辨率为\mathrm{1}/\mathit{\tau }；分段有限角度重建子图像表示患者状态的时长t_{sub}=\mathit{\tau }/{\mathit{N}}_{\mathit{s}\mathit{e}\mathit{g}}，对应的时间分辨率为N_{seg}/\mathit{\tau }。显然，这些子图像的时间分辨率更高。

图1B展示了传统图像域DNN的解决思路，它们以全角度重建图像为模型输入，直接学习有伪影图像到无伪影图像的映射关系。不同于这类思路，本文将子图像沿通道方向堆叠成张量并作为模型输入，用“更丰富的原始信息”和“更高的时间分辨率”这两大特性提高模型输入的质量，如图1C所示。此外，在模型\mathrm{ℛ} 中，本文还引入了退化编码器解决传统图像域DNN对CBCT运动退化建模不足的问题，还设计了伪影一致性损失提升模型整体性能，详见1.2 SBP-MAC模型。

1.2 SBP-MAC模型

图2A展示了模型输入数据的准备，即，从原始投影出发，用分段有限角度重建得到若干子图像，并沿通道方向堆叠成一个张量，见1.1。图2B展示了本文所提SBP-MAC模型的构成，包括生成器^［32］、退化编码器^［31］和空间自适应调制模块^［33］（SAM）。图2C展示了模型的输出结果，包括预测的有运动伪影图像和预测的无运动伪影图像。

本文所提模型中的生成器采用残差学习^［34］的方式进行训练，将张量沿通道方向叠加后与生成器的初始输出相加，得到最终的预测结果，如图2B的白色方框所示。这种残差学习将“直接预测图像内容”的困难任务转换为“仅需预测输入图像与目标图像的残差”这类相对简单的任务，从而可以降低生成器的训练难度。整体模型 \mathrm{ℛ} 的具体训练过程可描述如下：

当双向开关S_{\mathrm{1}}向上拨动、S_{\mathrm{2}}向下拨动时，模型的任务是从张量出发，预测无运动伪影图像，如以下公式所示：

y^{\text{'}}=R\left(T,M_{\mathrm{0}}\right)（1）

其中，y^{\text{'}} 表示预测的无运动伪影图像，\mathrm{ℛ} 表示所提SBP-MAC模型，T 表示张量，M_{\mathrm{0}} 表示退化特征。此时，退化编码器提取张量中空间变化的退化特征，并输入先验上采样模块，得到适应生成器各级跳跃连接大小的退化特征 M_i；然后，将 M_i 输入SAM模块，预测相同大小的调制参数后对原始跳跃连接特征 f_i 进行空间自适应调制，得到调制后的跳跃连接特征 F_i（具体的调制过程详见1.3 关键模块网络结构）；最后，由于有退化编码器提取张量中空间变化的退化特征，并且SAM模块能将这些特征整合到生成器的跳跃连接中，所以模型推理时能够感知到不同张量的运动信息，并根据提取到的退化特征，自适应校正不同运动波形导致的伪影。

当双向开关 S_{\mathrm{1}} 向下拨动、S_{\mathrm{2}} 向上拨动时，模型的任务是从张量出发，预测有运动伪影图像，如以下公式所示：

x^{\text{'}}=R\left(T,\mathrm{0}\right)（2）

其中，x^{\text{'}} 表示预测的有运动伪影图像。此时，M_{\mathrm{0}} 被双向开关 S_{\mathrm{1}} 置零，张量中空间变化的运动信息不再整合到生成器的跳跃连接中，模型感知不到这些运动信息，没有进行伪影校正；生成器仅仅完成张量融合这类相对简单的任务，在伪影一致性损失下，模型预测的有运动伪影图像与未校正图像的伪影表征一致。

1.3 关键模块网络结构

图3A表示退化编码器，由下采样卷积模块DownConv、卷积模块ConvBlock和残差模块ResnetBlock构成。其输入为张量T\in R^{C\times H\times W}，输出为退化特征 M_{\mathrm{0}}\in R^{\mathrm{2048}\times \frac{H}{\mathrm{32}}\times \frac{W}{\mathrm{32}}}。

图3B表示SAM，由插值模块Interpolate、多层感知机MLP、批量归一化模块BatchNorm构成。向SAM输入适应生成器各级跳跃连接f_i大小的退化特征M_i\in R^{\frac{\mathrm{512}}{{\mathrm{2}}^{i-\mathrm{1}}}\times \frac{H}{{\mathrm{2}}^{\mathrm{4}-i}}\times \frac{W}{{\mathrm{2}}^{\mathrm{4}-i}}}，通过插值、两级的多层感知机来预测相同大小的调制参数{\gamma }_i\in R^{\frac{\mathrm{512}}{{\mathrm{2}}^{i-\mathrm{1}}}\times \frac{H}{{\mathrm{2}}^{\mathrm{4}-i}}\times \frac{W}{{\mathrm{2}}^{\mathrm{4}-i}}}和 {\beta }_i\in R^{\frac{\mathrm{512}}{{\mathrm{2}}^{i-\mathrm{1}}}\times \frac{H}{{\mathrm{2}}^{\mathrm{4}-i}}\times \frac{W}{{\mathrm{2}}^{\mathrm{4}-i}}}，通过逐像素相加、相乘的方式，得到调制后的跳跃连接特征F_i\in R^{\frac{\mathrm{512}}{{\mathrm{2}}^{i-\mathrm{1}}}\times \frac{H}{{\mathrm{2}}^{\mathrm{4}-i}}\times \frac{W}{{\mathrm{2}}^{\mathrm{4}-i}}}，如以下公式所示：

F_i=(\mathrm{1}+{\gamma }_i)\odot f_i\oplus {\beta }_i（3）

图3C表示生成器，为编码-解码架构，其编码器由图像恢复模块NAFBlock^［32］和下采样卷积模块DownConv构成，解码器由图像恢复模块NAFBlock、PixelShuffle和上采样卷积模块UpConv构成。生成器对输入张量进行编码，用SAM调制f_i得到F_i后输入解码器参与解码；此外，生成器还将输入的张量沿通道方向叠加，用残差学习^［34］降低训练难度。

1.4 损失函数

我们用伪影校正损失 L_{MAC} 指导模型预测无伪影图像，如以下公式所示：

L_{MAC}=L_{MAE}\left(y^{\text{'}},y\right)+{\alpha }_{\mathrm{1}}{L}_{Contextual}\left(y^{\text{'}},y,l\right)+{\alpha }_{\mathrm{2}}{L}_{Gradient}\left(y^{\text{'}},y\right)（4）

其中，L_{MAE}为主要损失，我们将其权重系数固定为“1”，L_{Contextual}和L_{Gradient}为附加损失，我们用权重系数{\alpha }_{\mathrm{1}}和{\alpha }_{\mathrm{2}}调节其对L_{MAC}的影响；y^{\text{'}} 表示预测无运动伪影图像，y 表示标签无运动伪影图像。L_{MAE} 为平均绝对误差损失，它对处理预测图像异常值的鲁棒性较好；L_{Contextual} 为纹理损失^［35］，它能够更好地保留图像中的细节和纹理，提升视觉质量；L_{Gradient} 为图像梯度损失，它能更好地保持图像中牙齿的轮廓。这三个损失的定义如公式（5~7）所示：

L_{MAE}\left(y^{\text{'}},y\right)={‖y^{\text{'}}-y‖}_{\mathrm{1}},\#（5）

L_{Contextual}\left(y^{\text{'}},y,l\right)=-log\left(CX({\psi }^l\left(y^{\text{'}}\right),{\psi }^l\left(y\right))\right)（6）

L_{Gradient}\left(y^{\text{'}},y\right)={‖G_H\left(y^{\text{'}}\right)-G_H\left(y\right)‖}_{\mathrm{1}}+{‖G_W\left(y^{\text{'}}\right)-G_W\left(y\right)‖}_{\mathrm{1}}（7）

其中，\psi 表示VGG19预训练网络^［36］，{\psi }^l\left(y^{\text{'}}\right) 和 {\psi }^l\left(y^{\text{'}}\right) 表示VGG19网络第l层对图像 y^{\text{'}} 和 y^{\text{'}} 所提取的特征，CX 表示上下文相似度算子，G_H 和 G_W 分别表示高度方向、宽度方向上的图像梯度算子。

我们用伪影一致性损失 L_{Consistency} 指导模型预测有伪影图像，它的本质是感知损失，能帮助模型生成更接近标签的有伪影图像。其定义如以下公式所示：

L_{Consistency}=L_{perceptual}(x^{\text{'}},x,l)={‖{\psi }^l(x^{\text{'}}),{\psi }^l(x)‖}_{\mathrm{1}}（8）

其中，x^{\text{'}} 表示预测有运动伪影图像，x 表示标签有运动伪影图像。

最后，模型总损失 L_{total} 由 L_{MAC} 和 L_{Consistency} 构成，通过超参数 \lambda 调和，如以下公式所示：

L_{total}=\lambda L_{MAC}+\left(\mathrm{1}-\lambda \right)L_{Consistency}.（9）

1.5 运动伪影仿真

临床上，真实的有运动-无运动配对数据非常稀缺，不足以训练深度学习模型。为了解决这一问题，我们在大量的无运动数据基础上进行运动仿真，构建仿真的有运动-无运动配对数据集。因为人体头颅不易发生形变，且牙托限制了上下颌的咬动，所以我们将患者头部运动建模为6个自由度的全局刚性运动。因为患者运动多为随机抖动，且实际的牙科CBCT扫描过程中有牙托等固定装置的辅助，所以我们将运动波形建模为随机曲线，并确保运动以x轴平移、y轴平移、z轴旋转为主，其他轴向的运动幅度相应减小。

实际上，患者头部的全局刚性运动可以转化为系统机架的运动，所以我们将患者运动波形转化为机架的运动波形，并用机架运动矩阵修改原始系统矩阵，以此来仿真患者运动。用修改前后的系统矩阵对原始投影重建，可以获得配对的图像数据和张量数据。具体的系统矩阵修改过程如公式（10~12）所示：

P_i^{*}=P_i\cdot M{o}_i=K\cdot [R_i\left|T_i]\cdot \left(\begin{array}{c}{r}_i{t}_i\\ \mathrm{0}\mathrm{1}\end{array}\right)\right.（10）

t_i={\left[t_{xi},t_{yi},t_{zi}\right]}^T（11）

\begin{array}{c}\begin{array}{l}{r}_i=\left[\begin{array}{ccc}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& cos{\theta }_{xi}& -sin{\theta }_{xi}\\ \mathrm{0}& sin{\theta }_{xi}& cos{\theta }_{xi}\end{array}\right]\left[\begin{array}{ccc}cos{\theta }_{yi}& \mathrm{0}& sin{\theta }_{yi}\\ \mathrm{0}& \mathrm{1}& \mathrm{0}\\ -sin{\theta }_{yi}& \mathrm{0}& cos{\theta }_{yi}\end{array}\right]\\ \left[\begin{array}{ccc}cos{\theta }_{zi}& -sin{\theta }_{zi}& \mathrm{0}\\ sin{\theta }_{zi}& cos{\theta }_{zi}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{1}\end{array}\right]\end{array}\end{array}（12）

其中，P_i表示原始系统矩阵，P_i^{*} 表示修改后的系统矩阵，M{o}_i 表示运动矩阵。K表示P_i的内参矩阵，R_i表示机架固有的旋转，T_i表示机架固有的平移，[R_i|T_i]构成P_i的外参矩阵。t_i 表示平移运动矩阵，控制3个坐标的平移，r_i 表示旋转运动矩阵，控制3个坐标的旋转。所有公式的右下角标i表示投影序号。

图4A表示CBCT扫描示意图，x,y,z轴表示世界坐标系，x^{\text{'}},y^{\text{'}},z^{\text{'}}轴表示旋转坐标系，u,v轴表示探测器坐标系，灰色圆柱表示被扫描对象，其在扫描过程中发生了图4B、C所示的运动。其中，沿x和y轴平移的波形t_x、t_y幅度控制在1mm，沿z轴平移的波形t_z幅度控制在0.1mm；沿x和y轴旋转的波形{\theta }_x、{\theta }_y幅度控制在0.1°，沿z轴旋转的波形{\theta }_z幅度控制在1°。每套数据的运动波形均随机产生，并通过公式（10）转换到CBCT机架运动上。

1.6 实验设计

1.6.1 数据准备及实验设置

本文的所有数据均由菲森引擎（CBCT WuKong Matrix5000，Fussen）重建得到。源到探测器的距离（DSD）为635 mm，源到旋转中心的距离（DSO）为425 mm，投影角度为360°，投影角度个数 N_{angle} 为450，投影分段个数 N_{seg} 为9，探测器大小为768×768，探测器通道个数 N_{channel} 为768，探测器单元大小为0.2 mm×0.2 mm；重建体积大小为512×512×320，体素大小为0.3125 mm×0.3125 mm×0.3125 mm。训练数据包括9000套仿真的有运动张量数据和有运动-无运动配对图像数据；测试数据包括3200套仿真的有运动张量数据和有运动-无运动配对图像数据，9600套真实临床的有运动张量数据和有运动图像数据。其中，对于本文所提模型，采用有运动张量数据作为输入，无运动图像和有运动图像作为标签；对于其余4个对比方法，采用有运动CT图像作为输入，无运动CT图像作为标签。

本文所提模型的训练次数为200，一次训练的迭代个数为500；初始学习率设置为I_r=\mathrm{1}\times {\mathrm{10}}^{-\mathrm{4}}，每训练40次，学习率降低到原来的一半；模型训练采用Adam优化器，其中\left({\beta }_{\mathrm{1}},{\beta }_{\mathrm{2}}\right)=(\mathrm{0.9,0.999})。

本文在所提模型的训练过程中使用的计算机配置具体如下：CPU型号为Inter Xeon E5-2683（2.10 GHz），GPU型号为NVIDIA GPU RTX A6000（48 GB）。Python 版本为 3.9.18，CUDA 版本为 12.5，Pytorch版本为1.12.0。其他依赖包均在以上环境下进行安装配置。

1.6.2 对比方法及评价指标

本文采用了四种图像域恢复网络作为对比方法，分别是U-Net^［37］、REDCNN^［38］、HINet^［39］、Restormer^［40］，他们的模型参数设置均与本文所提模型相同。本文在仿真测试集上的定量指标为：结构相似度（SSIM）、峰值信噪比（PSNR）、均方根误差（RMSE）。

SSIM用于衡量预测图像与真值图像在亮度、对比度和结构上的相似性，其越高表示图像质量越好，如以下公式所示：

SSIM\left(I_{pred},I_{gt}\right)=\frac{\left(\mathrm{2}{\mu }_{pred}{\mu }_{gt}+C_{\mathrm{1}}\right)\left(\mathrm{2}{\sigma }_{I_{pred}{I}_{gt}}+C_{\mathrm{2}}\right)}{\left({\mu }_{I_{pred}}^{\mathrm{2}}+{\mu }_{I_{gt}}^{\mathrm{2}}+C_{\mathrm{1}}\right)\left({\sigma }_{I_{pred}}^{\mathrm{2}}+{\sigma }_{I_{gt}}^{\mathrm{2}}+C_{\mathrm{2}}\right)}（13）

其中，I_{pred} 表示模型预测的图像，I_{gt} 表示真值图像，{\mu }_{pred},{\mu }_{gt} 分别表示他们对应的均值，{\sigma }_{pred},{\sigma }_{gt} 分别表示他们对应的标准差，{\sigma }_{I_{pred}{I}_{gt}} 表示I_{pred} 和 I_{gt} 之间的协方差，C_{\mathrm{1}}=(K_{\mathrm{1}}{L)}^{\mathrm{2}}表示亮度相似度计算中的常数，C_{\mathrm{2}}=(K_{\mathrm{2}}{L)}^{\mathrm{2}}表示对比度相似度计算中的常数，且 K_{\mathrm{1}}=\mathrm{0.01}，K_{\mathrm{2}}=\mathrm{0.03}，L=\mathrm{255}。

PSNR用于衡量图像峰值信号与噪声的比值，其越高表示图像质量越好，如以下公式所示：

PSNR\left(I_{pred},I_{gt}\right)=\mathrm{10}lo{g}_{\mathrm{10}}\frac{MA{X}_{I_{pred}}^{\mathrm{2}}}{\frac{\mathrm{1}}{mn}{\sum }_{i=\mathrm{1}}^m{\sum }_{j=\mathrm{1}}^n{\left[I_{pred}\left(i,j\right)-I_{gt}\left(i,j\right)\right]}^{\mathrm{2}}}（14）

其中，m,n 表示图像两个维度的大小，i,j 表示图像像素的索引，MA{X}_{I_{pred}} 表示预测图像的最大像素值。

RMSE用于衡量预测图像与真值图像之间的均方根差异，其越低表示图像质量越好，如以下公式所示：

RMSE\left(I_{pred},I_{gt}\right)=\sqrt[]{\frac{\mathrm{1}}{mn}{\sum }_{i=\mathrm{1}}^m{\sum }_{j=\mathrm{1}}^n{\left[I_{pred}\left(i,j\right)-I_{gt}\left(i,j\right)\right]}^{\mathrm{2}}}（15）

此外，五位CT领域的专家在真实测试集上对本文所提方法及其他对比方法的伪影校正性能进行了评分。我们采用Mann-Whitney方法检验本文方法和对比方法在专家评分上的显著性差异。当P<0.05时，意味着本文方法的校正性能与对比方法相比有显著性差异。专家评分被分为五个等级：“5分”表示图像中无运动伪影，解剖结构恢复较好；“4分”表示图像中无明显条纹伪影，只在牙齿或骨组织附近有少量模糊；“3分”表示图像中存在少量条纹伪影，牙齿和骨组织存在模糊；“2分”表示图像中有显著条纹伪影，牙齿和骨组织结构比较模糊；“1分”表示图像中有严重运动伪影，牙齿和骨组织结构严重失真。

2 结果

2.1 仿真数据实验结果

图5~7分别展示了不同方法在仿真数据横断面、冠状面、矢状面上的伪影校正效果。其中，第1行为未校正结果，第2行为真值，第3~6行为对比方法结果，第7行为本方法结果；感兴趣区域ROI已放大展示，其关键细节由黄色箭头指示，图像显示窗均为［-1000，2400］HU。观察结果可发现，对比方法虽然可以去除一定程度的牙齿双轮廓，但它们的牙齿边缘清晰度不高，牙根管结构恢复程度有限；本文所提SBP-MAC模型在各个解剖面中，残留的牙齿双轮廓最少，牙齿边缘清晰度最高，牙根管结构恢复程度最高。观察结果可发现，本文所提的模型相比未校正的情况，3项指标均有显著的改善，PSNR提高了8.28%，RMSE降低了23.84%，SSIM提高了2.29%。相比其他对比方法，本方法所有定量结果最优，最优定量结果由宋体加粗标出（表1）。

2.2 真实数据实验结果

图8~10分别展示了不同方法在真实数据横断面、冠状面、矢状面上的伪影校正效果。其中，第1行为未校正的结果，第2~5行为对比方法结果，第6行为本方法结果；感兴趣区域ROI已放大展示，其关键细节由黄色箭头指示，图像显示窗均为［-1000，2400］HU。观察结果发现，对比方法虽然可以去除轻微程度的牙齿双轮廓，但它们无法消除较为明显的运动拖影，无法较好保持牙齿边缘清晰度，牙根管结构恢复程度有限；本文所提SBP-MAC模型在各个解剖面中，残留的牙齿双轮廓最少，牙齿边缘清晰度最高，牙根管结构恢复程度最高，还能消除明显的运动拖影。表2展示了五位相关领域专家在真实临床数据上对五种伪影校正方法的图像质量主观评分统计情况。其中，本文方法在真实临床数据中获得的评分最高（4.2421±0.3548），最高评分由宋体加粗标出；从U检验结果中可以发现，本文所提方法与其他对比方法的图像质量评分差异均具有统计学意义（P<0.05）。

2.3 伪影一致性损失有效性验证

图11展示了退化编码器不指导模型自适应伪影校正时，模型\mathrm{ℛ}的输出结果。回顾公式（2）可知，此时的输出结果应当逼近标签运动伪影图像。观察标签运动伪影图像 x 和生成运动伪影图像 x^{\text{'}}，以及他们分别与标签无伪影图像 y 之间的残差热图 x-y 和 x^{\text{'}}-y，可发现他们的伪影表征趋近一致；观察 x^{\text{'}} 和x 之间的残差热图 x^{\text{'}}-x，也可发现他们在牙齿、骨小梁附近的伪影都趋近一致。这些结果表明，在没有退化编码器起作用时，生成器确实完成了张量简单融合，伪影一致性损失有效。

2.4 超参数 \mathit{\lambda } 的必要性及有效性验证

表3展示了模型总损失L_{total}中调和超参数 \lambda不同取值对模型整体性能的影响，此时固定{\alpha }_{\mathrm{1}}=\mathrm{0.1}，{\alpha }_{\mathrm{2}}=\mathrm{1}，使得加权后的两个伪影校正附加损失保持在同一水平。结果表明，当\lambda =\mathrm{0}时，伪影校正损失L_{MAC}不参与模型训练，模型没有显著的伪影校正性能；当\lambda =\mathrm{1} 时，伪影一致性损失L_{Consistency}没有参与模型训练，模型有显著的伪影校正性能；当\lambda =\mathrm{0.8}时，超参数调制两种损失共同训练模型，此时伪影校正性能最佳，PSNR达到最高31.9936，RMSE达到最低6.5657，SSIM达到最高0.9547，最优定量结果由宋体加粗标出。以上实验结果说明了超参数 \lambda 对于提升模型性能的必要性和有效性。

2.5 各个关键设计的必要性及有效性验证

表4展示了不同关键设计对模型整体性能的影响，“--”表示不考虑该因素，“√”表示考虑该因素。结果表明，基础模型确实有一定程度的伪影校正性能，但逐步添加本文所提的3个关键设计后，各项指标均同步变好；在3个关键设计共同作用于基础模型时，各项定量指标达到最好，PSNR达到最高31.9936，RMSE达到最低6.5657，SSIM达到最高0.9547，最优定量结果由宋体加粗标出。

2.6 {\mathit{\alpha }}_{\mathrm{1}}和{\mathit{\alpha }}_{\mathrm{2}}对模型性能影响的探究性实验

表5展示了公式（4）中不同权重系数对模型性能的影响。结果表明，当用{\alpha }_{\mathrm{1}}=\mathrm{0.1}和{\alpha }_{\mathrm{2}}=\mathrm{1}加权L_{MAC}中的两个附加损失时，模型的各项定量指标达到最好，PSNR达到最高31.9936，RMSE达到最低6.5657，SSIM达到最高0.9547，最优定量结果由宋体加粗标出。

2.7 {\mathit{N}}_{\mathit{s}\mathit{e}\mathit{g}}对模型性能影响的探究性实验

表6和图12分别展示了不同分段个数对模型性能的影响。结果显示，随着分段个数的增加，模型各项指标逐渐变优，当N_seg=9时达到最优；若继续增加分段个数，定量指标有小幅度下降，伪影校正效果略低于N_seg=9。

3 讨论

传统的图像域DNN存在以下问题，一是没有分析CT图像重建链条中其他中间数据对模型性能的增益可能，二是试图用固定的网络权重处理复杂多变的退化过程，缺少具体的CBCT运动退化建模。为了解决这两个问题，本文分析了分段有限角度重建子图像的特性并作为模型的新输入，借鉴了自然图像去模糊中的可学习退化编码^［31］，提出一种新的基于分段反投影张量退化特征编码的牙科CBCT运动伪影校正模型（SBP-MAC）。

该模型由生成器和退化编码器构成，实现了从张量域到图像域的伪影校正，相比传统的图像域DNN校正模型有以下创新点。第一，将子图像堆叠成张量并作为模型的新输入。这能充分开发子图像的优点，即，更丰富的原始信息和更高的时间分辨率，直接从输入端提升现有图像域DNN的伪影校正性能。第二，引入退化编码器指导模型自适应伪影校正。这个可学习的编码器能提取张量中空间变化的运动信息，自适应调制生成器的各级跳跃连接特征，从而指导模型校正不同运动波形导致的伪影。第三，设计伪影一致性损失提升模型整体性能。规定在没有编码器起作用时，生成器仅输出有运动伪影图像，其与未校正图像的伪影表征一致。该损失一方面可以降低生成器的学习难度，另一方面可以使编码器专注于张量的特征提取。

为了验证本文所提模型对牙科CBCT运动伪影校正的有效性及鲁棒性，本文分别在仿真数据和真实数据上与4个图像域DNN模型进行了对比实验。不同模型的伪影校正效果以横断面、冠状面、矢状面展示。观察3个解剖面的ROI区域及箭头所指的关键解剖结构，可以发现无论是仿真数据还是真实数据，大部分案例中：UNet和REDCNN残留的牙齿双轮廓最多，牙齿边缘清晰度最低，牙根管结构恢复程度不高，整体伪影校正效果一般；HINet虽然能去除更多牙齿双轮廓，但牙齿边缘清晰度没有得到较好的保持，牙根管结构恢复程度较低；Restormer残留的牙齿双轮廓相对较少，牙齿边缘清晰度保持得较好，但对牙根管结构的恢复程度一般。本文所提SBP-MAC模型残留的牙齿双轮廓最少，牙齿边缘清晰度最高，牙根管结构恢复程度最高。针对一些牙齿粘连，牙齿边缘开裂和牙齿拖尾，本文所提模型的校正效果显著。

此外，仿真数据实验结果中的定量指标表明，本所提模型的所有指标结果均优于对比模型，PSNR提高了8.28%，RMSE降低了23.84%，SSIM提高了2.29%。真实数据实验结果中的专家评分表明，本文所提模型的评分最高，且与对比模型的评分差异均具有统计学意义（P<0.05）。

在“伪影一致性损失的有效性验证”实验中，本文检验了无退化编码器作用的条件下模型的输出，结果表明了该损失的有效性。在“超参数的必要性和有效性”实验中，结果表明了超参数对模型性能的必要性和有效性，找出了最优超参数。在“各个关键设计的必要性和有效性”实验中，结果验证了所有关键设计对基础模型性能的提升，证实了张量、退化编码器、伪影一致性损失的必要性和有效性。此外，我们探究了两个附加损失的权重因子对模型性能的影响，结果表明当附加损失维持相当水平时模型性能达到最优；我们还探究了不同N_{seg}对模型性能的影响，发现当N_{seg}大于9时，模型性能开始出现小幅下降。同时，本研究还存在一些不足。例如，模型的自适应伪影校正性能需要进一步验证，在未来工作中将针对同一患者的不同运动波形做相关方面的检验。

综上所述，本文所提SBP-MAC模型创新了传统图像域DNN的输入策略，解决了缺少CBCT运动退化建模的问题，可以有效解决患者在扫描过程中发生躯体运动导致的重建图像质量下降问题。

参考文献

[1] Ng SY. Comparison between conventional dental radiography and CBCT[M]// Cone Beam CT in Dentistry. Cham: Springer, 2023: 271-330.

[2] Schulze RKW, Drage NA. Cone-beam computed tomography and its applications in dental and maxillofacial radiology[J]. Clin Radiol, 2020, 75(9): 647-57.

[3] Nemtoi A, Czink C, Haba D, et al. Cone beam CT: a current overview of devices[J]. Dentomaxillofac Radiol, 2013, 42(8): 20120443.

[4] Spin-Neto R, Wenzel A. Patient movement and motion artefacts in cone beam computed tomography of the dentomaxillofacial region: a systematic literature review[J]. Oral Surg Oral Med Oral Pathol Oral Radiol, 2016, 121(4): 425-33.

[5] Keriş EY. Effect of patient anxiety on image motion artefacts in CBCT[J]. BMC Oral Health, 2017, 17(1): 73.

[6] Donaldson K, O'Connor S, Heath N. Dental cone beam CT image quality possibly reduced by patient movement[J]. Dentomaxillofac Radiol, 2013, 42(2): 91866873.

[7] Spin-Neto R, Hauge Matzen L, Hermann L, et al. Head motion and perception of discomfort by young children during simulated CBCT examinations[J]. Dentomaxillofac Radiol, 2021, 50(3): 20200445.

[8] Spin-Neto R, Costa C, Salgado DM, et al. Patient movement characteristics and the impact on CBCT image quality and interpretability[J]. Dentomaxillofac Radiol, 2018, 47(1): 20170216.

[9] Spin-Neto R, Matzen L H, Schropp L, et al. Factors affecting patient movement and re-exposure in cone beam computed tomography examination [J]. Oral Surg Oral Med Oral Pathol Oral Radiol, 2015, 119(5): 572-8.

[10] Nardi C, Taliani GG, Castellani A, et al. Repetition of examination due to motion artifacts in horizontal cone beam CT: comparison among three different kinds of head support[J]. J Int Soc Prev Community Dent, 2017, 7(4): 208-13.

[11] Yildizer Keris E, Demirel O, Ozdede M. Evaluation of motion artifacts in cone-beam computed tomography with three different patient positioning[J]. Oral Radiol, 2021, 37(2): 276-81.

[12] Jacobson MW, Stayman JW. Compensating for head motion in slowly-rotating cone beam CT systems with optimization transfer based motion estimation[C]//2008 IEEE Nuclear Science Symposium Conference Record. Dresden, Germany. IEEE, 2008: 5240-5.

[13] Kyme AZ, Fulton RR. Motion estimation and correction in SPECT, PET and CT[J]. Phys Med Biol, 2021, 66(18): 18TR02.

[14] Ouadah S, Jacobson M, Stayman JW, et al. Correction of patient motion in cone-beam CT using 3D-2D registration[J]. Phys Med Biol, 2017, 62(23): 8813-31.

[15] Sun T, Jacobs R, Pauwels R, et al. A motion correction approach for oral and maxillofacial cone-beam CT imaging[J]. Phys Med Biol, 2021, 66(12): 125008.

[16] Berger M, Müller K, Aichert A, et al. Marker-free motion correction in weight-bearing cone-beam CT of the knee joint[J]. Med Phys, 2016, 43(3): 1235-48.

[17] Niebler S, Schömer E, Tjaden H, et al. Projection-based improvement of 3D reconstructions from motion-impaired dental cone beam CT data[J]. Med Phys, 2019, 46(10): 4470-80.

[18] Hernandez D, Eldib ME, Hegazy MAA, et al. A head motion estimation algorithm for motion artifact correction in dental CT imaging[J]. Phys Med Biol, 2018, 63(6): 065014.

[19] Berger M, Xia Y, Aichinger W, et al. Motion compensation for cone-beam CT using Fourier consistency conditions[J]. Phys Med Biol, 2017, 62(17): 7181-215.

[20] Aichert A, Berger M, Wang J, et al. Epipolar consistency in transmission imaging[J]. IEEE Trans Med Imaging, 2015, 34(11): 2205-19.

[21] A'lvarez-Borrego J. Fast autofocus algorithm for automated microscopes[J]. Opt Eng, 2005, 44(6): 063601.

[22] Wicklein J, Kyriakou Y, Kalender WA, et al. An online motion- and misalignment-correction method for medical flat-detector CT[C]//SPIE Proceedings", "Medical Imaging 2013: Physics of Medical Imaging. Lake Buena Vista (Orlando Area), Florida, USA. SPIE, 2013: 466-72.

[23] Sisniega A, Stayman JW, Yorkston J, et al. Motion compensation in extremity cone-beam CT using a penalized image sharpness criterion[J]. Phys Med Biol, 2017, 62(9): 3712-34.

[24] Zhang Y, Zhang LY. A rigid motion artifact reduction method for CT based on blind deconvolution[J]. Algorithms, 2019, 12(8): 155.

[25] Gao C, Feng AQ, Liu XT, et al. A fully differentiable framework for 2D/3D registration and the projective spatial transformers[J]. IEEE Trans Med Imag, 2024, 43(1): 275-85.

[26] Ali ASRA, Landi C, Sarti C, et al. Non-iterative compensation for patient motion in dental CBCT imaging[C]//2023 IEEE Nuclear Science Symposium, Medical Imaging Conference and International Symposium on Room-Temperature Semiconductor Detectors (NSS MIC RTSD). Vancouver, BC, Canada. IEEE, 2023: 1.

[27] Preuhs A, Manhart M, Roser P, et al. Appearance learning for image-based motion estimation in tomography[J]. IEEE Trans Med Imaging, 2020, 39(11): 3667-78.

[28] Thies M, Wagner F, Maul N, et al. A gradient-based approach to fast and accurate head motion compensation in cone-beam CT[J]. IEEE Trans Med Imaging, 2024, doi: 10.1109/TMI.2024.3474250

[29] Amirian M, Montoya-Zegarra JA, Herzig I, et al. Mitigation of motion-induced artifacts in cone beam computed tomography using deep convolutional neural networks[J]. Med Phys, 2023, 50(10): 6228-42.

[30] Ko Y, Moon S, Baek J, et al. Rigid and non-rigid motion artifact reduction in X-ray CT using attention module[J]. Med Image Anal, 2021, 67: 101883.

[31] Li DS, Zhang Y, Cheung KC, et al. Learning degradation representations for Image deblurring[M]//Lecture Notes in Computer Science. Cham: Springer Nature Switzerland, 2022: 736-53.

[32] Chen LY, Chu XJ, Zhang XY, et al. Simple baselines for Image restoration[M]//Lecture Notes in Computer Science. Cham: Springer Nature Switzerland, 2022: 17-33.

[33] Park T, Liu MY, Wang TC, et al. Semantic image synthesis with spatially-adaptive normalization[C]//2019 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR). Long Beach, CA, USA. IEEE, 2019: 2332-41.

[34] He KM, Zhang XY, Ren SQ, et al. Deep residual learning for image recognition[C]//2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR). Las Vegas, NV, USA. IEEE, 2016: 770-8.

[35] Mechrez R, Talmi I, Shama F, et al. Maintaining natural image statistics with the contextual loss[C]//Asian Conference on Computer Vision. Cham: Springer, 2019: 427-443.

[36] Simonyan K, Zisserman A. Very deep convolutional networks for large-scale image recognition[EB/OL]. 2014: arXiv: 1409.1556.

[37] Ronneberger O, Fischer P, Brox T. U-net: Convolutional networks for biomedical image segmentation [C]//2015 Medical Image Computing and Computer Assisted Intervention (MICCAI). Munich, Germany. Springer International Publishing, 2015: 234-41.

[38] Chen H, Zhang Y, Kalra MK, et al. Low-dose CT with a residual encoder-decoder convolutional neural network[J]. IEEE Trans Med Imaging, 2017, 36(12): 2524-35.

[39] Chen LY, Lu X, Zhang J, et al. HINet: half instance normalization network for image restoration[C]//2021 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops (CVPRW). Nashville, TN, USA. IEEE, 2021: 182-92.

[40] Zamir SW, Arora A, Khan S, et al. Restormer: efficient transformer for high-resolution image restoration[C]//2022 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR). New Orleans, LA, USA. IEEE, 2022: 5718-29.

基金资助

RIGHTS & PERMISSIONS

版权所有©《南方医科大学学报》编辑部2021