简介：一套基于MATLAB实现的ADMM（交替方向乘子法）算法代码包，主要面向机器学习、图像处理、信号处理等领域的算法研究者和工程师，用于求解含等式或不等式约束的凸优化问题。压缩包共383个文件，内容约17.44MB，涵盖171个m源文件、70个p加密文件、22个dat数据、22个mat数据，以及面向Windows/Mac/Linux的mexw64、mexmaci64、mexa64预编译文件，并附带pdf、md、txt等格式的说明资料，目录按源码、数据、文档分类存放，便于对照学习。目前已有四千余人浏览学习，是实践ADMM的成熟参考资料。代码主体ADMM.m完整实现了问题定义、参数设置、子问题分解、交替更新与乘子更新等核心步骤，配套辅助函数用于特定子问题求解；大量dat/mat数据文件可供稀疏回归、低秩矩阵恢复等经典场景直接测试，mex文件则提供预编译加速便于性能比较。通过阅读代码与pdf文档，可以深入理解ADMM的收敛条件、惩罚参数rho的影响及加速策略，适合希望快速上手算法复现、二次开发乃至分布式并行化的中高级用户。1. ADMM算法：为什么分布式优化和约束问题都绕不开它ADMM算法全称是交替方向乘子法（Alternating Direction Method of Multipliers），这几年在压缩感知、低秩矩阵恢复、分布式机器学习、图像去模糊这些场景里几乎成了标配解法。你只要手上有一个带约束的优化问题，目标函数还能拆成两块，那ADMM基本就是最省心的那一个工具。它的聪明之处在于：不跟复杂的原问题硬碰硬，而是把大问题拆成小步迭代，每一步只需要处理一个简单子问题，MATLAB里十几行循环就能搭出主干。这份资源给的就是一套可直接运行的MATLAB实现，能帮你绕过从论文公式到可运行代码之间的那段空档期。适合正在做调参、算法对比或毕设复现的从业者。2. ADMM的数学原理与MATLAB落地：从对偶上升、增广拉格朗日到三段迭代2.1 ADMM的核心思想：对偶上升、增广拉格朗日与交替极小化要理解ADMM，得先看它从哪里来。最原始的拉格朗日乘子法，把约束罚进目标函数，然后同时对原变量和对偶变量做梯度上升，这就是对偶上升法的思路。问题是对偶上升要求目标函数有比较强的光滑性，而且收敛条件偏紧。后来有人加了二次罚项，得到增广拉格朗日函数，这样对偶上升的收敛范围宽了很多。但仍有一个麻烦：增广拉格朗日里的二次罚项让子问题不再可分，遇到大规模问题就没法并行。ADMM做的事情其实很直接：假设你要解决的问题是min f(x) + g(z)，约束是 Ax + Bz = c它不直接对x和z做联合优化，而是把x和z分开轮流更新，每一轮只解一个子问题。这种交替极小化的好处是，两个子问题往往都有闭式解，或者能用很便宜的数值方法求解。配合对偶变量的上升更新，最终收敛到原问题的最优解。MATLAB里实现这套逻辑，核心就是下面三段式子：x_new = argmin_x ( f(x) + (rho/2) * norm(A*x + B*z - c + u, 2)^2 ); z_new = argmin_z ( g(z) + (rho/2) * norm(A*x_new + B*z - c + u, 2)^2 ); u_new = u + A*x_new + B*z_new - c;第一个式子更新原变量x，第二个更新分裂变量z，第三个是对偶变量u的上升步。rho是惩罚参数，控制收敛速度和数值稳定性。代码看起来简单，但每一步的求解才是真正的工作量。2.2 MATLAB实现的最小骨架：单约束问题三段循环在一个标准的单约束问题里，ADMM的MATLAB实现可以缩到最小骨架。下面这段代码展示的是最简单形态：function [x, z, u, history] = admm_simple(A, B, c, rho, max_iter, tol) % A, B: 约束矩阵 % c: 约束右端项 % rho: 惩罚参数（ADMM的核心超参） % max_iter: 最大迭代次数 % tol: 相对终止容差 m = size(A, 1); n = size(A, 2); p = size(B, 2); x = zeros(n, 1); z = zeros(p, 1); u = zeros(m, 1); history = zeros(max_iter, 2); for k = 1:max_iter % 更新 x：这里假设 f(x) = ||x||_2^2 即平方正则 % 闭式解来自求解线性方程组 x = (A'*A + rho*eye(n)) \ (A'*(c - B*z - u) * rho); % 更新 z：这里假设 g(z) = indicator_{z>=0} % 投影到非负象限 z = max(0, c - A*x - u); % 对偶更新 u = u + (A*x + B*z - c); % 记录主残差和对偶残差 history(k, 1) = norm(A*x + B*z - c, 2); history(k, 2) = norm(rho * A' * B * (z - z_old), 2); end end注意上面这段代码里，x的闭式解写成了一个线性方程组的求解。实际使用时要先算(A'A + rhoI)的分解，不要每次迭代都重新求逆。另外z的更新用了非负投影，这是g(z)为指示函数时的特例，换成别的正则项就要改成对应的近端算子。对偶残差的定义在不同论文里写法略有差异，我在这里用了带rho和A'B的形式。2.3 为什么用ADMM：与你见过的梯度下降和近端梯度的区别很多人在接触ADMM之前，已经用过梯度下降或者近端梯度法。梯度下降只适合光滑目标函数，近端梯度法能处理一部分非光滑项，但它们都要求目标函数整体结构比较规矩。ADMM的出发点不同：不试图一次处理整个目标函数，而是把目标拆成f和g两块，各自有各自的解法，只要能算近端算子，就能塞进ADMM框架。举个例子，LASSO问题里f是平滑的二次损失，g是L1范数。用近端梯度做，L1近端算子本身很简单，但梯度步的步长选择会影响收敛；而ADMM里g对应的z-update就是一个软阈值操作，没有步长参数，调的东西只有rho。这一点在实际使用中的意义很大——少一个超参，调参的麻烦就少一截。还有一个实际场景上的区别：分布式或联邦式的优化环境里，数据可能分布在不同的节点上。梯度下降需要汇总梯度，ADMM可以通过交替更新的方式让每个节点独立计算自己的子问题，最后只同步对偶变量。这一点在集中式的MATLAB实现里体现不明显，但理解了这个背景，你就知道ADMM的变量分裂结构不是偶然。3. 以LASSO为例的完整MATLAB实现：目标函数拆解与主循环写法3.1 LASSO问题的ADMM重写：目标函数拆解与变量分裂LASSO问题的标准形式是：min 0.5 * ||A*x - b||_2^2 + lambda * ||x||_1要让ADMM上场，得先做一步变量分裂，把L1范数“隔离”到z上。常见做法是引入z，并加约束x - z = 0，于是问题变成：min f(x) + g(z)，约束 x - z = 0其中f(x) = 0.5 * ||A*x - b||_2^2，g(z) = lambda * ||z||_1。这样拆完之后，x-update变成了一个带二次正则的最小二乘问题，它的解就是线性方程组的解，而z-update是对L1范数的近端算子，即软阈值操作。拆解的关键是变量分裂这一步的选择。不同的问题，分裂方式可以完全不同。比如有的场景是拆出两块目标函数，有的场景是拆出约束矩阵A和B。LASSO这种拆法是最直观的：把难啃的L1正则项放到一边，把光滑的二次损失放到另一边。3.2 主迭代循环：x-update、z-update与对偶更新下面给出一个完整的MATLAB函数，输入设计矩阵A、观测向量b、正则系数lambda和ADMM参数，输出解和收敛历史：function [x, history] = admm_lasso(A, b, lambda, rho, alpha, max_iter, tol) % admm_lasso - 用ADMM求解LASSO问题 % 输入： % A - m*n 的测量矩阵 % b - m*1 的观测向量 % lambda - L1正则系数（非负） % rho - ADMM惩罚参数（正数） % alpha - 松弛因子，通常取1.0或1.5，>1会加速但可能不稳 % max_iter- 最大迭代次数 % tol - 相对残差容差 % 输出： % x - n*1 的最优解 % history - 迭代历史（主残差、对偶残差、目标值） [m, n] = size(A); x = zeros(n, 1); z = zeros(n, 1); u = zeros(n, 1); % 预分解：x-update需要反复求解(A'A + rho*I)这个系统 % 先算Cholesky因子，后面每次迭代只做回代，省掉重复分解 [L, p] = chol(A'*A + rho*eye(n), 'lower'); if p ~= 0 error('矩阵(A''A + rho*I)必须正定。'); end Atb = A'*b; history = zeros(max_iter, 3); for k = 1:max_iter % x-update: 求解 (A'A + rho*I) * x = A'b + rho*(z - u) % 这一项来自对f(x)求导并令导数为零 x = L \ (L' \ (Atb + rho*(z - u))); % z-update: 软阈值 % 对每个分量执行 shrink(z + u, lambda/rho) zold = z; z = max(0, 1 - lambda/(rho*abs(x + u))) .* (x + u); % u-update: 对偶变量沿残差方向上升 u = u + alpha * (x - z); % 收敛历史记录 history(k, 1) = norm(x - z, 2); history(k, 2) = norm(rho * (z - zold), 2); history(k, 3) = 0.5*norm(A*x - b, 2)^2 + lambda*norm(x, 1); % 终止判定：主残差和原始残差都小于阈值才停 if history(k,1) < tol && history(k,2) < tol break; end end history = history(1:k, :); end逻辑说明：这段代码的循环体完全对应ADMM的三个基本步骤。x-update使用预分解的Cholesky因子做快速回代，z-update实现的是软阈值近端算子，u-update则是对偶上升。整体上这个写法与标准ADMM论文里的形式一致，只是加入了alpha松弛因子，这让代码可以在alpha=1.5时获得比alpha=1更快的收敛，但代价是可能出现轻微震荡。参数说明：lambda控制稀疏度，lambda越大解越稀疏，lambda = norm(A'*b, Inf)时解全为0，这是一个有用的上界；rho是ADMM的收敛速度调节器，太大导致收敛快但数值精度下降，太小导致迭代次数激增；alpha是松弛因子，通常取1.0，加速到1.5之前要先确认残差曲线没有明显的震荡。3.3 收敛参数与输入输出的定义边界输入参数这一块有几个边界条件值得说明。首先是lambda的范围，lambda必须是大于等于0的实数，如果输入负值，软阈值操作里的(1 - lambda/(rho*abs(x+u)))可能为负，导致z的符号被反向修正，整个算法的收敛无从谈起。其次是rho，rho不能为0，否则增广拉格朗日的二次罚项消失，ADMM退化回原始拉格朗日法，收敛条件变得极其苛刻。另外要留意的是max_iter和tol的组合。新手容易把tol设得太小，比如1e-10，但数值精度和机器epsilon有限，这种设定往往白跑几千轮还停不下来。合理的做法是先跑到固定迭代次数，画残差曲线，再决定tol的等级。上面的代码里，如果矩阵A的规模较大，最好把Cholesky分解换成基于QR或SVD的预处理版本，否则每次调用函数时的一次分解开销依然不小。常见做法是在函数外部完成分解，把因子传入迭代主体。4. 参数设置与收敛性判定：惩罚参数rho、松弛因子alpha与终止条件4.1 惩罚参数rho：收敛速度与数值稳定性的天平rho是ADMM里最重要的超参数，它的作用在增广拉格朗日函数里看得最清楚：二次罚项的系数。rho越大，每一步迭代对约束违反的惩罚越重，迭代会更快地让Ax + Bz - c趋近于零，但代价是原变量x的更新步长变小，可能导致收敛精度下降甚至在最优解附近来回震荡。rho越小，约束违反在早期被容忍，迭代过程平滑但很慢。实际调rho的方法，我一般习惯先在log空间里试几个值，比如0.01、0.1、1、10、100，各跑200轮，对比主残差收敛曲线，然后选曲线下降最快且没有明显震荡的那个量级。注意一个问题：rho的选择和问题本身的尺度强相关。如果A是0到1之间的数据，rho=1通常合理；如果A的数值在10^3量级，rho可能需要调到10^3以上才有反应。若rho的取值与A的尺度严重不匹配，会出现一种极端的现象：主残差前期快速下降，后期长期不动，而z迭代还在缓慢漂移。这时候强行加大rho并不会改善，反而是先把数据标准化更有效。4.2 松弛因子alpha：加速还是添乱alpha是ADMM的松弛因子，标准的ADMM里alpha=1，就是普通的三段迭代。当alpha>1时，对偶更新的方向更激进，可以在某些问题里把收敛速度提升不少，特别是对于二次损失加L1正则这类结构的问题。但alpha不是越大越好。实践中alpha超过1.8以后，残差曲线往往出现锯齿状波动，甚至发散。原因在于对偶变量u更新里的(1 - alpha)项会引入额外的振荡模式。如果是第一次调这个参数，建议只试1.0和1.5两个值，对比总的迭代次数，不要贪心。很多现成实现里，alpha是写死为1.0的。如果你手里的资源代码里有alpha这个参数，说明作者考虑到加速场景，这是一个加分项。但用的时候要清楚：alpha只在目标函数光滑性较好的那段收敛前期有明显收益，到了末段收敛区域，它的作用会变小。4.3 终止条件怎么设：双重残差判定与自适应rho终止条件是ADMM实现里最容易糊弄的部分。只看主残差norm(x - z)小于阈值就停，是对偶残差可能还很大，算法实际上并没有收敛到鞍点。正确的做法是同时检查主残差和对偶残差，两个都小于各自阈值才算完成。阈值的设置也有讲究。相对残差比绝对残差更适合跨问题使用：把主残差除以max(norm(Ax), norm(Bz), norm(c))，得到一个尺度无关的指标。下面给出一个可用的判定片段：eps_pri = tol * max([norm(A*x,2), norm(B*z,2), norm(c,2)]) + tol_abs; eps_dual = tol * norm(A' * u, 2) + tol_abs; stop = (r_pri < eps_pri) && (r_dual < eps_dual);这里tol_abs是绝对容差，通常取1e-6量级。tol是相对容差，工程上取1e-4到1e-6之间比较常见。如果你发现算法迟迟不终止，但残差曲线已经走平，先检查是不是tol设得太小，再看残差是压在哪个水平走平的——这往往能直接告诉你rho应该增大还是减小。还有一种进阶做法是自适应rho，即每隔若干轮根据主残差与对偶残差的比值调整rho的值。常见策略是当主残差大于对偶残差若干倍时加大rho，反之减小rho。自适应rho能省去手动搜索的麻烦，代价是额外几个参数：调整间隔、放大倍数、缩小倍数。下面给一个常规实现片段：if mod(k, 10) == 0 if r_pri > 10 * r_dual rho = rho * 2; % 注意：改变了rho之后，之前预分解的Cholesky因子失效，需要重新分解 elseif r_dual > 10 * r_pri rho = rho / 2; end end这里有个容易翻车的细节：rho一变，x-update里的(A'A + rho*I)整个矩阵都变了，必须重新做Cholesky分解，否则用的是旧因子，算出来的x方向不对。每10轮重分解一次，代价通常可以承受。5. ADMM避坑指南：参数失误、残差误判与实现细节5.1 迭代发散或出现NaN：rho与步长设置不当现象：运行代码没几轮，变量就变成了NaN或者Inf，程序直接报错。原因：最常见的是rho设得过大，导致x-update的线性系统数值上接近奇异；或者是数据矩阵A本身的条件数太大，比如A里的元素跨越了10^6的量级，这时候(A'A + rho*I)的条件数依然很糟糕，解出来的x方向有巨大误差。另一个常见原因是软阈值操作里分母上出现了abs(x+u) = 0的情况，在MATLAB里虽然不会直接报错，但结果会变成NaN。解决：先把rho往小调，比如从10000调到1或者0.1。同时检查A的数值范围，如果A的元素量级差异太大，先对列做标准化，再进ADMM。对于软阈值的数值保护，可以在abs(x+u)下面加上一个小常数eps=1e-12：max(0, 1 - lambda/(rho*(abs(x+u)+eps)))。这一步是血泪经验，MATLAB里的0/0结果是NaN，NaN会顺着迭代一路传染。5.2 只盯着主残差：对偶残差滞后导致的误判现象：main残差已经到了1e-6，看起来一切正常，但最终的解和真实最优解相比有明显偏差。原因：主残差衡量的是x和z的一致性，对偶残差衡量的是对偶变量的收敛程度。如果只以主残差作为终止条件，可能在对偶变量还在漂移的时候就停了。对偶残差反映了当前迭代点到最优鞍点的距离，两个残差都小才是真正的收敛。解决：终止条件必须写成同时满足主残差和对偶残差的阈值。另外，在观察收敛曲线时，把两条残差曲线画在同一张图上，如果发现主残差已经下降但对偶残差还高出一个量级，那就说明当前的rho太小，对偶更新太慢，把rho调大一点就能拉平两条曲线。5.3 大矩阵场景每次重复求逆：缓存分解结果现象：矩阵A的规模在几千乘几千，代码每跑一轮迭代都要花很长时间，但看CPU占用率却不高。原因：典型的实现问题，每轮迭代都执行(A'A + rhoI) \ something。MATLAB里反斜杠运算符会先判断矩阵结构并做分解，反复调用意味着反复做LU或Cholesky分解，而这一步的复杂度是O(n^3)，完全淹没了回代的O(n^2)开销。解决：进入循环之前先做一次性分解。如果矩阵对称正定，用Cholesky；如果不保证正定，用LU。上面的LASSO代码已经展示了这个技巧。如果rho在迭代中会自适应变化，分解必须放在rho更新之后重新做。这里有一个取舍：自适应rho的更新间隔不能太短，否则频繁重构因子反而拖慢整体速度。5.4 可变rho的震荡：更新间隔与边界限制现象：使用自适应rho之后，残差曲线不是平滑下降，而是明显震荡，甚至发散。原因：rho的更新过于频繁或者调整倍数太大。比如每轮迭代都检查残差比并翻倍rho，会让增广拉格朗日的权重大幅摆动，收敛过程完全被rho的调整节奏带着走。解决：把rho的调整间隔加大，至少10轮一次。另外，调整倍数不要超过2倍。同时给rho加一个上下界，比如rho_min = 1e-4，rho_max = 1e4，防止极端情况下rho连续翻倍导致数值溢出。还有一个细节：改变rho之后，对偶变量u也应该按比例缩放，不然旧的u在新尺度下失去了参考意义。标准做法是u = u * (rho_old / rho_new)。5.5 不可分目标直接套用：用变量分裂改写现象：把某个带有非光滑耦合项的目标函数直接套进ADMM框架，结果迭代完全不收敛，或者z-update根本算不出来。原因：ADMM要求目标函数能拆成两个独立部分f(x)和g(z)，如果目标函数里有交叉项x.*z这种不可分项，把它硬塞进f或者g都是不对的，因为子问题不再可分离。解决：遇到这种情况，先做变量分裂，把交叉项用额外约束隔离。比如原问题里有x1*x2这种项，可以引入变量z并加约束x1 = z，让x2和z相互作用。具体做法因问题而异，总体上是要保证每个子问题都能解或者能近似解。如果实在无法拆开，ADMM不是你该用的工具，考虑换成近端梯度或者原始-对偶算法可能更合适。6. 验证ADMM正确性的实用技巧：从已知解到KKT残差6.1 构造已知解的反向验证拿到代码之后，第一步不应该是直接跑真实数据，而是构造一个已知的稀疏解。先生成随机的高斯测量矩阵A，用稀疏位置随机生成x_true，计算b = A*x_true，然后跑ADMM。因为b里没有噪声，lambda设一个很小的值，理论上最优解就是x_true本身。跑完之后比较恢复出的x和x_true，看最大误差在什么量级。6.2 用KKT条件残差做定量检查比解误差更严格的是KKT条件残差。ADMM的最优点必须满足：0属于A'(Ax - b) + lambda * sign(x)。这个包含集关系不好直接判断，但可以计算每个分量上的次梯度大小来近似。对于非零分量，残差应该接近0；对于零分量，次梯度的绝对值不应该超过lambda。这种检查能暴露主残差达标但解并不正确的问题。6.3 与现成求解器交叉对比如果环境里有优化工具箱，可以直接用quadprog或者CVX验证ADMM的解。LASSO可以用CVX一行代码求解，把两者的目标函数值对比，差值在1e-4以内就算合格。方法设计基础验证脚本：%% 基准测试：已知稀疏解 rng(2024); % 固定随机种子，保证可复现 m = 100; n = 300; A = randn(m, n); % 高斯测量矩阵 x_true = zeros(n, 1); support = randperm(n, 10); % 随机选择10个支撑位置 x_true(support) = randn(10, 1); b = A * x_true; % 无噪声观测 lambda = 0.01; rho = 0.5; % 实际应用中需要与A的尺度匹配 [x_est, history] = admm_lasso(A, b, lambda, rho, 1.0, 2000, 1e-6); recovery_err = norm(x_est - x_true, 2) / norm(x_true, 2); fprintf('恢复误差: %.4e\n', recovery_err);这种验证脚本的价值在于：一旦改动rho、终止条件或预处理方式，可以马上看出真实影响，而不会被真实数据的噪声干扰判断。从那以后，我拿到任何一份优化算法代码，第一步永远是先造已知解，跑通反向验证，再换真实数据，已经习惯到改不掉了。希望帮到你。本文还有配套的精品资源，点击获取
简介：一套基于MATLAB实现的ADMM（交替方向乘子法）算法代码包，主要面向机器学习、图像处理、信号处理等领域的算法研究者和工程师，用于求解含等式或不等式约束的凸优化问题。压缩包共383个文件，内容约17.44MB，涵盖171个m源文件、70个p加密文件、22个dat数据、22个mat数据，以及面向Windows/Mac/Linux的mexw64、mexmaci64、mexa64预编译文件，并附带pdf、md、txt等格式的说明资料，目录按源码、数据、文档分类存放，便于对照学习。目前已有四千余人浏览学习，是实践ADMM的成熟参考资料。代码主体ADMM.m完整实现了问题定义、参数设置、子问题分解、交替更新与乘子更新等核心步骤，配套辅助函数用于特定子问题求解；大量dat/mat数据文件可供稀疏回归、低秩矩阵恢复等经典场景直接测试，mex文件则提供预编译加速便于性能比较。通过阅读代码与pdf文档，可以深入理解ADMM的收敛条件、惩罚参数rho的影响及加速策略，适合希望快速上手算法复现、二次开发乃至分布式并行化的中高级用户。1. ADMM算法：为什么分布式优化和约束问题都绕不开它ADMM算法全称是交替方向乘子法（Alternating Direction Method of Multipliers），这几年在压缩感知、低秩矩阵恢复、分布式机器学习、图像去模糊这些场景里几乎成了标配解法。你只要手上有一个带约束的优化问题，目标函数还能拆成两块，那ADMM基本就是最省心的那一个工具。它的聪明之处在于：不跟复杂的原问题硬碰硬，而是把大问题拆成小步迭代，每一步只需要处理一个简单子问题，MATLAB里十几行循环就能搭出主干。这份资源给的就是一套可直接运行的MATLAB实现，能帮你绕过从论文公式到可运行代码之间的那段空档期。适合正在做调参、算法对比或毕设复现的从业者。2. ADMM的数学原理与MATLAB落地：从对偶上升、增广拉格朗日到三段迭代2.1 ADMM的核心思想：对偶上升、增广拉格朗日与交替极小化要理解ADMM，得先看它从哪里来。最原始的拉格朗日乘子法，把约束罚进目标函数，然后同时对原变量和对偶变量做梯度上升，这就是对偶上升法的思路。问题是对偶上升要求目标函数有比较强的光滑性，而且收敛条件偏紧。后来有人加了二次罚项，得到增广拉格朗日函数，这样对偶上升的收敛范围宽了很多。但仍有一个麻烦：增广拉格朗日里的二次罚项让子问题不再可分，遇到大规模问题就没法并行。ADMM做的事情其实很直接：假设你要解决的问题是min f(x) + g(z)，约束是 Ax + Bz = c它不直接对x和z做联合优化，而是把x和z分开轮流更新，每一轮只解一个子问题。这种交替极小化的好处是，两个子问题往往都有闭式解，或者能用很便宜的数值方法求解。配合对偶变量的上升更新，最终收敛到原问题的最优解。MATLAB里实现这套逻辑，核心就是下面三段式子：x_new = argmin_x ( f(x) + (rho/2) * norm(A*x + B*z - c + u, 2)^2 ); z_new = argmin_z ( g(z) + (rho/2) * norm(A*x_new + B*z - c + u, 2)^2 ); u_new = u + A*x_new + B*z_new - c;第一个式子更新原变量x，第二个更新分裂变量z，第三个是对偶变量u的上升步。rho是惩罚参数，控制收敛速度和数值稳定性。代码看起来简单，但每一步的求解才是真正的工作量。2.2 MATLAB实现的最小骨架：单约束问题三段循环在一个标准的单约束问题里，ADMM的MATLAB实现可以缩到最小骨架。下面这段代码展示的是最简单形态：function [x, z, u, history] = admm_simple(A, B, c, rho, max_iter, tol) % A, B: 约束矩阵 % c: 约束右端项 % rho: 惩罚参数（ADMM的核心超参） % max_iter: 最大迭代次数 % tol: 相对终止容差 m = size(A, 1); n = size(A, 2); p = size(B, 2); x = zeros(n, 1); z = zeros(p, 1); u = zeros(m, 1); history = zeros(max_iter, 2); for k = 1:max_iter % 更新 x：这里假设 f(x) = ||x||_2^2 即平方正则 % 闭式解来自求解线性方程组 x = (A'*A + rho*eye(n)) \ (A'*(c - B*z - u) * rho); % 更新 z：这里假设 g(z) = indicator_{z>=0} % 投影到非负象限 z = max(0, c - A*x - u); % 对偶更新 u = u + (A*x + B*z - c); % 记录主残差和对偶残差 history(k, 1) = norm(A*x + B*z - c, 2); history(k, 2) = norm(rho * A' * B * (z - z_old), 2); end end注意上面这段代码里，x的闭式解写成了一个线性方程组的求解。实际使用时要先算(A'A + rhoI)的分解，不要每次迭代都重新求逆。另外z的更新用了非负投影，这是g(z)为指示函数时的特例，换成别的正则项就要改成对应的近端算子。对偶残差的定义在不同论文里写法略有差异，我在这里用了带rho和A'B的形式。2.3 为什么用ADMM：与你见过的梯度下降和近端梯度的区别很多人在接触ADMM之前，已经用过梯度下降或者近端梯度法。梯度下降只适合光滑目标函数，近端梯度法能处理一部分非光滑项，但它们都要求目标函数整体结构比较规矩。ADMM的出发点不同：不试图一次处理整个目标函数，而是把目标拆成f和g两块，各自有各自的解法，只要能算近端算子，就能塞进ADMM框架。举个例子，LASSO问题里f是平滑的二次损失，g是L1范数。用近端梯度做，L1近端算子本身很简单，但梯度步的步长选择会影响收敛；而ADMM里g对应的z-update就是一个软阈值操作，没有步长参数，调的东西只有rho。这一点在实际使用中的意义很大——少一个超参，调参的麻烦就少一截。还有一个实际场景上的区别：分布式或联邦式的优化环境里，数据可能分布在不同的节点上。梯度下降需要汇总梯度，ADMM可以通过交替更新的方式让每个节点独立计算自己的子问题，最后只同步对偶变量。这一点在集中式的MATLAB实现里体现不明显，但理解了这个背景，你就知道ADMM的变量分裂结构不是偶然。3. 以LASSO为例的完整MATLAB实现：目标函数拆解与主循环写法3.1 LASSO问题的ADMM重写：目标函数拆解与变量分裂LASSO问题的标准形式是：min 0.5 * ||A*x - b||_2^2 + lambda * ||x||_1要让ADMM上场，得先做一步变量分裂，把L1范数“隔离”到z上。常见做法是引入z，并加约束x - z = 0，于是问题变成：min f(x) + g(z)，约束 x - z = 0其中f(x) = 0.5 * ||A*x - b||_2^2，g(z) = lambda * ||z||_1。这样拆完之后，x-update变成了一个带二次正则的最小二乘问题，它的解就是线性方程组的解，而z-update是对L1范数的近端算子，即软阈值操作。拆解的关键是变量分裂这一步的选择。不同的问题，分裂方式可以完全不同。比如有的场景是拆出两块目标函数，有的场景是拆出约束矩阵A和B。LASSO这种拆法是最直观的：把难啃的L1正则项放到一边，把光滑的二次损失放到另一边。3.2 主迭代循环：x-update、z-update与对偶更新下面给出一个完整的MATLAB函数，输入设计矩阵A、观测向量b、正则系数lambda和ADMM参数，输出解和收敛历史：function [x, history] = admm_lasso(A, b, lambda, rho, alpha, max_iter, tol) % admm_lasso - 用ADMM求解LASSO问题 % 输入： % A - m*n 的测量矩阵 % b - m*1 的观测向量 % lambda - L1正则系数（非负） % rho - ADMM惩罚参数（正数） % alpha - 松弛因子，通常取1.0或1.5，>1会加速但可能不稳 % max_iter- 最大迭代次数 % tol - 相对残差容差 % 输出： % x - n*1 的最优解 % history - 迭代历史（主残差、对偶残差、目标值） [m, n] = size(A); x = zeros(n, 1); z = zeros(n, 1); u = zeros(n, 1); % 预分解：x-update需要反复求解(A'A + rho*I)这个系统 % 先算Cholesky因子，后面每次迭代只做回代，省掉重复分解 [L, p] = chol(A'*A + rho*eye(n), 'lower'); if p ~= 0 error('矩阵(A''A + rho*I)必须正定。'); end Atb = A'*b; history = zeros(max_iter, 3); for k = 1:max_iter % x-update: 求解 (A'A + rho*I) * x = A'b + rho*(z - u) % 这一项来自对f(x)求导并令导数为零 x = L \ (L' \ (Atb + rho*(z - u))); % z-update: 软阈值 % 对每个分量执行 shrink(z + u, lambda/rho) zold = z; z = max(0, 1 - lambda/(rho*abs(x + u))) .* (x + u); % u-update: 对偶变量沿残差方向上升 u = u + alpha * (x - z); % 收敛历史记录 history(k, 1) = norm(x - z, 2); history(k, 2) = norm(rho * (z - zold), 2); history(k, 3) = 0.5*norm(A*x - b, 2)^2 + lambda*norm(x, 1); % 终止判定：主残差和原始残差都小于阈值才停 if history(k,1) < tol && history(k,2) < tol break; end end history = history(1:k, :); end逻辑说明：这段代码的循环体完全对应ADMM的三个基本步骤。x-update使用预分解的Cholesky因子做快速回代，z-update实现的是软阈值近端算子，u-update则是对偶上升。整体上这个写法与标准ADMM论文里的形式一致，只是加入了alpha松弛因子，这让代码可以在alpha=1.5时获得比alpha=1更快的收敛，但代价是可能出现轻微震荡。参数说明：lambda控制稀疏度，lambda越大解越稀疏，lambda = norm(A'*b, Inf)时解全为0，这是一个有用的上界；rho是ADMM的收敛速度调节器，太大导致收敛快但数值精度下降，太小导致迭代次数激增；alpha是松弛因子，通常取1.0，加速到1.5之前要先确认残差曲线没有明显的震荡。3.3 收敛参数与输入输出的定义边界输入参数这一块有几个边界条件值得说明。首先是lambda的范围，lambda必须是大于等于0的实数，如果输入负值，软阈值操作里的(1 - lambda/(rho*abs(x+u)))可能为负，导致z的符号被反向修正，整个算法的收敛无从谈起。其次是rho，rho不能为0，否则增广拉格朗日的二次罚项消失，ADMM退化回原始拉格朗日法，收敛条件变得极其苛刻。另外要留意的是max_iter和tol的组合。新手容易把tol设得太小，比如1e-10，但数值精度和机器epsilon有限，这种设定往往白跑几千轮还停不下来。合理的做法是先跑到固定迭代次数，画残差曲线，再决定tol的等级。上面的代码里，如果矩阵A的规模较大，最好把Cholesky分解换成基于QR或SVD的预处理版本，否则每次调用函数时的一次分解开销依然不小。常见做法是在函数外部完成分解，把因子传入迭代主体。4. 参数设置与收敛性判定：惩罚参数rho、松弛因子alpha与终止条件4.1 惩罚参数rho：收敛速度与数值稳定性的天平rho是ADMM里最重要的超参数，它的作用在增广拉格朗日函数里看得最清楚：二次罚项的系数。rho越大，每一步迭代对约束违反的惩罚越重，迭代会更快地让Ax + Bz - c趋近于零，但代价是原变量x的更新步长变小，可能导致收敛精度下降甚至在最优解附近来回震荡。rho越小，约束违反在早期被容忍，迭代过程平滑但很慢。实际调rho的方法，我一般习惯先在log空间里试几个值，比如0.01、0.1、1、10、100，各跑200轮，对比主残差收敛曲线，然后选曲线下降最快且没有明显震荡的那个量级。注意一个问题：rho的选择和问题本身的尺度强相关。如果A是0到1之间的数据，rho=1通常合理；如果A的数值在10^3量级，rho可能需要调到10^3以上才有反应。若rho的取值与A的尺度严重不匹配，会出现一种极端的现象：主残差前期快速下降，后期长期不动，而z迭代还在缓慢漂移。这时候强行加大rho并不会改善，反而是先把数据标准化更有效。4.2 松弛因子alpha：加速还是添乱alpha是ADMM的松弛因子，标准的ADMM里alpha=1，就是普通的三段迭代。当alpha>1时，对偶更新的方向更激进，可以在某些问题里把收敛速度提升不少，特别是对于二次损失加L1正则这类结构的问题。但alpha不是越大越好。实践中alpha超过1.8以后，残差曲线往往出现锯齿状波动，甚至发散。原因在于对偶变量u更新里的(1 - alpha)项会引入额外的振荡模式。如果是第一次调这个参数，建议只试1.0和1.5两个值，对比总的迭代次数，不要贪心。很多现成实现里，alpha是写死为1.0的。如果你手里的资源代码里有alpha这个参数，说明作者考虑到加速场景，这是一个加分项。但用的时候要清楚：alpha只在目标函数光滑性较好的那段收敛前期有明显收益，到了末段收敛区域，它的作用会变小。4.3 终止条件怎么设：双重残差判定与自适应rho终止条件是ADMM实现里最容易糊弄的部分。只看主残差norm(x - z)小于阈值就停，是对偶残差可能还很大，算法实际上并没有收敛到鞍点。正确的做法是同时检查主残差和对偶残差，两个都小于各自阈值才算完成。阈值的设置也有讲究。相对残差比绝对残差更适合跨问题使用：把主残差除以max(norm(Ax), norm(Bz), norm(c))，得到一个尺度无关的指标。下面给出一个可用的判定片段：eps_pri = tol * max([norm(A*x,2), norm(B*z,2), norm(c,2)]) + tol_abs; eps_dual = tol * norm(A' * u, 2) + tol_abs; stop = (r_pri < eps_pri) && (r_dual < eps_dual);这里tol_abs是绝对容差，通常取1e-6量级。tol是相对容差，工程上取1e-4到1e-6之间比较常见。如果你发现算法迟迟不终止，但残差曲线已经走平，先检查是不是tol设得太小，再看残差是压在哪个水平走平的——这往往能直接告诉你rho应该增大还是减小。还有一种进阶做法是自适应rho，即每隔若干轮根据主残差与对偶残差的比值调整rho的值。常见策略是当主残差大于对偶残差若干倍时加大rho，反之减小rho。自适应rho能省去手动搜索的麻烦，代价是额外几个参数：调整间隔、放大倍数、缩小倍数。下面给一个常规实现片段：if mod(k, 10) == 0 if r_pri > 10 * r_dual rho = rho * 2; % 注意：改变了rho之后，之前预分解的Cholesky因子失效，需要重新分解 elseif r_dual > 10 * r_pri rho = rho / 2; end end这里有个容易翻车的细节：rho一变，x-update里的(A'A + rho*I)整个矩阵都变了，必须重新做Cholesky分解，否则用的是旧因子，算出来的x方向不对。每10轮重分解一次，代价通常可以承受。5. ADMM避坑指南：参数失误、残差误判与实现细节5.1 迭代发散或出现NaN：rho与步长设置不当现象：运行代码没几轮，变量就变成了NaN或者Inf，程序直接报错。原因：最常见的是rho设得过大，导致x-update的线性系统数值上接近奇异；或者是数据矩阵A本身的条件数太大，比如A里的元素跨越了10^6的量级，这时候(A'A + rho*I)的条件数依然很糟糕，解出来的x方向有巨大误差。另一个常见原因是软阈值操作里分母上出现了abs(x+u) = 0的情况，在MATLAB里虽然不会直接报错，但结果会变成NaN。解决：先把rho往小调，比如从10000调到1或者0.1。同时检查A的数值范围，如果A的元素量级差异太大，先对列做标准化，再进ADMM。对于软阈值的数值保护，可以在abs(x+u)下面加上一个小常数eps=1e-12：max(0, 1 - lambda/(rho*(abs(x+u)+eps)))。这一步是血泪经验，MATLAB里的0/0结果是NaN，NaN会顺着迭代一路传染。5.2 只盯着主残差：对偶残差滞后导致的误判现象：main残差已经到了1e-6，看起来一切正常，但最终的解和真实最优解相比有明显偏差。原因：主残差衡量的是x和z的一致性，对偶残差衡量的是对偶变量的收敛程度。如果只以主残差作为终止条件，可能在对偶变量还在漂移的时候就停了。对偶残差反映了当前迭代点到最优鞍点的距离，两个残差都小才是真正的收敛。解决：终止条件必须写成同时满足主残差和对偶残差的阈值。另外，在观察收敛曲线时，把两条残差曲线画在同一张图上，如果发现主残差已经下降但对偶残差还高出一个量级，那就说明当前的rho太小，对偶更新太慢，把rho调大一点就能拉平两条曲线。5.3 大矩阵场景每次重复求逆：缓存分解结果现象：矩阵A的规模在几千乘几千，代码每跑一轮迭代都要花很长时间，但看CPU占用率却不高。原因：典型的实现问题，每轮迭代都执行(A'A + rhoI) \ something。MATLAB里反斜杠运算符会先判断矩阵结构并做分解，反复调用意味着反复做LU或Cholesky分解，而这一步的复杂度是O(n^3)，完全淹没了回代的O(n^2)开销。解决：进入循环之前先做一次性分解。如果矩阵对称正定，用Cholesky；如果不保证正定，用LU。上面的LASSO代码已经展示了这个技巧。如果rho在迭代中会自适应变化，分解必须放在rho更新之后重新做。这里有一个取舍：自适应rho的更新间隔不能太短，否则频繁重构因子反而拖慢整体速度。5.4 可变rho的震荡：更新间隔与边界限制现象：使用自适应rho之后，残差曲线不是平滑下降，而是明显震荡，甚至发散。原因：rho的更新过于频繁或者调整倍数太大。比如每轮迭代都检查残差比并翻倍rho，会让增广拉格朗日的权重大幅摆动，收敛过程完全被rho的调整节奏带着走。解决：把rho的调整间隔加大，至少10轮一次。另外，调整倍数不要超过2倍。同时给rho加一个上下界，比如rho_min = 1e-4，rho_max = 1e4，防止极端情况下rho连续翻倍导致数值溢出。还有一个细节：改变rho之后，对偶变量u也应该按比例缩放，不然旧的u在新尺度下失去了参考意义。标准做法是u = u * (rho_old / rho_new)。5.5 不可分目标直接套用：用变量分裂改写现象：把某个带有非光滑耦合项的目标函数直接套进ADMM框架，结果迭代完全不收敛，或者z-update根本算不出来。原因：ADMM要求目标函数能拆成两个独立部分f(x)和g(z)，如果目标函数里有交叉项x.*z这种不可分项，把它硬塞进f或者g都是不对的，因为子问题不再可分离。解决：遇到这种情况，先做变量分裂，把交叉项用额外约束隔离。比如原问题里有x1*x2这种项，可以引入变量z并加约束x1 = z，让x2和z相互作用。具体做法因问题而异，总体上是要保证每个子问题都能解或者能近似解。如果实在无法拆开，ADMM不是你该用的工具，考虑换成近端梯度或者原始-对偶算法可能更合适。6. 验证ADMM正确性的实用技巧：从已知解到KKT残差6.1 构造已知解的反向验证拿到代码之后，第一步不应该是直接跑真实数据，而是构造一个已知的稀疏解。先生成随机的高斯测量矩阵A，用稀疏位置随机生成x_true，计算b = A*x_true，然后跑ADMM。因为b里没有噪声，lambda设一个很小的值，理论上最优解就是x_true本身。跑完之后比较恢复出的x和x_true，看最大误差在什么量级。6.2 用KKT条件残差做定量检查比解误差更严格的是KKT条件残差。ADMM的最优点必须满足：0属于A'(Ax - b) + lambda * sign(x)。这个包含集关系不好直接判断，但可以计算每个分量上的次梯度大小来近似。对于非零分量，残差应该接近0；对于零分量，次梯度的绝对值不应该超过lambda。这种检查能暴露主残差达标但解并不正确的问题。6.3 与现成求解器交叉对比如果环境里有优化工具箱，可以直接用quadprog或者CVX验证ADMM的解。LASSO可以用CVX一行代码求解，把两者的目标函数值对比，差值在1e-4以内就算合格。方法设计基础验证脚本：%% 基准测试：已知稀疏解 rng(2024); % 固定随机种子，保证可复现 m = 100; n = 300; A = randn(m, n); % 高斯测量矩阵 x_true = zeros(n, 1); support = randperm(n, 10); % 随机选择10个支撑位置 x_true(support) = randn(10, 1); b = A * x_true; % 无噪声观测 lambda = 0.01; rho = 0.5; % 实际应用中需要与A的尺度匹配 [x_est, history] = admm_lasso(A, b, lambda, rho, 1.0, 2000, 1e-6); recovery_err = norm(x_est - x_true, 2) / norm(x_true, 2); fprintf('恢复误差: %.4e\n', recovery_err);这种验证脚本的价值在于：一旦改动rho、终止条件或预处理方式，可以马上看出真实影响，而不会被真实数据的噪声干扰判断。从那以后，我拿到任何一份优化算法代码，第一步永远是先造已知解，跑通反向验证，再换真实数据，已经习惯到改不掉了。希望帮到你。本文还有配套的精品资源，点击获取