SVD、LU 和 QR 等矩阵分解背后的理论很重要,但手动进行这些计算对于任何大尺寸的矩阵都不切实际。幸运的是,Python 的科学计算库,特别是 NumPy 和 SciPy,提供了这些算法的高效且数值稳定的实现。我们将主要使用 scipy.linalg 中的函数,该模块通常比 numpy.linalg 提供更多高级功能和优化,尽管 NumPy 也提供了基础版本。请确保您已安装 NumPy 和 SciPy 并导入它们:import numpy as np from scipy import linalg奇异值分解 (SVD)SVD 将矩阵 $A$ 分解为 $U \Sigma V^T$。scipy.linalg.svd 函数处理此分解。# 定义一个示例矩阵 A A = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]]) print("原始矩阵 A ({}):".format(A.shape)) print(A) # 执行 SVD # full_matrices=False 计算精简 SVD(对于非方阵更高效) U, s, Vh = linalg.svd(A, full_matrices=False) print("\nU 矩阵 ({}):".format(U.shape)) print(U) # s 包含奇异值(一维数组) print("\n奇异值 (s):", s.shape) print(s) # Vh 是 V 的转置(V.T) print("\nVh 矩阵 (V 转置) ({}):".format(Vh.shape)) print(Vh)理解输出:U:左奇异向量矩阵。其形状取决于 full_matrices。对于 $m \times n$ 矩阵且 $m \ge n$,当 full_matrices=False 时,$U$ 的形状将是 $m \times n$。s:一个包含按降序排列的奇异值($\sigma_1, \sigma_2, ..., \sigma_k$)的一维 NumPy 数组。它不是对角矩阵 $\Sigma$。Vh:右奇异向量矩阵,已转置($V^T$)。其形状取决于 full_matrices。对于 $m \ge n$,当 full_matrices=False 时,$Vh$ 的形状将是 $n \times n$。重构原始矩阵(验证):为了验证分解 $A = U \Sigma V^T$,您需要从奇异值 s 构建对角矩阵 $\Sigma$。# 构建 Sigma 矩阵(对角矩阵) # 需要匹配维度以进行乘法 Sigma = np.zeros(A.shape) # 以 A 的形状创建零矩阵 # 用奇异值填充对角线 # min(A.shape) 给出奇异值的数量 Sigma[:A.shape[1], :A.shape[1]] = np.diag(s) # 此处使用 Vh 的形状 (n x n) # 重构 A A_reconstructed = U @ Sigma @ Vh # 使用矩阵乘法运算符 @ print("\nSigma 矩阵 ({}):".format(Sigma.shape)) print(Sigma) print("\n重构矩阵 (U @ Sigma @ Vh):") print(A_reconstructed) # 检查是否与原始矩阵接近(由于浮点精度) print("\n重构矩阵是否接近原始矩阵?", np.allclose(A, A_reconstructed))np.allclose() 函数适用于比较浮点数组,同时考虑小的数值误差。当 $m \neq n$ 时,通常建议使用 full_matrices=False(精简 SVD),因为它避免了在 $U$ 或 $V$ 中计算不必要的零向量。LU 分解LU 分解将一个方阵 $A$ 分解为 $L$(下三角矩阵)和 $U$(上三角矩阵),使得 $A = LU$。通常,为了稳定性需要进行主元选择(行交换),从而得到 $PA = LU$,其中 $P$ 是一个置换矩阵。scipy.linalg.lu 处理此过程。# 定义一个方阵 B B = np.array([[2, 5, 8], [4, 2, 1], [9, 3, 7]]) print("原始矩阵 B ({}):".format(B.shape)) print(B) # 执行带置换的 LU 分解 P_mat, L, U = linalg.lu(B) # 注意:P_mat 就是置换矩阵本身 print("\n置换矩阵 P ({}):".format(P_mat.shape)) print(P_mat) print("\n下三角矩阵 L ({}):".format(L.shape)) print(L) print("\n上三角矩阵 U ({}):".format(U.shape)) print(U) # 验证:P @ B 应接近 L @ U print("\nP @ B 是否接近 L @ U?", np.allclose(P_mat @ B, L @ U))理解输出:P_mat:置换矩阵 $P$。将其应用于 $A$(作为 $P @ A$)表示高斯消元过程中执行的行交换。L:对角线上为一的下三角矩阵。U:上三角矩阵。LU 分解是线性方程组($Ax=b$)高效直接求解器的一个重要方法,尽管像 scipy.linalg.solve 这样的函数通常会将其抽象化。QR 分解QR 分解将任意 $m \times n$ 矩阵 $A$ 分解为 $Q$(一个正交矩阵,$Q^T Q = I$)和 $R$(一个上三角矩阵),使得 $A = QR$。scipy.linalg.qr 是应使用的函数。# 再次使用原始矩阵 A # A = np.array([[1, 2, 3], # [4, 5, 6], # [7, 8, 9], # [10, 11, 12]]) print("原始矩阵 A ({}):".format(A.shape)) print(A) # 执行 QR 分解 # mode='economic' 对于非方阵更高效 Q, R = linalg.qr(A, mode='economic') print("\n正交矩阵 Q ({}):".format(Q.shape)) print(Q) print("\n上三角矩阵 R ({}):".format(R.shape)) print(R) # 验证 Q 的正交性:Q.T @ Q 应接近单位矩阵 I = np.eye(Q.shape[1]) # 适当大小的单位矩阵 print("\nQ.T @ Q 是否接近单位矩阵?", np.allclose(Q.T @ Q, I)) # 验证分解:A 应接近 Q @ R print("\nA 是否接近 Q @ R?", np.allclose(A, Q @ R))理解输出:Q:正交(对于复数矩阵是酉)矩阵。当 mode='economic' 时,如果 $m \ge n$,则 $Q$ 是 $m \times n$;如果 $m < n$,则 $Q$ 是 $m \times m$。其列构成 $A$ 的列空间的正交规范基。R:上三角矩阵。当 mode='economic' 时,如果 $m \ge n$,则 $R$ 是 $n \times n$;如果 $m < n$,则 $R$ 是 $m \times n$。QR 分解在用于解决最小二乘问题、特征值计算(QR 算法)和正交化过程的算法中具有重要作用。这些 SciPy 函数提供了应用矩阵分解的工具。它们处理了数值复杂性,使您能够专注于这些分解如何用于解决机器学习中的问题,例如降维(SVD)、求解线性系统(LU)或寻找最小二乘解(QR)。下一节将提供一个侧重于 SVD 的实际示例。