随机化Kaczmarz算法在求解欠定线性方程组中的扩展与应用

字数 4997 2025-12-16 02:51:02

随机化Kaczmarz算法在求解欠定线性方程组中的扩展与应用

题目描述
我们考虑一个线性方程组 \(A\mathbf{x} = \mathbf{b}\)，其中 \(A \in \mathbb{R}^{m \times n}\)， \(\mathbf{b} \in \mathbb{R}^{m}\)，且 \(m < n\)，即方程组是欠定的（未知数个数多于方程个数）。通常情况下，欠定方程组有无穷多解。我们的目标是找到其中一个解，特别是希望找到具有某些期望性质（如较小范数）的解。经典Kaczmarz迭代算法通常用于超定或方阵系统，它通过逐行投影来逼近解。但对于欠定系统，标准Kaczmarz方法可能会收敛缓慢或不收敛到有意义的解。本题将探讨如何将随机化Kaczmarz算法（Randomized Kaczmarz, RK）扩展应用于欠定系统，以高效地找到一个解，并特别关注其收敛到最小范数解的性质。

解题过程循序渐进讲解

步骤1: 理解欠定系统与最小范数解
对于欠定线性方程组 \(A\mathbf{x} = \mathbf{b}\)，假设解存在（即 \(\mathbf{b} \in \text{range}(A)\)），其解集构成一个仿射空间。在所有解中，通常对最小欧几里得范数解（即 \(\min \|\mathbf{x}\|_2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\)）感兴趣，因为它通常对应于最“简单”或最稳定的解。这个最小范数解可以通过伪逆给出：\(\mathbf{x}^* = A^\dagger \mathbf{b}\)，其中 \(A^\dagger\) 是 \(A\) 的Moore-Penrose伪逆。

步骤2: 回顾标准随机化Kaczmarz算法（用于超定/方阵系统）
标准RK算法适用于 \(A\mathbf{x} = \mathbf{b}\)（通常 \(m \ge n\)），迭代格式为：

随机选取第 \(i\) 行，选取概率与行范数平方 \(\|A_i\|_2^2\) 成正比（或其他概率分布）。
更新解估计：\(\mathbf{x}_{k+1} = \mathbf{x}_k + \frac{b_i - A_i^T \mathbf{x}_k}{\|A_i\|_2^2} A_i\)，其中 \(A_i\) 是 \(A\) 的第 \(i\) 行（视为列向量）。
这个更新是当前估计向超平面 \(A_i^T \mathbf{x} = b_i\) 的正交投影。算法在一致相容系统中收敛到唯一解。

步骤3: 欠定系统的直接应用问题
如果直接将标准RK应用于欠定系统（\(m < n\)），迭代会在各行对应的超平面上投影，但由于解空间是无限维流形，迭代点可能会在解空间中“游荡”，而不收敛到一个特定点。实际上，如果初始猜测 \(\mathbf{x}_0 \in \text{row space}(A)\)，标准RK会收敛到解，但可能不是最小范数解。更常见的是，我们希望对欠定系统找到一个特定解，如最小范数解。

步骤4: 应用于增广系统——转化为方阵系统
一个常见技巧是将欠定系统转化为一个方阵系统，使RK可应用。考虑增广系统：

\[\begin{bmatrix} I & A^T \\ A & 0 \end{bmatrix} \begin{bmatrix} \mathbf{y} \\ \mathbf{x} \end{bmatrix} = \begin{bmatrix} \mathbf{0} \\ \mathbf{b} \end{bmatrix}. \]

这个系统的解满足 \(A\mathbf{x} = \mathbf{b}\) 且 \(\mathbf{y} = -A^T \mathbf{x}\)。但更直接的方法是利用对偶性：最小范数解 \(\mathbf{x}^*\) 位于 \(A\) 的行空间，即 \(\mathbf{x}^* = A^T \mathbf{z}\) 对某个 \(\mathbf{z} \in \mathbb{R}^m\) 成立。代入原系统得 \(A A^T \mathbf{z} = \mathbf{b}\)。这是一个 \(m \times m\) 的方阵系统（假设 \(A\) 行满秩，即秩为 \(m\)），可以对其应用标准RK。

步骤5: 针对 \(AA^T \mathbf{z} = \mathbf{b}\) 应用随机化Kaczmarz
系统 \(AA^T \mathbf{z} = \mathbf{b}\) 是方阵且对称正定（若 \(A\) 行满秩）。将 \(AA^T\) 视为矩阵，其行是 \(A\) 的行的线性组合。但更有效的是：注意到 \(AA^T\) 的第 \(i\) 行是 \(A_i^T A\)，其中 \(A_i\) 是 \(A\) 的第 \(i\) 行。对系统 \(AA^T \mathbf{z} = \mathbf{b}\) 应用RK的投影步为：
选取行索引 \(i\)，更新 \(\mathbf{z}_{k+1} = \mathbf{z}_k + \frac{b_i - (AA^T)_i^T \mathbf{z}_k}{\|(AA^T)_i\|_2^2} (AA^T)_i\)。
计算 \((AA^T)_i = A_i^T A\) 和其范数较繁琐。但我们可以简化：观察到对偶系统与原系统的联系。

步骤6: 随机化Kaczmarz对偶形式（RK Dual）
定义对偶变量 \(\mathbf{z} \in \mathbb{R}^m\)，并令原变量 \(\mathbf{x} = A^T \mathbf{z}\)。则 \(A\mathbf{x} = AA^T \mathbf{z} = \mathbf{b}\)。对 \(AA^T \mathbf{z} = \mathbf{b}\) 应用RK，选取行 \(i\) 的概率与 \(\|(AA^T)_i\|_2^2\) 成正比。但我们可以导出直接在原空间更新的等价形式：
RK更新为：\(\mathbf{z}_{k+1} = \mathbf{z}_k + \frac{b_i - A_i^T A \mathbf{z}_k}{\|A_i\|_2^2 \|A^T\|_F^2?} (A_i^T A)^T\)？这里需要仔细推导。

实际上，更常用的方法是采用随机化坐标下降（Randomized Coordinate Descent, RCD）视角。问题 \(\min_{\mathbf{x}} \frac{1}{2} \|\mathbf{x}\|_2^2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\) 的拉格朗日函数为 \(L(\mathbf{x}, \mathbf{z}) = \frac{1}{2} \|\mathbf{x}\|_2^2 + \mathbf{z}^T(\mathbf{b} - A\mathbf{x})\)。最优性条件为 \(\mathbf{x} = A^T \mathbf{z}\) 和 \(A\mathbf{x} = \mathbf{b}\)。将其结合为 \(AA^T \mathbf{z} = \mathbf{b}\)。对该系统应用随机化坐标下降（即对 \(\mathbf{z}\) 的每个分量更新）等价于对原问题应用随机化Kaczmarz的一种形式。

步骤7: 实际算法——两阶段更新
一个实用算法是交替更新原始变量和对偶变量：

初始化 \(\mathbf{z}_0 = \mathbf{0}\)，则 \(\mathbf{x}_0 = A^T \mathbf{z}_0 = \mathbf{0}\)。
在第 \(k\) 次迭代，随机选取行索引 \(i \in \{1,\dots, m\}\)，概率与 \(\|A_i\|_2^2\) 成正比（或均匀选取）。
更新对偶变量：\(z_i^{(k+1)} = z_i^{(k)} + \frac{b_i - A_i^T \mathbf{x}_k}{\|A_i\|_2^2}\)，其余 \(z_j^{(k+1)} = z_j^{(k)}\) 对 \(j \neq i\)。
更新原始变量：\(\mathbf{x}_{k+1} = A^T \mathbf{z}_{k+1}\)。
注意：由于 \(\mathbf{x}_{k+1} = A^T \mathbf{z}_{k+1} = \mathbf{x}_k + (z_i^{(k+1)} - z_i^{(k)}) A_i = \mathbf{x}_k + \frac{b_i - A_i^T \mathbf{x}_k}{\|A_i\|_2^2} A_i\)。这正是标准RK更新！因此，如果我们从 \(\mathbf{x}_0 = \mathbf{0}\) 开始，并保持 \(\mathbf{x}_k = A^T \mathbf{z}_k\)，那么对偶变量的坐标下降等价于原始变量的标准RK迭代。但关键点是：由于初始 \(\mathbf{x}_0 = \mathbf{0} \in \text{row space}(A)\)，且每次更新都在行空间中（因为加上 \(A_i\) 的倍数），因此所有迭代点 \(\mathbf{x}_k\) 都保持在行空间中。对于欠定系统，行空间中的解正是最小范数解（因为解空间是行空间与零空间直和，最小范数解是行空间分量）。

步骤8: 收敛性与解释
上述算法与标准RK形式完全相同，只是从零初始点开始。收敛性分析：对于欠定一致系统（解存在），从 \(\mathbf{x}_0 = \mathbf{0}\) 开始的标准RK迭代收敛到最小范数解 \(\mathbf{x}^* = A^\dagger \mathbf{b}\)。收敛速度取决于矩阵的条件数（在行空间上）。随机选取行可加速收敛。该算法可视为求解 \(\min \|\mathbf{x}\|_2^2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\) 的随机坐标上升法应用于对偶问题。

步骤9: 算法步骤总结

输入：\(A \in \mathbb{R}^{m \times n} (m < n)\), \(\mathbf{b} \in \mathbb{R}^m\)，确保解存在（否则算法会收敛到最小二乘解？这里假设一致）。
初始化：\(\mathbf{x}_0 = \mathbf{0} \in \mathbb{R}^n\)。
预计算行范数平方 \(w_i = \|A_i\|_2^2\) 对 \(i=1,\dots, m\)。
对于迭代 \(k=0,1,\dots\) 直到收敛：
a. 随机选取行索引 \(i\)，概率与 \(w_i\) 成正比（或使用均匀分布等其他分布）。
b. 计算 \(\alpha = \frac{b_i - A_i^T \mathbf{x}_k}{w_i}\)。
c. 更新 \(\mathbf{x}_{k+1} = \mathbf{x}_k + \alpha A_i\)。
输出：近似解 \(\mathbf{x}_k\)，它近似于 \(A\mathbf{x} = \mathbf{b}\) 的最小范数解。

步骤10: 扩展与注意事项

如果系统不一致（无解），该算法收敛到最小范数最小二乘解，但需要修改步长或使用其他变体。
加速技巧：使用随机化松弛（如随机化平均块Kaczmarz）或自适应采样（如基于残差的概率分布）。
与正则化结合：在噪声数据中，可引入Tikhonov正则化，相应修改步长公式。
实现时，注意数值稳定性，特别是当行范数很小时需处理除零问题。
该算法提供了一种简单迭代方法，无需矩阵求逆或分解，即可求解大规模欠定系统的最小范数解。

随机化Kaczmarz算法在求解欠定线性方程组中的扩展与应用题目描述我们考虑一个线性方程组 \(A\mathbf{x} = \mathbf{b}\)，其中 \(A \in \mathbb{R}^{m \times n}\)， \(\mathbf{b} \in \mathbb{R}^{m}\)，且 \(m < n\)，即方程组是欠定的（未知数个数多于方程个数）。通常情况下，欠定方程组有无穷多解。我们的目标是找到其中一个解，特别是希望找到具有某些期望性质（如较小范数）的解。经典Kaczmarz迭代算法通常用于超定或方阵系统，它通过逐行投影来逼近解。但对于欠定系统，标准Kaczmarz方法可能会收敛缓慢或不收敛到有意义的解。本题将探讨如何将随机化Kaczmarz算法（Randomized Kaczmarz, RK）扩展应用于欠定系统，以高效地找到一个解，并特别关注其收敛到最小范数解的性质。解题过程循序渐进讲解步骤1: 理解欠定系统与最小范数解对于欠定线性方程组 \(A\mathbf{x} = \mathbf{b}\)，假设解存在（即 \(\mathbf{b} \in \text{range}(A)\)），其解集构成一个仿射空间。在所有解中，通常对最小欧几里得范数解（即 \(\min \|\mathbf{x}\|_ 2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\)）感兴趣，因为它通常对应于最“简单”或最稳定的解。这个最小范数解可以通过伪逆给出：\(\mathbf{x}^* = A^\dagger \mathbf{b}\)，其中 \(A^\dagger\) 是 \(A\) 的Moore-Penrose伪逆。步骤2: 回顾标准随机化Kaczmarz算法（用于超定/方阵系统）标准RK算法适用于 \(A\mathbf{x} = \mathbf{b}\)（通常 \(m \ge n\)），迭代格式为：随机选取第 \(i\) 行，选取概率与行范数平方 \(\|A_ i\|_ 2^2\) 成正比（或其他概率分布）。更新解估计：\(\mathbf{x}_ {k+1} = \mathbf{x}_ k + \frac{b_ i - A_ i^T \mathbf{x}_ k}{\|A_ i\|_ 2^2} A_ i\)，其中 \(A_ i\) 是 \(A\) 的第 \(i\) 行（视为列向量）。这个更新是当前估计向超平面 \(A_ i^T \mathbf{x} = b_ i\) 的正交投影。算法在一致相容系统中收敛到唯一解。步骤3: 欠定系统的直接应用问题如果直接将标准RK应用于欠定系统（\(m < n\)），迭代会在各行对应的超平面上投影，但由于解空间是无限维流形，迭代点可能会在解空间中“游荡”，而不收敛到一个特定点。实际上，如果初始猜测 \(\mathbf{x}_ 0 \in \text{row space}(A)\)，标准RK会收敛到解，但可能不是最小范数解。更常见的是，我们希望对欠定系统找到一个特定解，如最小范数解。步骤4: 应用于增广系统——转化为方阵系统一个常见技巧是将欠定系统转化为一个方阵系统，使RK可应用。考虑增广系统： \[ \begin{bmatrix} I & A^T \\ A & 0 \end{bmatrix} \begin{bmatrix} \mathbf{y} \\ \mathbf{x} \end{bmatrix} = \begin{bmatrix} \mathbf{0} \\ \mathbf{b} \end{bmatrix}. \] 这个系统的解满足 \(A\mathbf{x} = \mathbf{b}\) 且 \(\mathbf{y} = -A^T \mathbf{x}\)。但更直接的方法是利用对偶性：最小范数解 \(\mathbf{x}^ \) 位于 \(A\) 的行空间，即 \(\mathbf{x}^ = A^T \mathbf{z}\) 对某个 \(\mathbf{z} \in \mathbb{R}^m\) 成立。代入原系统得 \(A A^T \mathbf{z} = \mathbf{b}\)。这是一个 \(m \times m\) 的方阵系统（假设 \(A\) 行满秩，即秩为 \(m\)），可以对其应用标准RK。步骤5: 针对 \(AA^T \mathbf{z} = \mathbf{b}\) 应用随机化Kaczmarz 系统 \(AA^T \mathbf{z} = \mathbf{b}\) 是方阵且对称正定（若 \(A\) 行满秩）。将 \(AA^T\) 视为矩阵，其行是 \(A\) 的行的线性组合。但更有效的是：注意到 \(AA^T\) 的第 \(i\) 行是 \(A_ i^T A\)，其中 \(A_ i\) 是 \(A\) 的第 \(i\) 行。对系统 \(AA^T \mathbf{z} = \mathbf{b}\) 应用RK的投影步为：选取行索引 \(i\)，更新 \(\mathbf{z}_ {k+1} = \mathbf{z}_ k + \frac{b_ i - (AA^T)_ i^T \mathbf{z}_ k}{\|(AA^T)_ i\|_ 2^2} (AA^T)_ i\)。计算 \((AA^T)_ i = A_ i^T A\) 和其范数较繁琐。但我们可以简化：观察到对偶系统与原系统的联系。步骤6: 随机化Kaczmarz对偶形式（RK Dual）定义对偶变量 \(\mathbf{z} \in \mathbb{R}^m\)，并令原变量 \(\mathbf{x} = A^T \mathbf{z}\)。则 \(A\mathbf{x} = AA^T \mathbf{z} = \mathbf{b}\)。对 \(AA^T \mathbf{z} = \mathbf{b}\) 应用RK，选取行 \(i\) 的概率与 \(\|(AA^T)_ i\| 2^2\) 成正比。但我们可以导出直接在原空间更新的等价形式： RK更新为：\(\mathbf{z} {k+1} = \mathbf{z}_ k + \frac{b_ i - A_ i^T A \mathbf{z}_ k}{\|A_ i\|_ 2^2 \|A^T\|_ F^2?} (A_ i^T A)^T\)？这里需要仔细推导。实际上，更常用的方法是采用随机化坐标下降（Randomized Coordinate Descent, RCD）视角。问题 \(\min_ {\mathbf{x}} \frac{1}{2} \|\mathbf{x}\|_ 2^2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\) 的拉格朗日函数为 \(L(\mathbf{x}, \mathbf{z}) = \frac{1}{2} \|\mathbf{x}\|_ 2^2 + \mathbf{z}^T(\mathbf{b} - A\mathbf{x})\)。最优性条件为 \(\mathbf{x} = A^T \mathbf{z}\) 和 \(A\mathbf{x} = \mathbf{b}\)。将其结合为 \(AA^T \mathbf{z} = \mathbf{b}\)。对该系统应用随机化坐标下降（即对 \(\mathbf{z}\) 的每个分量更新）等价于对原问题应用随机化Kaczmarz的一种形式。步骤7: 实际算法——两阶段更新一个实用算法是交替更新原始变量和对偶变量：初始化 \(\mathbf{z}_ 0 = \mathbf{0}\)，则 \(\mathbf{x}_ 0 = A^T \mathbf{z}_ 0 = \mathbf{0}\)。在第 \(k\) 次迭代，随机选取行索引 \(i \in \{1,\dots, m\}\)，概率与 \(\|A_ i\|_ 2^2\) 成正比（或均匀选取）。更新对偶变量：\(z_ i^{(k+1)} = z_ i^{(k)} + \frac{b_ i - A_ i^T \mathbf{x}_ k}{\|A_ i\|_ 2^2}\)，其余 \(z_ j^{(k+1)} = z_ j^{(k)}\) 对 \(j \neq i\)。更新原始变量：\(\mathbf{x} {k+1} = A^T \mathbf{z} {k+1}\)。注意：由于 \(\mathbf{x} {k+1} = A^T \mathbf{z} {k+1} = \mathbf{x}_ k + (z_ i^{(k+1)} - z_ i^{(k)}) A_ i = \mathbf{x}_ k + \frac{b_ i - A_ i^T \mathbf{x}_ k}{\|A_ i\|_ 2^2} A_ i\)。这正是标准RK更新！因此，如果我们从 \(\mathbf{x}_ 0 = \mathbf{0}\) 开始，并保持 \(\mathbf{x}_ k = A^T \mathbf{z}_ k\)，那么对偶变量的坐标下降等价于原始变量的标准RK迭代。但关键点是：由于初始 \(\mathbf{x}_ 0 = \mathbf{0} \in \text{row space}(A)\)，且每次更新都在行空间中（因为加上 \(A_ i\) 的倍数），因此所有迭代点 \(\mathbf{x}_ k\) 都保持在行空间中。对于欠定系统，行空间中的解正是最小范数解（因为解空间是行空间与零空间直和，最小范数解是行空间分量）。步骤8: 收敛性与解释上述算法与标准RK形式完全相同，只是从零初始点开始。收敛性分析：对于欠定一致系统（解存在），从 \(\mathbf{x}_ 0 = \mathbf{0}\) 开始的标准RK迭代收敛到最小范数解 \(\mathbf{x}^* = A^\dagger \mathbf{b}\)。收敛速度取决于矩阵的条件数（在行空间上）。随机选取行可加速收敛。该算法可视为求解 \(\min \|\mathbf{x}\|_ 2^2 \text{ s.t. } A\mathbf{x} = \mathbf{b}\) 的随机坐标上升法应用于对偶问题。步骤9: 算法步骤总结输入：\(A \in \mathbb{R}^{m \times n} (m < n)\), \(\mathbf{b} \in \mathbb{R}^m\)，确保解存在（否则算法会收敛到最小二乘解？这里假设一致）。初始化：\(\mathbf{x}_ 0 = \mathbf{0} \in \mathbb{R}^n\)。预计算行范数平方 \(w_ i = \|A_ i\|_ 2^2\) 对 \(i=1,\dots, m\)。对于迭代 \(k=0,1,\dots\) 直到收敛： a. 随机选取行索引 \(i\)，概率与 \(w_ i\) 成正比（或使用均匀分布等其他分布）。 b. 计算 \( \alpha = \frac{b_ i - A_ i^T \mathbf{x} k}{w_ i} \)。 c. 更新 \(\mathbf{x} {k+1} = \mathbf{x}_ k + \alpha A_ i\)。输出：近似解 \(\mathbf{x}_ k\)，它近似于 \(A\mathbf{x} = \mathbf{b}\) 的最小范数解。步骤10: 扩展与注意事项如果系统不一致（无解），该算法收敛到最小范数最小二乘解，但需要修改步长或使用其他变体。加速技巧：使用随机化松弛（如随机化平均块Kaczmarz）或自适应采样（如基于残差的概率分布）。与正则化结合：在噪声数据中，可引入Tikhonov正则化，相应修改步长公式。实现时，注意数值稳定性，特别是当行范数很小时需处理除零问题。该算法提供了一种简单迭代方法，无需矩阵求逆或分解，即可求解大规模欠定系统的最小范数解。