随机化Nesterov加速梯度下降在求解大规模稀疏线性方程组中的应用
我将为你讲解这个算法。首先，我们来理解问题的背景和描述。
问题描述
考虑求解大规模稀疏线性方程组：
\[ Ax = b \]
其中 \(A \in \mathbb{R}^{n \times n}\) 是一个大规模、稀疏、对称正定矩阵，\(b \in \mathbb{R}^n\) 是已知向量，需要求解 \(x \in \mathbb{R}^n\)。
- 大规模：\(n\) 非常大（例如百万甚至十亿级），使得直接法（如LU分解）在存储和计算上不可行。
- 稀疏：矩阵 \(A\) 绝大部分元素为零，这允许我们使用特殊的数据结构（如CSR、CSC格式）高效存储和进行矩阵-向量乘法。
- 对称正定：所有特征值为正，这保证了优化问题的强凸性，从而可以使用梯度类方法高效求解。
由于问题规模巨大，我们通常采用迭代法。将线性方程组求解转化为等价的凸优化问题：
\[ \min_{x \in \mathbb{R}^n} f(x) = \frac{1}{2} x^T A x - b^T x \]
其梯度为：
\[ \nabla f(x) = A x - b \]
梯度为零的点即为线性方程组的解。因此，我们可以用梯度下降法求解该优化问题。
但标准梯度下降法在病态条件数（即最大特征值与最小特征值之比 \(\kappa(A) = \frac{\lambda_{\max}}{\lambda_{\min}}\) 很大）时收敛很慢，收敛速度与条件数 \(\kappa(A)\) 相关，迭代复杂度为 \(O(\kappa \log(1/\varepsilon))\) 达到 \(\varepsilon\) 精度。
为了解决这个问题，随机化Nesterov加速梯度下降结合了：
- Nesterov加速技巧：通过在迭代中引入“动量项”，将梯度下降的收敛速度从 \(O(1/k)\) 提升到 \(O(1/k^2)\)（对于光滑强凸问题）。
- 随机化技术：每次迭代只使用部分梯度信息（例如通过随机坐标下降或随机梯度估计），大幅降低单次迭代的计算成本，特别适合大规模稀疏问题。
循序渐进讲解
第一步：回顾标准梯度下降法
对于问题 \(\min f(x)\)，梯度下降的迭代格式为：
\[ x_{k+1} = x_k - \alpha_k \nabla f(x_k) \]
其中 \(\alpha_k > 0\) 是步长（学习率）。
对于二次函数 \(f(x) = \frac{1}{2} x^T A x - b^T x\)，最优固定步长为 \(\alpha = \frac{2}{\lambda_{\min} + \lambda_{\max}}\)，收敛速度线性，收敛率 \(\rho = \frac{\kappa - 1}{\kappa + 1} \approx 1 - \frac{2}{\kappa}\)，需要 \(O(\kappa \log(1/\varepsilon))\) 次迭代。
第二步：引入Nesterov加速梯度法（确定型）
Nesterov加速梯度法（简称NAG）通过引入辅助序列 \(\{y_k\}\) 和动量系数，实现更快的收敛。对于强凸且梯度Lipschitz连续的函数，NAG的迭代格式为：
初始化 \(x_0 = y_0\)，令 \(q = \frac{\lambda_{\min}}{\lambda_{\max}} = \frac{1}{\kappa}\)。
对于 \(k = 0, 1, 2, \dots\)：
\[\begin{aligned} x_{k+1} &= y_k - \frac{1}{L} \nabla f(y_k) \\ y_{k+1} &= x_{k+1} + \frac{\sqrt{L} - \sqrt{q}}{\sqrt{L} + \sqrt{q}} (x_{k+1} - x_k) \end{aligned} \]
其中 \(L = \lambda_{\max}\) 是梯度的Lipschitz常数，\(q = \lambda_{\min}\) 是强凸系数。
该方法的收敛速度达到 \(O\left( \left(1 - \sqrt{\frac{1}{\kappa}}\right)^k \right)\)，即 \(O\left( \exp\left(-k / \sqrt{\kappa}\right) \right)\)，迭代复杂度为 \(O(\sqrt{\kappa} \log(1/\varepsilon))\)，比标准梯度下降的 \(O(\kappa \log(1/\varepsilon))\) 快得多。
第三步：结合随机化——随机坐标下降思想
在每次迭代中，计算全梯度 \(\nabla f(y_k) = A y_k - b\) 需要一次矩阵-向量乘法，复杂度 \(O(\text{nnz})\)，其中 nnz 是矩阵 \(A\) 的非零元个数。当 \(A\) 非常大时，即使单次矩阵-向量乘法也可能很昂贵。
随机坐标下降 每次只更新一个坐标，只计算梯度的一个分量。对于 \(f(x) = \frac{1}{2} x^T A x - b^T x\)，梯度的第 \(i\) 个分量为：
\[[\nabla f(x)]_i = (A x)_i - b_i = \sum_{j=1}^n A_{ij} x_j - b_i \]
如果 \(A\) 稀疏，且我们只更新第 \(i\) 个坐标，则计算这个分量只需要访问矩阵 \(A\) 的第 \(i\) 行，复杂度为 \(O(\text{该行的非零元个数})\)，通常远低于 \(O(\text{nnz})\)。
关键问题：如何将坐标下降与Nesterov加速结合？
第四步：随机化Nesterov加速梯度下降算法
一种有效方法是将Nesterov加速框架中的全梯度替换为随机梯度估计，并调整步长和动量参数。
算法步骤（以随机坐标下降为例）：
- 输入：对称正定矩阵 \(A\)，向量 \(b\)，初始猜测 \(x_0\)，最大迭代次数 \(K\)，目标精度 \(\varepsilon\)。
- 参数设置：
- 令 \(L_{\max} = \lambda_{\max}(A)\)（或上界），\(L_{\min} = \lambda_{\min}(A)\)（或下界）。
- 设 \(q = L_{\min} / n\)（这里除 \(n\) 是因为随机坐标下降的预期步长缩放，与随机梯度方差有关）。
- 令 \(L = L_{\max}\)。
- 设置动量参数 \(\beta = \frac{\sqrt{L} - \sqrt{q}}{\sqrt{L} + \sqrt{q}}\)。
- 初始化：\(x_0 = y_0\)。
- 迭代：对于 \(k = 0, 1, \dots, K-1\)：
a. 随机均匀选取坐标索引 \(i_k \in \{1, 2, \dots, n\}\)。
b. 计算随机梯度估计 \(g_k\)，它是全梯度 \(\nabla f(y_k)\) 在第 \(i_k\) 坐标方向上的近似：
\[ g_k = n \cdot ( (A y_k)_{i_k} - b_{i_k} ) \cdot e_{i_k} \]
其中 $e_{i_k}$ 是第 $i_k$ 个标准基向量。乘以 $n$ 是为了保证 $\mathbb{E}[g_k] = \nabla f(y_k)$（无偏估计）。
c. 更新 \(x_{k+1} = y_k - \frac{1}{L} g_k\)。
d. 更新 \(y_{k+1} = x_{k+1} + \beta (x_{k+1} - x_k)\)。
5. 输出：近似解 \(x_K\)。
注意：在实现中，计算 \((A y_k)_{i_k}\) 只需要矩阵 \(A\) 的第 \(i_k\) 行与向量 \(y_k\) 的内积，复杂度很低。由于每次只更新一个坐标，向量 \(x_k, y_k\) 的更新也很廉价。
第五步：收敛性分析
对于强凸且梯度Lipschitz连续的函数，随机化Nesterov加速梯度下降（上述变体）的期望收敛速度为：
\[\mathbb{E}[f(x_K) - f^*] \leq O\left( \left(1 - \sqrt{\frac{q}{L}}\right)^K \right) \]
其中 \(f^*\) 是最优值，\(q \approx \frac{\lambda_{\min}}{n}\)，\(L = \lambda_{\max}\)。因此，条件数变为 \(\kappa_{\text{eff}} = \frac{L}{q} \approx n \kappa\)，迭代复杂度为 \(O(\sqrt{n \kappa} \log(1/\varepsilon))\)。
虽然有效条件数增加了 \(n\) 倍，但每次迭代的成本从 \(O(\text{nnz})\) 降到了 \(O(\text{平均每行非零元个数})\)。对于非常稀疏的矩阵，这通常意味着单次迭代成本降低为原来的 \(1/n\) 量级，总体加速明显。
第六步：算法扩展与实用技巧
- 随机化策略：除了随机坐标，还可以采用随机小批量（一组坐标）或随机行采样（对最小二乘问题），平衡方差与单次迭代成本。
- 自适应步长：实际中 \(L\) 和 \(q\) 可能未知，可以使用回溯线搜索或自适应估计 Lipschitz 常数。
- 方差缩减技术：为了进一步加速，可以结合SVRG（随机方差缩减梯度）或SAGA等方法，在随机梯度中引入控制变量以减少方差，从而改善常数因子。
- 预热启动：可以先运行少量标准梯度下降或随机梯度下降，得到一个较好的初始点，再启动加速阶段。
第七步：总结与适用场景
- 适用：大规模稀疏对称正定线性方程组，尤其是矩阵-向量乘法昂贵但单行访问廉价的问题。
- 优势：结合了Nesterov加速的快速收敛和随机化的低迭代成本，总体时间通常远少于标准梯度下降和确定型Nesterov加速。
- 局限：需要对称正定性，对非对称或不定矩阵不直接适用；随机化导致收敛路径有噪声，可能需要适当调整参数。
这个算法体现了现代大规模优化中将经典加速技巧与随机化采样相结合的核心思想，是求解超大规模稀疏线性系统的重要工具之一。