简介面向信号处理与图像优化算法学习的研究者或开发者这份压缩包围绕交替方向乘子法、交替最小化算法、组稀疏去噪以及MM优化策略整理了一套便于直接运行的MATLAB代码工程。工程共计35个文件整体大小只有43KB其中以31个MATLAB脚本文件为主体覆盖矩阵补全、鲁棒主成分分析、全变分去噪、非凸正则化、彩色图像去噪、一维信号去噪、二阶总变分、分裂布雷格曼解卷积以及组稀疏去噪等典型图像处理任务同时附带说明文档和测试数据便于理解数据格式与实验流程。每个算法都配有独立演示脚本能够帮助读者快速观察不同迭代方法的收敛特性与去噪效果既适合作为优化理论课程的配套实验也可以作为论文复现的起步代码。该工程已有188人学习下载适合具备一定信号处理基础、希望通过代码加深对交替迭代和组稀疏模型理解的研究者与工程人员。1. 组稀疏信号去噪与 ADM、AMA、MM从“能跑”到“可调”的优化组合组稀疏信号去噪配上交替方向法ADM、交替最小化法AMA和 Majorization-MinimizationMM是图像信号处理里很常见的一套优化组合先假设信号在分组结构下稀疏再通过分裂迭代把目标函数拆成若干能单独求解的子问题。第一次接触的人容易卡在两个地方一是分不清 ADM 和 AMA 到底差在哪二是拿到图像数据后 λ、ρ 不知道怎么给跑出来要么没效果、要么整片抹平。这篇文章就从这两个问题入手把问题建模、迭代推导、完整可运行的 Python 实现、参数设定和踩坑点都摊开。内容面向做图像去噪、压缩感知重构、逆问题求解的工程师和同学。这套方法不需要 GPU不依赖深度学习框架一个 numpy 文件就能跑通每步迭代都有明确的物理含义样本量不大时比端到端网络更容易定位问题、控制结果值得作为基线方案沉淀下来。2. 从 ADM 到 AMA先分清两种交替迭代框架再动手组稀疏去噪的优化目标长这样min_x 0.5 * ||y - x||² λ * Σ_g ||x_g||₂y 是含噪观测x 是待恢复信号x_g 是第 g 组的分量。前一项要求恢复结果贴近观测后一项迫使每个分组要么整组保留、要么整组归零这就是组稀疏的含义。直接对这个目标做梯度下降很别扭因为组范数在零点不可导而且保真项和正则项对 x 的诉求不一致。ADM 和 AMA 都是为解决这类问题设计的。2.1 ADM 变量分裂把组稀疏去噪拆成三个子问题ADM 的核心操作是变量分裂。引入辅助变量 z令 z x把原问题改写成带等式约束的形式min_x 0.5 * ||y - z||² λ * Σ_g ||x_g||₂约束 x z这样保真项只依赖 z组稀疏正则只依赖 x两者通过约束耦合。对等式约束做增广拉格朗日处理得到三个更新步骤第一步更新 x只处理组范数项和二次惩罚项结果是组软阈值算子有闭式解。第二步更新 z只处理保真项和二次惩罚项也是解析解。第三步更新乘子 u把 x 与 z 的偏差反馈回去。# ADM 三步更新伪代码突出变量分裂结构 for k in range(max_iter): # 1) x-update组软阈值 x group_soft_threshold(z - u, groups, lam / rho) # 2) z-update保真项 二次惩罚的解析解 z (y rho * (x u)) / (1.0 rho) # 3) u-update乘子累积 u u (x - z)变量分裂的意义在于原本“既要贴近观测、又要分组稀疏”的矛盾被拆成了两个独立子问题各自都有高效解法不再需要手动权衡梯度步长。rho 是惩罚参数控制 x 和 z 之间的一致性强弱u 是乘子等价于在迭代中不断纠偏所以最终收敛时 x 与 z 的残差会趋近于零。2.2 AMA 与 ADM 的差别一个精确更新、一个近端步怎么选AMA 的全称是 Alternating Minimization Algorithm由 Tseng 在 1997 年前后提出和 ADM 同属分裂型迭代框架但在子问题的处理策略上有关键差别。ADM 要求两个子问题都尽量求到近似精确解AMA 则允许其中一个变量只做一步近端梯度更新不要求精确求解。这个差别在组稀疏去噪场景里很实际当正则项只是组范数时x 子问题有闭式解ADM 顺手但如果正则项是“组范数 全变分”这类组合x 子问题没有闭式解精确求解代价很高这时候用 AMA 对 x 做一步近端梯度步整体收敛反而是更划算的选择。对比项ADMAMA子问题求解两个子问题都尽量精确求解一个精确求解另一个只做近端梯度步收敛前提增广项 rho 适当残差单调下降梯度项满足 Lipschitz 条件或使用回溯步长适用场景两个子问题都有快速求解器其中一个子问题昂贵、无闭式解实现复杂度低稍高需要处理梯度步长选择原则可以简单记组范数本身有闭式 prox优先 ADM正则项叠加导致 prox 不再闭式换 AMA。后文第 4 章的落地代码以 ADM 为主第 3 章处理不光滑保真项时再引入 MM。2.3 ADM 求解组稀疏问题组软阈值算子的闭合解与实现组软阈值是组稀疏去噪里最核心的算子x 子问题的解。设 v 是输入向量对每一组下标 g先算整组的二范数然后按比例收缩def group_soft_threshold(v, groups, thresh): out v.copy() for idx in groups: norm_g np.linalg.norm(v[idx]) if norm_g thresh: # 整组按比例收缩方向保持 out[idx] v[idx] * (1.0 - thresh / norm_g) else: # 整组置零 out[idx] 0.0 return outthresh 对应 ADM 里的 lam / rho。当某组二范数小于阈值时整组归零这是“组”稀疏区别于逐元素 L1 稀疏的地方组内元素要么一起活要么一起死不会出现半个组保留、半个组被削平的情况。threshold 越大保留下来的组越少输出越平滑。调用时要注意 v 是 z - u也就是 x 子问题的输入不要传 y 本身否则会把乘子信息丢掉。3. Majorization-Minimization当保真项不光滑算法怎么继续下降ADM 在保真项是二次函数比如高斯噪声下的最小二乘时非常顺手。但图像信号处理里经常遇到非二次保真项椒盐噪声对应的 L1 保真项、低光图像泊松噪声对应的 KL 散度、带遮挡时的截断残差。这些函数不可导或没有解析解直接套 ADM 的 z 更新步骤会卡住。这时需要 Majorization-MinimizationMM也就是标题里的核心框架之一。3.1 MM 的两个条件和单调下降保证MM 的思路是原目标函数 f(x) 不好最小化就构造一个更容易处理的上界函数每次极小化这个上界函数来逼近原问题。上界函数需要满足两个条件一是处处不小于原函数二是在当前迭代点处取等号。满足这两个条件后每次迭代得到的新点不会让原目标函数值上升这是 MM 单调收敛的保证。用公式表达在第 k 轮找到 Q(x | x_k)满足 Q(x | x_k) ≥ f(x) 且 Q(x_k | x_k) f(x_k)然后求 x_{k1} argmin Q(x | x_k)。下一轮重新构造上界继续迭代。这个框架也叫序列凸近似在压缩感知、图像恢复的文献里经常和 ADM 配合出现MM 负责生成代理问题ADM 负责解代理问题。3.2 用 MM 把 L1 保真项替换成二次上界公式与参数含义以 L1 保真项为例。原问题变成 min_z ||y - z||₁ λ * Σ_g ||z_g||₂。L1 项在零点不可导z 子问题没有闭式解。MM 的做法是对 |y_i - z_i| 构造一个二次上界。常见构造是|t| ≤ (t² / (2M)) M / 2其中 M 是正的曲率参数。取 t y_i - z_i就把 L1 保真项替换成了加权二次项常数项与原函数当前值对齐。M 越大上界越“宽松”每次逼近的步幅越小M 越小上界越紧但小到一定程度不等式会失效迭代发散。工程上通常在当前迭代点取对等条件得到逐点权重 w_i 1 / (2|y_i - z_i^k| eps)分母加小常数防止除零。这样既满足上界条件又能让参数 M 随迭代自动调整。3.3 外层 MM 内层组软阈值最小可运行骨架下面是一段最小骨架外层用 MM 更新 L1 保真项的权重内层用近端梯度迭代解代理子问题组稀疏部分仍然走组软阈值算子def mm_l1_group_denoise(y, groups, lam0.05, tau0.1, outer30, inner30, eps1e-8): y y.ravel().astype(np.float64) z y.copy() for k in range(outer): # MM基于当前点构造 L1 保真项的二次上界权重 w 1.0 / (2.0 * np.abs(y - z) eps) # 内层近端梯度解 加权重二次保真 组稀疏正则 for _ in range(inner): grad 2.0 * w * (z - y) z group_soft_threshold(z - tau * grad, groups, lam * tau) return z.reshape(y.shape)梯度表达式来自代理二次项 2 * w_i * (z_i - y_i)对 z_i 求导得到。tau 是近端梯度步长受代理二次项的曲率限制必须小于 1 / (2 * max(w))否则内层迭代发散。如果运行时报 NaN优先把 tau 缩小一个数量级或者在外层循环里对 tau 做回溯。这段代码把 MM 和组软阈值算子串在了一起实际使用中直接用其中的 z 更新逻辑替换 ADM 里对应步骤即可不需要额外框架。4. 用 Python 在本地跑通组稀疏图像去噪完整代码与参数这一章给完整可运行的图像去噪流程。从读图、构造分组到 ADM 主循环、重叠聚合、收敛判断和 PSNR 评测按顺序落地。4.1 预处理与组构造图像怎么切成重叠块图像去噪里“组”常见的定义是滑动窗口切出的图像块每个块展平后作为一组。块内部像素被要求联合稀疏这样纹理和边缘不会被单独像素的噪声抹掉。重叠滑动能避免块与块之间的可见边界但也意味着同一个像素属于多个组去噪后需要聚合。import numpy as np import imageio.v2 as imageio img imageio.imread(cameraman.png).astype(np.float64) / 255.0 rng np.random.default_rng(0) noisy img rng.normal(0, 0.05, img.shape) def build_overlap_groups(shape, block8, step4): groups [] h, w shape for r0 in range(0, h - block 1, step): for c0 in range(0, w - block 1, step): idx [] for r in range(r0, r0 block): for c in range(c0, c0 block): idx.append(r * w c) groups.append(np.array(idx, dtypenp.int64)) return groups groups build_overlap_groups(img.shape, block8, step4)图像先转 float64 并除以 255归一化到 [0,1]后续加噪声、设正则参数和算 PSNR 都在同一尺度下进行避免 uint8 取值范围带来的参数漂移。step 控制重叠程度step 越小组越多、去噪效果更连续但迭代耗时线性上升。block8、step4 是性能和效果比较平衡的起点。4.2 ADM 主循环x、z、u 三步更新的完整实现def group_soft_aggregate(v, groups, thresh): out np.zeros_like(v, dtypenp.float64) cnt np.zeros_like(v, dtypenp.float64) for idx in groups: norm_g np.linalg.norm(v[idx]) if norm_g thresh: v_g v[idx] * (1.0 - thresh / norm_g) else: v_g np.zeros_like(v[idx]) out[idx] v_g cnt[idx] 1.0 return out / np.maximum(cnt, 1.0) def adm_group_denoise_2d(noisy, groups, lam0.05, rho0.2, max_iter80, tol1e-5): y noisy.flatten() z y.copy() x np.zeros_like(y) u np.zeros_like(y) history [] for k in range(max_iter): x_prev x.copy() # x-update组软阈值 重叠组平均 x group_soft_aggregate(z - u, groups, lam / rho) # z-update保真项二次 rho 惩罚项解析解 z (y rho * (x u)) / (1.0 rho) # u-update乘子累积 u u (x - z) # 收敛判断用相对变化量 rel np.linalg.norm(x - x_prev) / (np.linalg.norm(x_prev) 1e-12) history.append(rel) if rel tol: break return z.reshape(noisy.shape), historygroup_soft_aggregate 把组软阈值和重叠聚合合并到一步对每个组做阈值收缩后把结果累加到输出数组同时计数最后除以每个像素被覆盖的次数。这是重叠组稀疏去噪的工程近似写法严格意义上重叠组的 prox 没有闭式解但这种“收缩 平均”的做法在图像去噪里非常常用效果稳定。z 更新来自两个二次项之和的极小化分母 1.0 rho 对应保真项系数 1 和惩罚项系数 rho 的合并。4.3 组软阈值与重叠聚合让去噪图不留块状痕迹组软阈值算子本身不区分“这个块是纹理区还是平坦区”只看整组二范数。如果分组不重叠每个块独立处理后再拼回去块与块之间可能出现亮度不连续视觉上就是棋盘格。重叠聚合通过让每个像素被多个块共同决定来消除这个问题。块越多、重叠越大聚合结果越平滑但过大的重叠会让有效正则作用被平均稀释所以 block 和 step 通常保持 2:1 到 4:1 的关系。另外组软阈值对整组做比例收缩组内元素的相对关系保留。这意味着强边缘所在的小块不会被压成纯色块而是整组缩到一个较小的尺度边缘位置仍然保留。这个特性是组稀疏优于逐元素 L1 稀疏的关键也是去噪结果看起来更“像图”的原因。4.4 收敛判据与 PSNR 评测停止条件和客观指标一起看只用迭代次数退出不可靠。常见做法是同时看相对变化量和原始残差相对变化量反映 x 是否还在动原始残差 ||x - z||₂ 反映约束是否已经满足。工程上推荐输出每次迭代的 PSNR 曲线和残差曲线判断真正收敛的位置。def psnr(clean, rec, data_range1.0): mse np.mean((clean - rec) ** 2) if mse 1e-12: return 99.0 return 10.0 * np.log10(data_range ** 2 / mse) rec, history adm_group_denoise_2d(noisy, groups, lam0.05, rho0.2) print(PSNR:, psnr(img, rec)) print(last rel change:, history[-1])PSNR 计算先确保两个输入都是 [0,1] 范围内的浮点图像。data_range 填 1.0 而不是 255避免无量纲错误。如果 PSNR 低于带噪输入说明参数选择有问题优先检查 lam 是否过大或 rho 是否过小而不是去调迭代次数。5. 组稀疏去噪实装避坑五个绕不开的翻车点5.1 输出图出现棋盘格伪影现象去噪后的图能看出规则的方块边界特别是平坦区域像拼图没对齐。原因分组时 step 等于 block块与块不重叠每个块被独立阈值处理后再拼回去相邻块的噪声残留统计不一致形成可见边界。解决把 step 降到 block 的一半也就是至少 50% 重叠。block8 时 step 取 4如果还明显就取 2。同时确认聚合阶段除以了覆盖次数输出数组要用浮点累加而不是直接在原图上覆盖。5.2 lam 和 rho 调不出合适范围现象lam 调小输出和输入几乎一样噪声都在lam 调大图像整片发虚纹理全丢rho 调大收敛特别慢。原因lam 和 rho 都与图像数值尺度强相关直接套用网上代码里的数值经常对不上。图像尺度不一致同样的 lam 实际效果可能差一个量级。解决先把图像统一归一化到 [0,1]再按噪声水平给定参数。实际踩下来可以参考下面的起始范围噪声 Sigmalam 参考范围rho 参考范围0.010.005 - 0.020.1 - 0.30.050.02 - 0.080.2 - 0.50.100.08 - 0.200.3 - 0.8注意实际收缩阈值是 lam / rho所以 lam 和 rho 一起变大时收缩力度不一定增强可能只是收敛速度变化。调参时固定 rho等比扫描 lam比同时乱调两个参数高效得多。5.3 MM 迭代目标函数不降反升现象加入 MM 处理 L1 保真后目标函数曲线隔几步往上跳甚至出现 NaN。原因上界函数条件被破坏最常见的是权重 w 没有在当前迭代点重新计算用的是上一轮的值或者近端梯度步长 tau 超过 1 / (2 * max(w))导致内层迭代发散。解决外层每轮先重新算 w再做内层迭代。tau 设置保守一点比如取 0.5 / max(w)不要拍脑袋定死。内层迭代次数控制在 20 到 30 步不需要完全收敛外层 MM 框架会继续修正。5.4 停止条件太宽松导致提前退出现象程序跑了三四十轮就停了输出图看起来还有明显噪声但日志里相对变化量已经很小。原因相对变化量只反映 x 是否还在动。x 可能已经接近某个平坦区域但 z 还在慢慢调整乘子 u 也在漂移原始残差 ||x-z||₂ 仍然偏大。解决停止条件同时检查原始残差和对偶残差。原始残差 r_k ||x_k - z_k||₂对偶残差 s_k rho * ||z_k - z_{k-1}||₂。两个都小于阈值才算收敛同时设一个最小迭代次数比如 50 轮兜底避免刚进入稳定区就退出。5.5 整数图像直接开算导致参数尺度漂移现象同一组参数一张图效果好换一张图效果差很多有的图 PSNR 甚至为负。原因有的读图库返回 uint80-255有的返回 float0-1。代码里直接相加、减法和阈值操作在 uint8 上会溢出lam 和 rho 的绝对数值也无法跨尺度复用。解决所有图像数据在进入算法前统一转 float64如果最大值大于 1 就除以 255后续全程在 [0,1] 浮点域操作。加噪声、算 MSE、算 PSNR 都在这个域内完成。输出时先 clip 到 [0,1] 再乘 255 转 uint8。这个习惯能消除大部分跨环境复现不一致的问题。6. 把这套流程用在自己的数据上三步验证与热启动技巧6.1 仿真信号验证先确认算法没有逻辑错误不要一上来就跑真实图像。先构造一组已知真值的仿真信号确认算法恢复结果正确、收敛曲线正常。x_true np.zeros(512) x_true[10:40] np.random.default_rng(1).normal(size30) y_sim x_true 0.05 * np.random.default_rng(2).normal(size512) groups_sim [np.arange(i, i 8) for i in range(0, 512 - 8, 4)] x_rec, _ adm_group_denoise_2d(y_sim, groups_sim, lam0.05, rho0.2) print(np.linalg.norm(x_rec - x_true) / np.linalg.norm(x_true))这一步能快速暴露分组索引错误、阈值方向错误、聚合计数错误等低级问题。真实图像尺寸大、噪声复杂出错时不容易定位。6.2 真实图像验证PSNR、SSIM 和误差图一起看真实图像上除了 PSNR还要看 SSIM 和误差图。PSNR 对整体像素误差敏感SSIM 对局部结构保留更敏感。误差图能直接看出哪些区域被过度平滑、哪些区域噪声残留。常见做法是拿 cameraman 或 Set12 里的图分别在 sigma0.05 和 0.1 下跑和经典高斯滤波或 BM3D 做对比记录。一次结果不能说明问题至少跑三张不同内容风格的图结果一致才算可靠。6.3 视频去噪的热启动用上一帧结果省一轮迭代最后一个实用技巧是热启动。视频逐帧去噪时相邻帧内容变化不大完全没必要每帧从头开始迭代。把上一帧的 x、z、u 作为下一帧的初值传入通常能省 30% 到 50% 的迭代轮数。实现时只需要给主循环函数加几个可选初值参数迭代开始前判断是否传入即可。我调这类算法长期养成一个习惯跑完一组实验除了 PSNR一定把目标函数曲线和残差曲线打出来看。曾经因为只盯 PSNR错过了 rho 过大导致的振荡问题最后是靠这两条曲线才发现算法根本没进收敛区。希望帮到你少在这上面浪费一个下午。本文还有配套的精品资源点击获取