一、引言
医学图像到图像的翻译旨在将图像从一个域转移到另一个域，同时保持结构的完整性。其中一个应用是CBCT到CT的翻译，这对放射治疗应用至关重要。最近，更多的图像翻译方法是基于生成对抗网络(GAN)或变分自编码器(VAE)。CycleGAN 广泛用于医学图像翻译任务，它通过两个相互竞争的神经网络进行学习，间接捕获信息以提供关于源域和目标域分布的隐式表示。然而，这种隐式表示容易学习偏差，包括鉴别器过早收敛和模式崩溃，从而影响合成图像的质量。
最近，扩散模型被开发为一种新型的生成式AI方法，该方法可以生成具有优异Frechet Inception Distance (FID)的高质量图像。训练一个具有代表性的扩散模型去噪扩散概率模型(DDPM)，对不同程度的高斯噪声破坏的样本进行去噪。噪声损坏样本由称为前向扩散的马尔可夫链蒙特卡罗 (MCMC) 过程生成，通过随着高斯白噪声水平的增加而逐渐破坏图像。然后，基于损坏的图像，DDPM 通过学习基于 Langevin 动力学反转前向扩散过程（反向扩散），逐步去噪并将噪声转换为有意义的高质量图像。正向和反向扩散过程帮助DDPM学习图像的底层潜在空间和分布，并从学习的分布生成新的图像。
考虑到扩散模型可以在全面的真实性和质量方面生成优于基于 GAN 的模型的图像，因此将扩散模型应用于医学图像翻译的兴趣不断增长。通过DDPM，我们可能会对源图像域进行正向扩散，并使用学习到的反向扩散将噪声损坏的图像转换为目标图像域，实现直接域平移。这样的策略已经在自然图像翻译中取得了成功，例如猫狗图像翻译，它只需要保留一些与领域无关的特征（姿势、颜色等）。然而，扩散模型在医学图像到图像翻译任务中的应用目前是有限的，因为它对源图像的正向扩散过程导致结构细节损失，这在反向去噪过程中不能完全恢复。翻译(合成)图像的几何和结构完整性对于外科规划和放疗等医疗应用尤其重要，因为医学图像具有丰富的准确的解剖细节，对于这两种程序都是至关重要的。
为了保留结构信息，一些方法使用配对图像通过将源图像作为条件输入到反向扩散过程中来开发扩散模型，以约束目标域图像的解。然而，这种方法需要具有匹配解剖的精选源图像和目标图像对，由于扫描间解剖运动、成像成本、辐射剂量等，在医学成像中很难获得这些对。
对于扩散模型，正向扩散过程可以被认为是一个低通滤波器(在第三章节会进行解释)。相比之下，医学图像中的解剖结构轮廓大多嵌入高频信息，通过正向扩散过程进行过滤。对于相同的解剖结构，如果两个图像域/批次之间的差异大多在傅里叶域的中间频率上，我们可能会在条件目标域图像生成之前从源图像域中提取高频信息作为结构。我们可以使用仅通过目标域图像训练的扩散模型来填充中间频率空位。对于CBCT-toCT平移问题，我们定量和定性地分析了CBCT和CT图像在频域上的差异，发现它们之间的主要差异处于中间频率，如图2所示。基于这种频域观测，我们设计了一种由高通滤波器和低通滤波器(受控正向扩散过程)联合引导的图像平移扩散模型，其中低通滤波器捕获强度和整体语义信息，高通滤波器捕获解剖细节。生成的频率引导扩散模型 (FGDM) 以高频信息和低频信息为条件，为医学图像翻译生成中间频率信息。
与 GAN 或 VAE 相比，FGDM 提供了额外的好处，因为它仅使用来自目标域的数据来训练扩散模型以学习中间频率的分布。如果源域图像在两个频率端与目标域图像的末尾共享相似的信息，则得到的模型可用于来自不同源域的图像翻译（图 2）。对源域的松散限制允许自然零样本翻译 [14]，该模型可以应用于其他领域或分布中的源图像，而无需迁移学习或测试时间重新训练。据我们所知，我们是第一个在解剖水平上实现零样本医学图像翻译的扩散模型的人。本文的主要贡献总结如下：
1️⃣我们提出了一种频率引导的扩散模型，通过使用低频和高频信息来指导中间频率信号的生成，实现了保留准确解剖信息的医学图像翻译。
2️⃣我们管理和优化高通和低通滤波器以提取域不变特征信息。具体来说，我们两个滤波器的滤波阈值在测试阶段可以自由调整，以应对不同类型医学图像翻译任务的源域和目标域之间傅立叶空间中的不同域不变区域。
3️⃣我们的模型在CBCT到CT翻译任务中取得了最好的性能，超过了所有其他最先进的(SOTA)基准模型，包括基于gan的方法、基于VAE和基于扩散的方法。
二、方法
我们的模型的概述如图 4 所示。FGDM 可以将不可见的源域数据（CBCT）转换为目标域数据（CT）。对CBCT进行低通滤波和高通滤波，得到低频信息和高频信息作为扩散条件，通过反向扩散过程逐渐产生中频信息，得到最终的CT图像。
1.低通滤波器（Low Pass Filter）
在本节中，我们将展示前向扩散过程可以被视为低通滤波器（LPF）。首先，我们将 \(\textcolor{blue}{s_t(x, y)}\)、\(\textcolor{blue}{i(x, y)}\)和 \(\textcolor{blue}{z(x, y)}\) 分别设置为退化图像（在步骤 t 处）、原始图像和高斯白噪声。我们可以将前向扩散过程视为具有由递减序列 \(\textcolor{blue}{α_{1:t} \in (0, 1]^T}\)参数化的高斯转移的马尔可夫链。，其中前向过程的一个特殊属性可以从方程 2 导出为：
相应地，我们可以将 \(\textcolor{blue}{s_t}\) 表示为 \(\textcolor{blue}{s_0}\) 和噪声变量 \(\textcolor{blue}{z}\) 的线性组合：
在频域中，可以表示为：
其中\(\textcolor{blue}{S_t(u,v)}\)、\(\textcolor{blue}{I(u, v)}\) 和 \(\textcolor{blue}{Z(u, v)}\)分别是\(\textcolor{blue}{s_t(x, y)}\)、\(\textcolor{blue}{i(x, y)}\)和 \(\textcolor{blue}{z(x, y)}\)的傅里叶变换。高斯白噪声的自相关函数（ACF）可以定义为：
其中 \(\textcolor{blue}{τ_1}\) 和 \(\textcolor{blue}{τ_2}\) 表示滞后。在高斯白噪声的情况下，信号\(\textcolor{blue}{Z(u, v)}\)是均值为零且方差为常数 \(\textcolor{blue}{σ^2}\) 的随机变量：
因此，高斯白噪声的自相关函数为：
其中\(\textcolor{blue}{δ(τ_1, τ_2)}\)是狄拉克\(\textcolor{blue}{δ}\)函数。因此，高斯白噪声的功率在所有频率上都是相同的。
在一般图像中，功率谱密度（PSD）与空间频率的关系如下：
其中 \(\textcolor{blue}{k}\) 和 \(\textcolor{blue}{a}\) 是缩放/修改因子，表示 \(\textcolor{blue}{PSD(u, v)_I ∝ \frac{1} {∥(u,v)∥}}\) 。当高斯白噪声添加到图像中时，噪声的 PSD 会添加到图像的 PSD 中。如果\(PSD_I\) 是原始图像的功率，\(PSD_Z\) 是添加的高斯白噪声的功率，则步骤\(\textcolor{blue}{t}\)中的信噪比（SNR）可以定义为：
从方程19所示，\(SNR_t\)随着\(\textcolor{blue}{\alpha_t}\)增加。由于 \(\textcolor{blue}{\alpha_t}\) 是递减序列 (\(\textcolor{blue}{α_t ∝ \frac{1}{t}}\))，因此 \(SNR_t\) 和 \(\textcolor{blue}{t}\)之间的关系可以建立为：
从方程19所示，任何频率\(\textcolor{blue}{||(u，v)||}\)处的\(SNR_t\)也与该频率处原始信号的功率成正比，并且随着频率的增加而减小。
因此，随着前向扩散步长 \(\textcolor{blue}{t}\)的增加，\(SNR_t\)减小，并且较高频率的信号更早受到影响。将\(\textcolor{blue}{\phi}\)视为 SNR 阈值，低于该阈值信息将完全损坏，并且根据低通滤波器所需的频率阈值 \(\textcolor{blue}{ψ}\)，存在适当数量的步数 \(\textcolor{blue}{t}\)，超过该步数所有非低频信息都将丢失：
如果我们将这一步 \(\textcolor{blue}{t}\) 表示为\(\textcolor{blue}{\tilde{T}}\)，则前向扩散等于低通滤波器：
2.高通滤波器（High Pass Filter）
作者使用 Sobel 算子作为高通滤波器（HPF），它在技术上是一个离散微分算子，它计算图像强度函数的梯度以获得高频信息。具体来说，我们首先将图像分别与水平和垂直Sobel滤波器核\(\textcolor{blue}{K_x}\)和\(\textcolor{blue}{K_y}\)进行卷积，以获得滤波图像\(\textcolor{blue}{G_x}\)和\(\textcolor{blue}{G_y}\)。在图像中的每个点，可以使用以下方程组合所得的梯度近似值以给出梯度幅度。
然后通过阈值\(\textcolor{blue}{η}\)得到高频信息\(\textcolor{blue}{H}\)：
3.模型训练
作者选择使用混合扩散模型，去噪扩散GAN 作为主干网络，它采用一种特殊的 f 散度实例，称为软化反向 Kullback-Leibler 散度（SRKL）。它能够用不到 8 个步骤对图像进行采样，从而大大加快推理速度，而原始 DDPM 需要超过 1000 个步骤。它可以比原始DDPM快数百倍的推理速度实现接近原始DDPM的性能。训练是通过使用对抗性损失来匹配条件去噪模型 \(\textcolor{blue}{p_θ(s_{t−1}|s_t, H_η)}\) 和 \(\textcolor{blue}{q(s_{t−1}|s_t)}\) 来制定的，该对抗性损失可最小化每个去噪步骤的发散 Dadv：
SRKL 是一个与时间相关的判别器，它通过以下方程进行训练：
给定 SRKL，我们通过以下方式训练条件去噪扩散模型：
4.零样本图像翻译
在零样本图像转换中，模型在训练时仅在目标域数据上进行训练，因为来自源域的图像在训练期间不可用。该模型需要直接响应来自训练期间未观察到的领域的样本。零样本方法通常通过某种形式的辅助信息将观察到的类和未观察到的类关联起来，而在 CBCT 到 CT 图像转换任务中，我们发现它们的低频和高频信息大致相同，主要区别在于中频信息。两种模态/域之间的图像翻译任务的描述如下：
三、实验
1.数据集和预处理
头颈数据集：在头颈数据集上训练和测试了 FGDM。头颈部数据集包含患者特定的 CT 和 CBCT 图像，分为 70 个训练案例、9 个验证案例和 20 个测试案例，验证数据用于选择最佳性能检查点。 CBCT图像的体素间距为\(0.51×0.51×1.99mm^3\)，CT图像的体素间距为\(1.17×1.17×3.00mm^3\)。所有图像均以HU值[-1000,1000]截取，并线性归一化至[0,1]。我们使用带有 2D 切片的图像作为输入，所有切片大小调整为 192 × 192。为了测试，使用开源可变形图像配准工具箱将 CT 与 CBCT 对齐，并将配准的 CT 图像用作评估图像翻译性能的参考。
肺部数据集：为了进一步测试我们模型的零样本能力，我们将在头颈 CT 上训练的模型直接应用于肺部 CBCT 到 CT 的转换任务。我们总共收集了 18 名患者的肺部 CBCT 和 CT 图像，这些图像经历了与头颈部数据相同的预处理步骤，包括 HU 值截取、线性归一化和切片大小调整。与头颈部数据集类似，每个患者的肺部 CT 也变形配准到相应的 CBCT，作为评估参考。
处于危险中的器官 (OAR) 数据集：为了严格评估在具有挑战性和变化的条件下的性能，我们利用了美国医学物理学家协会 (AAPM) 胸部自动分割重大挑战的处于危险中的器官 (OAR) 数据集。该数据集来自 24 名患者，包含配对 CBCT 和解剖学 CT，包括食道、脊髓、心脏和肺。 OARs 数据集的预处理与头颈数据集和肺部数据集类似，其中涉及将 HU 值截取到 [-1000, 1000] 范围，然后线性归一化为 [0,1]。使用 2D 切片作为输入，所有切片大小调整为 192 × 192。在评估中包含该数据集有助于进一步评估 FGDM 在数据分布变化较多的场景中的稳健性和多功能性。
头部 MR T1 数据集：评估中使用的另一个数据集包括从两个不同存储库（中心）获取的头部 MR T1 图像：神经反馈 Skullstripped (NFBS)1 存储库和 MICCAI 2020 挑战：癌症传播的解剖脑屏障 (ABCs)2 存储库。选择该数据集是为了反映跨中心图像翻译的现实需求，其中来自不同机构的数据可能由于不同的成像协议而有所不同，并且需要协调一致以进行集体分析。 NFBS 数据集包含来自 125 名参与者的结构 T1 加权和匿名（污损）图像。该集合充当 FGDM 训练的源域数据。 ABCs 数据集包含 60 名诊断为胶质母细胞瘤和低级别胶质瘤的患者，作为测试的目标域数据。对于预处理，我们将图像强度剪裁在第 0 个百分位数和第 99 个百分位数之间，以消除异常值并标准化图像的动态范围。然后，我们将图像标准化为 [0,1] 并将所有切片的尺寸调整为 256 × 256。
2.实施细节
通过标准 Adam 优化器训练模型，在余弦退火学习率调度器下使用 10−4 的初始学习率，将最小学习率设置为 10−5。批量大小设置为 8，训练 200 个 epoch。在模型中，低频和高频的信号提取分别由前向扩散步数 T 和 Sobel 滤波阈值 η控制。为了训练一个对频率空间中源图像和目标图像之间的域不变区变化具有鲁棒性的模型，在训练过程中将 η 随机化为 1-25 的值，以提取不同级别的高频信号进行调节。对于低频信号，这种变化已经在逐步前向扩散过程中建立，可以通过改变测试源域图像上的前向扩散步数来直接进行测试时评估。在本研究中，在消融研究部分对这两个参数进行参数化研究后，作者选择 η 为 10，\(\tilde{T}\) 为 4 作为测试集。由于FGDM独特的设计，这种参数化研究是在测试阶段直接进行的，它对不同的数据输入具有高度的适应性和可调性，而不需要任何模型重新训练。所有网络都是使用 Pytorch 库实现的，在 NVIDIA Tesla V100 GPU 上运行了实验。
3.比较指标和方法
为了比较不同模型的图像翻译性能，作者测量了翻译图像的真实度和忠实度。为了真实起见，我们报告了翻译图像和目标参考图像之间广泛使用的 Frechet 初始距离分数 (FID) 。考虑到我们的目标参考图像是通过变形配准获得的，并且由于非变形引起的变化和变形配准误差而不完全是地面实况图像，我们通过将翻译图像与源 CBCT 图像和参考变形 CT 进行比较来评估忠实度图像，因为评估源 CBCT 图像结构保留的真实性至关重要。这种策略已被许多没有地面实况图像来计算忠实度的作品[11]、[12]、[39]使用。为了准确起见，我们报告了生成的图像与源图像和目标图像的 PSNR（峰值信噪比）和结构相似性指数度量（SSIM）。此外，为了进一步验证FGDM的性能，我们还计算了生成图像与目标图像在空间频域的MSE（均方误差）。为了与其他图像翻译方法进行全面比较，我们总共评估了三类广泛使用的图像翻译方法（基于 GAN、基于 VAE 和基于扩散模型）。其中，基于GAN和基于VAE的方法需要源域和目标域图像进行训练，而基于扩散模型的方法可以实现零样本，例如SDEdit