1. 从一张糊图说起：为什么成像逆问题需要生成先验

如果你拍过夜景、做过显微重建，或者处理过遥感图像，大概率遇到过这种场景：手里的观测数据少得可怜，噪声还大，但你要反推出一张干净、细节丰富的高质量图像。这就是典型的成像逆问题——已知退化过程（模糊核、降采样矩阵、编码掩膜等），要从观测值反解原始信号。

传统解法无非两类。一类是解析法，比如维纳滤波、Tikhonov正则化，思路是把逆问题写成最小二乘加一个手工设计的正则项，求解快但先验太弱，结果往往偏平滑、丢纹理。另一类是迭代法，比如TV正则、稀疏编码，靠人工设计的稀疏先验去约束解空间，效果比解析法好，但先验表达能力有限，遇到复杂自然图像就力不从心。

问题的根子在于：逆问题本身是病态的，观测数量远小于未知量，解不唯一。你必须注入额外的先验知识才能挑出那个“合理”的解。而手工设计的先验，无论是全变分、小波稀疏还是字典学习，本质上都是对自然图像统计规律的一种粗糙近似。它们能描述“图像应该是分段平滑的”或者“图像在小波域是稀疏的”，但描述不了“人脸应该有两只眼睛、一个鼻子”这种高层语义结构。

扩散模型的出现改变了这个局面。它在大规模图像数据上训练，学到的生成先验极其强大——你给它一个纯噪声，它能一步步去噪生成一张逼真的图。这意味着它隐式地掌握了自然图像的深层分布。那么一个自然的想法就冒出来了：能不能把这个生成先验拿来解逆问题？让扩散模型在去噪的过程中，同时受到观测数据的约束，最终生成一个既逼真又符合观测的解？

这就是**扩散后验采样（Diffusion Posterior Sampling，简称DPS）**要解决的核心问题。DPS的核心贡献在于，它给出了一种优雅的方式，把扩散模型的去噪过程和逆问题的数据保真项结合起来，不需要针对每个逆问题重新训练模型，也不需要知道退化算子的显式形式，就能实现高质量的图像重建。

这篇文章适合谁看？如果你正在做计算成像、图像恢复、医学影像重建，或者对扩散模型的应用感兴趣，想搞清楚“生成先验怎么跟逆问题结合”这件事，那接下来的内容应该对你有用。我会从设计思路、核心原理、实操细节到常见问题，把DPS这套方法拆开讲清楚。

2. DPS的整体设计思路：为什么是“后验采样”而不是“投影”

2.1 逆问题的贝叶斯视角

先把问题形式化。成像逆问题可以写成：

$$y = \mathcal{A}(x) + n$$

其中 $x$ 是真实图像，$\mathcal{A}$ 是已知的退化算子（比如模糊、降采样、掩膜编码），$n$ 是噪声，$y$ 是观测。我们的目标是从 $y$ 恢复 $x$。

从贝叶斯角度看，我们想要的是后验分布 $p(x|y)$。有了后验，就能取最大后验估计（MAP）或者后验均值作为重建结果。根据贝叶斯公式：

$$p(x|y) \propto p(y|x) \cdot p(x)$$

这里 $p(y|x)$ 是似然项，描述观测数据对解的约束；$p(x)$ 是先验项，描述我们对自然图像的认知。传统方法的问题就出在 $p(x)$ 上——手工设计的先验太弱。

扩散模型恰好提供了一个强大的 $p(x)$。它通过前向加噪过程把图像逐步破坏成噪声，再训练一个网络学习反向去噪过程。训练完成后，这个反向过程就隐式地定义了 $p(x)$ 的采样器。

2.2 为什么不能简单地在扩散采样中做投影

最朴素的想法是：在扩散模型的每一步去噪之后，把结果往观测数据上“投影”一下，比如做一次数据一致性修正。这种做法在早期工作中确实出现过，但问题很明显。

扩散采样是一个逐步去噪的过程，每一步的输出都带有随机性，而且中间步骤的图像并不是最终解，它只是去噪轨迹上的一个中间状态。如果你在中间步骤强行做投影，相当于用观测数据去“纠正”一个还没成型的中间结果，很容易破坏扩散模型本身的采样动力学，导致最终结果要么不真实，要么不满足观测。

更关键的是，投影操作通常假设退化算子是可逆的或者容易求伪逆的。对于复杂的退化（比如非线性编码、未知噪声水平），投影根本没法直接做。

2.3 DPS的核心洞察：用Tweedie公式搭桥

DPS的巧妙之处在于，它不去强行投影，而是利用扩散模型去噪器与后验得分之间的数学联系，把数据保真项自然地融入到采样过程中。

具体来说，扩散模型在每一步去噪时，实际上是在估计当前噪声水平下的干净图像。这个估计可以通过Tweedie公式与得分函数联系起来。而得分函数又恰好是后验分布梯度的组成部分。DPS利用这个关系，推导出一个近似的后验得分表达式，从而在每一步采样中同时考虑生成先验和观测数据。

这样做的好处是：不需要重新训练模型，不需要知道退化算子的显式形式（只需要能计算 $\mathcal{A}(\hat{x})$ 和它的梯度），而且适用于各种线性甚至非线性逆问题。一套框架，多个任务通用，这是DPS最吸引人的地方。

3. 核心原理拆解：DPS到底在每一步做了什么

3.1 扩散模型的前向与反向过程回顾

先快速过一遍扩散模型的基础。前向过程是一个马尔可夫链，逐步向图像加高斯噪声：

$$x_t = \sqrt{\bar{\alpha}_t} x_0 + \sqrt{1 - \bar{\alpha}_t} \epsilon, \quad \epsilon \sim \mathcal{N}(0, I)$$

其中 $\bar{\alpha}_t$ 是累积噪声系数，$t$ 从0到T。当T足够大时，$x_T$ 近似为标准高斯噪声。

反向过程则是从 $x_T$ 出发，逐步去噪得到 $x_0$。训练时，网络 $\epsilon_\theta(x_t, t)$ 学习预测当前噪声。采样时，利用这个预测噪声来估计干净图像：

$$\hat{x}_0 = \frac{x_t - \sqrt{1 - \bar{\alpha}t} \epsilon\theta(x_t, t)}{\sqrt{\bar{\alpha}_t}}$$

这个 $\hat{x}_0$ 就是扩散模型对当前步骤干净图像的估计，也是DPS中用来计算数据保真项的关键量。

3.2 后验得分的分解

我们真正想要的是后验分布 $p(x_t|y)$ 的得分函数 $\nabla_{x_t} \log p(x_t|y)$。根据贝叶斯公式：

$$\nabla_{x_t} \log p(x_t|y) = \nabla_{x_t} \log p(x_t) + \nabla_{x_t} \log p(y|x_t)$$

第一项是先验得分，扩散模型已经学到了。第二项是似然得分，描述观测数据对当前噪声图像的影响。

难点在于 $p(y|x_t)$ 的计算。因为 $x_t$ 是带噪的中间状态，而观测 $y$ 是针对干净图像 $x_0$ 的。DPS做了一个近似：用扩散模型估计的 $\hat{x}_0$ 来代替 $x_0$，然后假设

$$p(y|x_t) \approx p(y|\hat{x}_0)$$

在噪声水平较低时，这个近似是合理的。进一步假设观测噪声是高斯噪声，那么

$$\log p(y|\hat{x}_0) \propto -|y - \mathcal{A}(\hat{x}_0)|^2 / 2\sigma^2$$

于是似然得分就可以近似为：

$$\nabla_{x_t} \log p(y|x_t) \approx -\frac{1}{\sigma^2} \nabla_{x_t} |y - \mathcal{A}(\hat{x}_0)|^2$$

这个梯度可以通过链式法则计算，因为 $\hat{x}_0$ 是 $x_t$ 的函数。

3.3 采样更新公式

把先验得分和似然得分合起来，DPS的采样更新可以写成：

$$x_{t-1} = \text{DDIM_step}(x_t, \epsilon_\theta) - \zeta_t \nabla_{x_t} |y - \mathcal{A}(\hat{x}_0)|^2$$

其中 $\zeta_t$ 是一个步长系数，控制数据保真项的强度。DDIM_step是标准的扩散采样步骤，负责去噪；后面那一项则是数据一致性修正，把采样轨迹往满足观测的方向拉。

这个公式看起来简单，但里面有几个关键细节需要处理好，否则效果会大打折扣。

3.4 步长系数的选择

$\zeta_t$ 的选择直接影响重建质量。如果太大，数据保真项会主导采样过程，生成结果可能不真实；如果太小，观测约束不够，重建结果可能偏离真实解。

DPS原文给出的建议是：

$$\zeta_t = \frac{\zeta}{|\nabla_{x_t} |y - \mathcal{A}(\hat{x}_0)|^2|}$$

也就是用梯度的范数做归一化，让每一步的修正幅度相对稳定。$\zeta$ 是一个全局超参数，通常需要根据具体任务调。经验上，$\zeta$ 在0.1到1.0之间比较常见，噪声大的任务可以适当调大。

注意：这个归一化很重要。如果不做归一化，梯度范数在不同步骤可能差几个数量级，导致采样过程不稳定。

4. 实操过程：从零实现一个DPS采样器

4.1 环境准备与依赖

我平时用PyTorch做这类实验，扩散模型部分可以直接用预训练的DDPM或DDIM模型。以下是一个最小化的依赖列表：

pip install torch torchvision numpy scipy

如果你要用现成的扩散模型权重，可以加载OpenAI的guided-diffusion或者HuggingFace的diffusers库。我这里为了讲清楚原理，会写一个简化版的DDIM采样器，方便你理解每一步在做什么。

4.2 退化算子的定义

不同的逆问题对应不同的 $\mathcal{A}$。我整理了几个常见任务的算子实现方式：

以超分辨率为例，退化算子可以写成：

import torch import torch.nn.functional as F def downsample_operator(x, scale_factor=4): return F.avg_pool2d(x, kernel_size=scale_factor, stride=scale_factor)

对应的梯度计算需要用到自动微分，后面会讲。

4.3 核心采样循环

下面是DPS采样循环的伪代码，我加了详细注释：

def dps_sample(model, y, A, A_grad, num_steps=1000, zeta=0.3): """ model: 预训练扩散模型，输入(x_t, t)输出预测噪声 y: 观测数据 A: 退化算子 A_grad: 退化算子的梯度函数 num_steps: 采样步数 zeta: 数据保真项步长 """ # 初始化：从纯噪声开始 x = torch.randn_like(y) # 注意：这里需要根据任务调整尺寸 for t in reversed(range(num_steps)): t_tensor = torch.full((x.shape[0],), t, device=x.device) # 第一步：扩散模型预测噪声 with torch.no_grad(): eps_pred = model(x, t_tensor) # 第二步：估计干净图像 alpha_bar_t = get_alpha_bar(t) x0_hat = (x - torch.sqrt(1 - alpha_bar_t) * eps_pred) / torch.sqrt(alpha_bar_t) # 第三步：计算数据保真项梯度 x0_hat = x0_hat.detach().requires_grad_(True) residual = y - A(x0_hat) loss = torch.sum(residual ** 2) grad = torch.autograd.grad(loss, x0_hat)[0] # 第四步：DDIM去噪步骤 x = ddim_step(x, eps_pred, t) # 第五步：数据一致性修正 grad_norm = torch.norm(grad) + 1e-8 x = x - zeta * grad / grad_norm # 可选：对x做clip，防止数值溢出 x = torch.clamp(x, -1, 1) return x

这段代码里几个关键点值得展开说。

关于x0_hat的detach和requires_grad：我们需要对x0_hat求梯度，但不想让梯度回传到扩散模型参数上，所以先detach再开requires_grad。这样计算图只保留x0_hat到loss的部分，节省显存。

关于DDIM_step：这是标准的DDIM更新，根据当前x和预测噪声计算下一步的x。具体实现可以参考DDIM原文，核心是用确定性采样替代DDPM的随机采样，步数可以大幅减少。

关于梯度归一化：我用了L2范数归一化，这是DPS原文的做法。你也可以尝试其他归一化方式，比如用梯度的标准差，效果可能因任务而异。

4.4 参数选择与调优经验

跑通代码只是第一步，真正影响效果的是参数调优。我整理了几个关键参数的经验值：

我实测下来，zeta这个参数最敏感。在超分辨率任务上，zeta=0.3左右比较稳；在去模糊任务上，zeta=0.5效果更好。建议你先固定其他参数，单独扫zeta，画一条PSNR随zeta变化的曲线，找到峰值。

实操心得：如果重建结果出现明显的伪影或者不自然纹理，通常是zeta太小，数据保真不够；如果结果过于平滑、丢失细节，可能是zeta太大，生成先验被压制了。

5. 常见问题与排查技巧实录

5.1 重建结果全是噪声或完全偏离

这是最常见的问题，通常有几个原因。

原因一：退化算子实现有误。检查A(x)的输入输出尺寸是否匹配，卷积的padding是否正确，降采样的scale是否对得上。我踩过一次坑，avg_pool的kernel_size和stride设成不一样，导致尺寸对不上，梯度计算全乱。

原因二：alpha_bar的调度不对。扩散模型的噪声调度必须和训练时一致。如果你用的是预训练模型，一定要用模型对应的alpha_bar序列，不能自己随便设。

原因三：zeta太大导致发散。先把zeta调到0.01试试，如果结果正常了，说明是步长问题。然后再慢慢往上加。

5.2 重建结果模糊、细节丢失

这种情况通常是数据保真项太强，把生成先验压制了。可以尝试：

• 减小zeta
• 增加采样步数
• 检查噪声水平估计是否偏大（噪声水平越大，似然项权重越小，但DPS的近似可能不准）

另一个可能的原因是扩散模型本身的能力不够。如果你用的是在CIFAR-10上训练的模型去重建ImageNet图像，效果肯定好不了。生成先验的覆盖范围要和任务匹配。

5.3 采样速度太慢

DPS的主要计算开销在两块：扩散模型前向和梯度计算。扩散模型前向没法省，但梯度计算可以优化。

• 用DDIM替代DDPM，步数从1000降到100
• 梯度计算时只对x0_hat求导，不要对整个计算图求导
• 如果退化算子是线性的，可以解析计算梯度，不用自动微分

我实测过，在单张RTX 3090上，256x256图像、100步DDIM、超分辨率任务，一次重建大约需要15秒。如果降到50步，可以压到8秒左右，质量损失不大。

5.4 不同任务的效果差异

DPS不是万能的，不同逆问题上的表现差异很大。根据我的经验：

相位恢复这类非线性强的任务，DPS的近似假设不太成立，效果会打折扣。这时候可能需要考虑更复杂的后验采样方法，或者针对具体任务做改进。

5.5 常见问题速查表

6. 几个值得深挖的扩展方向

DPS打开了一扇门，但远不是终点。我在实际项目中遇到过几个值得继续探索的方向，这里分享出来供参考。

第一个方向是改进似然近似。DPS用 $\hat{x}_0$ 代替 $x_0$ 来计算似然，这在噪声水平低时还行，但在早期采样步骤（噪声很大）时误差明显。有些后续工作尝试用更精确的近似，比如考虑 $\hat{x}_0$ 的不确定性，或者用多个采样估计。如果你在做高噪声任务，这块值得花时间。

第二个方向是自适应步长。目前zeta是全局固定的，但不同步骤对数据保真项的需求不一样。早期步骤应该更依赖先验，后期步骤应该更依赖观测。设计一个随t变化的zeta调度，可能比固定值效果更好。我试过简单的线性衰减，有一定提升，但还没系统调优。

第三个方向是结合任务特定先验。DPS用的是通用图像先验，如果你做的是医学影像、遥感图像或者特定类型的自然图像，可以在扩散模型微调时注入领域知识，或者在采样时加入任务特定的约束。比如医学图像重建中，可以加入解剖结构的约束。

第四个方向是加速采样。100步DDIM还是慢，对于实时应用不够。可以尝试蒸馏、一致性模型或者少步采样技术。不过加速往往会损失一些质量，需要权衡。

第五个方向是理论分析。DPS的近似到底在什么条件下成立？误差上界是多少？这些问题目前还没有完整的理论回答。如果你偏理论方向，这是一个有价值的问题。

7. 我在实际项目中的几点体会

最后分享几个踩坑之后总结的经验，都是文档里不会写的。

关于噪声水平的估计。DPS假设观测噪声是已知的，但实际中往往不知道。我通常的做法是先用一个小网络估计噪声水平，或者把噪声水平也当作超参数，用验证集调。如果噪声水平设得偏大，重建结果会偏平滑；设得偏小，会保留更多噪声。这个参数对最终效果影响很大，值得多花时间。

关于初始化的选择。标准DPS从纯高斯噪声开始采样。但在某些任务中，从观测的某种粗估计开始（比如用传统方法先做一个初始重建，然后加噪到某个中间步骤）可以加速收敛，效果也更好。这个技巧在超分辨率任务上特别有用。

关于batch内的多样性。如果你一次重建多张图像，每张从不同的随机噪声开始，可以得到多个不同的解。这些解的均值通常比单个解更稳定，方差则可以用来评估不确定性。这个特性在医学影像中很有价值，医生可以看到重建的置信区间。

关于计算资源的规划。DPS的计算量大约是普通扩散采样的1.5到2倍，因为多了梯度计算。如果你要处理大批量数据，显存和时间的预算要提前算好。我的经验是，256x256图像、100步采样，batch size=4，大约需要8GB显存。512x512的话，batch size=1就需要10GB以上。

关于代码调试。DPS的调试比普通扩散模型麻烦，因为涉及梯度计算和算子实现。我的建议是先用一个极简的退化算子（比如恒等映射）验证采样器本身是否正确。如果恒等映射下能重建出原图，说明采样器没问题，再逐步增加算子复杂度。这个自顶向下的调试策略帮我省了很多时间。

这个方向后续还可以这样扩展：把DPS和变分推断结合起来，用变分分布去近似后验，可能比单点估计更稳定；或者把DPS用到视频重建中，利用时间维度的先验；再或者探索DPS在三维重建、计算摄影中的应用。生成先验和逆问题的结合才刚刚开始，后面还有很多有意思的问题值得做。