做贝叶斯推断和概率机器学习的同学大概率会在某个项目里遇到这种尴尬后验分布本身很干净但它的势能函数是非光滑的。最常见的就是 L1 正则也就是拉普拉斯先验再往深走还有全变差正则、带约束的模型、稀疏深度网络这类带非光滑惩罚项的问题。这时候标准的 Langevin 采样Langevin Monte Carlo / ULA会当场失效因为算法第一步要求你计算梯度的位置可能恰好落在导数不存在的点上。常见的绕路方案有三种把 |x| 改成 \sqrt{x^2\varepsilon} 之类的平滑近似、直接用次梯度硬算、或者引入辅助变量做增广。前两种的实现成本低但都会引入难以量化的偏差第三种数学上干净但会在变量维度上付出代价。这篇文章想谈的是另一条更优雅的技术路线用 Moreau--Yosida 正则化把非光滑势能磨平然后继续跑 Unadjusted Langevin SamplingULA。近年来这个方向在采样理论里非常活跃而最近一批论文又在思考一个更细粒度的问题当采样轨迹大部分时间都“躲开”了非光滑区域时传统的迭代复杂度上界是不是太悲观了于是就有了“Active-Trace Complexity Bounds”这类新度量。读完这篇文章你可以搞清楚三件事第一Moreau--Yosida 到底做了什么为什么它能替掉“求梯度”第二MYULA 算法长什么样怎么用 Python 在几十行内跑通第三所谓 active-trace 复杂度是在什么背景下提出来的它对工程判断有什么实际意义。整体不涉及太长篇的证明但会把需要的数学记号解释到能看懂论文的程度。1. 非光滑势能上的 Langevin 采样为什么是个真问题先回到问题本身。很多贝叶斯模型的后验可以写成[ \pi(x) \propto \exp(-U(x)), ]其中 (U:\mathbb{R}^d \to \mathbb{R}) 叫作势能函数。比如线性回归的贝叶斯版本中(U) 来自负对数似然加负对数先验如果先验是拉普拉斯分布那 (U) 里就有一个 (|x|_1) 项它在坐标轴处不可导。如果我们想用梯度类 MCMC 采样最自然的工具是 ULA。它的每一步只做两件事沿着负梯度走一小步再加一个高斯噪声。问题随之而来(U(x)|x|) 在 (x0) 处没有梯度。如果后验的支撑集足够复杂采样器会不断撞到这些不可导点。有人会反驳不可导点只是一个零测集随机游走“撞上”的概率也是零那不就行了吗问题在于离散化的 ULA 不是连续过程它每次从当前点起跳一旦当前迭代恰好落在这个集合附近或者势能本身不是凸的、存在折痕状结构算法行为就无法用光滑情形的那套分析来保证。工程上更麻烦真实问题里的不可导通常不是单个点而是一整片流形、一整条折线比如 (|x|_1) 的不可导区域是坐标平面。再看次梯度方案。把 (\nabla U) 换成某个次梯度理论上可以让迭代继续但采样理论对它的收敛性分析远不如光滑情形成熟而且当势能非凸时次梯度的选择还会影响算法的路径依赖。更关键的是次梯度在不可导点附近通常没有稳定的方向随机噪声与次梯度叠加后轨迹会变得非常“毛糙”。这就是 Moreau--Yosida 正则化登场的原因。它的思路不是绕开不可导点也不是给折线硬凑一个近似而是对势能的非光滑部分做一次“下半连续包络”变换得到一个新的、处处可微、且梯度有显式表达式的函数。更妙的是这个变换不会改变原问题的最优值结构也不会破坏凸性。我先把结论放在这里MYULA 本质上把“我需要知道非光滑势能某点的梯度”替换成了“我需要能求解该点的近端算子proximal operator”。在稀疏学习、约束估计、图像处理这类问题里prox 往往有闭式解这比求一个不存在的梯度靠谱得多。2. 三个核心概念ULA、Moreau--Yosida 与近端算子2.1 目标与记号设目标分布是[ \pi(x) \propto \exp(-U(x)), \quad Ugh, ]其中 (g) 是光滑部分具有 (L_g)-Lipschitz 梯度(h) 是非光滑但有良好结构的凸函数通常要求它的近端映射可以高效计算。我们的任务是产生近似服从 (\pi) 的样本。在经典的 ULA 分析中通常假设整个 (U) 是光滑的甚至要求强凸。引入分解 (Ugh) 的意义在于真正麻烦的只有 (h)而 (h) 往往是有特殊结构的。2.2 ULA从 Langevin 扩散到离散迭代连续时间的 Langevin 扩散满足随机微分方程[ dX_t -\nabla U(X_t),dt \sqrt{2},dB_t, ]其中 (B_t) 是标准布朗运动。在很宽的条件下这个扩散的平稳分布就是 (\pi)。数值上要用 Euler--Maruyama 离散化[ x_{k1} x_k - \gamma \nabla U(x_k) \sqrt{2\gamma},\xi