Matlab快速入门之线性代数：矩阵分解

本文属于Matlab快速入门之线性代数的第三篇，即矩阵分解，主要包括：Cholesky 分解、LU 分解、QR 分解、对分解使用多线程计算。讨论的这三种矩阵分解利用了三角形矩阵，其中对角线上下的所有元素都为零。涉及三角矩阵的线性方程组可以使用前代或回代方法轻松快捷地求解。

Cholesky 分解

Cholesky 分解将对称矩阵表示为三角矩阵与其转置的积A = R′R，其中，R 是上三角矩阵。

并非所有对称矩阵都可以通过这种方式进行分解；采用此类分解的矩阵被视为正定矩阵。这表明，A 的所有对角线元素都是正数，并且非对角线元素“不太大”。帕斯卡矩阵提供了有趣的示例。在本章中，示例矩阵 A 为 3×3 帕斯卡矩阵，暂时转换为 6×6，A 的元素为二项式系数。每个元素都是其北向和西向邻点之和。Cholesky 分解可以使用chol函数获得：

>> A = pascal(6)

A =

     1     1     1     1     1     1
     1     2     3     4     5     6
     1     3     6    10    15    21
     1     4    10    20    35    56
     1     5    15    35    70   126
     1     6    21    56   126   252

>> R = chol(A)

R =

     1     1     1     1     1     1
     0     1     2     3     4     5
     0     0     1     3     6    10
     0     0     0     1     4    10
     0     0     0     0     1     5
     0     0     0     0     0     1

这些元素同样为二项式系数。R'*R 等于 A 的情况说明了涉及二项式系数的积之和的单位矩阵。

注意
Cholesky 分解也适用于复矩阵。采用 Cholesky 分解的任何复矩阵都满足
A′ = A，并且被视为 Hermitian 正定矩阵。

通过 Cholesky 分解，可以将线性方程组Ax = b替换为R′Rx = b。

由于反斜杠运算符能识别三角形方程组，因此这可以在 MATLAB 环境中通过以下表达式快速进行求解：

x = R\(R'\b)

如果 A 为 n×n，则 chol(A) 的计算复杂度为 O(n³)，但后续的反斜杠解的复杂度仅为 O(n²)。

LU 分解

LU 分解（或高斯消去法）将任何方阵 A 都表示为下三角矩阵和上三角矩阵的置换之积：A = LU，

其中，L 是对角线元素为 1 的下三角矩阵的置换，U 是上三角矩阵。

出于理论和计算原因，必须进行置换。矩阵

[0  1
1  0]

在不交换其两行的情况下不能表示为三角矩阵的积。尽管矩阵

[ε  1
1   0]

可以表示为三角矩阵之积，但当ε很小时，因子中的元素也会很大并且会增大误差，因此即使置换并非完全必要，它们也是所希望的。部分主元消去法可确保 L 的元素的模以 1 为限，并且 U 的元素并不大于 A 的元素。

例如：

[L,U] = lu(B)

L =
    1.0000         0         0
    0.3750    0.5441    1.0000
    0.5000    1.0000         0

U =
    8.0000    1.0000    6.0000
         0    8.5000   -1.0000
         0         0    5.2941

通过对 A 执行 LU 分解，可以

A*x = b

使用以下表达式快速对线性方程组求解

x = U\(L\b)

行列式和逆矩阵是通过 LU 分解使用以下表达式进行计算的

det(A) = det(L)*det(U)

以及

inv(A) = inv(U)*inv(L)

也可以使用 det(A) = prod(diag(U)) 计算行列式，但行列式的符号可能会相反。

QR 分解

正交矩阵或包含正交列的矩阵为实矩阵，其列全部具有单位长度并且相互垂直。如果 Q 为正交矩阵，则

QTQ = I

其中 I 是单位矩阵。

最简单的正交矩阵为二维坐标旋转：

[ cos(θ)  sin(θ)
 -sin(θ)  cos(θ)]

对于复矩阵，对应的术语为单位。由于正交矩阵和单位矩阵会保留长度、保留角度并且不会增大误差，因此适用于数值计算。

正交或 QR 分解将任何矩形矩阵表示为正交或酉矩阵与上三角矩阵的积。此外，也可能涉及列置换。

A = QR

或

AP = QR

其中，Q 为正交或单位矩阵，R 为上三角矩阵，P 为置换向量。

QR 分解有四种变化形式 – 完全大小或合适大小，以及使用列置换或不使用列置换。

超定线性方程组涉及行数超过列数的矩形矩阵，也即 m×n 并且 m > n。完全大小的 QR 分解可生成一个方阵（m×m 正交矩阵 Q）和一个矩形 m×n 上三角矩阵 R：

>> C=gallery('uniformdata',[5 4], 0);
>> [Q,R] = qr(C)

Q =

    0.6191    0.1406   -0.1899   -0.5058    0.5522
    0.1506    0.4084    0.5034    0.5974    0.4475
    0.3954   -0.5564    0.6869   -0.1478   -0.2008
    0.3167    0.6676    0.1351   -0.1729   -0.6370
    0.5808   -0.2410   -0.4695    0.5792   -0.2207


R =

    1.5346    1.0663    1.2010    1.4036
         0    0.7245    0.3474   -0.0126
         0         0    0.9320    0.6596
         0         0         0    0.6648
         0         0         0         0

在许多情况下，Q 的最后 m – n 列是不需要的，因为这些列会与 R 底部的零相乘。因此，精简 QR 分解可生成一个矩形矩阵（具有正交列的 m×n Q）以及一个方阵 n×n 上三角矩阵 R。对于 5×4 示例，不会节省太多内存，但是对于更大的大量矩形矩阵，在时间和内存方面的节省可能会很重要：

>> [Q,R] = qr(C,0)	
Q =

    0.6191    0.1406   -0.1899   -0.5058
    0.1506    0.4084    0.5034    0.5974
    0.3954   -0.5564    0.6869   -0.1478
    0.3167    0.6676    0.1351   -0.1729
    0.5808   -0.2410   -0.4695    0.5792


R =

    1.5346    1.0663    1.2010    1.4036
         0    0.7245    0.3474   -0.0126
         0         0    0.9320    0.6596
         0         0         0    0.6648

与 LU 分解相比，QR 分解不需要进行任何消元或置换。但是，可选的列置换（因存在第三个输出参数而触发）对检测奇异性或秩亏是很有帮助的。在分解的每一步，未分解的剩余矩阵的列（范数最大）将用作该步骤的基础。这可以确保 R 的对角线元素以降序排列，并且各列之间的任何线性相关性肯定能够通过检查这些元素来显示。对于此处提供的小示例，C 的第二列的范数大于第一列的范数，因此这两列被交换：

>> [Q,R,P] = qr(C)

Q =
   -0.3522    0.8398   -0.4131
   -0.7044   -0.5285   -0.4739
   -0.6163    0.1241    0.7777

R =
  -11.3578   -8.2762
         0    7.2460
         0         0

P =
     0     1
     1     0

组合了合适大小和列置换后，第三个输出参数为置换向量而不是置换矩阵：

>> [Q,R,p] = qr(C,0)

Q =
   -0.3522    0.8398
   -0.7044   -0.5285
   -0.6163    0.1241

R =
  -11.3578   -8.2762
         0    7.2460


p =
     2     1

QR 分解可将超定线性方程组转换为等效的三角形方程组。表达式