行列式
行列式是方阵的一个属性
矩阵可以表示一组向量，方阵表示n个n维向量，用一个数字表示这些向量组的不同？
比方说在二维平面中，这里有三组二维向量，每组都有两个向量，那么每组向量的面积就可以表示它们的不同。当然这里说面积是针对二维平面来说的，在三维空间中，就是体积；在更高维度中，可能就是一个体，但这个体比较抽象
行列式描述的就是n个n维向量所组成的这样一个体的属性。我们从二维空间入手，来看一看行列式是怎么计算的
像这样的两个向量所组成的行列式，我们记作
，这里我们可以看到，这两个向量，在行列式中，我们是按行来排列的。但其实按行排列还是按列排列都是可以的，是无所谓的。对于行列式，也可以表示成
，而求它的值，实际上就是求上图中平行四边形的面积。
通过做绿色的辅助线，可得这个平行四边形的面积=(a + c)(b + d) - 2bc - cd - ab = ad - bc
则
对于二阶方阵来说，它的行列式的值为主对角线上的乘积减去非主对角线上的乘积。
反过来，则
，由此，我们可以看出
行列式表示向量组在空间中形成的有向体积
对于二维空间，我们往往是先从x轴看向y轴，则对于上面的u和v两个向量来说，从u到v的逆时针方向上，我们就认为它的行列式就是正的。而从u到v的顺时针方向上，我们就认为它的行列式是负的。
但是在三维或以上的空间中，体积的方向将变得极其复杂。简单说，在行列式中，向量排列的顺序是有意义的。在n阶行列式中，任意交换两行，则行列式的值取反。
行列式的基本性质
- 单位矩阵的行列式等于1，即det I = 1
- 交换行列式的两行，则行列式的值取反。
- 方阵的某一行乘以一个数k，则其对应的行列式也缩放了k倍，即
，比如
(a,b)这个向量扩展了两倍，则这个平行四边形的面积也扩展了两倍。这里需要注意的是，是其中的某一行，而不是所有行，这和矩阵的数量乘法是不一样的。
(其中A为方阵)；而正确的写法为
，，其中n是指A为n阶方阵。 - 方阵中的某一行加上一行数，则有：
，它其实表示的是一个向量相加，当然我们也可以直接代进公式，左边 = (a + a')d - (b + b')c = ad + a'd - bc - b'c。右边= ad - bc + a'd - b'c，说明左边=右边
这里的第4和第3条是行列式的加法和数量乘法，这是非常有意义的两种基本运算。
行列式和矩阵的逆
性质一：如果行列式的两行相同，则行列式的值为0。
证明： 假设方阵A交换两行后为A'
det(A) = -det(A')，因为两行相同，所以A = A'
det(A) = -det(A) 则det(A) = 0，得证
从直观上考虑：二维空间中，两个向量共线，面积为0；在三维空间中，两个向量共线，形成一个平面，体积为0；n维空间中，两个向量共线，则退化成一个n-1维的体，它的n维的体积为0。
性质二：如果行列式的一行是另一行的k倍，则行列式的值为0.
证明：
这一条跟上面一条在直观上理解是一样的，因为它们共线。
性质三：如果行列式有一行为0，则行列式的值为0。
证明：
在直观上理解，就是一个n维体中，有一个零向量，则这个n维体退化成了n-1维体，它的n维的体积为0
性质四：如果行列式的一行是其他行的线性组合，则行列式的值为0。
证明：
这是根据行列式的加法以及性质二来证明的。从n维空间的角度考虑，如果有两行线性相关，则无法构成n维空间的一组基，就算其他行都线性无关，也退化成n-1维空间，它的n维空间的体积为0.
行列式形成一个向量组，如果这组向量线性相关，则行列式的值为0，det(A)=0
矩阵不可逆
如果这组向量线性无关，则行列式的值不为0,det(A)≠0矩阵可逆
现在对于方阵A的等价命题又可以增加一条
- 对于方阵A,矩阵A可逆，A是非奇异矩阵
- 线性系统Ax=0(齐次线性方程组)只有唯一解，x=0
- 矩阵的行最简形式rref(A)=I
- A可以表示成一系列初等矩阵的乘积
- Ax=b只有唯一解
- 方阵A的列向量线性无关
- 方阵A的列向量可以生成n维空间
- 方阵A的列向量是n维空间的基
- 方阵A为满秩矩阵(秩=n)
- 方阵A的行秩为n
- 方阵A的列秩为n
- 方阵A的行空间为
- 方阵A的列空间为
- 方阵A的零空间为{O},维度为0
- det(A)≠0
计算行列式的算法
如果一个行列式的一行加(减)另一行的k倍，行列式的值不变。
证明：
这条性质是计算行列式的值的基础。
一个方阵行列式的值
- 等于其进行高斯消元法后的结果
上三角矩阵U - 等于其进行高斯-约旦消元法后的结果对角矩阵D(只有主对角线上的元素为非零元素，其他位置上的元素都为0)
但是这个高斯-约旦消元法跟之前讲线性系统时候的有所不同，它不能进行归一化操作。因为一旦进行归一化操作，行列式的值就发生了变化，这跟矩阵的性质是不一样的，归一化的本质就是某一行乘以了某个常数k，行列式的值就扩大了k倍。并且进行行置换和列置换需要改变行列式的正负号，因为任意交换两行，行列式的值取反。如果消元结果有零行，行列式的值为0。
对角矩阵的行列式
上三角矩阵的行列式
其实上三角矩阵的行列式跟对角矩阵的行列式的值是相同的，也是
，这是因为在进行约旦消元法的过程中，将主元列上面的数减去主元值乘以一个系数，全部化成0，也就变成了对角矩阵，而主元列的主元值都不变。
下三角矩阵的行列式
跟上三角矩阵的求法类似，进行高斯消元法就好，最后依然是一个对角矩阵，而主元列的主元值不变，依然为
所以上面的计算行列式的值不需要进行约旦消元法，只需要进行高斯消元法即可。
简单用2阶行列式进行验证
新增一个接口Formula
public interface Formula<T extends Number> { /** * 求行列式的值 * @return */ T result(Addr<T> addr); }
行列式实现类
@ToString public class Determinant<T extends Number> implements Formula<T> { //原始方阵 private Ranks<T> source; //高斯消元法后的方阵 private Ranks<T> matrix; @SuppressWarnings("unchecked") public Determinant(Ranks<T> matrix) throws CloneNotSupportedException { if (matrix.rowNum() != matrix.colNum()) { throw new IllegalArgumentException("矩阵必须为方阵"); } this.source = matrix; this.matrix = ((Matrix) matrix).clone(); } @Override public T result(Addr<T> addr) { T res = addr.one(); int times = forword(addr); for (int i = 0; i < matrix.rowNum(); i++) { res = addr.mul(res,matrix.get(new int[]{i,i})); } //交换次数为奇数，结果取负 if ((times & 1) == 1) { return addr.neg(res); } return res; } /** * 高斯消元，此处不可做归一化处理 * @param addr * @return */ @SuppressWarnings("unchecked") private int forword(Addr<T> addr) { int i = 0; int k = 0; //交换位置的次数 int times = 0; while (i < matrix.rowNum() && k < matrix.colNum()) { int maxRow = maxRow(i, k, matrix.rowNum(), addr); T[] temp = (T[]) new Number[matrix.colNum()]; if (maxRow != i) { times++; } //每一行无论是否有主元，都与主元最大值所在行进行交换 for (int j = 0; j < matrix.colNum(); j++) { temp[j] = matrix.get(new int[]{i, j}); matrix.set(new int[]{i, j}, matrix.get(new int[]{maxRow, j})); matrix.set(new int[]{maxRow, j}, temp[j]); } if (addr.equals(matrix.get(new int[]{i,k}),addr.zero())) { k++; }else { //将主元所在行以下的行的所有非主元数据全部变成0 for (int m = i + 1; m < matrix.rowNum(); m++) { T u1 = matrix.get(new int[]{m, k}); for (int l = 0; l < matrix.colNum(); l++) { matrix.set(new int[]{m, l}, addr.sub(matrix.get(new int[]{m, l}), addr.mul(matrix.get(new int[]{i, l}), addr.div(u1,matrix.get(new int[]{i,k}))))); } } i++; } } return times; } /** * 获取主元最大值所在的行 * @param addr * @return */ private int maxRow(int indexI,int indexJ,int n,Addr<T> addr) { T best = matrix.get(new int[]{indexI,indexJ}); int ret = indexI; for (int i = indexI + 1; i < n; i++) { if (addr.compair(matrix.get(new int[]{i,indexJ}),best) == 1) { best = matrix.get(new int[]{i,indexJ}); ret = i; } } return ret; } public static void main(String[] args) throws CloneNotSupportedException { Addr<Double> addr = new Addr<Double>() { @Override public Double zero() { return 0.0; } @Override public Double one() { return 1.0; } @Override public Double add(Double lhs, Double rhs) { return lhs + rhs; } @Override public Double sub(Double lhs, Double rhs) { return lhs - rhs; } @Override public Double mul(Double k, Double hs) { return k * hs; } @Override public Double square(Double hs) { return hs * hs; } @Override public double prescription(Double hs) { return Math.sqrt(hs); } @Override public double div(Double hs, double nhs) { return hs / nhs; } @Override public Double div(Double hs, Double nhs) { return hs / nhs; } @Override public double acos(double hs) { return Math.acos(hs); } @Override public double toDegree(double hs) { return Math.toDegrees(hs); } @Override public int compair(Double lhs, Double rhs) { if (lhs - rhs > 0) { return 1; }else if (lhs - rhs < 0) { return -1; }else { return 0; } } @Override public boolean equals(Double lhs, Double rhs) { double dis = 1e-8; if (Math.abs(lhs - rhs) < dis) { return true; } return false; } @Override public Double neg(Double hs) { return -hs; } }; Ranks<Double> matrix = new Matrix<>(new Double[][]{{7.0,3.0},{2.0,5.0}}); Formula<Double> determinant = new Determinant<>(matrix); System.out.println(determinant.result(addr)); System.out.println(determinant); } }
运行结果
29.000000000000004
Determinant(source=Matrix{values=[[7.0, 3.0],[2.0, 5.0]]}, matrix=Matrix{values=[[7.0, 3.0],[0.0, 4.142857142857143]]})
Python代码
from ._global import is_zero from .Matrix import Matrix class Determinant: def __init__(self, matrix): assert matrix.row_num() == matrix.col_num(), "矩阵必须为方阵" self._source = matrix self._matrix = [matrix.row_vector(i) for i in range(matrix.row_num())] def result(self): # 求行列式的值 res = 1 times = self._forword() for i in range(len(self._matrix)): res *= self._matrix[i][i] if times % 2 == 1: return -res return res def _forword(self): # 高斯消元，此处不可做归一化处理 i, k, times = 0, 0, 0 while i < len(self._matrix) and k < len(self._matrix[0]): # 看matrix[i][k]位置是否可以为主元 max_row = self._max_row(i, k, len(self._matrix)) if max_row != i: times += 1 self._matrix[i], self._matrix[max_row] = self._matrix[max_row], self._matrix[i] if is_zero(self._matrix[i][k]): k += 1 else: for j in range(i + 1, len(self._matrix)): self._matrix[j] = self._matrix[j] - self._matrix[j][k] * (self._matrix[i] / self._matrix[i][k]) i += 1 return times def _max_row(self, index_i, index_j, n): # 主元最大值所在的行 best, ret = self._matrix[index_i][index_j], index_i for i in range(index_i + 1, n): if self._matrix[i][index_j] > best: best, ret = self._matrix[i][index_j], i return ret def __repr__(self): return "Determinant({},{})".format(self._source, Matrix(self._matrix)) __str__ = __repr__
from playLA.Determinant import Determinant from playLA.Matrix import Matrix if __name__ == "__main__": matrix = Matrix([[7, 3], [2, 5]]) determinant = Determinant(matrix) print(determinant.result()) print(determinant)
运行结果
29.000000000000004
Determinant(Matrix([[7, 3], [2, 5]]),Matrix([[7, 3], [0.0, 4.142857142857143]]))
初等矩阵与行列式
之前我们讲的性质都可以从几何方面来进行理解，但是对于这个性质对于几何来讲就没什么用了，因为对于两个空间体体积的乘积，在几何上的意义就不那么明确了。所以我们要从纯代数的角度出发，来进行证明。
首先，如果A或者B中的某一行和其他行线性相关，则
,假设A中某一行和其他行线性相关，则
的结果矩阵中依然会保持线性相关性。
，所以在这种情况下，等式的两边成立。
现在我们来看A和B中的所有行都线性无关，从方阵的等价命题中
- 对于方阵A,矩阵A可逆，A是非奇异矩阵
- 线性系统Ax=0(齐次线性方程组)只有唯一解，x=0
- 矩阵的行最简形式rref(A)=I
- A可以表示成一系列初等矩阵的乘积
- Ax=b只有唯一解
- 方阵A的列向量线性无关
- 方阵A的列向量可以生成n维空间
- 方阵A的列向量是n维空间的基
- 方阵A为满秩矩阵(秩=n)
- 方阵A的行秩为n
- 方阵A的列秩为n
- 方阵A的行空间为
- 方阵A的列空间为
- 方阵A的零空间为{O},维度为0
- det(A)≠0
我们可以看到第4条，这表示A或者B都能表示为一系列初等矩阵的乘积。而初等矩阵是对单位矩阵进行初等变换(行操作)，则我们把A拆解成一系列初等矩阵
，进而我们需要研究一下初等矩阵的行列式
初等矩阵的变换一共有三种形式
- 如果E是单位矩阵的某一行乘以k,很明显 det(E) = k
- 如果E是单位矩阵的某两行交换位置，很明显 det(E) = -1
- 如果E是单位矩阵的某行加(减)另一行的k倍 det(E) = 1
现在我们需要先证明这个命题
现在我们依然分三种情况来证明
- 如果E是单位矩阵的某一行乘以k，方阵EB是B中的某一行乘以k，则
，而等式右边
，左右相等，在该种情况下得证。 - 如果E是单位矩阵的某两行交换位置，方阵EB是B中某两行交换位置，
，等式右边
，左右相等，在该种情况下得证。 - 如果E是单位矩阵的某行加(减)另一行的k倍，方阵EB是B中的某行加(减)另一行的k倍，
，等式右边
，左右相等，在该种情况下得证。
从而，下面的式子就可以进行转化
将A拆成一系列的初等矩阵的乘积，就有
而A本身的行列式为
最终得到
得证
我们将B替换成A的逆，则有
这里可以计算A的逆的行列式的值，而A的行列式的值由于在分母的位置，所以它不能为0，换句话说，如果A的行列式的值为0，A的逆的行列式的值不存在，进而A的逆不存在。这进而说明了方阵A的等价命题中A若存在逆，det(A)≠0。
矩阵的LDU分解，PLU分解，PLUP分解
- LDU分解
之前我们讲过矩阵的LU分解，将一个矩阵分解成了一个单位下三角矩阵L和一个上三角矩阵U
在不做归一化处理的高斯消元法中，将矩阵A化成一个上三角矩阵U,同时产生一个初等矩阵的逆矩阵L。现在如果我们想把U也化成一个单位上三角矩阵的话，L不变，等于我们进行高斯消元法的时候进行归一化处理，此时也相当于进行了一系列的初等矩阵的变换，我们同样要获取这些初等矩阵的逆矩阵D
我们可以看到D是一个对角矩阵。此时
- PLU分解
我们再来看一个矩阵的高斯消元的过程
此时它的这两步操作的初等矩阵的逆矩阵L为
，当我们要对第二行的主元进行消元的时候，我们发现，第二行第二列的元素为0，如果按照正常的高斯消元法，我们会替换第二行和第三行。
但同时会产生一个交换操作，之前的初等矩阵的逆矩阵L也要交换位置，此时这个矩阵将不再是单位下三角矩阵。在不动L的情况下，我们需要新产生一个初等矩阵的逆矩阵
，我们称这个矩阵为置换矩阵P，虽然有了矩阵P，但由于我们还是交换了位置，矩阵L也就变成了
最终得到
- PLUP分解
PLU分解并不能解决所有的矩阵的问题，比如
我们在进行高斯消元的过程中，第1列1下面两个0，这个没有问题，但是在对第二列进行消元的过程中，第二行和第三行的值都是0，即便进行行交换也没有用。则PLU分解的方法就失效了。为了继续高斯消元，需要交换矩阵的两列。但是初等矩阵中没有交换矩阵两列的初等变换，只能交换矩阵的两行。例如
要交换矩阵的两列，需要右乘以置换矩阵，比如
所以A需要去右乘一个置换矩阵，完成列交换，而这个置换矩阵我们称之为P'，再根据之前的PLU分解，最后可得
这里的A可以是任意矩阵。
行式就是列式
对于一个方阵，无论它其中的向量是行排列还是列排列，它们的行列式的值相等。
即为
之前我们讲解了任意A可以分解成PLUP'四个矩阵，其中P是一个行变换矩阵，P'是一个列变换矩阵，则
在这里，L、U的行列式的值和它们转置的行列式的值肯定是相等的，因为L是下三角矩阵，U是上三角矩阵，它们的行列式都等于对角矩阵的行列式的值，而L转置后是一个上三角矩阵，U转置后是一个下三角矩阵，它们的行列式同样等于对角矩阵的行列式的值。
P和P'是单位矩阵经过几次交换位置的初等变换的初等矩阵的结果，所以我们可以先证明这个命题
如果E表示单位矩阵两行交换位置det(E) = -1，而一个单位矩阵的两行交换了位置之后再进行转置，转置的结果还是两行交换了位置，它的行列式的值同样为-1，得证。
由于P和P'和它们的转置后的行交换的次数是一样的，所以P和P'的行列式的值和它们转置后的行列式的值是相等的，得证。
我们之前讲的所有性质，都是基于行的。换成列一样存在。
- 交换行列式的两列，则行列式的值取反。
- 方阵的某一列乘以一个数k，则其对应的行列式也缩放了k倍，即
- 方阵中的某一列加上一列数，则有：
- 如果行列式的两列相同，则行列式的值为0。
- 如果行列式的一列是另一列的k倍，则行列式的值为0。
- 如果行列式的一列是其他列的线形组合，则行列式的值为0。
- 如果一个方阵加(减)另一列的k倍，行列式的值不变。
特征值和特征向量
特征值和特征向量也是方阵的一个属性，描述的是方阵的"特征"，把方阵当作变换的时候的一个"特征"
我们之前已经知道矩阵看作变换，可以把空间中的一个向量(一个点)转换成另外一个向量(一个点)
图中就是把(1,4)这个蓝色的向量转换到了(-4,5)这个绿色的向量。当然它也是一个基变换，从(4,1)、(-2,1)组成的一组基到标准基的变换。
我们再来看一组转换
在这组转换里面，转换后的向量和转换前的向量在方向上并没有发生改变，只不过变成了原来的向量的某一个常数倍，满足
这里
称为A的特征值(eigenvalue)
称为A对应于的特征向量(eigenvector)
求解特征值和特征向量
这里需要注意的是，零向量肯定满足。是平凡解。既然是A的特征，零向量自然无法反映A的特征。如果是零向量就跟A无关了，A可以任意取。所以，特征向量不考虑零向量。但=0并不平凡，此时等式就变成了
，由于不能为零向量，该方程变成了一个齐次线性方程组。
根据方阵的等价命题中
- 对于方阵A,矩阵A可逆，A是非奇异矩阵
- 线性系统Ax=0(齐次线性方程组)只有唯一解，x=0
- 矩阵的行最简形式rref(A)=I
- A可以表示成一系列初等矩阵的乘积
- Ax=b只有唯一解
- 方阵A的列向量线性无关
- 方阵A的列向量可以生成n维空间
- 方阵A的列向量是n维空间的基
- 方阵A为满秩矩阵(秩=n)
- 方阵A的行秩为n
- 方阵A的列秩为n
- 方阵A的行空间为
- 方阵A的列空间为
- 方阵A的零空间为{O},维度为0
- det(A)≠0
的第2条，我们知道，A可逆的话，只能为零向量，所以A一定不可逆。
所以这里特征值可以为0
我们再对这个等式进行一系列的变形(由于矩阵不能与常数相减，所以将常数乘以一个单位矩阵，等式两边依然成立)
我们其实是希望
这个齐次线性方程组有非零解，根据方阵的等价命题中的2和15可知，其实就是
不可逆，必然的行列式等于0，即
，我们称该方程为特征方程
我们依然以
为例来说明该行列式的计算过程
这里A是二阶方阵，所以I为二阶单位矩阵，根据矩阵的加减法，代入得
由于我们知道了可以为2或者3，依次代入齐次线性方程组，就可以求出相应的特征向量
经过高斯消元，它们都有零行，u是一个二维向量，表示未知数的个数为2，而系数矩阵的非零行为1，小于未知数个数，表示它们都有无数解，这里的s表示可以取任意的实数。可以理解成(1,1)、(2,1)都是一维空间的基，它们所生成的一维空间的所有的向量都是该齐次线性方程组的解。
特征值和特征向量的相关概念
如果u是A对应于的一个特征向量，则ku也是一个特征向量。这里的k≠0。
这里的k其实就等于上面的s，它可以取任意值，这里我们也来证明一下。
它同样满足的基本式，所以得证
该命题同样也说明了有无数解
特征向量组成了这个方阵的零空间。当然这里要刨除零向量。
准确的说这个零空间
称为对应的特征空间
特征方程是一个关于的n次方程，有n个解(这个范围是一个复数范围，而不仅仅是实数内)
之前我们例子中的
它就是一个二次方程的两个解，我们称为简单特征值
但也有可能
这个二次方程的两个解相同，我们称为多重特征值，它的重数为2
这个二次方程的两个解为复数，我们称为复数特征值，不过复数特征值和多重特征值我们讨论的比较少，而应用最多的依然是简单特征值
特征值和特征向量的性质
之前我们已经介绍过，如果特征值为0，A一定不可逆，则现在我们可以为方阵A的等价命题再增加一条
- 对于方阵A,矩阵A可逆，A是非奇异矩阵
- 线性系统Ax=0(齐次线性方程组)只有唯一解，x=0
- 矩阵的行最简形式rref(A)=I
- A可以表示成一系列初等矩阵的乘积
- Ax=b只有唯一解
- 方阵A的列向量线性无关
- 方阵A的列向量可以生成n维空间
- 方阵A的列向量是n维空间的基
- 方阵A为满秩矩阵(秩=n)
- 方阵A的行秩为n
- 方阵A的列秩为n
- 方阵A的行空间为
- 方阵A的列空间为
- 方阵A的零空间为{O},维度为0
- det(A)≠0
- 0不是A的特征值
根据特征值和特征向量的定义和特征方程
如果A是一个对角矩阵，它的行列式为
那么的行列式为
由此可见，对角矩阵的特征值是其对角线上的元素
同理，如果A是一个上三角矩阵，它的行列式为
那么的行列式为
则上三角矩阵的特征值是其对角线上的元素
同理，如果A是一个下三角矩阵，它的行列式为
那么的行列式为
则下三角矩阵的特征值是其对角线上的元素
若是A的特征值，则
是
的特征值
这里我们可以采用数学归纳法来证明
当m=1时，成立
假设m=k时
成立
当m=k+1时
得证
但是这里的证明，m都是大于1的，它并不能说明m=-1时的情况，但是m=-1时，等式依然成立，我们需要有另外的证明
若是A的特征值，则
是
(这里叫A的逆)的特征值
证明:
得证
直观理解特征值和特征向量
我个人觉得最直观的理解可以参考https://m.bilibili.com/video/BV1ys411472E?p=14
即一个向量在经过线性变换后(被一个变换矩阵点乘)，它的方向不发生变化，而模的长度发生了伸缩，这个伸缩的倍数就是特征值，而以这个向量为基的一维空间的除了零向量以外的所有向量就是特征向量。现在我们来理解几个特例
投影变换
二维空间，把任何向量投影到(2,1)所在的直线，这个变换对应一个矩阵，虽然这个矩阵我们并不知道是什么，但是根据线性变换的定义
T(u + v) = T(u) + T(v)
T(cu) = cT(u) c
R
我们将任意两个向量，先相加再投影，和先投影再相加，所得到的结果一定是相同的。同一个向量乘以一个常数再投影，与先投影再乘以一个常数也是相同的。
使用代码来验证，向量类重写Clone()方法
@Getter @AllArgsConstructor public class Vector<T extends Number> implements Attribute<T>,Cloneable { private T[] values; @SuppressWarnings("unchecked") public static <T extends Number> Attribute<T> zero(int dim,Addr<T> addr) { T[] values = (T[])new Number[dim]; Arrays.fill(values,addr.zero()); return new Vector<>(values); } @Override public String toString() { return "Vector{" + "values=" + Arrays.toString(values) + '}'; } @Override @SuppressWarnings("unchecked") public Vector clone() throws CloneNotSupportedException { T[] values = Arrays.copyOf(this.values,this.len()); return new Vector(values); } @Override @SuppressWarnings("unchecked") public boolean equals(Attribute<T> another,Addr<T> addr) { if (this.len() != another.len()) { return false; } for (int i = 0; i < len(); i++) { if (!addr.equals(this.get(i),another.get(i))) { return false; } } return true; } @Override public int len() { return values.length; } @Override public T get(int index) { if (index >= len() || index < 0) { throw new IllegalArgumentException("索引超出范围"); } return values[index]; } @Override public void set(int index, T value) { if (index >= len() || index < 0) { throw new IllegalArgumentException("索引超出范围"); } this.values[index] = value; } @Override @SuppressWarnings("unchecked") public Attribute<T> add(Attribute<T> another,Addr<T> addr) { if (this.len() != another.len()) { throw new IllegalArgumentException("加法错误，两个向量的维度必须相等"); } T[] newValues = (T[])new Number[this.len()]; for (int i = 0; i < this.len(); i++) { newValues[i] = addr.add(this.get(i),another.get(i)); } this.values = newValues; return this; } @Override @SuppressWarnings("unchecked") public Attribute<T> sub(Attribute<T> another, Addr<T> addr) { if (this.len() != another.len()) { throw new IllegalArgumentException("减法错误，两个向量的维度必须相等"); } T[] newValues = (T[])new Number[this.len()]; for (int i = 0; i < this.len(); i++) { newValues[i] = addr.sub(this.get(i),another.get(i)); } this.values = newValues; return this; } @Override @SuppressWarnings("unchecked") public Attribute<T> div(T k, Addr<T> addr) { T[] newValues = (T[])new Number[this.len()]; for (int i = 0; i < this.len(); i++) { newValues[i] = addr.div(this.get(i),k); } this.values = newValues; return this; } @Override @SuppressWarnings("unchecked") public Attribute<T> mul(T k, Addr<T> addr) { T[] newValues = (T[])new Number[this.len()]; for (int i = 0; i < this.len(); i++) { newValues[i] = addr.mul(k,this.get(i)); } this.values = newValues; return this; } @Override public double norm(Addr<T> addr) { T total = addr.zero(); for (int i = 0; i < this.len(); i++) { total = addr.add(total,addr.square(this.get(i))); } return addr.prescription(total); } @Override public Attribute<Double> normalize(Addr<T> addr) { Double[] newValues = new Double[this.len()]; double norm = this.norm(addr); for (int i = 0; i < this.len(); i++) { newValues[i] = addr.div(this.get(i),norm); } return new Vector<>(newValues); } @Override public T dot(Attribute<T> another, Addr<T> addr) { if (this.len() != another.len()) { throw new IllegalArgumentException("点乘错误，两个向量的维度必须相等"); } T total = addr.zero(); for (int i = 0; i < this.len(); i++) { total = addr.add(total,addr.mul(this.get(i),another.get(i))); } return total; } @Override public double degree(Attribute<T> another, Addr<T> addr) { T dot = this.dot(another,addr); double selfNorm = this.norm(addr); double anotherNorm = another.norm(addr); double productNorm = selfNorm * anotherNorm; return addr.toDegree(addr.acos(addr.div(dot,productNorm))); } @Override public Attribute<Double> project(Attribute<T> another, Addr<T> addr) { T dot = this.dot(another,addr); double selfNorm = this.norm(addr); Attribute<Double> unitVector = this.normalize(addr); double d = addr.div(dot, selfNorm); return mul(d,unitVector); } private Attribute<Double> mul(double d,Attribute<Double> v) { Double[] values = new Double[v.len()]; for (int i = 0; i < v.len(); i++) { values[i] = d * v.get(i); } return new Vector<>(values); } @Override @SuppressWarnings("unchecked") public Attribute<Double> rectangular(Attribute<T> another, Addr<T> addr) { Attribute<Double> att = (Attribute<Double>) another; Addr<Double> add = (Addr<Double>) addr; return att.sub(this.project(another,addr),add); } @SuppressWarnings("unchecked") public static void main(String[] args) throws CloneNotSupportedException { Attribute<Integer> vector = new Vector<>(new Integer[]{5,2}); System.out.println(vector.toString()); System.out.println(vector.len()); System.out.println(vector.get(0)); Attribute<Integer> vector1 = new Vector<>(new Integer[]{6,3}); Addr<Integer> addr = new Addr<Integer>() { @Override public Integer zero() { return 0; } @Override public Integer one() { return 1; } @Override public Integer add(Integer lhs, Integer rhs) { return lhs + rhs; } @Override public Integer sub(Integer lhs, Integer rhs) { return lhs - rhs; } @Override public Integer mul(Integer k, Integer hs) { return k * hs; } @Override public Integer square(Integer hs) { return hs * hs; } @Override public double prescription(Integer hs) { return Math.sqrt(hs); } @Override public double div(Integer hs, double nhs) { return hs / nhs; } @Override public Integer div(Integer hs, Integer nhs) { return hs / nhs; } @Override public double acos(double hs) { return Math.acos(hs); } @Override public double toDegree(double hs) { return Math.toDegrees(hs); } @Override public int compair(Integer lhs, Integer rhs) { if (lhs - rhs > 0) { return 1; }else if (lhs - rhs < 0) { return -1; }else { return 0; } } @Override public boolean equals(Integer lhs, Integer rhs) { return lhs == rhs; } @Override public Integer neg(Integer hs) { return -hs; } }; Addr<Double> addr1 = new Addr<Double>() { @Override public Double zero() { return 0.0; } @Override public Double one() { return 1.0; } @Override public Double add(Double lhs, Double rhs) { return lhs + rhs; } @Override public Double sub(Double lhs, Double rhs) { return lhs - rhs; } @Override public Double mul(Double k, Double hs) { return k * hs; } @Override public Double square(Double hs) { return hs * hs; } @Override public double prescription(Double hs) { return Math.sqrt(hs); } @Override public double div(Double hs, double nhs) { return hs / nhs; } @Override public Double div(Double hs, Double nhs) { return hs / nhs; } @Override public double acos(double hs) { return Math.acos(hs); } @Override public double toDegree(double hs) { return Math.toDegrees(hs); } @Override public int compair(Double lhs, Double rhs) { if (lhs - rhs > 0) { return 1; }else if (lhs - rhs < 0) { return -1; }else { return 0; } } @Override public boolean equals(Double lhs, Double rhs) { double dis = 1e-6; if (Math.abs(lhs - rhs) < dis) { return true; } return false; } @Override public Double neg(Double hs) { return -hs; } }; //(5,2)*(6,3) System.out.println(vector.dot(vector1,addr)); //(5,2),(6,3)的夹角 System.out.println(vector.degree(vector1,addr)); //(6,3)在(5,2)上的映射 System.out.println(vector.project(vector1,addr)); //(5,2)+(6,3) vector.add(vector1, addr); System.out.println(vector); vector.mul(3,addr); System.out.println(vector); System.out.println(Vector.zero(4,addr)); System.out.println(vector.norm(addr)); Attribute<Double> unitVector = vector.normalize(addr); System.out.println(unitVector); System.out.println(unitVector.norm(addr1)); Attribute<Double> v1 = new Vector<>(new Double[]{5.0,2.0}); Attribute<Double> v2 = new Vector<>(new Double[]{6.0,3.0}); Attribute<Double> v3 = v1.rectangular(v2,addr1); //(5,2)的垂直向量 System.out.println(v3); //垂直的两个向量点乘为0 System.out.println(v1.dot(v3,addr1)); System.out.println(addr1.equals(addr1.zero(),v1.dot(v3,addr1))); Attribute<Double> line = new Vector<>(new Double[]{2.0,1.0}); Attribute<Double> v = new Vector<>(new Double[]{2.0,5.0}); Attribute<Double> u = new Vector<>(new Double[]{-4.0,-4.0}); Attribute<Double> vc = ((Vector) v).clone(); Attribute<Double> uc = ((Vector) u).clone(); //先相加，再投影 System.out.println(line.project(v.add(u,addr1),addr1)); //先投影，再相加 System.out.println(line.project(vc,addr1).add(line.project(u,addr1),addr1)); //先相乘，再投影 System.out.println(line.project(u.mul(3.0,addr1),addr1)); //先投影，再相乘 System.out.println(line.project(uc,addr1).mul(3.0,addr1)); } }
运行结果
Vector{values=[5, 2]}
2
5
36
4.763641690726143
Vector{values=[6.20689655172414, 2.4827586206896552]}
Vector{values=[11, 5]}
Vector{values=[33, 15]}
Vector{values=[0, 0, 0, 0]}
36.24913792078372
Vector{values=[0.9103664774626047, 0.413802944301184]}
0.9999999999999999
Vector{values=[-0.20689655172413968, 0.5172413793103448]}
-8.881784197001252E-15
true
Vector{values=[-1.2, -0.6]}
Vector{values=[-1.2000000000000002, -0.6000000000000001]}
Vector{values=[-14.399999999999999, -7.199999999999999]}
Vector{values=[-14.399999999999999, -7.199999999999999]}
从最后四行结果可以看出，它们满足线性变换的定义。
现在我们可以考虑一下，什么向量在经过投影之后的方向跟原向量还在一个直线上呢？
很明显，在(2,1)相同直线上的向量的投影跟原向量在同一直线上，此时它的伸缩长度不变，
；另外就是跟(2,1)直线垂直的向量的投影(此时为零向量)，跟原向量在同一直线上，由于是零向量，所以
。
这是一个初等矩阵，代表了交换两行。对于向量来说，就代表了交换两个维度的值。
从二维空间的角度来看，它代表着所有的向量关于(1,1)这个向量所在直线的翻转。那么在什么情况下，向量翻转后跟原向量还在同一个直线上呢，答案是很明显的。
要么，这个向量就在(1,1)这个向量所在的直线上，它翻转后的向量为原向量，伸缩长度不变，；要么跟(1,1)这个向量所在的直线垂直，它翻转后的向量为原向量的反方向，长度不变，所以
。我们用代数的方式来验证一下
所以这里特征值为正负1.
如果矩阵A含有两个不同的特征值，则它们对应的特征向量线性无关。
现在我们要证明u和v是线性无关的
反证法： 假设u和v是线性相关的，则
，则有
由于
，所以
由于，左右两边同时乘以
，得到
则
这里k≠0，
和是两个不同的特征值，所以它们相减也不为0，那剩下的就只有v为零向量了，这与特征向量不包含零向量相矛盾，所以u和v不可能线性相关，就只能是线性无关的了，得证。
"不简单"的特征值
当一个矩阵的特征值为不同的实数的时候，我们称之为简单特征值，否则为多重特征值和复数特征值。
在该图中，我们要将浅色的F变成深色的F，等于是逆时针旋转了90度，等于顺时针旋转了270度，之前我们讲过旋转矩阵
因为
为270度，则该旋转矩阵具体为
，通过直观的方式，我们看不到任何向量在旋转后还跟原向量在同一直线上，根本不存在特征值和特征向量。通过代数的方式，我们来看一下
则这个特征值为两个复数
，它所对应的几何含义已经消失了。
单位矩阵
这是一个二阶单位矩阵就是将二维平面中的原向量转化成它自己，向量根本不发生变化。那么在二维平面中所有的向量都是A的特征向量，而特征值为1。通过代数的形式来看为
我们将特征值代入得
由于这个齐次线性方程组没有主元列，都是自由列，表示u可以任意取，则二维平面上所有向量都是特征向量。
用一组基来表示，特征空间的基为
由这组基的线性组合所组成的向量都是A的特征向量。
由这里可以看出当特征值为多重特征值的时候，它的特征向量有可能将不再是一根直线了，它有可能是一个高维空间
但是不能保证多重特征值的特征向量一定是一个高维空间，例如
根据它是一个上三角矩阵，我们可知它的特征值就为3，这也是一个重数为2的多重特征值
代入
得
得到的(1,0)是特征空间的一组基，它是一维的
如果矩阵A的某个特征值的重数=k，则对应的特征空间的维度<=k.我们把特征值的重数定义为代数重数，把特征空间的维度定义为几何重数。
则几何重数不大于代数重数
若重数为1，则很简单，其特征空间维度为1.也就是说简单特征值的特征空间一定是一维的，而多重特征值的特征空间最多可以达到特征值的重数。
numpy中求特征值和特征向量
之前我们已经知道求解特征值和特征向量的方法
它实际上就是求解关于
的一元n次多项式方程，有n个解。由于一元n次多项式方程的底层算法较为复杂，它也不是一个线性代数的算法，它已经是数值分析领域的一个算法，所以这里直接使用numpy来求解矩阵的特征值和特征向量
numpy是Python所独有的，所以这里只有Python代码
import numpy as np from numpy.linalg import eig if __name__ == "__main__": A1 = np.array([[4, -2], [1, 1]]) eigenvalues1, eigenvectors1 = eig(A1) print(eigenvalues1) print(eigenvectors1) print() # 关于y=x翻转 A2 = np.array([[0, 1], [1, 0]]) eigenvalues2, eigenvectors2 = eig(A2) print(eigenvalues2) print(eigenvectors2) print() # 旋转270度 A3 = np.array([[0, -1], [1, 0]]) eigenvalues3, eigenvectors3 = eig(A3) print(eigenvalues3) print(eigenvectors3) print() # 单位矩阵 A4 = np.array([[1, 0], [0, 1]]) eigenvalues4, eigenvectors4 = eig(A4) print(eigenvalues4) print(eigenvectors4) print() # 几何重数为1 A5 = np.array([[3, 1], [0, 3]]) eigenvalues5, eigenvectors5 = eig(A5) print(eigenvalues5) print(eigenvectors5) print()
运行结果(带注释)
[3. 2.] //这里的3，2都是特征值
[[0.89442719 0.70710678]
[0.4472136 0.70710678]]
//这里看特征向量是需要看该矩阵的列向量的，对应着相同位置的特征值的特征空间的一组基
[ 1. -1.]
[[ 0.70710678 -0.70710678]
[ 0.70710678 0.70710678]]
[0.+1.j 0.-1.j] //这里是两个复数
[[0.70710678+0.j 0.70710678-0.j ]
[0. -0.70710678j 0. +0.70710678j]]
//这里特征向量也是两个复数
[1. 1.] //这里具有相同的特征值，重数为2
[[1. 0.]
[0. 1.]]
//特征向量的矩阵的列向量共同构成了二维空间的一组基
[3. 3.] //这里也具有相同的特征值，重数为2
[[ 1.00000000e+00 -1.00000000e+00]
[ 0.00000000e+00 6.66133815e-16]]
//这两个特征向量是线性相关的
矩阵相似型
如果矩阵A、B满足
，则称A和B相似
这其实是以不同的视角观察相同的内容，我们可以把P想象成一个坐标系，A是P坐标系下的变换，而B是标准坐标系下的变换。则A变换是在P坐标系下观察的B变换。它们本质是同一个变换，只是坐标系不同而已。
这里x是P坐标系下的一个向量
该向量x在P坐标系下经过A变换，变成了
。如果我们把x向量转成标准坐标系下的向量，则有
，我们将标准坐标系下的
进行B变换，变成了
，自然
，最后我们将标准坐标系变换后的转成P坐标系下的向量就有P的逆右乘，就有
，这个就等于。此时A、B其实是同一个变换，只不过是坐标系不同，则有A、B相似。
更多的时候，我们都是从标准坐标系下去观察的，所以将进行变形，先右乘P的逆，再左乘P，就有了
，而这样在已知P坐标系和P坐标系下的A变换的时候，我们就知道了向量在标准坐标系中的是如何变换的，而它的变换矩阵就是B。
因为A和B本质是同一个变换，则A和B的特征方程相同，特征值相同。
证明：
得证
矩阵对角化
这里的A是一个标准坐标系中的变换矩阵，而D是P坐标系中的对角矩阵，也是A变换在P坐标系下的变换。D的样式如下
如果A有n个线性无关的特征向量，则A和某个D相似。称A可以被对角化
就是A矩阵的对角化。而D的对角元素就是A的各个特征值，P是由A的特征向量排成列向量所组成的矩阵，也称为特征向量矩阵
这里的u就是A的特征向量
证明： 我们将等式两边同时右乘P,得到
AP = PD，则有
因为左边等于右边，得证
之前我们说过，如果矩阵A含有两个不同的特征值，则它们对应的特征向量线性无关。推广之
如果A有n个不同的特征值，则A可以被对角化。这里需要注意的是，如果A没有n个不相同的特征值，A不一定不能被对角化。这是因为如果A有多重特征值，有可能A的特征向量是特征空间的一组基，它们彼此之间依然是线性无关的。但如果A有多重特征值时，而A的特征向量是线性相关的，则A不能被对角化。
则如果A有n个线性无关的特征向量，A就可以被对角化。
由于A的对角化依然要用到特征值和特征向量，所以我们这里就不再使用Java来编写了，只使用Python，而我们之前自己实现的向量、矩阵、线性系统的代码却依然可以派上用场。
import numpy as np from numpy.linalg import eig, inv from playLA.Matrix import Matrix def diagonalize(A): # 矩阵的对角化 # 必须是一个二维数组 assert A.ndim == 2 assert A.shape[0] == A.shape[1] eigenvalues, eigenvectors = eig(A) # 特征向量 P = eigenvectors assert Matrix(P.tolist()).rank() == A.shape[0], "特征向量彼此必须线性无关" # 使用特征值创建一个对角矩阵 D = np.diag(eigenvalues) # 特征向量的逆 Pinv = inv(P) return P, D, Pinv if __name__ == "__main__": A1 = np.array([[4, -2], [1, 1]]) P1, D1, Pinv1 = diagonalize(A1) print(P1) print(D1) print(Pinv1) print(P1.dot(D1).dot(Pinv1)) A2 = np.array([[3, 1], [0, 3]]) P2, D2, Pinv2 = diagonalize(A2) print(P2) print(D2) print(Pinv2) print(P2.dot(D2).dot(Pinv2))
运行结果
[[0.89442719 0.70710678]
[0.4472136 0.70710678]]
[[3. 0.]
[0. 2.]]
[[ 2.23606798 -2.23606798]
[-1.41421356 2.82842712]]
[[ 4. -2.]
[ 1. 1.]]
Traceback (most recent call last):
File "/Users/admin/Downloads/PycharmProjects/matrix/main_numpy_diag.py", line 32, in <module>
P2, D2, Pinv2 = diagonalize(A2)
File "/Users/admin/Downloads/PycharmProjects/matrix/main_numpy_diag.py", line 13, in diagonalize
assert Matrix(P.tolist()).rank() == A.shape[0], "特征向量彼此必须线性无关"
AssertionError: 特征向量彼此必须线性无关
矩阵对角化的应用
求解矩阵的幂
由可知
对于一个对角矩阵而言，
则有
则有
这整个的过程，我们只进行了两次矩阵的乘法运算，其他的过程都是对一个数字进行幂运算。所以这种运算是比一个矩阵要进行多次乘法运算要简单很多的。这个计算的时间复杂度是要大大降低的。
大量问题的形式是：
这是一个动态系统
，它都是随着时间的推动，而不断发生变化的。比如经济系统、粒子系统、生态系统、随机系统。
由于D的对角数据都是A的各个特征值，所以特征值反应动态系统的各个分量的变化速率。
对称矩阵与矩阵的SVD分解
完美的对称矩阵
什么是对称矩阵，形如下面形式的矩阵就是一个对称矩阵，它同样也是一个方阵
我们可以看到主对角线之外的所有元素，它们都是根据对角线对称的
所以对称矩阵和它的转置矩阵相等
为什么说对称矩阵是完美的，它们有下面几点
- 对称矩阵的特征值一定是实数。
- 对称矩阵的多重特征值，其对应的特征空间的维度一定等于重数。即对称矩阵的几何重数等于代数重数。
- 对称矩阵一定有n个线性无关的特征向量。
- 对称矩阵一定可以被对角化。
正交对角化
之前我们在讲特征值和特征向量的时候讲过这样一个例子。很显然，在这个例子中的两个特征向量所构成的坐标系是很扭曲的。
很明显，对于坐标系来说，我们更喜欢正交坐标系，那什么样的特征向量会构成正交坐标系呢？
再之后我们又讲解了两个变换，一个是投影变换，一个是翻转变换，它们的特征向量都是正交的。
对称矩阵的所有不同的特征值对应的特征向量互相垂直。
现在我们来证明这个性质，假设矩阵A的两个特征向量
对应不同的特征值
要证明两个向量垂直，则只需要证明它们点乘的结果为0，即
我们这里的向量都是列向量，也可以看成是一个n行1列的矩阵，所以把
向量看成矩阵的时候，用
来表示，
是一个1行n列的矩阵。这里的
依然是一个n行1列的矩阵。所以有
。由于
是A的一个特征值，则
，再根据转置的性质
，由于A是一个对称矩阵，则
，再
也是A的一个特征值，则有
，由于本身就是向量，最后可得
。
由于不等于，则，得证
在有相同的特征值的情况下，几何重数等于代数重数，k重相同特征值的特征空间为k维。我们依然可以找到n个互相垂直的特征向量。
对称矩阵一定可以被对角化
，由于对称矩阵的特征向量是正交的，所以P是一个正交矩阵，我们又可以写成
，这里Q是一个标准正交矩阵，它的基的模为1。由于标准正交矩阵的逆就是它的转置，所以有
，这又叫做正交对角化。
则有对称矩阵一定可以被正交对角化。
如果A可以被正交对角化，则A一定是对称矩阵。
证明：
这里D是一个对角矩阵，所以D的转置一定为D，得证。
A是对称矩阵
A可以被正交对角化。这是一个等价命题。
奇异值
之前很长一段内容，都是基于方阵的，比如行列式、特征值、特征向量、矩阵相似型、对角化、对称矩阵、正交对角化。
然而大部分时候我们的真实数据都很难是一个方阵。那么对于非方阵，我们该怎么处理呢？
若A是一个m*n的矩阵，则
是一个n*n的方阵，且对称。
是一个n*m的矩阵点乘一个m*n的矩阵。
它的第i行第j列元素为：
第i行点乘A第j列，其实就是A第i列点乘A第j列
它的第j行第i列元素为：第j行点乘A第i列，其实就是A第j列点乘A第i列
由于向量的点乘满足乘法交换律，所以A第i列点乘A第j列=A第j列点乘A第i列，换句话说的第i行第j列元素=第j行第i列元素。则是对称矩阵。
可以被正交对角化，拥有n个实数特征值
(有可能有重复)，n个互相垂直的标准特征向量
。
我们来看一下这样一段式子的推导
A(非方阵)点乘(对称矩阵)的标准特征向量的模的平方，由于两个向量的点乘是一个标量，所以
；
向量可以看成是一个n行1列的矩阵，所以把向量看成矩阵的时候，用
来表示，
是一个1行n列的矩阵，所以有
；再根据转置的性质
；根据特征值和特征向量的定义
，这里
是的一个特征值，是对应的一个特征向量；本身就是，所以
；由于是标准特征向量，所以它的模等于1，则
。
最后可得
，由于这里是一个平方操作，由此可见的特征值一定>=0。则一定有
，我们将这个值称为奇异值(Singular Value),奇异值就是
的模。
奇异值的几何含义
之前我们一直在讲一个很奇怪的向量式子，A是原矩阵，是的特征值对应的标准特征向量，而的模是的特征值开方，即奇异值。那么对于A来说到底是什么呢？
其实所组成的集合
是A的列空间的一组正交基，不过此时
，依然是的特征值。
我们先来证明正交性：我们假设、
是这组集合中的两个向量
这个推导过程，前面已经说过很多次了，这里就不多说了，由于
、
是的标准特征向量，所以
=0，最终就等于0。这里它的正交性得证。
现在我们要证明是A的列空间的一组基。
由于是标准特征向量，则
是n维空间的一组基，而且是一组标准正交基。那么在该n维空间中，任意一个向量都可以表示为这组基的一个线性组合
，A的列空间是该n维空间的一个子空间。我们假设
是A的列空间的一个向量，它是n维空间的向量
转换到A的列空间的一个向量，则A就是这个转换矩阵。
这个时候，我们可以看到A的列空间的向量已经可以写成的一个线性组合了。但
并不一定是线性无关的，这是因为其中有一些可能为0。很明显只有的特征值为0的时候，才可能为0。如果我们把特征值为0的抛弃的话，那么剩下的一定是线性无关的，因为上面已经证明了是两两正交、两两垂直的。所以最终转化为
结论就是，如果A有r个不为零的奇异值：
是A的列空间的一组正交基。A的列空间的维度为r，rank(A) = r
是A的列空间的一组标准正交基，通常把奇异值按照从大到小排序
是A的列空间的一组标准正交基，这里
，我们从这里也可以看到奇异值
不能为0。
矩阵的SVD分解
矩阵的SVD分解 - Singular Value Decomposition，又叫做矩阵的奇异值分解
它和矩阵的LU分解或者矩阵的QR分解不同，这两种分解方式都是把矩阵分解成两部分(LU是分解成L单位下三角矩阵,U上三角矩阵；QR分解是分解成Q标准正交矩阵,R上三角矩阵)，SVD分解的本意就是直译——奇异值分解。它对任意形状的矩阵都适用，而LU或者QR分解都会有一定的限制。
矩阵的SVD分解是分解成这么几个部分。如果A是m*n的矩阵，这里U是m*m的方阵；
是m*n的矩阵，形状与A相同；V是n*n的方阵。
V是的特征向量矩阵进行标准化，V是一个n*n标准正交矩阵，它里面的列向量都是一个标准特征向量
这里u1、u2...um就是之前所说的其中
，虽然u是A的列空间的一组标准正交基，而A的列空间的秩为r，所以这里的u只有r个，而r是小于m且小于n的。所以u肯定不是整个m维空间的一组基，要凑满m个u，则需要使用格拉姆-施密特过程进行扩展，在已知有r个正交向量的基础上，去求后面的更高维度的正交向量。
是A的所有奇异值组成的对角矩阵，由于我们知道奇异值只有r个，所以在填充m*n的矩阵的过程中，我们只需要在其他部分填充0就可以了。
现在我们来证明
证明： 在等式的两边都右乘V，V是一个标准正交矩阵，所以
其实就是V的逆，则有
，我们假设
是V中的列向量，也是的标准特征向量,则AV可以进行如下转化
这里奇异值只有r个非零的
，如果为0，则后面的列向量就一定为零向量。
可进行如下变换
虽然U进行了格拉姆-施密特过程的扩展，但是只有r个非零的奇异值，所以最终结果，除前r个有值的列向量，后面依然为零向量。
由于AV和都可以化成相同的形式，所以得证。
如果给定一个任意矩阵A，对其进行SVD分解，需要经过以下步骤
- 先求解的特征值和特征向量
- 把所有非零的特征值开根后(奇异值)，在放入m*n的的对角矩阵中，不够的填充0。在这里我们通常把奇异值按照从大到小进行排序，而也叫做奇异值矩阵。
- 将求得的特征向量标准化(每一个特征向量除以它的模)后，按照奇异值大小的顺序(每一个特征值对应一个特征向量)得到n*n的标准正交矩阵V。
- U的求解比较麻烦一点，，其中是在V中排好顺序的标准特征向量，在经过格拉姆-施密特过程扩展后得到的m*m的标准正交矩阵U
这样我们就完成了整个矩阵A的SVD分解
因为这个过程比较复杂，我们就不从底层来实现它的算法了，直接使用scipy来写对一个矩阵的分解代码
import numpy as np from scipy.linalg import svd if __name__ == "__main__": A = np.array([[1, 2], [3, 4], [5, 6], [7, 8]]) U, s, VT = svd(A) print(U) # 这里的s只是按照从大到小排列的奇异值 print(s) print(VT) # 构建一个零矩阵，将奇异值填充进该零矩阵 Sigma = np.zeros(A.shape) for i in range(len(s)): Sigma[i][i] = s[i] print(Sigma) print(U.dot(Sigma).dot(VT))
运行结果
[[-0.15248323 -0.82264747 -0.39450102 -0.37995913]
[-0.34991837 -0.42137529 0.24279655 0.80065588]
[-0.54735351 -0.0201031 0.69790998 -0.46143436]
[-0.74478865 0.38116908 -0.5462055 0.04073761]]
[14.2690955 0.62682823]
[[-0.64142303 -0.7671874 ]
[ 0.7671874 -0.64142303]]
[[14.2690955 0. ]
[ 0. 0.62682823]
[ 0. 0. ]
[ 0. 0. ]]
[[1. 2.]
[3. 4.]
[5. 6.]
[7. 8.]]
SVD分解的应用
矩阵的SVD分解的应用非常广泛，在各个领域，只要是用线性的地方，几乎都离不开SVD分解的应用。
一般说我们可以将矩阵A看成是一个变换。如果A是m*n的矩阵，当然一般我们说矩阵是一个变换的时候，A是一个方阵的情况居多。
关于变换这里可以看一下直观的表述https://m.bilibili.com/video/BV16A411T7zX
由于A有n列，它只能乘以一个n维的向量，结果向量是一个m维的向量，即把一个n维的向量转换成一个m维的向量。
V是n维空间的一个标准正交基，则n维空间内的任何向量x都可以写成V的线性相关的形式
那么矩阵A点乘该向量x就有(这里V是一个标准正交矩阵，所以它的转置就等于它的逆)
该式子的解读为：原来的x在正交的坐标系V下，相应的坐标为k1、k2...kn。经过了A转换后，就变成了在正交坐标系U下看，相应的每一个维度的坐标值，不是k1、k2...kn了，而是要拉伸奇异值那么多倍。这里我们可以记住标准正交矩阵的变换的作用是旋转，而对角矩阵的变换是拉伸。
现在我们不把矩阵A看成一个变换了，而是把矩阵A看成一片数据。这个应用在机器学习领域和数据分析时更常使用。
A是一个m*n的二维数据
最后的u、v不再是一个向量，而是一个二维的矩阵，u是一个m*1的矩阵，v是一个n*1的矩阵，所以v的转置就是一个1*n的矩阵。每一个u乘以v的转置得到的结果都是一个m*n的矩阵。
由于我们在进行奇异值分解的时候已经将奇异值按照从大到小排列了，所以我们可以理解成这些奇异值是一些权值。所以这里是一系列m*n的矩阵相加，每一个m*n的矩阵有一个不同的权值，它们是从大到小来组合的。可能排在后面的奇异值很小，我们可以将这些很小的奇异值权值的矩阵去除，相应的A的改动并不是很大，这样一来，我们就对A这些数据进行了压缩、去噪、降维。
对于这个应用可以非常好的直接应用于图像领域，如果A是一个高为m，宽为n的图像的话，那么我们对它进行奇异值分解之后，再把奇异值比较小的项扔掉的话，剩下的这些项依然可以很好的表达这个图像的语义，所谓语义是指我们人的肉眼去看这个图像，我们还能看懂这个图像的意思是什么，但是我们表征这个图像的数据量，存储空间都变低了。这是SVD分解的一个非常好的应用。当我们扔的项越来越多，我们可以看到图像在一步一步的变模糊。
SVD分解还在很多领域有着重要的作用，比如定义伪逆。之前我们说过逆矩阵只有方阵才有，但是对于非方阵来说，我们可以定义一个伪逆，伪逆跟逆非常的像，在非方阵的计算上，我们可以用伪逆来代替逆。
除此之外，在推荐系统、自然语言处理、搜索引擎等跟计算机高度相关的领域，其中的很多相关理论都和SVD分解有关。
