CVXPY 最小长度最小二乘（Minimum-Length Least Squares）求解实战：基于 DQCP 的准凸规划指南

科学计算【免费下载链接】cvxpyA Python-embedded modeling language for convex optimization problems.项目地址https://gitcode.com/gh_mirrors/cv/cvxpy点击查看免费下载导读本文以 CVXPY 官方示例 doc/source/examples/dqcp/minimum_length_least_squares.rst 为骨架讲解如何在 CVXPY 中建模并求解最小长度最小二乘问题在保证均方误差MSE不超过给定阈值的前提下寻找长度最后一个非零分量的下标最短的向量解。该问题是典型的拟凸规划Quasiconvex Program, QCP可通过**有纪律的拟凸规划DQCP**框架声明式建模并借助二分法bisection求解。读完本文你将掌握cp.length拟凸原子、DQCP 建模语法、problem.solve(qcpTrue)求解入口以及 CVXPY 内部将 QCP 归约为参数化 DCP 问题并二分求解的完整原理并可在自己的稀疏解恢复类任务中直接复用。一、问题定义什么是最小长度最小二乘最小长度最小二乘问题的目标是找到一个长度最短的向量 x ∈ Rⁿ使其对给定最小二乘系统的均方误差足够小。其标准形式为minimize len(x) subject to (1/n)·‖Ax - b‖₂² ≤ ε其中变量x ∈ Rⁿ问题数据矩阵 A、观测向量 b、样本数 n、误差容忍度 ε目标函数len(x)即向量 x 的长度——最后一个非零分量的下标按 1 计数若 x 全为零向量则长度为 0约束归一化均方误差(1/n)·‖Ax - b‖₂² ≤ ε。直观理解给定存在噪声或冗余的线性方程组 Ax ≈ b我们希望用尽量稀疏有效分量集中在开头、尾部尽快归零的解逼近它这与稀疏信号恢复、特征选择等场景目标一致。从 CVXPY 源码看length原子定义于 cvxpy/atoms/length.py其数值实现使用容差ATOM_EVAL_TOL判定非零分量再取最大非零下标加 1仅接受向量输入对非向量会抛出ValueError(length can only be applied to vectors.)cvxpy/atoms/length.py数值计算中np.abs(x) ATOM_EVAL_TOL判定非零若无非零元素返回 0否则返回np.max(nz) 1cvxpy/atoms/length.py。二、为什么是 QCP/DQCP曲率分析关键问题是len(x)与 MSE 约束组合后问题本身不是凸问题但属于拟凸问题。len(x)既不凸也不凹源码中is_atom_convex与is_atom_concave均返回False但它满足**拟凸quasiconvex**性质is_atom_quasiconvex返回Truecvxpy/atoms/length.py。其子水平集{x : len(x) ≤ t}是凸集——当 t ≥ n 时退化为全空间否则等价于x[t:] 0显然为凸集MSE 约束(1/n)·‖Ax - b‖₂² ≤ ε是凸约束平方范数凸拟凸目标 凸约束构成的极小化问题即为拟凸规划QCP。CVXPY 通过**有纪律的拟凸规划DQCP**语法支持此类问题在length上还规定其单调性——对非负参数单调不减、对非正参数单调不增is_incr/is_decr见 cvxpy/atoms/length.py这是 DQCP 组合规则单调函数复合拟凸原子成立的前提。注意length的单调性依赖参数符号。示例 cvxpy/tests/test_dqcp.py 验证了length(abs(x))单调不减故为 DQCP而length(abs(x) - 1)不再单调is_dqcp()返回False。三、完整可运行示例建模与求解以下代码完整复刻官方示例问题数据构造与示例一致可直接在安装了 CVXPY 与 NumPy 的 Python 3 环境中运行。# 安装最新版 CVXPY如尚未安装 # !pip install --upgrade cvxpy import cvxpy as cp import numpy as np # 1. 构造问题数据 n 10 np.random.seed(1) A np.random.randn(n, n) x_star np.random.randn(n) # 真实解作为对照 b A x_star # 无噪声观测 epsilon 1e-2 # 均方误差容忍度 # 2. 声明变量与目标/约束表达式 x cp.Variable(n) mse cp.sum_squares(A x - b) / n problem cp.Problem(cp.Minimize(cp.length(x)), [mse epsilon]) # 3. 验证 DQCP 属性并求解 print(Is problem DQCP?: , problem.is_dqcp()) problem.solve(qcpTrue) print(Found a solution, with length: , problem.value)求解输出示例文档中的运行结果Is problem DQCP?: True Found a solution, with length: 8.0进一步查看均方误差与最优解print(MSE: , mse.value) print(x: , x.value) print(x_star: , x_star)示例运行输出MSE: 0.00926009328813662 x: [-2.58366030e-01 1.38434327e00 2.10714108e-01 9.44811159e-01 -1.14622208e00 1.51283929e-01 6.62931941e-01 -1.16358584e00 2.78132907e-13 -1.76314786e-13] x_star: [-0.44712856 1.2245077 0.40349164 0.59357852 -1.09491185 0.16938243 0.74055645 -0.9537006 -0.26621851 0.03261455]结果解读problem.value为 8.0即最优解 x 的最后一个非零分量下标为 8第 9、10 个分量在数值容差内为 0与官方示例结果一致mse.value ≈ 0.00926 ε 1e-2约束满足对照x_star可见解的前 8 个分量与真实解接近后两个分量被压缩至 1e-13 量级数值零实现了尾部截断式的稀疏化。官方示例被 cvxpy/tests/test_dqcp.py 的test_length_example完整收录注释Fix #1760作为回归测试断言problem.value与 8 数值接近进一步印证了上述输出具有可复现性。四、求解原理DQCP → 参数化 DCP 二分法problem.solve(qcpTrue)并非直接求解 QCP而是走一条归约 二分的流水线DQCP 归约Dqcp2Dcp将拟凸目标len(x)替换为带参数 t 的不等式len(x) ≤ t原问题变成一族参数化的 DCP 可行性问题{x : mse ≤ ε, len(x) ≤ t}若目标非负则 t 声明为非负 Parameter否则为一般 Parametercvxpy/reductions/dqcp2dcp/dqcp2dcp.py。其中惰性约束lazy constraints以可调用对象形式存在仅在 t 取具体值时求值子水平集改写拟凸原子length的子水平集按 cvxpy/reductions/dqcp2dcp/sets.py 的length_sub展开当 t 为 Parameter 时若 t 0 视为不可行、t ≥ n 视为无约束否则生成约束x[floor(t):] 0对 t 向下取整当 t 为常数时直接生成该等式约束。这就是len(x) ≤ t ⟺ 尾部 t 个分量全为 0的 DCP 化表达二分求解bisectcvxpy/reductions/solvers/bisection.py 的bisect函数在 t 上做二分先求解可行性问题判断原问题是否可行若未给定区间则通过_find_bisection_interval指数扩张寻找 [low, high]cvxpy/reductions/solvers/bisection.py随后在_bisect中反复取中点 t (lowhigh)/2 求解参数化 DCP直至区间宽度小于 eps默认 1e-6或达到最大迭代次数默认 100整数目标加速length属于整数取值原子与ceil、floor并列cvxpy/reductions/dqcp2dcp/tighten.py 为其提供收紧函数不可行侧向上取整np.ceil、可行侧向下取整np.floor使二分区间快速收敛到精确整数值这正是示例能精确得到 8.0 而非 7.9 之类近似值的原因。可选的verboseTrue参数会逐次打印二分迭代的下界、上界与查询点cvxpy/reductions/solvers/bisection.pycvxpy/tests/test_dqcp.py 还提供了verboseTrue下length原子不崩溃的回归测试。五、DQCP 建模规则速查与注意事项从示例及 CVXPY 源码可归纳出以下使用要点入口开关DQCP 问题的求解必须显式传qcpTrue即problem.solve(qcpTrue)否则 CVXPY 会按 DCP 路径处理并报错属性自检建模后先用problem.is_dqcp()、expr.is_quasiconvex()/expr.is_quasiconcave()验证避免进入求解阶段才发现语法错误。源码曲率测试见 cvxpy/tests/test_dqcp.pylength(x)曲率为 QUASICONVEX-length(x)为 QUASICONCAVE单调性组合规则DQCP 允许单调函数 ∘ 拟凸原子 ∘ DCP 表达式的复合但单调性必须成立见 cvxpy/reductions/dqcp2dcp/dqcp2dcp.py 的约束改写逻辑单调函数通过求逆作用到约束两侧拟凸原子替换为其子水平集/超水平集。破坏单调性如length(abs(x) - 1)会导致is_dqcp()为Falselength 的输入限制length只能作用于向量源码显式校验矩阵输入会抛出异常容差语义length的数值结果依赖ATOM_EVAL_TOL容差非零判定是按数值容差而非精确非零因此打印出的解尾部会出现 1e-13 量级的数值残差求解器选择二分内部使用常规凸求解器求解每个子问题默认选择可用求解器即可若需手动指定可传problem.solve(qcpTrue, solvercp.SCS)等。测试 cvxpy/tests/test_dqcp.py 展示了length加等式约束的小规模验证期望值为 2可作为最小复现样例。六、延伸应用从示例到实际场景最小长度最小二乘示例虽小其拟凸目标 凸约束的建模范式可推广到大量实际问题稀疏解恢复将目标换成len(x)或拟凸的sign、dist_ratio等原子配合 L2 拟合约束可在不做 L1 松弛的情况下直接求最短有效支撑解特征/模型选择在误差上限约束下最小化有效参数量得到尾部截断的稀疏系数向量通用 DQCP 工具箱CVXPY 在 cvxpy/reductions/dqcp2dcp/sets.py 中注册了multiply乘积、DivExpression比率、length、sign、dist_ratio、gen_lambda_max、condition_number等原子的子水平集/超水平集改写同类拟凸目标如条件数最小化、广义最大特征值均可套用本文的建模与求解流程。结语本文从 doc/source/examples/dqcp/minimum_length_least_squares.rst 出发完整给出了最小长度最小二乘问题的数学模型、CVXPY 可运行代码、DQCP 属性分析与结果解读并深入源码揭示了其底层实现length拟凸原子cvxpy/atoms/length.py→ DQCP 归约cvxpy/reductions/dqcp2dcp/dqcp2dcp.py→ 子水平集改写cvxpy/reductions/dqcp2dcp/sets.py→ 二分求解cvxpy/reductions/solvers/bisection.py与整数收紧cvxpy/reductions/dqcp2dcp/tighten.py。掌握这套拟凸建模 qcpTrue 二分归约的方法论后你可以在 CVXPY 中自信地处理更大规模的拟凸规划问题。赞分享科学计算【免费下载链接】cvxpyA Python-embedded modeling language for convex optimization problems.项目地址https://gitcode.com/gh_mirrors/cv/cvxpy点击查看免费下载相关推荐CVXPY 天线阵列波束宽度最小化基于二分法求解准凸优化问题的完整实战CVXPY 天线阵列波束宽度最小化基于二分法求解准凸优化问题的完整实战 导读 本文以 CVXPY 官方示例 ant_array_min_beamwidth.r科学计算使用CVXPY求解最小二乘问题原理与实践使用CVXPY求解最小二乘问题原理与实践 最小二乘问题概述 最小二乘法Least squares是统计学和机器学习中最基础的线性回归方法之一。该方法通过最科学计算使用 CVXPY 求解时钟网格尺寸优化最小化功耗的半定规划SDP实战指南使用 CVXPY 求解时钟网格尺寸优化最小化功耗的半定规划SDP实战指南 本文基于 CVXPY 官方示例 doc/source/examples/appl科学计算上一篇认知心理画像用「意识阶梯 ELM 双路径」精准校准信息传播策略的 AAS 技能实战指南下一篇xbar 插件升级指南从 BitBar 迁移到 xbarmetadata、shell、paramN、变量与键盘快捷键创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考