格拉斯曼流形上的核范数正则化最小二乘优化光看这个标题就知道不是个省油的灯——三个词拆开都认识合在一起就让人有点发怵格拉斯曼流形、核范数、最小二乘优化。我最初接触这篇工作的时候也犯嘀咕流形优化那套黎曼梯度、指数映射、切空间投影已经够折腾了再叠一个核范数这还能收敛吗直到把目标函数拆成两个子问题逐个看才发现事情远比想象中干净一个子问题是Procrustes问题闭式解直接写在脸上另一个子问题是奇异值软阈值一步到位连步长都不需要调。这篇是我们的系列第三篇早期工作技术路线很朴素但正是这种“朴素”让它成了后续很多工作的地基。如果你正在做子空间学习、低秩建模、或者视觉特征处理这篇文章值得你坐下来慢慢读我尽量把里面的几何直觉和实现细节都讲透。1. 这套组合拳要解决的问题子空间估计和低秩系数为什么必须一起优化先摆出最常见的场景。你有一堆观测数据排成矩阵 $X \in \mathbb{R}^{D \times n}$行是特征维度列是样本或者帧。无论是图像序列、基因表达谱、还是用户行为日志这类数据经常存在一个天然的低维结构列向量其实落在某个 $r$ 维子空间附近只是被噪声、缺失值、离群点糊住了。教科书做法是两步走先用主成分分析得到子空间 $U$把 $X$ 投影到 $U$ 上得到系数 $W$再做后续的分类、聚类、回归。听起来顺理成章但问题恰恰出在“先估计子空间”这件事上。PCA估计出的 $\hat{U}$ 本身就带噪声子空间方向只要偏一点点后续所有系数都会继承这个偏差而且是系统性偏差没办法用下游模型救回来。更麻烦的是当你希望系数矩阵额外满足低秩结构时——比如视频帧之间高度线性相关、多视角特征共享同一个低维表示——两步走的方案彻底失效因为你没有任何机制在估计子空间的同时逼着系数低秩。所以这篇工作想的很直接把子空间 $U$ 当作流形上的变量把系数矩阵 $W$ 当作另一个变量两者放进同一个目标函数里联合优化。目标函数长这样$$\min_{U \in \mathrm{Gr}(r,D),; W \in \mathbb{R}^{r \times n}} \frac{1}{2} | X - U W |F^2 \lambda | W |*$$其中 $\mathrm{Gr}(r,D)$ 是格拉斯曼流形$| W |_*$ 是核范数。第一项是重构误差第二项是低秩正则。$U$ 在流形上跑$W$ 在矩阵空间里跑中间靠 $X$ 耦合成一个整体。为什么这个耦合值得做因为子空间估计和低秩系数本质上是同一枚硬币的两面。如果 $U$ 估计得准$WU^T X$ 的奇异值谱会自然衰减核范数稍微推一把就能得到明显的低秩结果反过来如果 $W$ 的低秩约束反馈到 $U$ 的更新里$U$ 会更偏向那些“能承载低秩表达”的方向而不是被噪声方向带跑。这个思路在数学上还有个非常优雅的副产品两个子问题都有闭式解下面会详细展开。但先记住核心结论——联合优化不是锦上添花而是低秩结构存在的场景下两步法在原理上就输了一截。2. 必备数学背景格拉斯曼流形与核范数到底扮演什么角色2.1 格拉斯曼流形子空间本身是一个“点”大部分人对流形的第一印象是曲面上的点在流动但格拉斯曼流形 $\mathrm{Gr}(r,D)$ 的定义有点与众不同它上面的每个“点”是整个 $\mathbb{R}^D$ 里的一个 $r$ 维线性子空间。你没看错点是子空间而不是子空间里的某个向量。这个定义带来的第一个问题是一个子空间怎么用矩阵表示最自然的办法是找它的正交基 $U \in \mathbb{R}^{D \times r}$满足 $U^T U I_r$。但正交基不唯一——右乘任意正交矩阵 $Q \in O(r)$ 后$UQ$ 的列空间和 $U$ 完全一样。所以格拉斯曼流形本质上是Stiefel流形再取商$$\mathrm{Gr}(r,D) \mathrm{St}(D,r) / O(r)$$其中 $\mathrm{St}(D,r)$ 是所有 $D \times r$ 列正交矩阵的集合。商结构意味着优化时不能拿 $U$ 和 $UQ$ 当成两个不同的点它们必须被视为同一个对象。这也是为什么很多人直接用Stiefel流形上的优化器会踩坑自由度冗余会让迭代步态歪歪扭扭收敛变慢。格拉斯曼流形上怎么衡量两个子空间差多远用主角度principal angles。对 $U_1, U_2$先算 $U_1^T U_2$ 的奇异值 $\sigma_i$主角度就是 $\theta_i \arccos(\sigma_i)$。所有主角度合起来就是两个子空间之间的几何距离。这个度量在后面做收敛判定时非常重要千万不要用 $|U_1 - U_2|_F$ 去判断因为 $U_1$ 和 $U_2$ 各自右乘不同的正交矩阵都表示同一个子空间F范数会给出完全错误的“距离”。2.2 核范数向量稀疏的“矩阵版本”核范数 $|W|_* \sum_i \sigma_i(W)$ 是奇异值的 $L_1$ 范数。直接约束 $\mathrm{rank}(W) \le k$ 是非凸的、带组合爆炸的问题而核范数是它在单位球内的凸包络松驰后的优化问题保持凸性又能得到尽量低秩的解。核范数能发挥作用的关键是它的近端算子有闭式解。给定矩阵 $Y$考虑最小化$$\min_Z \frac{1}{2} |Y - Z|F^2 \lambda |Z|*$$答案是奇异值软阈值算子SVT$$\mathrm{SVT}_\lambda(Y) U_Y , \mathrm{diag}\left(\max(\sigma_i - \lambda, 0)\right) V_Y^T$$其中 $Y U_Y \Sigma V_Y^T$ 是奇异值分解。效果很直观小于 $\lambda$ 的奇异值直接清零大于 $\lambda$ 的每个减去 $\lambda$。这就是把 L1 软阈值推广到奇异值谱上——向量里的元素被单独软阈值化矩阵里的奇异值也被单独软阈值化逻辑完全一致。2.3 两者怎么搭起来流形管理的是“子空间的方向”核范数管理的是“子空间系数矩阵的谱”。一个是几何对象一个是谱对象看似八竿子打不着但在我们的目标函数里正好各管一摊$U$ 负责描述数据落在哪个子空间附近$W$ 负责描述在这个子空间里的表达要尽量低秩。$U$ 动不动用流形几何$W$ 动不动用谱分析然后通过最小二乘项耦合。这种分工带来的直接好处是两个变量的更新规则各自独立、各自闭式不用像一般非凸优化那样在切空间上小心翼翼地挪步子。3. 模型与优化U子问题是ProcrustesW子问题是软阈值3.1 固定 U 更新 W投影后再做奇异值软阈值先固定 $U$看关于 $W$ 的子问题。因为 $U^T U I$可以做一个非常舒服的正交分解$$| X - U W |_F^2 | U^T X - W |_F^2 | X |_F^2 - | U^T X |_F^2$$推导只需要两步展开$X - UW$ 先投影到 $U$ 的列空间得到 $U(U^TX - W)$再投影到正交补空间得到 $(I - UU^T)X$两个部分正交所以F范数平方直接相加。第二项和第三项都与 $W$ 无关于是整个子问题等价于$$\min_W \frac{1}{2} | U^T X - W |F^2 \lambda | W |*$$这就是前面说的近端问题答案标准得不能再标准$$W^{(k1)} \mathrm{SVT}_\lambda\left(U^{(k)T} X\right)$$注意这里的 $U^T X$ 天然是一个 $r \times n$ 的小矩阵。$r$ 通常远小于特征维度 $D$所以这个奇异值分解非常便宜每次迭代只求一个 $r \times n$ 矩阵的 SVD而不是在原始 $D \times n$ 尺度上操作。3.2 固定 W 更新 U一个标准的 Procrustes 问题再看固定 $W$ 时分关于 $U$ 的子问题。展开重构误差$$| X - U W |_F^2 | X |_F^2 - 2 , \mathrm{tr}(W^T U^T X) | W |_F^2$$因为 $U^T U I$ 把交叉项直接简化成 $W^T U^T U W W^T W$。常数项去掉后问题变成$$\max_{U \in \mathrm{St}(D,r)} \mathrm{tr}(U^T X W^T)$$令 $M X W^T$这是一个 $D \times r$ 矩阵。对 $M$ 做紧凑 SVD$$M P \Sigma Q^T$$其中 $P \in \mathbb{R}^{D \times r}$、$Q \in \mathbb{R}^{r \times r}$。于是限制在Stiefel流形上的最优解是$$U^{(k1)} P Q^T$$这叫做正交 Procrustes 解。验证逻辑很直接$\mathrm{tr}(U^T P \Sigma Q^T) \mathrm{tr}(Q^T U^T P \Sigma)$令 $Z U^T P$则 $|Z|_2 \le 1$迹最大化在 $Z Q$ 时取到——前提是 $\Sigma$ 的对角元非负。在格拉斯曼流形上$P$ 和 $P Q^T$ 的列空间完全相同它们代表同一个子空间。所以实际实现时我会直接取 $U P$少一次矩阵乘法还能避开零奇异值带来的符号不唯一问题。但如果你想保持某种规范比如让 $U$ 的某列指向特定方向$P Q^T$ 也有优势它在正交群 $O(r)$ 的表示下是“对齐”的。两种取法在数学上都正确后续有没有数值差异我会在第4节详细说。3.3 完整迭代与复杂度整个算法结构非常简洁每次迭代只有四步输入数据矩阵 X子空间维度 r正则系数 λ初始 U^(0)比如对随机高斯矩阵做 QR 迭代 1. Y : U^(k)T X 2. W^(k1) : SVT_λ(Y) 3. M : X W^(k1)T 4. 对 M 做紧凑 SVDM : P Σ Q^T令 U^(k1) : P或 P Q^T 直到收敛每次迭代的计算量由几个矩阵乘法主导$U^T X$ 是 $O(Dnr)$$X W^T$ 也是 $O(Dnr)$对 $M$ 的紧凑 SVD 是 $O(D r^2)$对 $Y$ 的软阈值 SVD 是 $O(r^2 n)$。总体复杂度大约 $O(Dnr)$当 $r \ll \min(D,n)$ 时这比直接在 $D \times n$ 矩阵上做全量 SVD 便宜太多。更关键的是这个算法不需要学习率。传统流形优化每个迭代步都要选步长、做回退retraction这套流程既慢又敏感。这里两个子问题各自全局最优迭代一次就完成一次“跃迁”而不是滑动一小步。3.4 算法为什么收敛以及为什么不保证全局最优目标函数对单个变量而言都是凸的固定 $U$ 时关于 $W$ 是凸的核范数正则的二阶凸函数固定 $W$ 时关于 $U$ 也是凸的虽然是凸集上的线性函数的最大化。交替最小化的经典性质就是目标值单调不增每次子问题更新都严格下降或持平。加上格拉斯曼流形是紧的所以迭代序列必有收敛子列。但别高兴太早整体目标函数对 $(U, W)$ 联合起来不是凸的所以算法保证的是收敛到稳定的驻点不保证全局最优。实践里初始化影响确实存在一个稳妥做法是拿 PCA 或随机初始化跑 5 到 10 次取目标函数值最小的一次。4. 数值实现里的关键细节谱投影、再正交化、收敛判断4.1 谱投影的数值稳定性为什么我建议取 P 而不是 P Q^T在固定 $W$ 更新 $U$ 时如果 $M X W^T$ 的奇异值出现为 0 的情况比如当前 $W$ 的秩已经小于 $r$那么 $P$ 和 $Q$ 的零奇异值部分是不唯一的。取 $U P Q^T$ 时这种不唯一性会被矩阵乘法放大到所有方向下一轮 $U^T X$ 可能有随机跳变。取 $U P$ 则不一样零奇异值对应的列虽然在分解时也不唯一但它们在子空间表示里根本不参与——列空间由非零奇异值对应的列张开。格拉斯曼流形上的等价关系恰好把这种不唯一性吸收掉了数值上更稳。所以代码里我会写成U, _, _ torch.linalg.svd(M, full_matricesFalse)直接拿左奇异向量矩阵当新的 $U$。另一个隐性坑是浮点误差让 $U^T U$ 渐渐偏离单位阵。如果采用 $U P$ 的更新方式SVD 出来的 $U$ 天然列正交不会漂移。但如果你的实现里加了其他操作比如初始化用 QR 但后续手工修正务必每隔几轮做一次 QR 或 SVD 再正交化否则 $U^T U I$ 的假设崩塌整套推导全废。4.2 零矩阵退化问题实际运行中会遇到一种让人抓狂的情况$\lambda$ 设太大$W$ 直接变成零矩阵。此时 $M X W^T 0$对零矩阵做 SVD$P$ 任意选择一列正交基$U$ 更新完全随机下一轮 $Y U^T X$ 又是满的$W$ 被软阈值压成零……目标值不降但 $U$ 到处乱跳。遇到这种情况要意识到算法在告诉你当前 $\lambda$ 太大了有效维数被压到 0。排查顺序是先查 $W$ 的奇异值谱是不是全为零再看 $U$ 的主角度变化是不是超过 $90^\circ$。不要盲目调迭代次数。4.3 收敛判定必须用主角度别用 F 范数很多人在流形优化里栽跟头都是因为收敛判据选错。直接比较 $|U^{(k1)} - U^{(k)}|_F$ 会犯两类错误一是 $U$ 右乘一个正交矩阵后表示同一个子空间但 F 范数差异很大二是符号翻转在 $r1$ 时特别常见比如 $U$ 和 $-U$ 是同一个一维子空间F 范数却等于 $2\sqrt{D}$。正确做法是用主角度最大值或 Frobenius 范数形式def subspace_distance(U1, U2): s torch.linalg.svdvals(U1.T U2).clamp(min0.0, max1.0) angles torch.arccos(s) return angles.max().item(), torch.norm(angles).item()注意s一定要clamp到 $[0,1]$因为浮点矩阵乘法可能让奇异值略微超过 1直接arccos会得到NaN这个坑我踩过不止一次。4.4 终止条件的选择我通常同时跟踪三个指标$U$ 的最大主角度变化小于 $10^{-4}$ 弧度$W$ 的相对变化 $|W^{(k1)} - W^{(k)}|_F / \max(1, |W^{(k)}|_F)$ 小于 $10^{-5}$目标函数相对下降小于 $10^{-6}$。三个条件不必全部满足取先到者。但从经验看$W$ 的相对变化是最敏感的早期信号$U$ 的主角度是最终是否稳定的硬指标。只看目标函数很容易被骗因为非凸问题里目标函数可能在一个很长的“高原”上缓慢蠕动人眼看不出变化但 $U$ 其实还在慢慢飘。5. 实验怎么设计λ怎么调一组可复现的经验5.1 模拟数据生成与指标设定要说明这套方法确实有效先造一组干净可控的数据。参考设定如下维度 $D 256$样本数 $n 2000$真实子空间维度 $r_{\text{true}} 6$生成随机列正交矩阵 $U_{\text{true}} \in \mathbb{R}^{256 \times 6}$生成 $W_{\text{true}} \in \mathbb{R}^{6 \times 2000}$让它的奇异值按指数衰减前 3 个奇异值大后面迅速变小模拟“有效秩其实只有 2~3”的低秩结构观测 $X U_{\text{true}} W_{\text{true}} \eta E$其中 $E \in \mathbb{R}^{256 \times 2000}$ 是标准高斯噪声$\eta$ 用来控制信噪比常用的取法让信噪比在 5 到 20 之间。评估指标要有两个层次。一个是重构误差 $|X - U W|_F / |X|_F$证明模型拟合数据另一个是子空间估计误差用主角度向量 $(\theta_1, \dots, \theta_6)$ 的 Frobenius 范数来度量这个指标比重构误差更能说明联合优化的价值。5.2 不做联合优化的代价同一份数据上对照组跑三个模型PCA 估计 $U$再固定 $U$ 做最小二乘求 $W$对应“两步走”标准流程$\lambda 0$本算法的 $\lambda 0$ 版本纯粹在格拉斯曼流形上做最小二乘没有低秩正则完整模型$\lambda$ 取适当值。模拟结果对照起来非常直观噪声中等时联合优化比 PCA 两步法的子空间误差低 20%~40%加入核范数后如果真实 $W_{\text{true}}$ 本身谱衰减明显重构误差还能进一步下降。原因也很清楚PCA 只按方差大小选方向不关心这些方向能不能承载低秩系数而联合优化把“方向是否正确”和“系数是否低秩”放在一起秤噪声方向即使方差大如果它不能和低秩系数兼容也会被优化过程慢慢丢掉。5.3 λ的调参策略核范数系数 $\lambda$ 是这个模型里最需要经验的一个量。我试了不同的数据后总结出一套流程先跑一次 $\lambda 0$ 的版本得到 $W_0$记下 $\sigma_{\max}(W_0)$也可以直接取 $U^T X$ 的最大奇异值。这个值大概是 $\lambda$ 的一个自然上界因为当 $\lambda \ge \sigma_{\max}$ 时 $W$ 会被整个压成零。然后在 $[0.01 \sigma_{\max}, 0.2 \sigma_{\max}]$ 范围内按对数网格选候选值。选择标准看两个信号一个是 $W$ 的有效秩奇异值大于 $10^{-3} \sigma_{\max}$ 的数目另一个是验证集重构误差。太小的 $\lambda$ 对谱几乎没有影响有效秩接近 $n$等于没做太大的 $\lambda$ 有效秩骤降重构误差曲线出现明显拐弯或回升。一个我一直推荐的技巧是画奇异值谱图。把 $\lambda 0$ 和候选 $\lambda$ 下的 $W$ 奇异值画在同一张图上理想的 $\lambda$ 会让谱图出现一个明显的“悬崖”——前面几个奇异值还在后面的几乎贴零线。没有悬崖说明数据本身谱衰减不够核范数帮不上什么忙这时候别强行加大 $\lambda$而是要考虑换模型。5.4 真实场景视频帧到底怎么用拿视频动态背景建模举例。假设有 $n$ 帧灰度图每帧 $w \times h$ 像素把每帧拉平成 $D w \times h$ 的列向量堆成矩阵 $X$。背景部分高度线性相关所以它的子空间维度 $r$ 很小$W$ 的奇异值衰减快。前景物体、光照扰动就是噪声项。套用这套模型时$D$ 可能高达几万甚至几十万像素直接对 $M$ 做全量 SVD 内存不够。两种解决办法一是先做一次 PCA 把维度压到比如 100 维再在压缩后的张量上跑流形优化二是用随机 SVD 只求前 $r$ 个奇异值分量。由于 $r$ 通常不到 10随机 SVD 的误差完全可以接受。6. 复现中踩过的坑以及从这个思路延伸出去的方向6.1 软阈值别写错奇异值软阈值最容易犯的错误是把阈值逻辑写成“小于 $\lambda$ 的奇异值直接置零大于 $\lambda$ 的保持不变”。这是硬阈值不是软阈值。核范数近端算子的正确逻辑是$$\sigma_i^{\text{new}} \max(\sigma_i - \lambda, 0)$$每个保留下来的奇异值都要减去 $\lambda$而不是原样保留。我曾经在这上面吃过亏硬阈值版本下目标函数根本不下降因为近端算子选错了目标迭代像无头苍蝇一样乱撞。检查方法很简单看 $W$ 的奇异值是不是精确落在“原奇异值减 $\lambda$”这条曲线上。6.2 别在 $U$ 上做“额外再正交化”如果是按 $U P$ 的方式更新每次 SVD 得到的左奇异向量自动满足列正交不需要再做 Gram–Schmidt。额外再正交化反而可能破坏格拉斯曼流形上的等价关系或引入不必要的数值漂移。如果你是看了其他流形优化的教程习惯在每次梯度更新后都做一次 QR 再正交化那要意识到那是针对 Stiefel 流形上小步梯度更新场景的。本文这种 Procrustes 一步更新自带正交性没必要画蛇添足。6.3 主角度计算的 clamp 细节前面提过acos之前要clamp到 $[0,1]$。这里有一个更隐蔽的问题当 $U_1^T U_2$ 的奇异值精确等于 1 时对应主角度为 0表示两个子空间在某几个方向上完全重合当奇异值等于 0 时主角度为 $\pi/2$表示两个子空间在某几个方向完全正交。但 SVD 每次返回的奇异值排序可能不一样如果子空间距离在迭代后期会震荡要先对奇异值排序再和上一轮比较别被顺序变化迷惑。6.4 后续可以怎么扩展这个模型其实是一个“母体结构”换不同场景能长出很多变体。常见扩展方向有三个一是缺失数据场景。把最小二乘项换成只在观测项上计算 $|P_\Omega(X - UW)|_F^2$此时 $W$ 子问题不再有闭式解要在近端梯度迭代里跑几步内层循环$U$ 子问题仍然是 Procrustes 类问题但不一定有分析解需要用加权版本或迭代 SVD。二是把核范数换成非凸惩罚比如 $\sum_i \sigma_i^p$$0 p 1$。非凸惩罚对奇异值的压低更狠低秩效果通常比核范数更好但代价是近端算子没有闭式解需要数值求解收敛分析也更复杂。三是加图正则化。如果样本之间有关系结构可以在目标函数里加 $\gamma \mathrm{tr}(W L W^T)$$L$ 是拉普拉斯矩阵。$W$ 子问题那一步近端算子依然可以闭式求解近端算子的线性项处理是标准操作整套算法框架一点都不受影响。我在实际项目里最常用的还是缺少数据那一版。视频背景建模时总有部分像素被前景遮挡把遮挡位置直接排除在观测集外剩下的矩阵照样能用这套交替优化框架求解只是 $W$ 的内层循环多花一点算力。最后再分享两个小建议。第一把算法套到自己的数据上前先画一下 $X$ 的奇异值谱如果谱从头到尾都是平的、根本没有衰减趋势那核范数大概帮不上什么忙这个模型的选择要从根上再想想。第二判断收敛或者比较不同初始化结果的时候一定要用主角度这类流形上的几何量别用 $U$ 的逐元素差——我见过太多人在这个细节上绕了远路还以为是算法不收敛。这个方向做下来最深的体会就是把几何结构和谱结构放在同一个模型里并不冲突只要子问题拆得干净数学会告诉你路怎么走。
格拉斯曼流形上的核范数正则化最小二乘优化光看这个标题就知道不是个省油的灯——三个词拆开都认识合在一起就让人有点发怵格拉斯曼流形、核范数、最小二乘优化。我最初接触这篇工作的时候也犯嘀咕流形优化那套黎曼梯度、指数映射、切空间投影已经够折腾了再叠一个核范数这还能收敛吗直到把目标函数拆成两个子问题逐个看才发现事情远比想象中干净一个子问题是Procrustes问题闭式解直接写在脸上另一个子问题是奇异值软阈值一步到位连步长都不需要调。这篇是我们的系列第三篇早期工作技术路线很朴素但正是这种“朴素”让它成了后续很多工作的地基。如果你正在做子空间学习、低秩建模、或者视觉特征处理这篇文章值得你坐下来慢慢读我尽量把里面的几何直觉和实现细节都讲透。1. 这套组合拳要解决的问题子空间估计和低秩系数为什么必须一起优化先摆出最常见的场景。你有一堆观测数据排成矩阵 $X \in \mathbb{R}^{D \times n}$行是特征维度列是样本或者帧。无论是图像序列、基因表达谱、还是用户行为日志这类数据经常存在一个天然的低维结构列向量其实落在某个 $r$ 维子空间附近只是被噪声、缺失值、离群点糊住了。教科书做法是两步走先用主成分分析得到子空间 $U$把 $X$ 投影到 $U$ 上得到系数 $W$再做后续的分类、聚类、回归。听起来顺理成章但问题恰恰出在“先估计子空间”这件事上。PCA估计出的 $\hat{U}$ 本身就带噪声子空间方向只要偏一点点后续所有系数都会继承这个偏差而且是系统性偏差没办法用下游模型救回来。更麻烦的是当你希望系数矩阵额外满足低秩结构时——比如视频帧之间高度线性相关、多视角特征共享同一个低维表示——两步走的方案彻底失效因为你没有任何机制在估计子空间的同时逼着系数低秩。所以这篇工作想的很直接把子空间 $U$ 当作流形上的变量把系数矩阵 $W$ 当作另一个变量两者放进同一个目标函数里联合优化。目标函数长这样$$\min_{U \in \mathrm{Gr}(r,D),; W \in \mathbb{R}^{r \times n}} \frac{1}{2} | X - U W |F^2 \lambda | W |*$$其中 $\mathrm{Gr}(r,D)$ 是格拉斯曼流形$| W |_*$ 是核范数。第一项是重构误差第二项是低秩正则。$U$ 在流形上跑$W$ 在矩阵空间里跑中间靠 $X$ 耦合成一个整体。为什么这个耦合值得做因为子空间估计和低秩系数本质上是同一枚硬币的两面。如果 $U$ 估计得准$WU^T X$ 的奇异值谱会自然衰减核范数稍微推一把就能得到明显的低秩结果反过来如果 $W$ 的低秩约束反馈到 $U$ 的更新里$U$ 会更偏向那些“能承载低秩表达”的方向而不是被噪声方向带跑。这个思路在数学上还有个非常优雅的副产品两个子问题都有闭式解下面会详细展开。但先记住核心结论——联合优化不是锦上添花而是低秩结构存在的场景下两步法在原理上就输了一截。2. 必备数学背景格拉斯曼流形与核范数到底扮演什么角色2.1 格拉斯曼流形子空间本身是一个“点”大部分人对流形的第一印象是曲面上的点在流动但格拉斯曼流形 $\mathrm{Gr}(r,D)$ 的定义有点与众不同它上面的每个“点”是整个 $\mathbb{R}^D$ 里的一个 $r$ 维线性子空间。你没看错点是子空间而不是子空间里的某个向量。这个定义带来的第一个问题是一个子空间怎么用矩阵表示最自然的办法是找它的正交基 $U \in \mathbb{R}^{D \times r}$满足 $U^T U I_r$。但正交基不唯一——右乘任意正交矩阵 $Q \in O(r)$ 后$UQ$ 的列空间和 $U$ 完全一样。所以格拉斯曼流形本质上是Stiefel流形再取商$$\mathrm{Gr}(r,D) \mathrm{St}(D,r) / O(r)$$其中 $\mathrm{St}(D,r)$ 是所有 $D \times r$ 列正交矩阵的集合。商结构意味着优化时不能拿 $U$ 和 $UQ$ 当成两个不同的点它们必须被视为同一个对象。这也是为什么很多人直接用Stiefel流形上的优化器会踩坑自由度冗余会让迭代步态歪歪扭扭收敛变慢。格拉斯曼流形上怎么衡量两个子空间差多远用主角度principal angles。对 $U_1, U_2$先算 $U_1^T U_2$ 的奇异值 $\sigma_i$主角度就是 $\theta_i \arccos(\sigma_i)$。所有主角度合起来就是两个子空间之间的几何距离。这个度量在后面做收敛判定时非常重要千万不要用 $|U_1 - U_2|_F$ 去判断因为 $U_1$ 和 $U_2$ 各自右乘不同的正交矩阵都表示同一个子空间F范数会给出完全错误的“距离”。2.2 核范数向量稀疏的“矩阵版本”核范数 $|W|_* \sum_i \sigma_i(W)$ 是奇异值的 $L_1$ 范数。直接约束 $\mathrm{rank}(W) \le k$ 是非凸的、带组合爆炸的问题而核范数是它在单位球内的凸包络松驰后的优化问题保持凸性又能得到尽量低秩的解。核范数能发挥作用的关键是它的近端算子有闭式解。给定矩阵 $Y$考虑最小化$$\min_Z \frac{1}{2} |Y - Z|F^2 \lambda |Z|*$$答案是奇异值软阈值算子SVT$$\mathrm{SVT}_\lambda(Y) U_Y , \mathrm{diag}\left(\max(\sigma_i - \lambda, 0)\right) V_Y^T$$其中 $Y U_Y \Sigma V_Y^T$ 是奇异值分解。效果很直观小于 $\lambda$ 的奇异值直接清零大于 $\lambda$ 的每个减去 $\lambda$。这就是把 L1 软阈值推广到奇异值谱上——向量里的元素被单独软阈值化矩阵里的奇异值也被单独软阈值化逻辑完全一致。2.3 两者怎么搭起来流形管理的是“子空间的方向”核范数管理的是“子空间系数矩阵的谱”。一个是几何对象一个是谱对象看似八竿子打不着但在我们的目标函数里正好各管一摊$U$ 负责描述数据落在哪个子空间附近$W$ 负责描述在这个子空间里的表达要尽量低秩。$U$ 动不动用流形几何$W$ 动不动用谱分析然后通过最小二乘项耦合。这种分工带来的直接好处是两个变量的更新规则各自独立、各自闭式不用像一般非凸优化那样在切空间上小心翼翼地挪步子。3. 模型与优化U子问题是ProcrustesW子问题是软阈值3.1 固定 U 更新 W投影后再做奇异值软阈值先固定 $U$看关于 $W$ 的子问题。因为 $U^T U I$可以做一个非常舒服的正交分解$$| X - U W |_F^2 | U^T X - W |_F^2 | X |_F^2 - | U^T X |_F^2$$推导只需要两步展开$X - UW$ 先投影到 $U$ 的列空间得到 $U(U^TX - W)$再投影到正交补空间得到 $(I - UU^T)X$两个部分正交所以F范数平方直接相加。第二项和第三项都与 $W$ 无关于是整个子问题等价于$$\min_W \frac{1}{2} | U^T X - W |F^2 \lambda | W |*$$这就是前面说的近端问题答案标准得不能再标准$$W^{(k1)} \mathrm{SVT}_\lambda\left(U^{(k)T} X\right)$$注意这里的 $U^T X$ 天然是一个 $r \times n$ 的小矩阵。$r$ 通常远小于特征维度 $D$所以这个奇异值分解非常便宜每次迭代只求一个 $r \times n$ 矩阵的 SVD而不是在原始 $D \times n$ 尺度上操作。3.2 固定 W 更新 U一个标准的 Procrustes 问题再看固定 $W$ 时分关于 $U$ 的子问题。展开重构误差$$| X - U W |_F^2 | X |_F^2 - 2 , \mathrm{tr}(W^T U^T X) | W |_F^2$$因为 $U^T U I$ 把交叉项直接简化成 $W^T U^T U W W^T W$。常数项去掉后问题变成$$\max_{U \in \mathrm{St}(D,r)} \mathrm{tr}(U^T X W^T)$$令 $M X W^T$这是一个 $D \times r$ 矩阵。对 $M$ 做紧凑 SVD$$M P \Sigma Q^T$$其中 $P \in \mathbb{R}^{D \times r}$、$Q \in \mathbb{R}^{r \times r}$。于是限制在Stiefel流形上的最优解是$$U^{(k1)} P Q^T$$这叫做正交 Procrustes 解。验证逻辑很直接$\mathrm{tr}(U^T P \Sigma Q^T) \mathrm{tr}(Q^T U^T P \Sigma)$令 $Z U^T P$则 $|Z|_2 \le 1$迹最大化在 $Z Q$ 时取到——前提是 $\Sigma$ 的对角元非负。在格拉斯曼流形上$P$ 和 $P Q^T$ 的列空间完全相同它们代表同一个子空间。所以实际实现时我会直接取 $U P$少一次矩阵乘法还能避开零奇异值带来的符号不唯一问题。但如果你想保持某种规范比如让 $U$ 的某列指向特定方向$P Q^T$ 也有优势它在正交群 $O(r)$ 的表示下是“对齐”的。两种取法在数学上都正确后续有没有数值差异我会在第4节详细说。3.3 完整迭代与复杂度整个算法结构非常简洁每次迭代只有四步输入数据矩阵 X子空间维度 r正则系数 λ初始 U^(0)比如对随机高斯矩阵做 QR 迭代 1. Y : U^(k)T X 2. W^(k1) : SVT_λ(Y) 3. M : X W^(k1)T 4. 对 M 做紧凑 SVDM : P Σ Q^T令 U^(k1) : P或 P Q^T 直到收敛每次迭代的计算量由几个矩阵乘法主导$U^T X$ 是 $O(Dnr)$$X W^T$ 也是 $O(Dnr)$对 $M$ 的紧凑 SVD 是 $O(D r^2)$对 $Y$ 的软阈值 SVD 是 $O(r^2 n)$。总体复杂度大约 $O(Dnr)$当 $r \ll \min(D,n)$ 时这比直接在 $D \times n$ 矩阵上做全量 SVD 便宜太多。更关键的是这个算法不需要学习率。传统流形优化每个迭代步都要选步长、做回退retraction这套流程既慢又敏感。这里两个子问题各自全局最优迭代一次就完成一次“跃迁”而不是滑动一小步。3.4 算法为什么收敛以及为什么不保证全局最优目标函数对单个变量而言都是凸的固定 $U$ 时关于 $W$ 是凸的核范数正则的二阶凸函数固定 $W$ 时关于 $U$ 也是凸的虽然是凸集上的线性函数的最大化。交替最小化的经典性质就是目标值单调不增每次子问题更新都严格下降或持平。加上格拉斯曼流形是紧的所以迭代序列必有收敛子列。但别高兴太早整体目标函数对 $(U, W)$ 联合起来不是凸的所以算法保证的是收敛到稳定的驻点不保证全局最优。实践里初始化影响确实存在一个稳妥做法是拿 PCA 或随机初始化跑 5 到 10 次取目标函数值最小的一次。4. 数值实现里的关键细节谱投影、再正交化、收敛判断4.1 谱投影的数值稳定性为什么我建议取 P 而不是 P Q^T在固定 $W$ 更新 $U$ 时如果 $M X W^T$ 的奇异值出现为 0 的情况比如当前 $W$ 的秩已经小于 $r$那么 $P$ 和 $Q$ 的零奇异值部分是不唯一的。取 $U P Q^T$ 时这种不唯一性会被矩阵乘法放大到所有方向下一轮 $U^T X$ 可能有随机跳变。取 $U P$ 则不一样零奇异值对应的列虽然在分解时也不唯一但它们在子空间表示里根本不参与——列空间由非零奇异值对应的列张开。格拉斯曼流形上的等价关系恰好把这种不唯一性吸收掉了数值上更稳。所以代码里我会写成U, _, _ torch.linalg.svd(M, full_matricesFalse)直接拿左奇异向量矩阵当新的 $U$。另一个隐性坑是浮点误差让 $U^T U$ 渐渐偏离单位阵。如果采用 $U P$ 的更新方式SVD 出来的 $U$ 天然列正交不会漂移。但如果你的实现里加了其他操作比如初始化用 QR 但后续手工修正务必每隔几轮做一次 QR 或 SVD 再正交化否则 $U^T U I$ 的假设崩塌整套推导全废。4.2 零矩阵退化问题实际运行中会遇到一种让人抓狂的情况$\lambda$ 设太大$W$ 直接变成零矩阵。此时 $M X W^T 0$对零矩阵做 SVD$P$ 任意选择一列正交基$U$ 更新完全随机下一轮 $Y U^T X$ 又是满的$W$ 被软阈值压成零……目标值不降但 $U$ 到处乱跳。遇到这种情况要意识到算法在告诉你当前 $\lambda$ 太大了有效维数被压到 0。排查顺序是先查 $W$ 的奇异值谱是不是全为零再看 $U$ 的主角度变化是不是超过 $90^\circ$。不要盲目调迭代次数。4.3 收敛判定必须用主角度别用 F 范数很多人在流形优化里栽跟头都是因为收敛判据选错。直接比较 $|U^{(k1)} - U^{(k)}|_F$ 会犯两类错误一是 $U$ 右乘一个正交矩阵后表示同一个子空间但 F 范数差异很大二是符号翻转在 $r1$ 时特别常见比如 $U$ 和 $-U$ 是同一个一维子空间F 范数却等于 $2\sqrt{D}$。正确做法是用主角度最大值或 Frobenius 范数形式def subspace_distance(U1, U2): s torch.linalg.svdvals(U1.T U2).clamp(min0.0, max1.0) angles torch.arccos(s) return angles.max().item(), torch.norm(angles).item()注意s一定要clamp到 $[0,1]$因为浮点矩阵乘法可能让奇异值略微超过 1直接arccos会得到NaN这个坑我踩过不止一次。4.4 终止条件的选择我通常同时跟踪三个指标$U$ 的最大主角度变化小于 $10^{-4}$ 弧度$W$ 的相对变化 $|W^{(k1)} - W^{(k)}|_F / \max(1, |W^{(k)}|_F)$ 小于 $10^{-5}$目标函数相对下降小于 $10^{-6}$。三个条件不必全部满足取先到者。但从经验看$W$ 的相对变化是最敏感的早期信号$U$ 的主角度是最终是否稳定的硬指标。只看目标函数很容易被骗因为非凸问题里目标函数可能在一个很长的“高原”上缓慢蠕动人眼看不出变化但 $U$ 其实还在慢慢飘。5. 实验怎么设计λ怎么调一组可复现的经验5.1 模拟数据生成与指标设定要说明这套方法确实有效先造一组干净可控的数据。参考设定如下维度 $D 256$样本数 $n 2000$真实子空间维度 $r_{\text{true}} 6$生成随机列正交矩阵 $U_{\text{true}} \in \mathbb{R}^{256 \times 6}$生成 $W_{\text{true}} \in \mathbb{R}^{6 \times 2000}$让它的奇异值按指数衰减前 3 个奇异值大后面迅速变小模拟“有效秩其实只有 2~3”的低秩结构观测 $X U_{\text{true}} W_{\text{true}} \eta E$其中 $E \in \mathbb{R}^{256 \times 2000}$ 是标准高斯噪声$\eta$ 用来控制信噪比常用的取法让信噪比在 5 到 20 之间。评估指标要有两个层次。一个是重构误差 $|X - U W|_F / |X|_F$证明模型拟合数据另一个是子空间估计误差用主角度向量 $(\theta_1, \dots, \theta_6)$ 的 Frobenius 范数来度量这个指标比重构误差更能说明联合优化的价值。5.2 不做联合优化的代价同一份数据上对照组跑三个模型PCA 估计 $U$再固定 $U$ 做最小二乘求 $W$对应“两步走”标准流程$\lambda 0$本算法的 $\lambda 0$ 版本纯粹在格拉斯曼流形上做最小二乘没有低秩正则完整模型$\lambda$ 取适当值。模拟结果对照起来非常直观噪声中等时联合优化比 PCA 两步法的子空间误差低 20%~40%加入核范数后如果真实 $W_{\text{true}}$ 本身谱衰减明显重构误差还能进一步下降。原因也很清楚PCA 只按方差大小选方向不关心这些方向能不能承载低秩系数而联合优化把“方向是否正确”和“系数是否低秩”放在一起秤噪声方向即使方差大如果它不能和低秩系数兼容也会被优化过程慢慢丢掉。5.3 λ的调参策略核范数系数 $\lambda$ 是这个模型里最需要经验的一个量。我试了不同的数据后总结出一套流程先跑一次 $\lambda 0$ 的版本得到 $W_0$记下 $\sigma_{\max}(W_0)$也可以直接取 $U^T X$ 的最大奇异值。这个值大概是 $\lambda$ 的一个自然上界因为当 $\lambda \ge \sigma_{\max}$ 时 $W$ 会被整个压成零。然后在 $[0.01 \sigma_{\max}, 0.2 \sigma_{\max}]$ 范围内按对数网格选候选值。选择标准看两个信号一个是 $W$ 的有效秩奇异值大于 $10^{-3} \sigma_{\max}$ 的数目另一个是验证集重构误差。太小的 $\lambda$ 对谱几乎没有影响有效秩接近 $n$等于没做太大的 $\lambda$ 有效秩骤降重构误差曲线出现明显拐弯或回升。一个我一直推荐的技巧是画奇异值谱图。把 $\lambda 0$ 和候选 $\lambda$ 下的 $W$ 奇异值画在同一张图上理想的 $\lambda$ 会让谱图出现一个明显的“悬崖”——前面几个奇异值还在后面的几乎贴零线。没有悬崖说明数据本身谱衰减不够核范数帮不上什么忙这时候别强行加大 $\lambda$而是要考虑换模型。5.4 真实场景视频帧到底怎么用拿视频动态背景建模举例。假设有 $n$ 帧灰度图每帧 $w \times h$ 像素把每帧拉平成 $D w \times h$ 的列向量堆成矩阵 $X$。背景部分高度线性相关所以它的子空间维度 $r$ 很小$W$ 的奇异值衰减快。前景物体、光照扰动就是噪声项。套用这套模型时$D$ 可能高达几万甚至几十万像素直接对 $M$ 做全量 SVD 内存不够。两种解决办法一是先做一次 PCA 把维度压到比如 100 维再在压缩后的张量上跑流形优化二是用随机 SVD 只求前 $r$ 个奇异值分量。由于 $r$ 通常不到 10随机 SVD 的误差完全可以接受。6. 复现中踩过的坑以及从这个思路延伸出去的方向6.1 软阈值别写错奇异值软阈值最容易犯的错误是把阈值逻辑写成“小于 $\lambda$ 的奇异值直接置零大于 $\lambda$ 的保持不变”。这是硬阈值不是软阈值。核范数近端算子的正确逻辑是$$\sigma_i^{\text{new}} \max(\sigma_i - \lambda, 0)$$每个保留下来的奇异值都要减去 $\lambda$而不是原样保留。我曾经在这上面吃过亏硬阈值版本下目标函数根本不下降因为近端算子选错了目标迭代像无头苍蝇一样乱撞。检查方法很简单看 $W$ 的奇异值是不是精确落在“原奇异值减 $\lambda$”这条曲线上。6.2 别在 $U$ 上做“额外再正交化”如果是按 $U P$ 的方式更新每次 SVD 得到的左奇异向量自动满足列正交不需要再做 Gram–Schmidt。额外再正交化反而可能破坏格拉斯曼流形上的等价关系或引入不必要的数值漂移。如果你是看了其他流形优化的教程习惯在每次梯度更新后都做一次 QR 再正交化那要意识到那是针对 Stiefel 流形上小步梯度更新场景的。本文这种 Procrustes 一步更新自带正交性没必要画蛇添足。6.3 主角度计算的 clamp 细节前面提过acos之前要clamp到 $[0,1]$。这里有一个更隐蔽的问题当 $U_1^T U_2$ 的奇异值精确等于 1 时对应主角度为 0表示两个子空间在某几个方向上完全重合当奇异值等于 0 时主角度为 $\pi/2$表示两个子空间在某几个方向完全正交。但 SVD 每次返回的奇异值排序可能不一样如果子空间距离在迭代后期会震荡要先对奇异值排序再和上一轮比较别被顺序变化迷惑。6.4 后续可以怎么扩展这个模型其实是一个“母体结构”换不同场景能长出很多变体。常见扩展方向有三个一是缺失数据场景。把最小二乘项换成只在观测项上计算 $|P_\Omega(X - UW)|_F^2$此时 $W$ 子问题不再有闭式解要在近端梯度迭代里跑几步内层循环$U$ 子问题仍然是 Procrustes 类问题但不一定有分析解需要用加权版本或迭代 SVD。二是把核范数换成非凸惩罚比如 $\sum_i \sigma_i^p$$0 p 1$。非凸惩罚对奇异值的压低更狠低秩效果通常比核范数更好但代价是近端算子没有闭式解需要数值求解收敛分析也更复杂。三是加图正则化。如果样本之间有关系结构可以在目标函数里加 $\gamma \mathrm{tr}(W L W^T)$$L$ 是拉普拉斯矩阵。$W$ 子问题那一步近端算子依然可以闭式求解近端算子的线性项处理是标准操作整套算法框架一点都不受影响。我在实际项目里最常用的还是缺少数据那一版。视频背景建模时总有部分像素被前景遮挡把遮挡位置直接排除在观测集外剩下的矩阵照样能用这套交替优化框架求解只是 $W$ 的内层循环多花一点算力。最后再分享两个小建议。第一把算法套到自己的数据上前先画一下 $X$ 的奇异值谱如果谱从头到尾都是平的、根本没有衰减趋势那核范数大概帮不上什么忙这个模型的选择要从根上再想想。第二判断收敛或者比较不同初始化结果的时候一定要用主角度这类流形上的几何量别用 $U$ 的逐元素差——我见过太多人在这个细节上绕了远路还以为是算法不收敛。这个方向做下来最深的体会就是把几何结构和谱结构放在同一个模型里并不冲突只要子问题拆得干净数学会告诉你路怎么走。
平方复杂度无法被常规优化消除｜TaoToken 统一 Key 通道实测)
