DPS算法解析：用扩散模型先验求解成像逆问题

1. 项目背景与核心思路

1.1 为什么成像逆问题需要新的求解思路

搞计算成像的人应该都有同感：我们手上的相机、显微镜、CT设备，本质上都在做一件事——把真实场景经过一套物理系统变成测量值。超分辨率、去模糊、去噪、CT重建、相位恢复、压缩感知，统统可以写进同一个数学框架：观测模型 y = A(x) + n。其中 x 是我们要恢复的原始图像，A 是前向成像算子，n 是噪声。

难点在于，A 这个算子往往是病态的。要么是欠定（测量数不够），要么是条件数极差（微小噪声被放大），要么是信息本身就丢了（比如超分里高频细节直接不存在）。传统方法一般加一个正则项，比如 Tikhonov、全变分（TV）、稀疏先验，把解空间往“看起来合理”的方向掰。但这类先验是手设计的，表征能力有限，遇到复杂结构（纹理、边缘、语义内容）经常力不从心。

所以要换思路：能不能让机器自己从海量自然图像里学出一个先验，再拿这个先验去约束逆问题的求解？这正是生成先验（Generative Prior）的出发点。扩散模型在后训练采样这一块，天然适合干这个活。

1.2 扩散模型凭什么能当先验用

先快速过一遍扩散模型的基本逻辑。前向过程就是不断给图像加高斯噪声，加到最后变成纯噪声。反向过程就是学一个神经网络去一步步去噪，从纯噪声里把图像“雕”出来。

关键在于：这个去噪网络在训练过程中，被迫学会了大量关于自然图像结构的统计知识——哪些像素组合是常见的，哪些边缘形状是合理的，哪些纹理模式会同时出现。所以它本质上是一个极其强大的隐式先验，比手设计的全变分不知道高到哪里去了。

把扩散模型用于逆问题求解，最直接的想法就是：在反向采样的每一步，除了让图像符合去噪网络的先验，还强迫它符合观测方程 y = A(x)。这样采样出来的结果，既真实又跟测量数据一致。

这里有一个技术难点：标准扩散采样过程中，我们手里只有带噪的 x_t，没有干净图像 x_0。观测模型 y = A(x_0) + n 是定义在 x_0 上的，怎么在每一步把观测约束塞进去？DPS 给出的答案是：先用当前 x_t 估计一个 x_0（叫“干净图像预测”），然后在这个估计上施加观测约束，再把梯度传回 x_t。讲白了就是把“后验采样”的梯度用近似方式算出来，每一步都往测量一致的方向推一把。

一句话总结：DPS 的核心贡献是把扩散先验和后验观测约束拆开，用一步 x_0 预测构建近似梯度，从而把任意线性（甚至部分非线性）成像逆问题塞进扩散采样的迭代框架里。

2. 数学原理与算法拆解

2.1 从贝叶斯后验到采样过程

逆问题的贝叶斯视角很优雅。我们想要求解的是后验分布 p(x|y)，根据贝叶斯公式：

p(x|y) ∝ p(y|x) · p(x)

其中 p(y|x) 是似然项，编码了测量模型和噪声分布；p(x) 是先验项，告诉算法什么样的图像是“自然的”。传统变分法直接算这个后验的最大值（MAP），但 MAP 估计容易丢失多样性，而且对多模态后验无能为力。更好的做法是——直接从这个后验里采样，得到多个合理解。

扩散模型天生是采样器。标准扩散采样是从 p(x_T)（纯噪声）开始，逐步采样 p(x_{t-1}|x_t)。如果能把观测信息 y 注入每一步的条件分布 p(x_{t-1}|x_t, y)，那最终采样出来的 x_0 就会同时满足先验和观测约束。

用贝叶斯公式展开：

p(x_{t-1}|x_t, y) ∝ p(x_{t-1}|x_t) · p(y|x_{t-1})

这里 p(y|x_{t-1}) 是似然项，但它跟模型观测 y = A(x_0) + n 之间还隔着一个距离。我们需要把 x_{t-1} 映射到观测空间。

2.2 DPS 的近似梯度推导

DPS 的关键技巧是对数似然梯度的近似。考虑我们要在反向采样过程中对 p(x_t | y) 取梯度，标准做法是：

∇_{x_t} log p(y|x_t) = ∇_{x_t} log E_{p(x_0|x_t)} [ p(y|x_0) ]

这个期望没法解析算。DPS 做了两件事：

第一，用 Tweedie 公式把条件期望换成点估计。扩散模型的去噪网络 ε_θ(x_t, t) 实际上给出了 E[x_0 | x_t] 的估计：

x̂_0 = (x_t - √(1-ᾱ_t) · ε_θ(x_t, t)) / √ᾱ_t

这一步不是近似，是 Tweedie 公式的精确结果。

第二，用点估计 x̂_0 代替对 x_0 的积分，把梯度近似为：

∇_{x_t} log p(y|x_t) ≈ ∇_{x_t} log p(y|x̂_0(x_t))

如果噪声是高斯的，似然项 p(y|x_0) ∝ exp(-||y - A(x_0)||² / 2σ²)，那么：

∇_{x_t} log p(y|x̂_0(x_t)) ≈ - (1/σ²) · A^T (A(x̂_0) - y) · (∂x̂_0/∂x_t)

这里 A^T 是 A 的伴随算子。这个近似相当于：先看当前估计的 x̂_0 跟观测 y 差多远，然后把差值通过伴随算子投影回图像空间，再沿着 x̂_0 对 x_t 的依赖关系反向传播。

有人会问：这个近似靠谱吗？毕竟点估计代替期望是有误差的。原文实验证明了在多个任务上这个近似足够好用，因为扩散模型的去噪网络对 x_0 的估计相当准确，尤其在采样后期。我自己的实践体会是，这个近似在早期采样步（图像还很模糊时）会有一定偏差，但迭代修正会逐步拉回来。

2.3 算法流程总览

把整个算法串起来看，DPS 的推理流程非常清晰：

1. 从高斯噪声采样 x_T。
2. 对 t = T, T-1, ..., 1： a. 用去噪网络估计噪声 ε_θ(x_t, t)，并根据 Tweedie 公式计算 x̂_0。 b. 计算观测残差 r = y - A(x̂_0)。 c. 把残差通过伴随算子投影回图像域，得到梯度项 g = A^T(r)。 d. 把 g 经过 ∂x̂_0/∂x_t 的反传，乘上步长 ξ_t，加到采样过程中。 e. 执行一步标准 DDPM/DDIM 采样更新。
3. 输出 x_0。

注意第 d 步的实际实现不需要显式求雅可比矩阵，直接用自动微分工具（PyTorch 的 autograd）就能搞定。这也是 DPS 能快速推广的关键——换一个前向算子 A，只需要定义 A 和 A^T，其他逻辑基本不动。

算法伪代码可以写成这样：

# 核心采样循环（伪代码） for t in reversed(range(1, T)): t_tensor = torch.full((batch_size,), t, device=device) # 1. 预测噪声并计算 x0 估计 eps_pred = model(x_t, t_tensor) x0_hat = (x_t - torch.sqrt(1 - alpha_bar[t]) * eps_pred) / torch.sqrt(alpha_bar[t]) # 2. 计算观测梯度 grad_term = A_T(A(x0_hat) - y) # 伴随算子 # 通过 x0_hat 对 x_t 的依赖反传 grad_x0 = torch.autograd.grad(x0_hat, x_t, grad_outputs=grad_term)[0] # 3. 更新采样 x_0 = x_t - beta_t / torch.sqrt(1 - alpha_bar[t]) * grad_x0 # 简化形式 # 再加噪声/DDIM 步骤

3. 实操配置与关键参数

3.1 前向算子与伴随算子的正确写法

实操中第一个坑就是算子 A 和伴随算子 A^T 的实现。很多逆问题算法只需要 A 和 A^T 成对出现，但写反了或者写错了，结果会非常离谱。

以超分辨率为例。如果我们想把 256×256 图像降采样到 64×64，A 就是一个降采样算子。在 PyTorch 里最省事的是用F.interpolate：

def A_sr(x): # x: (B, C, H, W) return F.interpolate(x, scale_factor=0.25, mode='bicubic', align_corners=False) def AT_sr(y): # 上采样回原图尺寸，注意这里只是近似伴随 return F.interpolate(y, scale_factor=4, mode='bicubic', align_corners=False)

严格来说，bicubic 的降采样和上采样并不互为精确伴随，但在 DPS 框架里这个近似也能收敛，只是质量略有损失。如果追求严格，应该预计算降采样矩阵，用矩阵转置实现精确伴随。

另一个典型场景是修复（inpainting）。A 是一个掩膜算子，只保留可见像素：

def A_inpaint(x, mask): return x * mask # mask 中 1 表示可见 def AT_inpaint(y, mask): return y * mask # 伴随算子还是掩膜操作

修复任务的伴随算子特别简单，但要注意：残差 r = y - A(x̂_0) 只在掩膜为 1 的地方非零，所以梯度只回传到掩膜可见区域，其余区域完全靠先验生成。这就是为什么修复看起来“凭空造内容”——其实是在靠扩散模型的先验填空。

3.2 关键超参数调节指南

DPS 里最核心的超参数是梯度步长 ζ（论文里通常写作 ζ_t，或者代码里叫grad_scale）。它的作用是控制观测约束的强度。取值太小，结果跟测量数据对不上；取值太大，采样过程不稳定，甚至发散。

我的实践经验：用 DDPM 采样时，梯度应该乘以系数，形式为：

x_t = x_t - ζ_t · ∇_{x_t} ||y - A(x̂_0)||²

ζ_t 需要随噪声水平缩放。论文里的建议是 ζ_t = ζ / sqrt(||A^T(y - A(x̂_0))||)，归一化到相对稳定的尺度。实际代码里，通常设置一个基础步长比如 0.1 到 1.0，然后除以梯度范数。这个归一化技巧非常重要，因为不同任务、不同测量维度下，残差量级差异巨大。

我调参的心得：

• 先固定其他条件，只扫描 ζ。从 0.1 开始，看重建结果的 PSNR 和视觉质量。
• 如果结果噪声大、跟观测差异大，增大 ζ。
• 如果中间过程出现 NaN 或剧烈跳变，减小 ζ，同时检查学习率是否压过了 DDPM 自身的采样步。
• 配合 DDIM 加速采样时，ζ 可以稍微调大一点，因为 DDIM 本身更稳定。

第二个重要参数是噪声控制。DPS 假设观测噪声 y = A(x) + n，n 的方差 σ² 也影响梯度权重。在实现中，||y - A(x̂_0)||² / (2σ²)分母上的 σ² 越小，约束越强。如果实际噪声是 σ=0.01，但代码里写死了 σ=0.05，那结果会欠拟合观测。建议把 σ 设为可配置项，按实际数据调整。

提示：大量 DPS 复现效果不理想，八成是 ζ 没归一化，或者 σ 跟实际噪声不匹配。先检查这两个参数，再怀疑网络权重。这是我在多个任务上排查过的经验总结。

3.3 采样步数与调度策略

DPS 原文用的是完整 1000 步 DDPM 采样，效果最好但耗时最长。实跑的时候可以搭配 DDIM 加速：

• 完整 1000 步：质量最高，适合离线处理和高要求场景。每张图在单卡上大概要几十秒到几分钟，视分辨率和网络规模而定。
• DDIM 200 步：质量损失很小，速度提升约 5 倍。大多数情况下推荐先用这个配置跑通。
• DDIM 50 步：速度快，但 DPS 的梯度近似误差会被放大，容易出现细节崩坏。

还有一个技巧：可以把梯度项只加在采样过程的前 60% 步（噪声大时先让先验生效，后段再用观测约束精细调整）。我实测这个策略在部分任务上能改善细节，但论文没有采用，属于非标准技巧。新手不建议改，先用默认配置跑通，再慢慢探索变体。

4. 实战案例：超分辨率与修复

4.1 超分辨率实测：256 到 64 再回来

我拿一个模拟项目做过验证，取 256×256 的自然图像，用 4 倍双三次降采样得到 64×64 的观测，再跑 DPS 恢复。整个过程分成三步：

第一步，准备预训练扩散模型。直接加载社区公开的预训练权重，不需要针对逆问题重新训练。这是 DPS 最爽的地方——一个通用的预训练模型就能处理一大类任务，完全不需要为每个任务单独训练。

第二步，定义前向算子。用F.interpolate做 bicubic 降采样，伴随算子用上采样。这里要特别注意，bicubic 不是严格线性的吗？它不是严格线性算子，因为align_corners参数的引入让映射关系不那么直白。如果要严格，可以采用矩阵乘法形式的预计算插值权重。

第三步，跑采样。关键代码片段：

# 核心梯度更新模块 def dps_grad_step(x_t, t, y, A, AT, model, zeta, sigma_y): x_t.requires_grad_(True) # 预测噪声和 x0 eps = model(x_t, t) x0_hat = (x_t - torch.sqrt(1 - alpha_bar[t]) * eps) / torch.sqrt(alpha_bar[t]) # 观测残差 residual = A(x0_hat) - y grad_term = 2 * AT(residual) / (sigma_y ** 2) # 反传 grad = torch.autograd.grad(x0_hat, x_t, grad_outputs=grad_term)[0] return grad, x0_hat

第四步，评估。超分任务上，用 64×64 低分辨率输入恢复到 256×256，相比裸 bicubic 插值，PSNR 直接提升 4~6 dB 是常态。在纹理丰富区域（比如树叶、毛发），DPS 生成的细节真实感远超传统插值。而且多次采样能产生不同的细节版本，这对艺术创作和增强类应用特别有价值。

一个明显的短板是速度。每张 256 图像跑 200 步 DDIM，在中等显卡上大约需要 40 秒到 1 分钟。如果想用到视频或者交互式场景，需要配合加速采样和模型蒸馏，那是另一条线了。

4.2 修复任务：掩膜与内容生成

修复任务是 DPS 最直观、最惊艳的应用场景。给定一张图和一块破损区域（掩膜），模型需要在看不到任何内容的区域“想象”出合理的内容。

我的实测：在 128×128 图像中间挖掉一块 32×32 的矩形区域，DPS 能生成跟周围环境协调一致的纹理。关键点在于——掩膜区域完全没有观测约束，算法输出的内容完全由扩散先验决定。如果先验足够强（模型训练数据丰富），生成结果会自然融入周围。

修复任务调参有两个注意点。第一，掩膜的过渡边缘要平滑。直接二值化掩膜会在边缘处产生接缝，建议对掩膜做少量高斯模糊，让观测约束在边缘处渐变。这个技巧在代码里就一行：

mask_smooth = F.gaussian_blur(mask.float(), kernel_size=5, sigma=2.0)

第二，如果破损区域面积很大（超过图像面积 30%），单次重建可能出现内容不连续。这时候可以增大 ζ，或者考虑分块修复策略——把大区域拆成小块逐块修复，每块保留边缘重叠区做约束融合。

4.3 其他算子的扩展思路

DPS 框架对逆问题的类型没有太多限制。除了超分和修复，还可以做：

• 图像去模糊：A 是模糊核的卷积运算，A^T 是相关运算。在 PyTorch 里用F.conv2d和F.conv_transpose2d实现。关键是把模糊核参数绑定，确保 A 和 A^T 共享同一组核参数。
• CT 重建：A 是 Radon 变换（投影），A^T 是滤波反投影算子。可以利用已有的torch-radon类库实现。DPS 在稀疏视角 CT 上的效果比传统 FBP 好得多，伪影大幅减少。
• 相位恢复：测量是强度信息，丢失了相位，前向算子是非线性的。DPS 的梯度近似框架依然能工作，只是 A 的定义变成非线性映射，autograd 自动计算雅可比向量积，实现上几乎无感。

非线性算子的支持是 DPS 区别于很多传统方法的重要优势。传统 MRI 重建、CT 重建等方法对算子线性性质有依赖，而 DPS 只需要前向算子可微即可。这一点让它的应用面宽了很多。

5. 常见问题与性能优化

5.1 踩坑实录：重影、模糊和发散

我在多次复现和调优中遇到了几类高频问题，这里直接整理成速查表：

最坑的是“发散但不报错”——输出结果看起来正常，但 PSNR 极低，反复检查代码没发现问题。后来发现是torch.autograd.grad时把x_t的requires_grad设成了 True，但每一步循环重新创建了计算图，导致梯度内存累积。解决办法是每次迭代后手动释放计算图，确保 x_t 是叶子变量。

另一个容易踩的坑是alpha_bar的数值精度。直接按公式计算np.cumprod(1 - beta)会随着步数增加产生累积误差，尤其在使用 FP16 精度时。推荐的做法是预计算好alpha_bar存成 float32 张量，采样时再转到对应设备。

5.2 显存与速度优化

DPS 比普通扩散采样多了一次梯度反传，显存开销大约增加 50%。我在 8GB 显存的卡上跑 256×256 图像、200 步，勉强能过。几个优化手段：

• 梯度检查点：对去噪网络和 x0 预测部分使用torch.utils.checkpoint，用少量计算换显存。
• 混合精度：网络推理用 FP16，梯度反传也走 FP16，但不建议整个采样循环都用 FP16，会出现数值不稳。
• 减少批次大小：逆问题重建大多可以逐张处理，batch size 用 1 就够了。

速度优化方面，推荐把 DDIM 步数从 1000 降到 200，配合eta=0的确定性采样。DPS 在 DDIM 上的表现比 DDPM 略差一点，但速度差 5 倍，综合收益是正的。如果对质量要求极高，再切回 DDPM。

5.3 与 Condition 类扩散模型的对比思考

DPS 解决问题的路径是“采样时注入观测梯度”，它的优势是没有改变扩散模型本身的训练目标，预训练权重直接复用，所以零成本适配新任务。另一类做法是在训练阶段就把条件信息（低分辨率、掩膜等）融入模型，比如超分模型、修复模型的专用架构。专用模型质量好、速度快，但每换一个逆问题就要重新训练或微调。

在实际工程选型中，我的建议是：

• 如果只有一两个固定逆问题、对时延敏感，选专用条件扩散模型。
• 如果有很多不同类型的逆问题、需要快速验证效果，DPS 是正确的起点。
• 研究阶段先用 DPS 验证先验有效性，再做专用模型优化，是很高效的路径。

个人体验最深刻的一点：DPS 把“算法归算法，先验归先验”解耦得干净。换任务时只要改算子定义，其他一切不动。这种模块化程度在同代方法里相当少见，也让它成为我日常工作流里一个非常趁手的底座。后续如果要做模型蒸馏或采样加速，在这个底座上扩展也很顺手。

DPS 之后还有一系列改进，比如通过更好的 x0 估计器降低梯度噪声，或者用块坐标下降处理高维问题。不过那些都是锦上添花，先把基础版跑通、把参数吃透，才是最务实的路径。以上是我在实际项目里反复折腾后的经验沉淀，希望对做计算成像的朋友有实际帮助。