总结机器学习中的线性代数和矩阵论基础

• 线性代数一直迷的不行……so
• 这个博客主要参考MIT18.06和MIT18.065两门课的内容，构建并梳理各个知识点的联系。一定要学过一遍线性代数再看，这里着重记录“感觉”，信息密度大，多停多想，不要贪快。遇到不知道且没有解释的的名词可以查下百科。对于没有证明（或许提供了思路）的定理可以当作练习
• MIT18.06讲了线性代数，18.065讲了线性代数在数据科学中的诸多用法。本文计算机视觉和机器学习有关的内容还参考了很多Prince的Computer vision: models, learning, and inference的附录C，鼓励大家看原书，内容非常清楚。本文同时穿插了一些信号、图像方面的知识点
• 这里主要介绍为应用服务的理论工具。部分应用需要少量机器学习基础。如果对某部分应用不感兴趣，可以连带相关章节略去不看
• 作者水平非常有限，如有错误或任何批评，恳请读者一定要指出来

---

---

线性变换

（逻辑上应该先讲线性空间，但是会不明意义，所以先从几何角度建立线性变换的感觉，强烈建议去看3B1B 线性代数的本质理解线性变换，无前置基础）
这里给出$T$满足线性变换的条件
$$\begin{array}{rl}T(\mathbf{v}+\mathbf{w})& =T(\mathbf{v})+T(\mathbf{w})\\ T(c\mathbf{v})& =cT(\mathbf{v})\end{array}$$

• 只需要确定$T$对所有基向量的影响，就能完全掌握这个变换
• 两个定理（对几何直觉构建帮助不大，可以先略去）
  • 对于同一个元素，在基$\alpha$下坐标为$x$，在基$\beta$下坐标为$x^{\prime }$，如果$\beta =\alpha P$，那么$x=P{x}^{\prime }$.
    证明：因为$\alpha x=\beta x^{\prime }$. $P$也叫过渡矩阵或基变换矩阵
  • 如果$\beta =\alpha P$，线性变换$T$在两组基下的矩阵依次为$A,B$，即$T(\alpha )=\alpha A,T(\beta )=\beta B$，则$B=P^{-1}AP$.
    证明：$\beta B=T(\beta )=T(\alpha P)=T(\alpha )P=\alpha AP=\beta P^{-1}\alpha P$，因为$\beta$满秩，所以$B=P^{-1}AP$
    即$A$和$B$是两组基下的同一个变换，它们一定相似，且$P$就是相似变换矩阵
• $n$维线性空间和${\mathbb{R}}^n$都同构，即维数相等的线性空间都同构，线性空间的结构完全由维数决定

---

线性空间

数域$F$中对加法和数乘封闭的非空集合

• Ax=b的解空间不构成线性子空间，注意$b$的存在，这种情况不满足加法和乘法的封闭法则

子空间的加和交

对于两个子空间$S$和$U$
$$dim(S)+dim(U)=dim(S+U)+dim(S\cap U)$$

矩阵空间

以矩阵为元素的线性空间

4个基本子空间

1. 列空间Column Space $C(A)$: 所有$Ax$
2. 行空间Raw Space $C(A^T)$
3. 零空间Null space of a matrix $N(A)$: $Ax=0$ for all $x$
4. 左零空间Left Null space $N(A^T)$: $x^{T}A=0$ for all $x$

关系如图所示（图片摘自网络），这张图非常重要，理解了这张图就理解了很多抽象的概念。

注意$Ax=0$就已经揭示了$C(A^T)$和$N(A)$是正交的。

A=CR分解，行秩等于列秩

对于$A\in {\mathbb{R}}^{m\times n}$，如果列秩为$r$，$r$维的列空间中的无关列向量拼成矩阵$C=[c_{v_1},\cdots ,c_{v_r}]$，$A$的每一列可以由$C$中的列线性表达，则有$A=CR$，其中$C\in {\mathbb{R}}^{m\times r},R\in {\mathbb{R}}^{r\times n}$. 注意该式也可看作$A$的每一行由$R$中的$r$行进行线性组合得到，所以行空间维度不超过$r$，也即行秩不超过列秩；类似方法可证列秩不超过行秩。即证明$行秩=列秩=r$。所以$A=CR$分解中$C$是列无关的，$R$也是行无关的

• 从$C$和$R$两个角度看$A$的构成，也说明了秩的乘法性质： r a n k [ A B ] ⩽ min ⁡ { r a n k [ A ] , r a n k [ B ] } rank[AB]\leqslant \min\{rank[A], rank[B] \} rank[AB]⩽min{rank[A],rank[B]}
• $Ax=b$有解即$b$在$A$的列空间中

几个性质

• 当 $B$ 可逆时，$N(BA)=N(A)$
• 行空间和列空间维度都为$r$
• 行空间与零空间正交，所以零空间$N(A)$维度为$n-r$，也即线性方程组$Ax=0$解的构成维度；同理左零空间$N(A^T)$与列空间正交，维度为$m-r$

---

子空间投影

• 对于向量$a$和$b$，$b$在$a$上投影向量p为
  $$a\frac{a^{T}b}{a^{T}a}=\frac{a{a}^T}{a^{T}a}b=Pb\text{\hspace{1em}(1)}$$
  其中$P$看作是投影矩阵，可以自行证明。

最小二乘问题

• 求解$Ax=b$时，如果无解，$b$不在$C(A)$当中，那么改求$A\hat{x}=p$，其中$p$是$b$在$C(A)$上的投影。注意到$b-A\hat{x}$垂直于$C(A)$，在Left Null Space。所以$A^T(b-A\hat{x})=0$，即$A^{T}A\hat{x}=A^{T}b$
• 所以$\hat{x}=(A^{T}A{)}^{-1}{A}^{T}b$，$p=A\hat{x}=A(A^{T}A{)}^{-1}{A}^{T}b=Pb$，注意该式和式(1)的相似性. $P=A(A^{T}A{)}^{-1}{A}^T$称之为投影矩阵. $\hat{x}$也称之为$A$的左逆
• 注意$P$的两个性质，这两个性质也是判定投影矩阵的充分条件：$P^T=P,P^2=P$，（注意$P^2=P$不足以判定投影矩阵，例如矩阵$(1,1;0,0)$）. 因为$P^2=P$，特征值只能为0或1. 此外，$I-P$也是投影矩阵
• 因为$b=p+e$. ($e$是投影误差)，而$p=Pb$，所以$e=(I-P)b$. 且有 $(I-P{)}^2=I-P,(I-P{)}^T=I-P$. 注意和$P$的类似性
• $A^{T}A$ 可逆当且仅当$A$的所有列独立. 一个可靠的证明方法是两者的Null Space相同. 这里似乎$r(A^{T}A)=r(A)$是成立的

分析角度看最小二乘问题

$$Ax=b$$
以最小二乘形式看待
$$\begin{array}{rl}\hat{x}& =\underset{x}{\mathrm{argmin}}[(Ax-b{)}^T(Ax-b)]\\ & =\underset{x}{\mathrm{argmin}}[x^T{A}^{T}Ax-2x^T{A}^{T}b]\end{array}$$
标量对$x$求导得
$$x=(A^{T}A{)}^{-1}{A}^{T}b$$
结果和投影方法是一样的。
其中$A$的行数如果比$x$少，则奇异

应用：线性回归

这一节摘自PRML P143，为了与PRML配图一致，符号标注不符合本文的惯用标注。

几何解释

记训练集标注为$\mathbf{t}=({\mathbf{t}}_1,...,{\mathbf{t}}_{\mathbf{N}}{)}^{\mathbf{T}}$，并构成标注空间${\mathbb{R}}^N$，$\mathcal{S}$是能在训练集的标注空间中用广义线性回归张成的超平面，也即训练集的列空间

• 这里线性回归的基可以是带核$\phi (X)$的，实际上带核的仍然是张成超平面，而不是曲面，超平面的第$i$个基由${\phi }_i(X)$决定，${\phi }_i$表示第$i$个特征，$X$表示所以的N个数据。
  这样线性回归是求了标注空间中训练集所在位置在超平面上的投影，垂直距离即为误差向量，范数为最小二乘的结果。（图片摘自Bishop的Pattern Recognition and Machine Learning）
  
• 注意误差向量垂直于$\mathcal{S}$，所以如果截距项全1向量在$S$当中，那么必有误差向量自身之和为0. 也即线性回归拟合均值一定等于标注均值

多重共线性缺陷

之前只知道多重共线性不好，到底哪里不好一直说不清楚。这里把它讲清楚。
多重共线性的灾难在于参数值爆炸。
我们记（可选经过核变换后）训练集为$\Phi \in {\mathbb{R}}^{N\times M}$，其中$M$是特征维度。$rank[\Phi ]<M$时，即产生了多重共线性问题，也即特征之间线性相关。因为$rank[\Phi ]=rank[{\Phi }^T\Phi ]=rank[\Phi {\Phi }^T]$，（注：方法为证明$\Phi x=0$与${\Phi }^T\Phi x=0$同解）。所以如何判断$rank[\Phi ]$与$M$的关系，只需要计算${\Phi }^T\Phi$是否奇异。

• 此外，如果${\Phi }^T\Phi$接近奇异，即行列式很小，那么线性回归的参数闭式解$({\Phi }^T\Phi {)}^{-1}{\Phi }^T\mathbf{t}$会非常大。关于行列式和逆的定义参考下一节。
• 从几何角度解释，即两个基向量方向非常近，那么为了表达出与这两个基向量几乎垂直的方向上的位置，这两个向量需要不断抵消， 系数会增长非常快！

这里如果无法求逆，则可以求伪逆，后文伪逆部分会继续讨论最小二乘求逆问题。

---

行列式和逆

最基本的性质

1. $det[I]=1$
2. 交换两行，行列式取反
3. a) 某一行乘$k$倍，行列式乘$k$倍.
   b) $$det\left[\begin{array}{c}⋮\\ a_{i_1}+a_{i_2}\\ ⋮\end{array}\right]=det\left[\begin{array}{c}⋮\\ a_{i_1}\\ ⋮\end{array}\right]+det\left[\begin{array}{c}⋮\\ a_{i_2}\\ ⋮\end{array}\right]$$

用这三条性质能推出剩下一大堆东西，包括基本算法和其他性质，都能推出来，列几个常用的性质

• $det[A^T]=det[A]$
• $det[AB]=det[A]det[B]$，所以$det[A^{-1}]=1/det[A]$
• 利用性质2和3可以证明$$det\left[\begin{array}{c}⋮\\ a_{i_1}\\ ⋮\\ a_{i_2}\\ ⋮\end{array}\right]=det\left[\begin{array}{c}⋮\\ a_{i_1}\\ ⋮\\ a_{i_1}+a_{i_2}\\ ⋮\end{array}\right]$$
  由此进一步可证不满秩矩阵的行列式为0。

行列式表示矩阵组成的体积

可以证明体积也满足上述三个重要性质。于是得证。（3b不是很好证）
建议参考3B1B——线性代数的本质

• 三阶行列式即混合积：$a\times (b\cdot c)=det[[a,b,c]]$

行列式算法

对于方阵$A\in {\mathbb{R}}^{n\times n}$，一种写法：
$$det[A]=\sum _{n!\text{ terms}}\pm a_{1\alpha }{a}_{2\beta }\cdots a_{n\omega }$$
其中
$$(\alpha ,\beta ,\gamma ,\cdots ,\omega )=\text{perm of }(1,2,\cdots ,n)$$
另一种写法：
$$det[A]=a_{11}{C}_{11}+a_{12}{C}_{12}+\cdots +a_{1n}{C}_{1n}$$
这里只以第一行为例，也可以取其他行。$C_{ij}$是第$i$行第$j$列的代数余子式cofactor，注意有正有负。

逆矩阵

$$A{C}^T=det[A]I$$
注意$C$是代数余子式形成的矩阵，$C^T$称之为伴随矩阵Adjugate Matrix. 正确性可以自行验证.

• 伴随矩阵和逆矩阵只差一个系数。然而，伴随矩阵对不可逆的矩阵也有定义，并且不需要用到除法（摘自百度百科）

克拉默法则

$$Ax=b$$
有唯一解时
$$x=A^{-1}b=\frac{1}{det[A]}{C}^{T}b$$
注意$C^{T}b$中代数余子式的含义，可以得到
$$\begin{gathered} x_i=\frac{det[B_i]}{det[A]} \\ {B}_i=\left[\begin{array}{c}{A}_{\cdot 1},\cdots ,A_{\cdot (i-1)},b,A_{\cdot (i+1)},\cdots ,A_{\cdot n}\end{array}\right] \end{gathered}$$
几何解释参考3B1B——线性代数的本质

特殊等式

• $det[I_n+u{v}^T]=det[1+u^{T}v]$。它的扩展形式为$det[I_n+A{B}^T]=det[I_m+A^{T}B]$，其中$A,B\in {\mathbb{R}}^{n\times m}$. 只需证明左右特征值连乘相等即可。注意$A^{T}B$与$A{B}^T$特征值相等（Sylvester定理）

---

正交矩阵

$Q$是方阵，且每一列的范数为1，不同列正交。故有$Q^{T}Q=I,Q^{-1}=Q^T$
例子：哈达玛矩阵，参考博客

• 实际上每一行范数也为1，也与其他行正交

• 复空间的正交矩阵叫做酉矩阵unitary，即满足${\bar{A}}^{T}A=I$
• 正交矩阵特征值绝对值都为1，因为${\overline{(Qx)}}^{T}Qx=\bar{\lambda }\lambda {\bar{x}}^{T}x$，${\bar{x}}^{T}x>0$，所以$\Vert \lambda \Vert =1$. 注意旋转矩阵这一类，考虑其旋转和特征向量的物理意义，除了$I$外，无法找出实特征向量和特征值

• 注意：如果$Q^{T}Q=I$，当$Q$不为方阵的时候，不一定能推出$Q{Q}^T=I$

旋转矩阵与正交变换

正交矩阵行列式绝对值为1，所以线性变换体积不变。左乘的线性变换效果是可看作是一种反演（改变手性）+旋转，这个在SVD中会有感受。行列式为1时，为旋转矩阵，为-1则包含了反演
正交矩阵的线性变换叫正交变换，正交变换有以下好性质：

• 距离不变
• 夹角不变
• 上述两条也能推出内积不变，即$x_1^T{x}_2=(A{x}_1{)}^T(A{x}_2)$

反射矩阵

即正交对称阵$r^T=r$.

• 注意这意味着特征值都为实数±1，且特征向量都存在。这是一个不带旋转的变换，如果带旋转的话不可能特征值都为±1

A=QR与Gram-schmitt正交化

A=QR

$Q$是正交矩阵，其中$R$是上三角矩阵.
这个上三角矩阵成立的原因可以思考一下Gram-schmitt正交化的过程，对$Q$的列向量$q_1,q_2,\cdots ,q_n$一个一个进行正交。每一次正交只依赖于前面的列向量

数值计算：Gram-schmitt中的过小分母问题

• 在Gram-schmitt时，会有分母$q_i^T{q}_i$的出现，如果两个向量离得很近，这个值有很小，会在数值计算中出现问题。一种解决方法是：假定上一个被正交化的是$q_i$，在计算下一个要被正交化的向量$q_j$时，从没有正交化的向量当中把$a_{\cdot }-\frac{q_i{q}_i^T}{q_i^T{q}_i}{a}_{\cdot }$都计算一遍，选范数最大的。注意这一步计算不引入额外开销，因为这个计算是必须要经历的。

应用：信号处理中的变换

信号处理中有一大类变换是线性变换，其中的部分变换是正交变换，可以看作是对原始信号$x$施加了一组新的标准正交基，例如$y={\bar{A}}^{T}x$，新的基为正交矩阵$A$中的各列，而且$A$和$A^T$互为正反变换。实际例子包括：

• 傅里叶变换，记$\omega =e^{2\pi /n}$，则
  $$A=\frac{1}{\sqrt{n}}\left[\begin{array}{ccccc}1& 1& 1& \cdots & 1\\ 1& \omega & {\omega }^2& \cdots & {\omega }^{n-1}\\ 1& {\omega }^2& {\omega }^4& \cdots & {\omega }^{2(n-1)}\\ & & & ⋮& \\ 1& {\omega }^{n-1}& {\omega }^{2(n-1)}& \cdots & {\omega }^{(n-1)(n-1)}\end{array}\right]$$
  注意$A$是复数矩阵，虽然对称，但不是共轭对称，所以不是埃尔米特矩阵，但它是酉矩阵，反变换为${\bar{A}}^T$
• 小波变换。参考我的另一篇博客：小波变换——公式整理和简单介绍

---

特征值与特征向量

对于方阵$A\in {\mathbb{R}}^{n\times n}$
$$Av=\lambda v$$
$v$为特征向量，$\lambda$为对应特征值

• $\sum \lambda =tr[A]$，$tr$表示矩阵的迹trace：对角线元素之和
• $\prod \lambda =det[A]$
• 对称矩阵特征值一定为实数，反对称矩阵则为零或纯虚数；其他矩阵介于对称与反对称之间
• $A$和$A^T$的特征值相同，但特征向量不一定，因为Left Null Space not equals to Null Space.
• Sylvester定理：$AB$和$BA$的非零特征值相等
  
  （摘自程云鹏《矩阵论》，笔记来自zhenliang师兄。感谢！）

特征值分解

$$A=V\Lambda V^{-1}$$
其中，$V$是特征向量列拼接组成的矩阵，$\Lambda$是对应特征值。这里需要假设$n$个线性无关的特征向量存在，即$A$可逆。

• $A^k=V{\Lambda }^k{V}^{-1}$，所以$A^k$特征向量不变，特征值乘$k$次。
• $A^{k}u=(V{\Lambda }^k)(V^{-1}u)=({\lambda }_1{v}_1,\cdots ,{\lambda }_n{v}_n)(c_1,\cdots ,c_n{)}^T={\sum }_{i=1}^n{c}_i{\lambda }_i^k{v}^i$
• 若$A$的所有特征值都不同，那么一定有$n$个线性无关特征向量（可对角化）；但特征值相同时（代数重数大于1），需要再研究，线性无关的特征向量可能不够，即几何重数低于代数重数
• 特征向量组唯一当且仅当所有的特征值都是唯一的
• 给定特征值和特征向量唯一确定一个矩阵，可以从几何角度感知，参考3B1B线性代数本质

实对称矩阵的特征值分解

谱定理Spectral Theorem

如果$S$是实对称矩阵，$S$一定可以分解成$n$个线性无关的正交实特征向量和实特征值（不一定都不同）
$$S=Q\Lambda Q^T$$

正规矩阵Normal Matrix

实际上当$A$满足$A{\bar{A}}^T={\bar{A}}^{T}A$时，特征向量都正交，这样的矩阵叫正规矩阵normal matrix，其满足
$$A=Q\Lambda Q^T$$
正规矩阵的典型特例：

• 实对称矩阵，$\lambda \in \mathbb{R}$
• 正交矩阵，$|\lambda |=1$
• $A={\bar{A}}^T$，其中$\bar{A}$是$A$的共轭矩阵。满足上述条件的矩阵也叫埃尔米特矩阵Hermitian Matrix或自共轭矩阵
• $A$是反对称矩阵，特征向量也正交，注意特征值一定为0或纯虚数

主轴定理

$Q$是$A$的特征向量组成的正交矩阵，作用在单位圆上后，可以将$A$看作是沿$v_i$方向扩展${\lambda }_i$倍。如图所示（摘自GoodFellow《深度学习》）

解释：$Au=Q\Lambda Q^{T}u$，从右往左看，这里正交矩阵$Q^T$对单位圆$u$的旋转作用看不出来，$\Lambda$把坐标轴拉长，$Q$再把拉长后的椭圆旋转

• $\Lambda$通常按降序排列
• 矩阵奇异当且仅当有0特征值
• 所有特征值为正，则正定；不为负，则半正定；不为正，则半负定；为负，则负定。感受一下正负定对线性变换的影响

实对称矩阵的谱分解

实对称矩阵的特征分解：
$$S=(Q\Lambda )Q^T=\sum _i{\lambda }_i{q}_i{q}_i^T$$
注意谱是对称的，并且全是秩一矩阵

且有
$$S{q}_i=(\sum _i{\lambda }_i{q}_i{q}_i^T)q_i=\lambda q_i$$

应用：瑞利商Rayleigh Quotient

$$R(S,x)=\frac{x^{T}Sx}{x^{T}x}$$
其中$S$是实对称矩阵。

• 易证${\lambda }_{\min }⩽R(S,x)⩽{\lambda }_{\max }$，$x$分别指向最大和最小特征值对应的特征向量时，取最大和最小值。

• 当$x$指向其他特征向量方向时，$x$位于鞍点saddle point，可以把$x+\Delta x$用特征向量表示，然后发现函数增减和$\Delta x$方向有关。不过鞍点可以是maxmin问题的解，例如当$S$是对角阵时，对于$i$维子空间$V^i$，$$\underset{V^i}{\max }\underset{x\in V^i}{\min }R(S,x)={\lambda }_i$$
  其中${\lambda }_i$是第$i$大的特征值。

• 广义瑞利商：
  $$R(S,B,x)=\frac{x^{T}Sx}{x^{T}Bx}=\frac{{\tilde{x}}^T(B^{-1/2}S{B}^{-1/2})\tilde{x}}{{\tilde{x}}^T\tilde{x}}$$
  仍然能化成一般形式

• 很多聚类或降维任务，包括线性判别分析、主成分分析等，经过推导，往往能推导出瑞利商的形式

应用：马尔可夫矩阵Markov Matrix（随机矩阵Stochastic Matrix）与转移过程

马尔可夫矩阵$A$：

1. 每列的和为1
2. 所有元素大于0

物理意义是列到行概率转移。

• 上述定义1保证了1是特征值，原因在于$A-I$奇异. 并且其对应的特征向量每个值都不小于0（总之不能异号），根据马尔可夫矩阵的物理意义可以理解。特征值为1的特征向量也是稳态$A^{\infty }{u}_0$的组成部分
• 两条定义共同保证了特征值不大于1. 注意马尔可夫矩阵在无限次放后的稳态性，$A^n=V{\Lambda }^k{V}^{-1}$，特征值必须收敛

数值计算：快速特征值计算方法

记$A_0=A$，做QR分解得到$A_0=Q_0{R}_0$，取$A_1=R_0{Q}_0$. 则$A_1=R_0{A}_0{R}^{-1}$，和$A_0$相似，所以特征值相同，同样方法构造$A_2,A_3,\cdots$，下三角部分会越来越小（证明略，可能比较复杂），从而对角线越来越接近特征值。

改进

引入平移矩阵，得到$A_0-sI$，进行QR分解，取$A_1=R_0{Q}_0+sI=R_0(A_0-sI)R_0^{-1}+sI=R_0{A}_0{R}_0^{-1}$，会发现这样这样的平移不改变特征值。这样做的好处是得到的特征值收敛更快。

这一节内容参考了博客：MIT 18.065—机器学习中的矩阵方法12 计算特征值和奇异值

这样的计算方法在SVD分解中也有应用。MATLAB对这一类问题就用的是这种方法。这里Mark一下，不再展开。

相似矩阵

$$A\sim B\iff A=MB{M}^{-1}$$
注意相似的传递性，这有点像是抽象代数里的群作用在集合上形成的“轨道”，相似的矩阵在一个轨道里

• 相似矩阵特征值相同。对于$Ax=\lambda x$，有$M^{-1}AM{M}^{-1}x=\lambda M^{-1}x$，进而$B(M^{-1}x)=\lambda M^{-1}x$，这说明$A$的特征值$\lambda$也是$B$的特征值。特征向量不同。
• 上述命题反过来不成立，如果特征值有重根，特征向量个数有可能不同。当特征向量不够的时候，甚至无法对角化。例如
  $$\left[\begin{array}{cc}4& 1\\ 0& 4\end{array}\right]$$
  这也叫Jordan标准型Jordan form
• 每个矩阵都像相似于一个Jordan矩阵，Jordan矩阵由Jordan块组成。而Jordan块数则等于不同的特征向量数。这里仅mark一下，不再展开，等用到时再详细了解吧。
• $AB$和$BA$特征值相同，证明方法是$B(AB)B^{-1}=BA$，两者相似，这里假定$B$可逆. SVD和这里有一定的相通之处，即$A^{T}A,A{A}^T$有相同非零特征值。

---

正定矩阵

一种特殊的实对称矩阵。几个等价定义：

• $\forall x\ne 0,x^{T}Sx>0$，即二次型quadratic form（相对于线型linear form）大于0
• 特征值全为正，即$S=Q\Lambda Q^T$中$\Lambda$对角线元素全正。所以迹恒正
• 主元pivots全为正。主元实际上和LU分解，和二次型配方都是一一对应的
• 顺序主子式leading determinants全为正。注意主元和顺序主子式的关系，前$k$主元的乘积是第$k$顺序主子式，可以用高斯消元法解释，高斯消元法不会影响行列式。如果有一个主元为负，那么该主子式为负
• $S=A^{T}A$， $A$列无关.
  如果$A$不满秩，亦能半正定positive semi definite.考虑定义1和$N(A)$易证. 所以协方差总是半正定，就像是方差总是正数

性质：

• 如果$A,B$都正定，那么$A+B$正定，用定义1秒证
• 正定矩阵可逆，且逆矩阵正定

说明：

• 二次型$x^{T}Sx$也叫做能量Energy，许多机器学习和深度学习的问题就是在找参数$x$来最小化能量
• 二阶偏导数的Hessian矩阵半正定时，原函数是凸函数convex function

应用：图的拉普拉斯矩阵Laplacian Matrix

对于$v$个点的图简无向单simple graph，拉普拉斯矩阵$L\in {\mathbb{R}}^{v\times v}$存储了图的信息
$$L=D-A$$
其中$D\in {\mathbb{R}}^{v\times v}$是对角矩阵，表示每个顶点的度；$A$是图的邻接矩阵，仅由1或0组成。

• $L$一定是半正定矩阵。
  证明方法为$L=M^{T}M$，其中$M\in {\mathbb{R}}^{e\times v}$是关联矩阵Incidence matrix，$e$是图的边数，$M_{i\cdot }$是第$i$条边（连接$p,q$两点，且$p>q$），$M_{ip}=1,M_{iq}=-1$，行中其他元素为0.
  另外注意$rank[L]=v-1$，所以只能半正定. $L$有特征值0，对应特征向量为${\text{1}}_n$
• 拉普拉斯矩阵还有一些性质，和图的部分数量有关，另外和谱聚类等方法有联系，可以自行查阅相关资料，如ESL第14章。

扩展：隔行扫描定理Cauchy Interlacing Theorem

该定理揭示了秩1正定矩阵所带来的变化会如何影响特征值。
对于对称矩阵$S$，特征值为${\lambda }_1⩾{\lambda }_2⩾\cdots$，对于$S+\theta u{u}^T$得到特征值${\mu }_1⩾{\mu }_2⩾\cdots$，其中$\theta$是一个系数，$u$是单位向量。$u{u}^T$一定正定，并可以看作是positive change，会导致特征值上升。但是特征值的上升有上界，实际上有
$${\lambda }_1+\theta ⩾{\mu }_1⩾{\lambda }_1⩾{\mu }_2⩾\lambda 2⩾\cdots \text{\hspace{1em}(3)}$$
如何证明？

Weyl不等式

（该定理不再证明）
对于对称矩阵$S,T$，降特征值降序排列，有${\lambda }_{i+j-1}^{S+T}⩽{\lambda }_i^S+{\lambda }_j^T$.

• 当$j=1$时，${\lambda }_i^{S+T}⩽{\lambda }_i^S+{\lambda }_{\max }^T$.
  可以发现式(3)中${\lambda }_1+\theta ⩾{\mu }_1$
  该定理给出了对称矩阵带来变化的特征值变化上界
• 当$j=2$时，${\lambda }_{i+1}^{S+T}⩽{\lambda }_i^S+{\lambda }_2^T$
  式(3)中秩1矩阵的${\lambda }_2^T=0$，所以序关系得证

---

奇异值分解Singular Value Decomposition(SVD)

奇异值分解非常重要，在后文应用非常多。
对任意$m\times n$矩阵$A$，总能
$$A=U\Sigma V^T$$
其中，$U$是$m$阶正交矩阵，$\Sigma$是$m\times n$对角矩阵且${\sigma }_1⩾{\sigma }_2⩾\cdots ⩾0$，$V$是$n$阶正交矩阵。奇异值分解并不唯一，往往让$\Sigma$对角元素沿着从大到小的顺序排，在某些情况下唯一。参考：可逆矩阵的奇异值分解唯一吗？存在奇异值分解唯一的矩阵吗？

• 性质：$A^{T}A=V({\Sigma }^T\Sigma )V^T$，$A{A}^T=U(\Sigma {\Sigma }^T)U^T$
• 如果$A$是正定，则SVD即是正交特征分解

代数解释

如果$rank[A]=r$，奇异值分解让$A$把行空间$C(A^T)$的一组$r$个单位正交基同构映射到列空间$C(A)$中去得到另一组$r$个单位正交基，注意$C(A)$和$C(A^T)$维度相同，即
$$A[v_1,v_2,\cdots ,v_r]=[u_1,u_2,\cdots ,u_r]\left[\begin{array}{cccccc}{\sigma }_1& & & & & \\ & \ddots & & & & \\ & & {\sigma }_r& & & \\ & & & 0& & \\ & & & & \ddots & \\ & & & & & 0\end{array}\right]\text{\hspace{1em}(2)}$$
也可写成
$$AV=U\Sigma$$
其中$V$中多余的列如$v_{r+1},\cdots ,v_n$从$A$的零空间$N(A)$中拿即可，其维度正好是$n-r$，而且所有列都正交；$U$中多余的列类似，从左零空间$N(A^T)$中拿. 对于$i>r$，${\Sigma }_{ii}=0$，$r$列之后等式两边全是0

• 如何找出$V$的前$r$列?
  从$A^{T}A$的非零特征值对应的特征向量中拿，$A^{T}A=V({\Sigma }^T\Sigma )V^T$，$V$也即实对称矩阵$A^{T}A$的特征向量，从中也能看出$\Sigma$；

  • ${\sigma }_i^2{v}_i=A^T(A{v}_i)$，所以$v_i$一定在行空间$C(A^T)$中

• 如何找出$U$的前$r$列？
  令$u_i=\frac{A{v}_i}{{\sigma }_i}$

  • 可以看出$u_i$一定在$A$的列空间里
  • 可以证明这$r$列是单位正交的，不妨记前$r$列的各个矩阵为$V_r\in {\mathbb{R}}^{n\times r},U_r\in {\mathbb{R}}^{m\times r},{\Sigma }_r\in {\mathbb{R}}^{r\times r}$，则$$\begin{gathered} U_r=A{V}_r{\Sigma }_r^{-1} \\ {U}_r^T{U}_r=(A{V}_r{\Sigma }_r^{-1}{)}^{T}A{V}_r{\Sigma }_r^{-1}={\Sigma }_r^{-1}{V}_r^T(V{\Sigma }^T\Sigma V^T)V_r{\Sigma }_r^{-1}=({\Sigma }_r^{-1}{V}_r^{T}V{\Sigma }^T)({\Sigma }_r^{-1}{V}_r^{T}V{\Sigma }^T{)}^T \end{gathered}$$.
    其中$V_r^{T}V=I_{r\times n}$，${\Sigma }_r^{-1}{I}_{r\times n}={\Sigma }_r^{-1}$, ${\Sigma }_r^{-1}{\Sigma }^T=I_{r\times r}$，所以$U_r^T{U}_r=I_{r\times r}$. 即证明了这$r$列一定单位正交，式(2)是一定可以构造出来的

• 注意$V$和$U$一起构成了四个基本子空间，SVD实现了四个基本子空间的分解，

• $U:[C(A),N(A^T)]$: ， $V：[C(A^T),N(A)]$

• 从矩阵分解的角度揭示了矩阵从四个记本子空间中如何生成的，这也是SVD各种秒用的原理之一

几何解释

从线性变换的角度看奇异值分解揭示了矩阵的组成方法， 如下图所示（摘自Prince的Computer vision: models, learning, and inference）：

• $V^T$ 旋转（+反演，反射与否可用行列式判断，下同）
• $L$ 伸缩每个维度（假定$L$对角线都是正值，也即把反演放到旋转中去）
• $U$ 旋转（+反演）

• 只改变某部分，对最终矩阵$A_3$ 的影响
  

奇异值$\sigma$

• 个数：$min(m,n)$
• 非0奇异值的数量是矩阵的秩
• 最大奇异值与最小奇异值之比反映了可逆性的度量，称为条件数。这个概念和后文的数值秩有联系
• 根据行列式的定义，知，奇异值的乘积为行列式（因为只有$L$拉伸了矩阵，$U$和$V$不改变面积）

逆

$$A^{-1}=V{L}^{-1}{U}^T$$
从几何角度和数值角度都很直观

极分解Polar Decomposition生成SVD

$$A=QS$$
其中$Q$正交（酉矩阵），$S$半正定矩阵（半正定埃尔米特矩阵）。
所以
$$A=QV\Lambda {\bar{V}}^T$$
$QV=U$正交，$V$正交，这也即SVD分解

• 关于求$S$，${\bar{A}}^{T}A={\bar{S}}^T{\bar{Q}}^{T}QS=S^2$，用矩阵开方的知识，正规矩阵都可以开方，$S=\sqrt{{\bar{A}}^{T}A}=\sqrt{V{\Lambda }^2{\bar{V}}^T}=V\Lambda {\bar{V}}^T$

---

矩阵范数

这里插入矩阵范数，以便应用部分使用

范数定义

满足这三条的叫做范数

• $\Vert A\Vert ⩾0$
• $\Vert cA\Vert =|c|\Vert A\Vert$
• $\Vert A+B\Vert ⩽\Vert A\Vert +\Vert B\Vert$

如果定义$\Vert A\Vert =({\sum }_{i,j}{a}_{i,j}^p{)}^{1/p}$，当$p<1$时，会违背第三条性质，所以不是范数；$p⩾1$是可以的。
接下来给出3个重要范数，这些范数都可以用奇异值表示。

2范数（谱范数）

$$\Vert A{\Vert }_2=\underset{x}{\max }\frac{\Vert Ax{\Vert }_2}{\Vert x{\Vert }_2}={\sigma }_1$$
即最大的奇异值，不妨假定$\Vert x{\Vert }_2=1$，证明方法为记$x={\sum }_i{\lambda }_i{v}_i$，则${\sum }_i{\lambda }_i^2=1$，根据SVD分解性质有$Ax={\sum }_i{\lambda }_i{\sigma }_i{u}_i$，所以$\Vert Ax{\Vert }_2=\sqrt{{\sum }_i{\lambda }_i^2{\sigma }_i^2}⩽{\sigma }_1$

F范数Frobenius Norm

$\Vert A{\Vert }_F$是所有元素平方和再开方，$\Vert A{\Vert }_F=\sqrt{tr[A^{T}A]|}=\sqrt{tr[V{\Sigma }^T\Sigma V^T]}$，根据迹的轮换对称性 $tr[AB]=tr[BA]$，所以$$\Vert A{\Vert }_F=\sqrt{tr[V^{T}V{\Sigma }^T\Sigma ]}=\sqrt{\sum _{i=1}^r{\sigma }^2}$$

核范数Nuclear Norm

$$\Vert A{\Vert }_N=\sum _{i=1}^r{\sigma }_i$$

这三种范数在正交变换下不变

$QA=(QU)\Sigma V^T$，这是$QA$的SVD分解，如果范数只由$\Sigma$决定，那么$Q$不影响范数

向量特例

这三个范数相等，所以经常省略下标，简写为$\Vert x\Vert$

• 基追踪问题Basis pursuit
  $$\begin{gathered} min\Vert x{\Vert }_{p}s.t.Ax=b \\ p=1or2 \end{gathered}$$
  这个问题似乎不太好解

应用：SVD矩阵低秩压缩

这里矩阵可以是一幅图像或DataFrame型数据。如果$rank[A]=r$，那么
$$A=U\Sigma V^T=\sum _{i=1}^r{\sigma }_i{u}_i{v}_i^T$$
有点像谱分解，上述$\sigma$按绝对值从大到小排列，越大的$\sigma$占据了数据的主要成分。

• 如果$r\ll \min (m,n)$，则该矩阵可以很好的无损压缩。
• 注意${\Vert u_i{v}_i^T\Vert }_F=tr[v_i{u}_i^T{u}_i{v}_i^T]=tr[v_i{v}_i^T]=tr[v_i^T{v}_i]=1$. 所以每一个谱的数据量主要由$\sigma$衡量。可以只保留最大的$k$个$\sigma$，对$A$进行有损压缩得到$A_k=U_k{\Sigma }_k{V}_k^T$

Eckart-Young定理

如果$A_k=U_k{\Sigma }_k{V}_k^T$，对任意秩$k$矩阵，
$$\Vert A-B\Vert ⩾\Vert A-A_k\Vert$$
$A_k$是秩$k$矩阵中最接近$A$的。该定理不再证明。它说明了SVD矩阵低秩压缩和PCA的最优性。

数值秩Numerical Rank

对于奇异值很小的部分，数值秩采取软处理。定义容忍度$ϵ\in (0,1)$，$ran{k}_{ϵ}[A]=r$含义是：
$$\begin{array}{r}{\sigma }_{k+1}⩽ϵ{\sigma }_1(x)\\ {\sigma }_k>ϵ{\sigma }_1(x)\end{array}$$

• $ran{k}_0[x]=rank[x]$
• ${\sigma }_{k+1}(x)=\Vert A-A_k{\Vert }_2$，根据Eckart-Young定理，${\sigma }_{k+1}$度量了最优拟合差异的2范数

数值秩类似于模糊数学的概念，如果一个矩阵满秩，但是很多奇异值很小，导致矩阵不容易求逆。在数值秩看来这就是一个低秩矩阵。典型例子：

• 希尔伯特矩阵$H$，其中$H_{ij}=1/(i+j-1)$，该矩阵满秩，但数值秩很低。对于1000维希尔伯特矩阵，$rank[H_{1000}]=1000,ran{k}_{10^{-15}}[H_{1000}]=28$
• 范德蒙矩阵。注意范德蒙矩阵也是多项式线性回归中扩展后的数据矩阵

其他判断低数值秩矩阵的方法包括Sylvester等式，不再展开细讲。

应用：主成分分析Principal Components Analysis(PCA)

假设存在数据集$X=[x_1,\cdots ,x_n]\in {\mathbb{R}}^{m\times n}$，$m$是特征数，$n$是样本数，希望通过线性变换$W\in {\mathbb{R}}^{k\times m}$把特征降低到$k$维，左乘得到降维结果$WX$，并希望降维再升维度的数据能通过线性变换尽可能代表原始数据，距离最小。
因为降维后的矩阵秩不超过$k$，如果再升维回去，秩仍然不超过$k$，根据Eckart-Young定理，$X_k$是最接近$X$的矩阵。
其中$X_k=U_k{\Sigma }_k{V}_k^T$通过SVD分解得到。
什么样的线性变换能在降维和升维后得到$X_k$？

• 降维：实际上可以取$U_k^T$左乘，$U_k^{T}X$把数据集降维到${\mathbb{R}}^{k\times n}$，即$u_i^{T}x$是降维后第$i$维的结果
• 升维：取矩阵$U_k$作为升维矩阵，$$U_k(U_k^{T}X)=U_k(U_k^{T}U)\Sigma V^T=U_k{I}_{k\times m}\Sigma V^T=U_k[{\Sigma }_k,{\text{0}}_{k\times (n-k)}]V^T=U_k{\Sigma }_k{V}_k^T=X_k$$.

所以这样的变换是成立的，即对$X{X}^T$进行特征分解，取前$k$大特征值对应的特征向量$u_1,\cdots ,u_r$作为降维的投影向量。
如果变换之前，先对$X$减去各个特征均值，得到$\hat{X}=X-{\mu }_m{\text{1}}_{1\times n}$，即$\hat{X}{\text{1}}_n=0$，再进行降维，这就是PCA，其中${\text{1}}_{1\times n}$和${\text{1}}_n$表示全1行向量和列向量。注意$\hat{X}{\hat{X}}^T/(m-1)$是特征间的协方差矩阵。

应用：多维缩放Multiple Dimensional Scaling

对于$n$个坐标未知的样本点，给定两两之间的距离平方，构成距离矩阵$D\in {\mathbb{R}}^{n\times n}$，其中$D_{ij}=\Vert x_i-x_j{\Vert }^2$. 求出这$n$个样本点的坐标。往往我们要求坐标空间的维度$m<n$.

• 这个问题在降维中有实际意义。参考低维嵌入和流形学习。考虑高维空间中的点，这些点可能分布在低维流形manifold上，如图所示（摘自周志华《机器学习》），通过计算相邻样本点的距离，构成图，计算任意两点间的测地线距离（最短路径长度）当作这两点的距离，就获得了距离矩阵，但是不知道低维样本点坐标。

已知$D_{ij}=x_i^T{x}_i+x_j^T{x}_j-2x_i^T{x}_j$，思路是先求出$X^{T}X$. 假设${\sum }_{i=1}^n{x}_i=0$，则
$$\sum _{i=1}^n\sum _{j=1}^n{D}_{ij}=2n\cdot tr[X^{T}X]$$
可以求出$tr[X^{T}X]$. 列求和则得到
$$\sum _{j=1}^n{D}_{ij}=n{x}_i^T{x}_i+tr[X^{T}X]$$
又可以求出任何$x_i^T{x}_i$，从而任何$x_i^T{x}_j$都可以求出。写成矩阵结果为
$$X_i^T{X}_j=-\frac{1}{2}(D-\frac{1}{n}D{\text{1}}_{n\times n}-\frac{1}{n}{\text{1}}_{n\times n}D+\frac{1}{n^2}tr[D^{T}D])$$
其中${\text{1}}_{n\times n}$表示全1矩阵。
如果$D$是满足三角不等式的合理距离矩阵，那么$X^{T}X$是半正定对称阵。
从而特征值分解（也即SVD分解）得到$X^{T}X=V\Lambda V^T$，如果我们希望降到$m$维，根据前文Eckart-Young定理，取特征值最大的前$m$个特征值和特征向量，得到$\hat{X}={\Lambda }_m^{1/2}{V}_m^T$.
（注意找低维近似的时候，没有办法做矩阵开方）

应用：主方向/最小值方向问题

主方向问题：
$$\underset{b}{\mathrm{argmax}}\Vert Ab\Vert s.t.\Vert b\Vert =1$$
最小值方向问题则是求$\arg \min \Vert Ab\Vert$。

• 作为工具，该问题和最小二乘的“地位”类似。可以用于求解几何模型中针孔摄像机的外在参数，双视图匹配等问题。下一节正交Procrustes问题也在这些类似问题上求参数时会用到。详见Computer vision: models, learning, and inference的14-16章。
• 根据限制，$b$必须在单位圆上，问题变为找一个方向，对应变换后的椭圆长轴。
  根据上述SVD分解，待求方向经过$V^T$变换后，为坐标系主轴方向$X=(1,0,...,0)$
  把它变回去，该方向为$(V^T{)}^{-1}X=VX$，即$V$的第一列。
• 同理，最小值方向为$V$的最后一列

应用：正交Procrustes问题

Procrustes是希腊语Προκρούστης.
寻找正交矩阵$Q$，使
$$\hat{Q}=\underset{Q}{\mathrm{argmin}}{\Vert QA-B\Vert }_F$$

首先
$$\begin{array}{rl}\hat{Q}& =\underset{Q}{\mathrm{argmin}}tr\left[(QA-B{)}^T(QA-B)\right]\\ & =\underset{Q}{\mathrm{argmax}}tr\left[Q^{T}B{A}^T\right]\end{array}$$
计算SVD分解$B{A}^T=U\Sigma V^T$，则
$$\begin{array}{rl}\hat{Q}& =\underset{Q}{\mathrm{argmax}}tr\left[Q^{T}U\Sigma V^T\right]\\ & =\underset{Q}{\mathrm{argmax}}tr\left[V^T{Q}^{T}U\Sigma \right]\end{array}$$
注意到
$$tr\left[V^T{Q}^{T}U\Sigma \right]=tr\left[Z\Sigma \right]=\sum _{i=1}^I{z}_{ii}{l}_{ii}$$
其中$Z=V^T{Q}^{T}U$，$Z$是三个正交矩阵乘积，所以也正交，对角线上每个数值小于等于1，所以选择$Z=I$来最大化上述目标，全局解为
$$\hat{Q}=U{V}^T$$

• 本问题的一个特殊形式是给定方阵$B$，寻找最接近的正交矩阵$Q$，即最优化
  $$\hat{Q}=\underset{Q}{\mathrm{argmin}}\Vert Q-B{\Vert }_F$$
  用上述方法，得到解为$\hat{Q}=U{V}^T$，其中$U\Sigma V^T=B$

应用：压缩感知Compressed Sensing

压缩感知的目的是通过部分信息恢复全部信息。思想是把全信息空间转换成稀疏表达空间（例如空域转频域、把读者喜好的书籍条目转为读者喜好的书籍分类等）。这种表达理应是稀疏的，所以当只观察到部分数据时，可以以此去寻找稀疏表达空间中尽量稀疏的解。
压缩感知可以用于矩阵补全，进一步可用于填补数据缺失值。对于矩阵$A$，部分数据缺失，压缩感知的处理方式是填补缺失值并使$rank[A]$最小。这是NP难问题，注意到$rank[A]$在集合 { X ∈ R m × n : ∥ X ∥ F 2 ⩽ 1 } \{X\in \mathbb R^{m\times n}:\|X\|^2_F \leqslant 1\} {X∈Rm×n:∥X∥F2​⩽1}的凸包是$X$的核范数，所以可通过最小化核范数的形式近似求解原问题。可用半正定规划求解。秩和核范数的关系有点像是0范数和1范数的关系。
更多讨论参考周志华《机器学习》11.6节。

---

求逆技巧

左右逆

对于$A\in {\mathbb{R}}^{m\times n}$

• 左逆: $(A^{T}A{)}^{-1}{A}^T\in {\mathbb{R}}^{n\times m}$，这需要$A^{T}A$可逆，也即$rank[A]=n$
• 右逆: $A^T(A{A}^T{)}^{-1}\in {\mathbb{R}}^{n\times m}$，需要$rank[A^T]=rank[A]=m$

伪逆

考虑行空间和列空间，从行空间中拿元素$x$，$Ax$一定在列空间中。而且$x$和$Ax$能构成一一对应，两个空间维度相同。注意该定理和SVD代数解释中的原理联系。证明：

• 对于$x_1,x_2\in C(A^T)$，$A{x}_1\ne A{x}_2$，否则$A(x_1-x_2)=0$，$x_1-x_2$不可能在$N(A)$当中，所以单射！
• 另外只取$x\in C(A^T),Ax$就能张成列空间，证明：
  对于$x^{\prime }\in {\mathbb{R}}^n$，$A{x}^{\prime }$是表示各种列向量的组合情况，必然能张成列空间。另一方面，必有分解$x^{\prime }=x+x_N$，其中$x\in C(A^T),x_N\in N(A)$，$A{x}^{\prime }=Ax+A{x}_N=Ax$，也即只需要$Ax$就能张成列空间。满射！
• 所以双射，即一一对应成立。

此时限制在这两空间上的逆就是伪逆，即对于$Ax=y$，伪逆$A^{+}$使$x=A^{+}y=(A^{+}A)x$

应用：无需$A^{T}A$可逆的最小二乘问题

先来介绍下SVD伪逆

SVD伪逆

伪逆有很多求法，一种求法是SVD分解：
$A=U\Sigma V^T$，$U,V$可逆， ${\Sigma }^{+}$对角线元素为$1/{\sigma }_1,1/{\sigma }_2,\cdots ,1/{\sigma }_r,0,\cdots ,0$，其余元素为0，注意$\Sigma \in {\mathbb{R}}^{m\times n},{\Sigma }^{+}\in {\mathbb{R}}^{n\times m}$，$A^{+}=V{\Sigma }^{+}{U}^T$

• $AV=U\Sigma =[u_1{\sigma }_1,\cdots ,u_r{\sigma }_r,0,\cdots ,0]$，$A$把$v_{r+1},\cdots ,v_n$变换成0，这是必然的，因为$v_{r+1},\cdots ,v_n$是在$N(A)$中。有趣的是$A^{+}U=V{\Sigma }^{+}=[v_1/{\sigma }_1,\cdots ,u_r/{\sigma }_r,0,\cdots ,0]$，这说明$A^{+}$把$u_{r+1},\cdots ,u_m$变换为0，进一步说明:
  • 性质1： $A^{+}$把$N(A^T)$变换成0

对于$Ax=b$，有最优解$\hat{x}=(A^{T}A{)}^{-1}{A}^{T}b$，这需要$A^{T}A$可逆。用伪逆可以消除该限制，而且思考过程非常漂亮！
重新审视最小二乘问题，仍用投影的思路，但是换一个思考方式。记$b=b_C+b_N$，其中$b_C\in C(A),b_N\in N(A^T)$，根据投影，$Ax=b$与$Ax=b_C$同解。
如果取$x=A^{+}{b}_C$，则$Ax=A{A}^{+}{b}_C$，因为$b_C\in C(A)$，所以$A{A}^{+}{b}_C=b_C$，说明$x=A^{+}{b}_C$是$Ax=b_C$的解，也即$Ax=b$的解。
根据性质1，$A^{+}{b}_C=A^{+}b$，所以$x=A^{+}b=V{\Sigma }^{+}Ub$是$Ax=b$的一个最小二乘可行解

• 讨论：
  1. 如果$A^{T}A$可逆，$rank[A]=n$，意味着$m⩾n$，则$x=A^{+}b=(A^{T}A{)}^{-1}{A}^{T}b$
     此时$A^{+}=(A^{T}A{)}^{-1}{A}^T$；进一步如果$A$可逆，则$A^{+}=A^{-1}$
  2. 如果$A^{T}A$不可逆，则$x=A^{+}b$是可行解中的一个，是唯一的$C(A^T)$中的解

注意当$\sigma$从0扰动到非0时，伪逆变化非常大。一种处理方法是把$\sigma$小于某个阈值的都当作0. 另一种方法是引入L2正则化

L2正则化

对于$Ax=b$问题，其中$A\in {\mathbb{R}}^{m\times n}$，如果$rank[A]<n$或者$A$接近不列满秩，$A^{T}A$行列式很小，求逆后结果爆炸，则可以改求
$$\underset{x}{\min }\Vert Ax-b{\Vert }_2^2+{\delta }^2\Vert x{\Vert }_2^2(\delta >0)$$
上式对$x$进行正则化regularization，这是统计学家的观点。但从线性代数的角度可以写成成矩阵形式，即最小二乘
$$A^{\ast }x=\left[\begin{array}{c}A\\ \delta I_{n\times n}\end{array}\right]x=\left[\begin{array}{c}b\\ 0_{n\times 1}\end{array}\right]=b^{\ast }$$
这里$A^{\ast T}{A}^{\ast }$的可逆性很好

• $rank[A^{\ast }]=n$，从而$A^{\ast T}{A}^{\ast }$可逆
• 另一种更推荐的解释方法：$A^{\ast T}{A}^{\ast }=A^{T}A+{\delta }^2{I}_{n\times n}$，这里$A^{T}A$半正定，$\delta$越大，$A^{\ast T}{A}^{\ast }$越正定，越远离不满秩。

从而$x=(A^{\ast T}{A}^{\ast }{)}^{-1}{A}^{\ast T}{b}^{\ast }=(A^{\ast T}{A}^{\ast }{)}^{-1}{A}^{T}b=(A^{T}A+{\delta }^2{I}_{n\times n}{)}^{-1}{A}^{T}b$

逆关系

Sherman-Morrison恒等式

对于$U,V\in {\mathbb{R}}^{n\times k}$，$rank[U]=rank[V]=k$
$$(I_n-U{V}^T{)}^{-1}=I_n+U(I_k-V^{T}U{)}^{-1}{V}^T$$

只需证明
$$\begin{array}{rl}{I}_n& =I_n-U{V}^T+(I_n-U{V}^T)U(I_k-V^{T}U{)}^{-1}{V}^T\\ & =I_n-U{V}^T+U(I_k-V^{T}U{)}^{-1}{V}^T-U{V}^{T}U(I_k-V^{T}U{)}^{-1}{V}^T\\ & =I_n-U{V}^T+U(I_k-V^{T}U)(I_k-V^{T}U{)}^{-1}{V}^T=I_n\end{array}$$

• 秩1矩阵$u,v$
  $$(I-u{v}^T{)}^{-1}=I+\frac{u{v}^T}{1-v^{T}u}$$

Sherman–Morrison–Woodbury恒等式

该式也叫矩阵求逆引理、Woodbury矩阵恒等式。
令Sherman-Morrison恒等式中$U=AX,V^T=C^{-1}Y$，其中$A\in {\mathbb{R}}^{n\times n},X\in {\mathbb{R}}^{n\times k},Y\in {\mathbb{R}}^{k\times n},C\in {\mathbb{R}}^{k\times k}$，推导几步，易得
$$(A^{-1}+X{C}^{-1}Y{)}^{-1}=A-AX(C+YAX{)}^{-1}YA$$

• 该式的一个用处是降低求逆维度，可以把左侧对$n$维求逆改为右侧对$k$维求逆，如果$n>k$，则有助于求逆。

如果$C=I_k$，易得
$$(A^{-1}+XY{)}^{-1}=A-AX(I_k+YAX{)}^{-1}YA$$

• 该公式揭示了如果对矩阵$A^{-1}$有秩$k$的修改$XY$，例如$XY$是修改的$CR$分解，那么其逆矩阵也会有秩$k$的修改过程。

如果$Y=X^T$，易得
$$(A^{-1}+X{C}^{-1}{X}^T{)}^{-1}=A-AX(X^{T}AX+C{)}^{-1}{X}^{T}A$$

应用：因子分析Factor Analysis的概率密度函数

对于$n$维高斯分布，如果维度过大（例如把整个图像像素建模成高斯分布），协方差参数太多，达到$\mathcal{O}(n^2)$量级，可以用因子分析模型进行有损压缩，把参数降低到线性。该模型形式化为
$$\begin{array}{rl}P(x|h)& =\mathcal{N}(x|\mu +\Phi h,\Sigma )\\ P(h)& =\mathcal{N}(h|0,I)\end{array}$$
其中观测数据$x$的维度为$n$，隐变量$h$的维度为$k(n>k)$，参数包括：$\Phi \in {\mathbb{R}}^{n\times k}$，$\Sigma$ 是对角阵.

• $x=\mu +\Phi h+ϵ$，$x$可由隐变量$h$解释
• 当$\Sigma$是单位阵的常数倍时，该模型对球形协方差建模，称为概率主成分分析probabilistic principal component analysis或概率PCA.
• 因子分析模型和PCA类似，可以用于人脸形状建模等问题。

可以求出$P(x)=\mathcal{N}(x|\mu ,\Phi {\Phi }^T+\Sigma )$ ，仍然是高斯分布，且协方差参数进行了压缩。如何求解呢？

求解方法1：高斯线性模型

注意到这是高斯线性模型（linear Gaussian model)！直接套结论。结果为$\Phi {\Phi }^T+\Sigma$. 可以参考PRML 2.3.3节，用了高斯分布逆辨识，比较tricky.
这个过程中需要用舒尔补恒等式Schur complement identity对分块矩阵求逆，而舒尔补恒等式则和Sherman-Morrison-Woodbury恒等式有很大联系。

求解方法2：配平方法：边缘分布推导

直接暴力求解，需要耐心
$$\begin{array}{rl}P(x)& =\int \mathcal{N}(x|\mu +\Phi h,\Sigma )\mathcal{N}(h|0,I)dh\\ =& C_1\int exp\left\{-\frac{1}{2}(x-\mu -\Phi h{)}^T{\Sigma }^{-1}(x-\mu -\Phi h)\right\}exp\left\{-\frac{1}{2}{h}^{T}h\right\}dh\\ =& C_2\int exp\left\{-\frac{1}{2}[(x-\mu {)}^T{\Sigma }^{-1}(x-\mu )-2(x-\mu {)}^T{\Sigma }^{-1}\Phi h+h^T({\Phi }^T{\Sigma }^{-1}\Phi +I)h]\right\}dh\\ =& C_3\exp \left\{-\frac{1}{2}[(x-\mu {)}^T[{\Sigma }^{-1}-{\Sigma }^{-1}\Phi ({\Phi }^T{\Sigma }^{-1}\Phi +I{)}^{-1}{\Phi }^T{\Sigma }^{-1}](x-\mu )]\right\}\\ & \times \int exp\left\{-\frac{1}{2}[h-({\Phi }^T{\Sigma }^{-1}\Phi +I{)}^{-1}{\Phi }^T{\Sigma }^{-1}(x-\mu ){]}^T({\Phi }^T{\Sigma }^{-1}\Phi +I)[h-({\Phi }^T{\Sigma }^{-1}\Phi +I{)}^{-1}{\Phi }^T{\Sigma }^{-1}(x-\mu )]\right\}dh\end{array}$$
其中$C_1,C_2,C_3$表示$x,h$无关的系数。注意到把$x$有关而$h$无关的项提到积分外，积分内是关于$h$的一个高斯分布。积分后得到
$$\begin{array}{r}{C}_4\exp \left\{-\frac{1}{2}[(x-\mu {)}^T[{\Sigma }^{-1}-{\Sigma }^{-1}\Phi ({\Phi }^T{\Sigma }^{-1}\Phi +I{)}^{-1}{\Phi }^T{\Sigma }^{-1}](x-\mu )]\right\}\end{array}$$
这是一个高斯分布，均值为$\mu$，方差为
$$[{\Sigma }^{-1}-{\Sigma }^{-1}\Phi ({\Phi }^T{\Sigma }^{-1}\Phi +I{)}^{-1}{\Phi }^T{\Sigma }^{-1}{]}^{-1}$$
这里套用Sherman–Morrison–Woodbury恒等式，上式化为$\Phi {\Phi }^T+\Sigma$

应用：卡尔曼滤波

卡尔曼滤波在跟踪问题中作为时序模型常见。卡尔曼滤波中含有两个高斯线性模型，和因子分析模型类似，Sherman–Morrison–Woodbury恒等式可以化简卡尔曼滤波的表达式。更多内容可以参考Prince的Computer vision: models, learning, and inference第19章或Bishop的Pattern Recognition and Machine Learning第13.3节

---

分解方法Factorizations汇总

$$\begin{array}{rlr}A=& CR& \text{ 列/行向量无关}\\ A=& LU& \text{ elimination，下上三角分解}\\ A=& QR& \text{Gram-Schmitt, Q正交}\\ A=& V\Lambda V^{-1}& 一般特征分解（\Lambda 为Jordan标准型）\\ S=& Q\Lambda Q^T& \text{S对称, 特征分解，Q正交}\\ A=& U\Sigma V^T& \text{SVD}\\ A=& QS& \text{Polar Decomposition}\end{array}$$

注意每种分解前后的自由度相等，可以自行验证。
对于矩阵$A\in {\mathbb{R}}^{m\times n},rank[A]=r$，其自由度为$(m+n-r)r$，可以用SVD分解或者从几何角度找子空间中的向量组来证明。

A=LU的解释与高斯消元法

$$\begin{gathered} LU=[l_1,l_2,\cdots ]\left[\begin{array}{c}{u}_1\\ u_2\\ ⋮\end{array}\right] \\ =l_1{u}_1+l_2{u}_2+\cdots \end{gathered}$$
$L$是单位下三角矩阵（对角线为1）。$l_i{u}_i$维度与$A$同，可以看作是$l_1{u}_1$剥离了$A$的第一行第一列，……这样，$l_i{u}_i$只需面对已经消除了$i-1$行和$i-1$列的矩阵，所以$l_i$的前$i-1$个元素都是0，$u_i$同理。这样即实现了$LU$分解.

• 当A的所有顺序主子式都不为0时，主元不为0，矩阵$A$可以进行LU分解，且是唯一分解（摘自百度百科）
  
  上图也可以从高斯消元的角度看待，$U$是消元结果，$L$是初等矩阵的拼接

---

矩阵微分

一般情况

矩阵微分不难理解，但是很多资料疏于整理，导致这套技术不统一也难成体系，学习者不易掌握一套通用方法。
矩阵求导有很多种形式：

• 标量对矩阵（向量）求（偏）导：这是机器学习最常用的形式，即Loss对参数求导
• 矩阵（向量）对标量求导
• 向量对向量求导

推荐两个文章：

• 矩阵求导术（上）：这篇博客讲授Jacobian辨识，并将其作为一种较为通用的标量对矩阵求导的工具。可以处理点乘、矩阵乘法等多种运算。
• Matrix Calculus：师兄的文章，除了标量对矩阵求导外，对各类矩阵求导形式进行了精炼的推导和总结。适合自学和速查

杂知识：

• 标量对向量求导后，再求二阶导的Jacobian矩阵即为Hessian矩阵
• 对于，$f(A)=\ln det[A]$，$\nabla f=(A^{-1}{)}^T$，可以自行验证

应用：神经网络反向传播向量化

利用该方法可以推导出向量化的反向传播，见神经网络反向传播向量化（CS231n A1 Q4）

应用：动量Momentum优化二次型

对于函数$f(x)=\frac{1}{2}{x}^{T}Sx$，$S$是正定矩阵。
采用梯度下降的线搜索会以zigzag的形式逼近极小值。（如果每次线搜索都搜到最好的位置，那么zigzag每两个相邻迭代的方向是垂直的）
更高效的方法是采用momentum优化器，迭代方式为
$$\begin{gathered} x_{k+1}=x_k-\alpha z_k \\ {z}_{k+1}=\nabla f_{k+1}+\beta z_k \end{gathered}$$
用向量化写法可形式化为
$$\left[\begin{array}{cc}I& 0\\ -S& I\end{array}\right]{\left[\begin{array}{c}x\\ z\end{array}\right]}_{k+1}=\left[\begin{array}{cc}I& -\alpha I\\ 0& \beta I\end{array}\right]{\left[\begin{array}{c}x\\ z\end{array}\right]}_k$$

记$S$特征值$\lambda$，对应特征向量$q$。注意$S$是正定矩阵，只考虑一个特征分量，如果$x_k=c_{k}q,z_k=d_{k}q$，那么$S{x}_k=c_k\lambda q$，易证
$$\begin{array}{rl}\left[\begin{array}{cc}1& 0\\ -\lambda & 1\end{array}\right]\left[\begin{array}{c}{c}_{k+1}\\ d_{k+1}\end{array}\right]& =\left[\begin{array}{cc}1& -\alpha \\ 0& \beta \end{array}\right]\left[\begin{array}{c}{c}_k\\ d_k\end{array}\right]\\ \left[\begin{array}{c}{c}_{k+1}\\ d_{k+1}\end{array}\right]& =\left[\begin{array}{cc}1& -\alpha \\ \lambda & \beta -\lambda \alpha \end{array}\right]\left[\begin{array}{c}{c}_k\\ d_k\end{array}\right]\\ & =R\left[\begin{array}{c}{c}_k\\ d_k\end{array}\right]\end{array}$$
$S$是正定的，记最小和最大特征值分别为$m,M(>0)$，$\frac{M}{m}$即条件数，为1时，$S$是单位阵；当条件数很大时，则不好用梯度下降优化。
这里希望找到最好的$\alpha ,\beta$，使对任意$\lambda$(满足$m⩽\lambda ⩽M$)，$c_k,d_k$减少最快，即$R$的所有特征值尽可能小。
可以求得最优解为
$$\begin{gathered} {\alpha }_{opt}={\left(\frac{2}{\sqrt{M}+\sqrt{m}}\right)}^2 \\ {\beta }_{opt}={\left(\frac{\sqrt{M}-\sqrt{m}}{\sqrt{M}+\sqrt{m}}\right)}^2 \end{gathered}$$
此时$$\text{R的特征值}<{\left(\frac{1-\sqrt{b}}{1+\sqrt{b}}\right)}^2$$
其中$b=\frac{m}{M}$。
（该式正确性存疑，望有凸优化的大佬告知）

逆矩阵

对于矩阵函数$A(t)$，易得
$$\begin{array}{rl}(A+dA{)}^{-1}-A^{-1}& =(A+dA{)}^{-1}[A-(A+dA)]A^{-1}\\ & =-A^{-1}dA{A}^{-1}\end{array}$$
所以
$$\frac{d{A}^{-1}}{dt}=-A^{-1}\frac{dA}{dt}{A}^{-1}$$

特征值

对于矩阵$A(t)$，转置后特征值不变，但特征向量则有可能会变，所以对于同一特征值，有列特征向量和行特征向量：
$$\begin{gathered} A(t)x(t)=\lambda (t)x(t) \\ y(t{)}^{T}A(t)=\lambda (t)y(t{)}^T \end{gathered}$$
对相同特征值对应的行列特征向量，引入限制
$$y(t{)}^{T}x(t)=1$$
则
$$y(t{)}^{T}A(t)x(t)=\lambda (t)$$
所以
$$\begin{array}{rl}\frac{d\lambda }{dt}& =\frac{d{y}^T}{dt}Ax+y^T\frac{dA}{dt}x+y^{T}A\frac{dx}{dt}\\ & =y^T\frac{dA}{dt}x+\frac{d{y}^T}{dt}\lambda x+\lambda y^T\frac{dx}{dt}\\ & =y^T\frac{dA}{dt}x+\lambda \frac{d(y^{T}x)}{dt}\\ & =y(t{)}^T\frac{dA(t)}{dt}x(t)\end{array}$$
该式揭示了特征值的变化如何对应于矩阵的变化

奇异值

对于$A$的同一奇异值$\sigma$对应的左右奇异向量$u,v$，有$u(t{)}^{T}A(t)v(t)=u(t{)}^T[\sigma (t)u(t)]=\sigma (t)$，所以
$$\begin{array}{rl}\frac{d\sigma }{dt}& =\frac{d{u}^T}{dt}Av+u^T\frac{dA}{dt}v+u^{T}A\frac{dv}{dt}\\ & =\sigma \frac{d{u}^T}{dt}u+u^T\frac{dA}{dt}v+\sigma v^T\frac{dv}{dt}\end{array}$$
注意到$\frac{d(u^{T}u)}{dt}=0=\frac{d{u}^T}{dt}u+u^T\frac{du}{dt}$，所以$\frac{d{u}^T}{dt}u=0$，同理$\frac{d{v}^T}{dt}=0$，所以
$$\begin{array}{rl}\frac{d\sigma }{dt}& =u^T(t)\frac{dA(t)}{dt}v(t)\end{array}$$
其中$u,v$是对应$\sigma$的左右奇异向量。

行列式

$$d|A|=|A|tr(A^{-1}dA)$$
从而
$$d\ln |A|=tr(A^{-1}dA)$$

---

特殊矩阵

循环矩阵Circulant Matrix

行或列逐级轮换的矩阵
$$C=\left[\begin{array}{ccccc}{c}_0& c_{n-1}& \cdots & c_2& c_1\\ c_1& c_0& c_{n-1}& & c_2\\ ⋮& c_1& c_0& \ddots & ⋮\\ c_{n-2}& & \ddots & \ddots & c_{n-1}\\ c_{n-1}& c_{n-2}& \cdots & c_1& c_0\end{array}\right]$$
$C$是一个正规矩阵。
一种特殊的循环矩阵是
$$P_4=\left[\begin{array}{cccc}& 1& & \\ & & 1& \\ & & & 1\\ 1& & & \end{array}\right]$$
这是一个正交矩阵。该矩阵作用的效果是轮换
$$P_4\left[\begin{array}{c}{x}_0\\ x_1\\ x_2\\ x_3\end{array}\right]=\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\\ x_0\end{array}\right]$$

• 如果$C_1,C_2$是两个循环矩阵，那么$C_1{C}_2$仍然是循环矩阵。
  注意$P^2$相当于$P$对自己的列轮换，所以
  $$P_4^2=\left[\begin{array}{cccc}& & 1& \\ & & & 1\\ 1& & & \\ & 1& & \end{array}\right]$$

• 任意循环矩阵$C\in {\mathbb{R}}^{n\times n}$，可以分解为$C={\sum }_{i=0}^{n-1}{c}_i{P}_n^i$

• 循环矩阵是一种特殊的常对角矩阵Toeplitz matrix

应用：离散傅里叶变换DFT

正交矩阵一节提到了傅里叶变换，这里我们从循环矩阵的角度来考察。
思考$P_n$的特征值和特征向量，因为$P_n$的效果是轮换，所以可以想到特征值等分复平面单位圆，一个特征向量的每个元素依特征值等差取在复平面单位圆上。形式化：记$\omega =e^{2\pi /n}$，则特征值为$0,\omega ,{\omega }^2,\cdots ,{\omega }^{n-1}$，对应特征向量的矩阵为
$$V=\frac{1}{\sqrt{n}}\left[\begin{array}{ccccc}1& 1& 1& \cdots & 1\\ 1& \omega & {\omega }^2& \cdots & {\omega }^{n-1}\\ 1& {\omega }^2& {\omega }^4& \cdots & {\omega }^{2(n-1)}\\ & & & ⋮& \\ 1& {\omega }^{n-1}& {\omega }^{2(n-1)}& \cdots & {\omega }^{(n-1)(n-1)}\end{array}\right]$$
这恰是离散傅里叶变换的基。
另外，因为任意循环矩阵$C={\sum }_{i=0}^{n-1}{c}_i{P}_n^i$，所以特征向量和$P$都是相同的。可以证明特征值则满足$$diag[\Lambda ]=V\left[\begin{array}{c}{c}_0\\ c_1\\ ⋮\\ c_{n-1}\end{array}\right]$$

克罗内克积Kronecker Product（直积）

$$A\otimes B=\left[\begin{array}{ccc}{a}_{11}B& \cdots & a_{1n}B\\ ⋮& \ddots & ⋮\\ a_{m1}B& \cdots & a_{mn}B\end{array}\right]$$
其中$A\in {\mathbb{R}}^{m\times n},B\in {\mathbb{R}}^{p\times q},A\otimes B\in {\mathbb{R}}^{mp\times nq}$
有性质$(A\otimes B{)}^T=A^T\otimes B^T$
克罗内克积在求二阶导中常见，例如矩阵求导术（下）

---

参考文献：
[1] MIT18.06 Linear Algebra
[2] MIT18.065 Matrix Methods in Data Analysis, Signal Processing, and Machine Learning
[3] 3Blue1Brown 线性代数的本质
[4] Prince S J D. Computer vision: models, learning, and inference. Cambridge University Press, 2012.
[5] Goodfellow I, Bengio Y, Courville A. Deep learning. MIT press, 2016.
[6] Gonzales, Rafael C., and Richard E. Woods. Digital image processing. Fourth Edition 2018.
[7] Bishop C . Pattern Recognition and Machine LearningStat Sci. 2006.
[8] 周志华. 机器学习. 清华大学出版社. 2016.
[9] 同济大学数学系. 线性代数第五版. 高等教育出版社. 2007.
[10] 程云鹏等. 矩阵论. 西北工业大学出版社. 2006.
[11] 居余马等. 线性代数第2版. 清华大学出版社. 2002.
[12] 维基百科，如Woodbury matrix identity等