随机化正交投影算法在求解病态线性最小二乘问题中的应用
题目描述:
考虑一个病态(ill-conditioned)的线性最小二乘问题:给定矩阵 \(A \in \mathbb{R}^{m \times n}\)(通常 \(m \gg n\))和向量 \(b \in \mathbb{R}^m\),求 \(x \in \mathbb{R}^n\) 使得残差 \(\|Ax - b\|_2\) 最小化。当 \(A\) 的条件数很大时,传统算法(如正规方程法、QR分解法)会因舍入误差导致数值不稳定。本题目要求讲解一种基于随机化正交投影的算法,其核心思想是:先用一个随机矩阵 \(\Omega \in \mathbb{R}^{n \times d}\)(\(d \ll n\))对 \(A\) 进行投影,得到一个低维子空间,在该子空间中求解一个更小规模且良态的最小二乘问题,从而快速稳定地获得原问题的近似解。
解题过程循序渐进讲解:
步骤1:问题形式化与动机
- 原最小二乘问题为:
\[ \min_{x \in \mathbb{R}^n} \|Ax - b\|_2^2 \]
其理论解为 \(x^* = A^\dagger b\)(\(A^\dagger\) 为 Moore-Penrose 伪逆)。
2. 当 \(A\) 病态时,直接计算 \(A^\dagger\) 或求解正规方程 \(A^T A x = A^T b\) 会放大误差。例如,若 \(A = U \Sigma V^T\) 是奇异值分解(SVD),小奇异值对应的分量在计算中会引入大的扰动。
3. 随机化正交投影算法的目标:通过随机投影将问题压缩到一个低维空间,在压缩后的空间中求解一个条件数更小的近似问题,从而兼顾效率与稳定性。
步骤2:算法核心思想——随机投影与子空间构造
- 选取一个随机矩阵 \(\Omega \in \mathbb{R}^{n \times d}\),其中 \(d\) 是目标子空间维度(\(d \ll n\)),通常 \(d\) 略大于数值秩(numerical rank)以保证精度。\(\Omega\) 的常见选择:
- 高斯随机矩阵(元素独立同分布服从 \(N(0,1)\))。
- 稀疏随机矩阵(如每行只有一个随机 ±1)。
- 随机哈达玛矩阵(快速变换)。
- 计算投影矩阵 \(Y = A \Omega \in \mathbb{R}^{m \times d}\)。\(Y\) 的列张成了 \(A\) 的列空间的一个随机子空间。
- 对 \(Y\) 进行经济型 QR 分解:\(Y = Q R\),其中 \(Q \in \mathbb{R}^{m \times d}\) 列正交(\(Q^T Q = I_d\)),\(R \in \mathbb{R}^{d \times d}\) 上三角。此时 \(Q\) 的列构成了 \(A\) 的列空间的一个近似基。
步骤3:问题压缩与近似解计算
- 将原问题投影到子空间 \(\text{span}(Q)\) 上。因为 \(A x \approx Q Q^T A x\),我们寻找形如 \(x = \Omega z\)(\(z \in \mathbb{R}^d\))的近似解,但更稳健的做法是直接求解压缩后的问题:
\[ \min_{z \in \mathbb{R}^d} \| A Q z - b \|_2^2 \]
这里 \(A Q \in \mathbb{R}^{m \times d}\) 是一个“瘦”矩阵(\(d \ll n\)),其条件数通常远小于 \(A\) 的条件数,因为随机投影保留了 \(A\) 的主要谱结构而抑制了小奇异值的负面影响。
2. 求解压缩问题:计算 \(B = A Q \in \mathbb{R}^{m \times d}\),然后解最小二乘问题 \(\min_z \| B z - b \|_2^2\)。由于 \(d\) 很小,可用稳定方法(如 QR 分解)高效求解:
- 对 \(B\) 做 QR 分解:\(B = U T\)(\(U \in \mathbb{R}^{m \times d}\) 列正交,\(T \in \mathbb{R}^{d \times d}\) 上三角)。
- 计算 \(\tilde{z} = T^{-1} U^T b\)。
- 得到原问题的近似解:\(\tilde{x} = Q \tilde{z}\)。
步骤4:算法流程总结
输入:\(A \in \mathbb{R}^{m \times n}\), \(b \in \mathbb{R}^m\), 目标维度 \(d\)(\(d \ll n\))。
输出:近似解 \(\tilde{x} \in \mathbb{R}^n\)。
- 生成随机矩阵 \(\Omega \in \mathbb{R}^{n \times d}\)(如高斯随机矩阵)。
- 计算 \(Y = A \Omega \in \mathbb{R}^{m \times d}\)。
- 对 \(Y\) 做经济型 QR 分解:\([Q, R] = \text{qr}(Y, 0)\)。
- 计算 \(B = A Q \in \mathbb{R}^{m \times d}\)。
- 对 \(B\) 做经济型 QR 分解:\([U, T] = \text{qr}(B, 0)\)。
- 计算 \(\tilde{z} = T^{-1} (U^T b)\)。
- 返回 \(\tilde{x} = Q \tilde{z}\)。
步骤5:误差与稳定性分析
- 理论上,该算法给出的近似解满足(以高概率):
\[ \| A \tilde{x} - b \|_2^2 \leq \min_{x} \| A x - b \|_2^2 + \epsilon \| b \|_2^2 \]
其中 \(\epsilon\) 与随机矩阵的分布、维度 \(d\) 有关。当 \(d\) 大于数值秩时,误差可控。
2. 稳定性提升原因:
- 随机投影相当于对 \(A\) 的奇异值进行软阈值过滤,减少了小奇异值在求解中的贡献。
- 压缩后的问题维度低,QR 分解等直接法的舍入误差影响小。
- 计算复杂度:
- 矩阵乘法 \(A \Omega\) 成本为 \(O(m n d)\),但若 \(A\) 稀疏可加速。
- QR 分解成本为 \(O(m d^2)\),远低于原问题的 \(O(m n^2)\)。
因此算法适用于大规模病态问题。
步骤6:扩展与变体
- 幂迭代技巧:若 \(A\) 的奇异值衰减缓慢,可先计算 \(Y = (A A^T)^q A \Omega\)(\(q=1,2\)),以增强大奇异值方向的权重,得到更精确的子空间。
- 自适应确定 \(d\):可基于 \(Y\) 的奇异值阈值自动选择 \(d\),避免手动设定。
- 与正则化结合:在压缩问题中可加入 Tikhonov 正则化,进一步抑制噪声。
总结:随机化正交投影算法通过随机压缩将原高维病态问题转化为低维良态问题,在保证精度的前提下显著提升计算效率和数值稳定性,特别适合大规模病态最小二乘问题。