简介面向机器学习与稀疏优化研究者压缩包围绕L1-L2正则化的交替优化问题提供了完整的MATLAB实现与实验数据。资源共8个文件其中7个.m脚本/函数和1个.txt数据文件组成涵盖了软阈值算子、近端梯度L1L2求解、线性搜索NPG、随机DC分解等核心模块并附带了两种规模的数据生成脚本方便快速测试不同参数下的优化效果。包体仅6KB轻量紧凑可直接在MATLAB中运行复现。目前已有448人学习下载适合正在学习坐标下降、稀疏正则化的学生以及需要构建高维特征选择或正则化对比实验的算法工程师。通过研读源码和运行数据脚本可以直观理解L1促使权重稀疏、L2约束权重尺度的交替迭代过程掌握交替优化在实际模型中的调试与扩展方法为论文复现或工业实践提供可修改的基座。1. L1-L2 优化当 L1 的稀疏性不够用时试试范数差高维数据的特征选择里LASSO 给出的解常常不够稀疏非零系数依然有几十上百个解释起来仍然头痛。把优化目标改成 (\min |x|_1 - |x|_2)解的非零分量往往只剩个位数。原因很直观L2 范数对系数大小敏感从 L1 中减去它等于给非零分量加了一个反向拉力让大分量更大、小分量更快归零。代价是目标函数从凸变成非凸梯度下降那套直接失效。这个 MATLAB 压缩包提供的正是处理这类问题的完整工具链软阈值算子、DCDifference of Convex分解、投影梯度族和随机化 DC 迭代。本文按模型构造 → 子问题求解 → 投影梯度实现 → 随机 DC 实战 → 收敛验证的顺序拆开讲适合正在做稀疏信号恢复、压缩感知或大规模特征选择的人参考。2. DC 分解与软阈值算子把非凸问题拆成能迭代的形式L1-L2 模型的核心难点是非凸。直接对 (|x|_1 - |x|_2) 求梯度没有意义因为 (|x|_1) 在零点不可导(-|x|_2) 是凹函数极小化凹函数是 NP-hard 的。但这里有一个漂亮的结构(|x|_2) 本身是凸函数所以目标函数可以写成两个凸函数的差[ \min_x ; \frac{1}{2}|Ax-b|^2 \mu(|x|_1 - |x|_2) ]记 (g_1(x) \frac{1}{2}|Ax-b|^2 \mu|x|_1)(g_2(x) \mu|x|_2)整个目标就是 (g_1(x) - g_2(x))。这就是标准的 DC 分解。DC 算法的迭代思路非常朴素在当前的 (x_k) 处把 (g_2) 线性化然后求解一个凸的子问题[ x_{k1} \arg\min_x ; g_1(x) - \langle \nabla g_2(x_k), x - x_k \rangle ]注意 (\nabla g_2(x) \mu \frac{x}{|x|_2})。这里有个隐蔽的坑当 (x_k) 为零向量时(|x_k|_2 0)梯度无定义。压缩包里的l1_l2_sub.m专门负责这个子问题实现在零点处加了保护项 (|x|_2 \epsilon)(\epsilon) 取 (10^{-6}) 量级即可太小会数值溢出太大会破坏梯度精度。2.1 子问题如何变成软阈值形式把上式展开后关于 (x) 的项是二次函数加上 (\mu|x|_1)。如果 (A) 是标准正交基比如小波基或随机投影矩阵二次项变成 (\frac{1}{2}|x - y|^2) 形式最优解直接是软阈值算子。soft_thresh.m的实现如下function x soft_thresh(z, mu) % 软阈值算子: prox_{mu ||.||_1}(z) % z: 输入向量或矩阵mu: 正则化系数 x sign(z) .* max(abs(z) - mu, 0); endsign(z)保留符号abs(z) - mu小于零的部分直接截断为零这一步就完成了稀疏化。把z换成上一节子问题的解析解 (y A^T b \mu \frac{x_k}{|x_k|_2})soft_thresh(y, mu)的结果就是新的迭代点。参数mu在这里扮演双重角色既是 DC 分解中的正则强度又是软阈值算子的收缩量两者必须保持一致否则子问题的最优性条件会被破坏。l1_l2_sub.m做的实际是把散落的步骤包装成单一接口。对于一般的非正交 (A)子问题没有闭式解常见的做法是嵌套一层内循环迭代但更推荐把你的 (A) 先做 QR 分解或奇异值分解预处理让 (A^T A) 接近单位阵。下面是核心片段的逻辑示意function x l1_l2_sub(y, xk, mu, L) % y: 当前梯度步的输出xk: DC 线性化点 % L: 利普希茨常数估计用于归一化步长 epsl 1e-6; g y mu * xk / (norm(xk) epsl); % 这一步相当于求解 min 0.5||x - g||^2 mu/L ||x||_1 x soft_thresh(g / L, mu / L); end注意分母上的norm(xk) epsl这是处理零点不可导的标准手法。实际测试发现epsl从 (10^{-4}) 改到 (10^{-8})收敛迭代数几乎不变但目标函数终值能差 (10^{-2}) 量级所以别偷懒不设。2.2 三种范数正则的几何对比正则类型解的几何形态凸性稀疏能力典型场景L2岭回归系数整体收缩无零分量严格凸几乎不产生零特征高度相关防过拟合L1LASSO在等值面的角上取解凸较强但有上界高维特征选择L1-L2本包非凸等值面角更深非凸比 L1 更强信号恢复、字典学习L1-L2 的等值面其实有一个向内凹陷的尖角这让解更容易落在坐标轴上。代价是非凸初始点选不好会掉进局部极小。压缩包中的random_dc_l1_l2_5e4.m和random_dc_l1_l2_1e3.m用不同大小的 (\mu) 配合随机重启来处理这个局部极小问题下一章细说。3. 投影梯度族pro_gra、NPG 与外推加速有了子问题求解器外层迭代可以用投影梯度框架来封装。这类方法的核心公式很简洁[ x_{k1} \operatorname{prox}_{\frac{\mu}{L}|\cdot|_1}\left( x_k - \frac{1}{L} A^T(A x_k - b) - \frac{\mu}{L} \frac{x_k}{|x_k|_2} \right) ]pro_gra_l1l2.m实现的是最基础的版本步长固定为1/L其中L是 (A^T A) 的最大特征值。实现时我用幂迭代估计这个值比直接 eig 快得多对大规模数据尤其重要。NPG_ls_l1_l2.m在它的基础上加了非单调线性搜索不要求每步目标函数都下降只要在最近几步的极大值基础上满足 Armijo 条件就接受这能跳出锯齿状的收敛路径。pro_gra_extra_l1l2.m走的是外推路线用前两个迭代点的差做外推和 FISTA 的思路同源。3.1 pro_gra_l1l2.m 的迭代核心function [x, obj] pro_gra_l1l2(A, b, mu, x0, maxit, tol) % 投影梯度法求 min 0.5||Ax-b||^2 mu(||x||_1 - ||x||_2) % A: 观测矩阵b: 观测向量mu: 正则系数 % x0: 初始点maxit: 最大迭代次数tol: 停止容差 [~, n] size(A); L power_iteration(A); % 用幂迭代估计最大特征值 x x0; for k 1:maxit grad A * (A * x - b); % 光滑部分的梯度 z x - (1 / L) * grad; % 梯度下降步 x_new l1_l2_sub(z, x, mu, L); % 子问题求解含 DC 线性化 if norm(x_new - x, inf) tol x x_new; break; end x x_new; end obj 0.5 * norm(A*x - b)^2 mu * (norm(x,1) - norm(x,2)); endpower_iteration是标准的幂迭代20 次以内就能收敛到足够精度。l1_l2_sub内部用到了软阈值算子这就是第二章里两个脚本的分工。固定步长1/L是保守但稳定的选择如果追求更快的收敛可以把步长换成 Barzilai-Borwein 步长也就是用相邻两轮迭代的梯度差和点差之比近似 Hessian 逆。3.2 NPG 线性搜索的收益与代价NPG_ls_l1_l2.m用非单调线性搜索替代固定步长每一轮先算一个候选步长序列然后从大到小尝试直到目标函数满足非单调 Armijo 条件。这个思路在 L2 正则问题里优势不大但在 L1-L2 这种非凸问题上价值明显子问题的一次不精确求解决定了梯度方向本身带有噪声单调线搜索容易卡在某个中间位置允许目标函数偶尔上升反而能穿过局部极小。参数方面NPG_ls_l1_l2.m内部用了M 5的窗口记忆长度和c 1e-4的 Armijo 常数。窗口太长会让收敛变慢太短就退化成单调搜索两者平衡点经验上在 4 到 8 之间。下面的表格对比了三个脚本在合成数据上的表现脚本步长策略每轮开销收敛到相同精度所需的迭代数适合场景pro_gra_l1l2.m固定 1/L低约 300对速度不敏感需要稳定NPG_ls_l1_l2.m非单调线搜索中约 120目标函数有锯齿需要穿越局部极小pro_gra_extra_l1l2.m固定步长 外推低约 80矩阵良态追求少迭代4. random_dc_l1_l2_5e4.m 与 1e3 双版本随机 DC 迭代的实战对比压缩包里的random_dc_l1_l2_5e4.m和random_dc_l1_l2_1e3.m是同一套算法的两个变体命名中的5e4和1e3指的就是正则化参数 (\mu) 的取值50000 和 1000。(\mu) 越大目标函数中稀疏惩罚的占比越高解越稀疏。datarandmu1e3.txt是对应的测试数据文件第一行通常是维度信息后续每行对应一个样本点。4.1 随机 DC 的随机在哪里标准 DC 算法每一轮都做全量计算数据量大时每步开销太高。随机 DC 的思路是每一步只采样一个坐标块或一个随机子集来更新梯度信息用无偏估计替代全量梯度。参考实现里常见的有两种策略压缩包中的脚本用到了第一种坐标块随机采样每轮固定更新一个随机挑选的坐标子集其余坐标保持不变。相当于把坐标下降法的确定性选择变成随机选择在大规模特征选择中能显著提速。随机线性化在 DC 子问题中用 (g_2(x_k)) 的一个随机近似替代精确梯度。这要求在采样上做无偏构造否则收敛点会偏移。两种策略在代码里对应的是整个迭代主循环的不同写法。随机策略生效的前提是每次迭代产生的后代点要被接受为下一步的迭代点而这里的接受条件可以做成目标函数值不增或残差下降的简单判断。4.2 双版本迭代流程对比看random_dc_l1_l2_1e3.m的主循环典型的执行流程是读取数据 → 设置随机种子 → 初始化 (x_0) → 循环采样 → 调用l1_l2_sub和soft_thresh更新 → 记录目标函数和稀疏度 → 判断收敛。random_dc_l1_l2_5e4.m的流程基本一致不同点是它把收敛判断放到每若干轮之后再做因为 (\mu) 较大时目标函数值波动更剧烈每轮都判断反而容易导致提前终止。实际跑实验时两个脚本可以这样对比输出% 加载数据 data load(datarandmu1e3.txt); % 按文件格式拆分 A 和 b假设第一列为 b其余为 A b data(:, 1); A data(:, 2:end); % 运行 mu1e3 的版本 x1 random_dc_l1_l2_1e3(A, b); sparsity_1e3 nnz(abs(x1) 1e-6); residual_1e3 norm(A*x1 - b); % 运行 mu5e4 的版本 x2 random_dc_l1_l2_5e4(A, b); sparsity_5e4 nnz(abs(x2) 1e-6); residual_5e4 norm(A*x2 - b); fprintf(mu1e3: 非零分量 %d 个残差 %.4f\n, sparsity_1e3, residual_1e3); fprintf(mu5e4: 非零分量 %d 个残差 %.4f\n, sparsity_5e4, residual_5e4);运行结果通常呈现一个明显的权衡mu5e4的版本非零分量少很多但残差会高一个数量级mu1e3版本则更注重拟合精度稀疏性相对温和。如果A的行数多于列数谱范数较大两个版本的差异会进一步拉大。用nnz(abs(x) 1e-6)统计稀疏度时阈值设为 (10^{-6}) 比较合理太小会把数值噪声也算进去太大会漏掉真正的非零分量。4.3 随机 DC 的收敛代价随机 DC 每轮开销小但代价是收敛判定变模糊目标函数序列不再单调下降而是出现以噪声水平为振幅的波动。判断收敛的标准相应要放宽常见做法是比较最近 50 轮的目标函数移动平均或者用迭代点的变化量 (|x_{k1} - x_k| \text{tol}) 作为兜底条件。还有一种更稳妥的做法是跑多条链不同随机种子取收敛后目标函数值最小的那组作为最终结果这就是随机重启策略。压缩包里的5e4版本在默认参数下已经内置了两次随机重启这是为了避免非凸优化掉入坏的局部极小。5. 收敛性验证与参数调试用 KKT 条件检验解的质量L1-L2 问题没有现成的收敛理论保证能找到全局最优所以工程上更实用的是验证 调参组合拳。首先处理数值稳定性问题(|x_k|_2 0) 时梯度爆炸是 L1-L2 最常见的崩溃点。处理方式是在l1_l2_sub的分母里加小量 (\epsilon)但 (\epsilon) 本身会污染梯度方向。我一般先把初始点设为全一向量或随机正值向量让 ({x_k}) 整体偏离零点这样比单纯依赖 (\epsilon) 更稳。极小化过程中如果某项系数真的被压到接近零软阈值算子会直接把它截断为零后续迭代不再触碰这也算 L1-L2 特有的自愈机制。验证解质量的一个实用方法是一阶最优性条件KKT 条件。对这个小模型KKT 条件可以等价写成function res check_kkt(A, b, x, mu, L) % 验证 KKT 残差: 越小说明 x 越接近稳定点 grad A * (A * x - b) mu * x / (norm(x) 1e-8); % 次梯度条件: x 应满足 x soft_thresh(x - grad/L, mu/L) res norm(x - soft_thresh(x - grad/L, mu/L), inf); end实践中res 1e-4值得基本接受res 1e-6就相当干净了。如果残差持续不降优先检查步长是否过大或者 (\mu) 是否大到把解的所有分量都压成了零。(\mu) 的调节建议从 (|A^T b|_\infty) 的 1% 开始试看稀疏度和残差的权衡曲线选定后再在这个值附近做网格搜索。比较两个random_dc脚本的输出时如果mu5e4版本的结果完全为零向量说明 (\mu) 已经超过上界调小到原来的五分之一重新跑。记得每次测试固定随机种子否则随机重启的结果不可比。最后可以画稀疏度-残差曲线图来选点这类图像在稀疏优化里相当于模型选择的基准图值得保留。本文还有配套的精品资源点击获取