.jpg)
矩阵分解
矩阵分解（Matrix Factorizations）就是将一个矩阵用两个以上的矩阵相乘的等式来表达。而矩阵乘法涉及到数据的合成(即将两个或多个线性变换的效果组合成一个矩阵)，因此可以说，矩阵分解是对数据的一种分析。
| 矩阵分解 |
LU分解 |
| 分解形式 |
A=LU
(L代表下三角矩阵，U代表上三角矩阵) |
| 目的 |
提高计算效率 |
| 前提 |
(1)矩阵是方阵（LU分解主要是针对方阵）；
(2)矩阵是可逆的，也就是该矩阵是满秩矩阵，每一行都是独立向量；
(3)消元过程中没有0主元出现，也就是消元过程中不能出现行交换的初等变换。 |
LU分解
分解形式
A=LU=⎣⎢⎢⎡1l21l31l4101l32l42001l430001⎦⎥⎥⎤⎣⎢⎢⎡u11000u12u2200u13u23u330u14u24u34u44⎦⎥⎥⎤
简述
LU分解常用于求解工业和商业问题中的序列方程。它是最常见的求解线性系统 Ax=b 的方法，主要思路是：把A 分解成一个下三角矩阵（Lower Triangular Matrix）和一个上三角矩阵（Upper Triangular Matrix），简称LU。分解后，等式Ax=b可以写成L(Ux)=b，这样我们令y=Ux，这样就可以通过分开求解两个等式来得到x：
Ly=bUx=y
即可以先通过Ly=b求出y，最后再用Ux=y，求出x。
- 本质上，LU分解是高斯消元法的一种表达方式，为了更好地理解LU分解，我们需要解释一下高斯消元法。
高斯消元法（Gaussian Elimination）
简述
为了求解线性系统，我们基本的策略是用相等的式子去替代要求解的式子，来让求解变得更容易。对于一个线性系统中有若干个未知数，我们需要将方程组中的一方程的未知数用含有另一未知数的代数式表示，并将其代入到另一方程中，这就消去了一未知数，得到一解；或将方程组中的一方程倍乘某个常数加到另外一方程中去，也可达到消去一未知数的目的，并最终得到线性系统的解。这一求解算法称之为高斯消元法（Gaussian Elimination）。
示例
下文将展示高斯消元法的求解过程，设需求解的线性方程组如下：
⎩⎪⎨⎪⎧x1−2x2+x3=02x2−8x3=85x1−5x3=10
将其表达为矩阵形式为：
⎣⎡105−2201−8−50810⎦⎤
- 为了将多余的未知数x1消除，我们需要将第三行的5消除，矩阵的第一行乘−5倍再和第三行相加，得到：
⎣⎡100−22101−8−100810⎦⎤
- 接着再进行对多余的未知数x2的消除，这里将第二行乘−5倍再与第三行相加，得到：
⎣⎡100−2201−83008−30⎦⎤
至此，我们得到了一个新的三角形矩阵，转换为线性方程组的形式为：
⎩⎪⎨⎪⎧x1−2x2+x3=0x2−4x3=4x3=−1
通过简单的回代，我们最终可以得到线性方程组的解：
⎩⎪⎨⎪⎧x1=1x2=0x3=−1⎣⎡10001000110−1⎦⎤
高斯消元法与LU分解之间的联系
那么为什么说LU分解其实就是高斯消元法的一种表达方式呢？这是因为我们可以将上述示例中每一步的操作用矩阵相乘的形式（行变换的本质就是左乘矩阵）表达出来。
例如上文示例中的第一步：“为了将多余的未知数x1消除，我们需要将第三行的5消除，矩阵的第一行乘−5倍再和第三行相加” ，可以表达为：
⎣⎡10−5010001⎦⎤⎣⎡105−2201−8−50810⎦⎤=⎣⎡100−22101−8−100810⎦⎤
同样地，第二步：“继续进行对多余的未知数x2的消除，这里将第二行乘−5倍再与第三行相加” 可以表达为：
⎣⎡10001−5001⎦⎤⎣⎡100−22101−8−100810⎦⎤=⎣⎡100−2201−83008−30⎦⎤
因此，我们可以用一个矩阵相乘的等式来记录整个高斯消元的过程：
⎣⎡10001−5001⎦⎤⎣⎡10−5010001⎦⎤⎣⎡105−2201−8−50810⎦⎤=⎣⎡100−2201−83008−30⎦⎤
将三个消元矩阵合并，用E表示则为：
EA=⎣⎡10−501−5001⎦⎤⎣⎡105−2201−8−50810⎦⎤=⎣⎡100−2201−83008−30⎦⎤=U
那么我们可以得到关于A的表达式：
A=⎣⎡105−2201−8−50810⎦⎤=E−1U=⎣⎡105015001⎦⎤⎣⎡100−2201−83008−30⎦⎤=LU
这样一来我们就得到了A的LU分解，即将矩阵A分解为一个下三角矩阵（L）和一个上三角矩阵（U）的乘积。
伪代码
U = A, L = I
for j = 1 : n - 1 do
for i = j + 1 : n do
l_ij = u_ij / u_jj
for k = j : n do
u_ik = u_ik - l_ij * u_jk
end for
end for
end for
- 此处是矩阵index从1开始计数。
- Credit to http://iacs-courses.seas.harvard.edu/courses/am205/slides/am205_lec07.pdf
算法效率
接下来我们来思考一下算法效率的问题，假如我们有一个100∗100的矩阵，并且各项均非0，那我们需要计算多少次才能得到矩阵U呢？
对于这个问题我们先从列的角度来考虑，第一列消元运算结束以后，矩阵变为：
⎣⎢⎢⎢⎡10⋮0⋯⋱⋮⋯⋯⋯⋱⋯⋯⋯⋮⋱⎦⎥⎥⎥⎤
这一步中，第一列的元素运算了100次，而第一行共有100个元素，于是仅第一行与第一列消元结束后，我们就计算了1002次。之后我们要研究剩下的99∗99的矩阵，以此类推可知，最后的计算量为∑k=1nk2。
接着，我们可以利用微积分计算得到，A的操作次数为31n3，b的操作次数为n2次。
结合上文伪代码， 因为有三重循环嵌套，因此LU分解的复杂度为O(n3)，其中加法操作为31n3，乘法操作为31n3，总共为32n3次，求解b和y各需要n2次
前提条件
(1)矩阵是方阵（LU分解主要是针对方阵）；
(2)矩阵是可逆的，也就是该矩阵是满秩矩阵，每一行都是独立向量；
(3)消元过程中没有0主元出现，也就是消元过程中不能出现行交换的初等变换。
LUD分解（LDU Decomposition）
上文中我们已经得到了A=LU，我们还可以继续分解下去得到A=LDU，D为对角矩阵(Diagonal Matrix)，矩阵的LDU分解是在LU分解之后，把U再次分解，目的是把U的对角线元素都化为1，上文示例中我们得到的A=LU如下
A=⎣⎡105−2201−8−50810⎦⎤=LU=⎣⎡105015001⎦⎤⎣⎡100−2201−83008−30⎦⎤
接下来，让U矩阵中位于对角线上的元素都为1，就可以提出一个对角矩阵D
LU=⎣⎡105015001⎦⎤⎣⎡100−2201−83008−30⎦⎤=⎣⎡105015001⎦⎤⎣⎡1000200030⎦⎤⎣⎡100−2101−4104−1⎦⎤=LDU
LUP分解（LU Decomposition with Partial Pivoting）
我们来考虑一个矩阵$ A=\left[\begin{matrix} 0 & 1 \ 1 & 1 \end{matrix}\right]，尽管矩阵A非奇异，并且很简单，但是使用LU分解的算法会失败，因为第一个主元A_{0,0} = 0$。在高斯消元法中，我们常常通过人为的行变换来消除主元为0时的影响。例如我们需要把矩阵A的第一行和第二行进行交换得到$ A=\left[\begin{matrix} 1 & 1 \ 0 & 1 \end{matrix}\right]，这样我们就能将消元进行下去。这样的操作也可以用矩阵相乘的方法来表达，那么在该示例中，即为PA=[0110][0111]=[1011]，这里的矩阵P称为置换矩阵（PermutationMatrix）。那么就有了LUP分解（LUdecompositionwithpartialpivoting）：PA = LU。这一算法的出现显著改善了LU$分解算法的稳定性。
伪代码
U = A, L = I, P = I
for j = 1 : n - 1 do
Select i (>= j) that maximizes |u_ij|
Interchange rows of U: u_(j,j:n) <-> u_(i,j:n)
Interchange rows of L: l_(j,1:j-1) <-> l_(i,1:j-1)
Interchange rows of P: p_(j,:) <-> p_(i,:)
for i = j + 1 : n do
l_ij = u_ij / u_jj
for k = j : n do
u_ik = u_ik - l_ij * u_jk
end for
end for
end for
