奇异值分解（SVD）与 主成分分析（PCA）

线性映射

线性映射同构于矩阵乘。线性空间的一组基$\alpha =({\alpha }_1,\cdots ,{\alpha }_n)$，一个点的坐标为$x=(x_1,\cdots ,x_n)$，那么点可以记做两者内积$\alpha \cdot x={\sum }_{i=1}^n{x}_i{\alpha }_i$。注意，基不一定是向量，可以是任何线性无关的对象（比如，三角函数）。

线性映射$\mathcal{A}:\mathbb{U}\to \mathbb{V}$，在空间$U$下的基$({\alpha }_1,\cdots ,{\alpha }_n)$下和空间$V$下的基$({\beta }_1,\cdots ,{\beta }_m)$的矩阵是$A$，在空间$U$下的基$({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)$下和空间$V$下的基$({\tilde{\beta }}_1,\cdots ,{\tilde{\beta }}_m)$的矩阵是$B$，即
$$\begin{gathered} \mathcal{A}({\alpha }_1,\cdots ,{\alpha }_n)=({\beta }_1,\cdots ,{\beta }_m)A \\ \mathcal{A}({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)=({\tilde{\beta }}_1,\cdots ,{\tilde{\beta }}_m)B \end{gathered}$$
假设$({\alpha }_1,\cdots ,{\alpha }_n)$到$({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)$的过渡矩阵为$P$，$({\beta }_1,\cdots ,{\beta }_m)$到$({\tilde{\beta }}_1,\cdots ,{\tilde{\beta }}_m)$的过渡矩阵为$Q$，即
$$\begin{gathered} ({\alpha }_1,\cdots ,{\alpha }_n)P=({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n) \\ ({\beta }_1,\cdots ,{\beta }_m)Q=({\tilde{\beta }}_1,\cdots ,{\tilde{\beta }}_m) \end{gathered}$$
那么$B=QAP$（相抵），相抵矩阵代表相同的线性映射。对于基$\alpha$下坐标为$x$的点$P_1\in \mathbb{U}$，映射到了$P_2=\mathcal{A}(\alpha \cdot x)=\mathcal{A}(\alpha )\cdot x=({\beta }_1,\cdots ,{\beta }_m)Ax$，在基$\beta$下点$P_2\in \mathbb{V}$的坐标为$y=Ax$

线性变换$\mathcal{A}:\mathbb{V}\to \mathbb{V}$，空间$V$的两组基$({\alpha }_1,\cdots ,{\alpha }_n),({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)$，若
$$\begin{gathered} \mathcal{A}({\alpha }_1,\cdots ,{\alpha }_n)=({\alpha }_1,\cdots ,{\alpha }_n)A \\ \mathcal{A}({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)=({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)B \end{gathered}$$
假设从基$\alpha$到基$\tilde{\alpha }$的过渡矩阵为可逆方阵$P$，即
$$({\alpha }_1,\cdots ,{\alpha }_n)P=({\tilde{\alpha }}_1,\cdots ,{\tilde{\alpha }}_n)$$
那么$B=P^{-1}AP$（相似），相似矩阵代表相同的线性变换。对于基$\alpha$下坐标为$x$的点$P_1\in \mathbb{V}$，映射到了$P_2=\mathcal{A}(\alpha \cdot x)=\mathcal{A}(\alpha )\cdot x=({\alpha }_1,\cdots ,{\alpha }_n)Ax$，在基$\alpha$下点$P_2\in \mathbb{V}$的坐标为$y=Ax$

SVD

对角化：对于$n$维方阵$A$，如果存在$n$个线性无关的特征向量$w_1,\cdots ,w_n$，以及对应的特征值${\lambda }_1\le \cdots \le {\lambda }_n$，那么可以表示为：$A=W\Sigma W^{-1}$，其中$W=[w_1,\cdots ,w_n]$，$\Sigma =diag({\lambda }_1,\cdots ,{\lambda }_n)$

一个实对称矩阵（$A=A^T\in R^{n\times n}$），它满足：

1. 不同特征值对应的特征向量彼此正交
2. 若特征值$\lambda$的重数为$r$，那么就存在$r$个线性无关的特征向量$\lambda$（几何重数=代数重数）

因此，一个$n$阶实对称方阵中一定可以找到$n$个单位正交特征向量！或者说，存在$W^{T}W=I$（酉矩阵，这里是共轭转置）

对于长矩阵$A\in R^{m\times n}$，有

特征值分解：

1. 计算协方差矩阵$C=\frac{1}{n-1}{A}^{T}A\in R^{n\times n}$（无偏估计，不除以$n-1$或者除以$n$都不影响特征值和特征向量）
2. 协方差矩阵$C$是$n$阶实对称方阵，因此可以对角化：$C=W\Sigma W^T$

奇异值分解（Singular Value Decomposition，SVD）：

1. 对于矩阵$A^{T}A\in R^{n\times n}$，计算$n$个单位正交的特征向量（右奇异向量$v_i$），按列组合成酉矩阵$V\in R^{n\times n}$，对应的特征值为${\lambda }_i$

2. 对于矩阵$A{A}^T\in R^{m\times m}$，计算$m$个单位正交的特征向量（左奇异向量$u_i$），按列组合成酉矩阵$U\in R^{m\times m}$，对应的特征值为${\tilde{\lambda }}_i$

3. 根据${\sigma }_i=A{v}_i/u_i$或者${\sigma }_i^2={\lambda }_i={\tilde{\lambda }}_i$计算奇异值，组合成$\Sigma \in R^{m\times n}$，它的对角线是对应的奇异值

4. 最终得到$A=U\Sigma V^T$（易知，$AV=U\Sigma$，$A^{T}A=V{\Sigma }^T{U}^T\cdot U\Sigma V^T=V{\Sigma }^2{V}^T$，$A{A}^T=U\Sigma V^T\cdot V{\Sigma }^T{U}^T=U{\Sigma }^2{U}^T$）

一般地，特征值分解以及奇异值分解都将$\Sigma$中的特征值或奇异值按照从大到小的顺序排列。并且，奇异值会快速衰减（前10%甚至1%的奇异值的加和，可以占全部奇异值之和的99%以上），可用于压缩数据。

用numpy计算，

import numpy as np

A = np.array([[-1, 1, 1],
              [-4, 3, 2],
              [1, 0, 3]])

#特征值、特征向量
eigenvalue, featurevector = np.linalg.eig(A)	#计算
index = list(reversed(np.argsort(eigenvalue)))	#从大到小排序
eigenvalue = eigenvalue[index]
featurevector = featurevector.T[index]
print("特征值：\n", eigenvalue)
print("特征向量：\n", featurevector)

#逆矩阵
det = np.linalg.det(W)
print("det:",det)
W_inv = np.linalg.inv(W) 
print("W_inv:\n",W_inv)

#特征值的对角阵
Sigma = np.diag(eigenvalue) 
print("Sigma:\n",Sigma)

# A2 = A
A2 = W@Sigma@W_inv
print("A2:\n",A2)

其实，更简单的

# 特征值分解
S,W = np.linalg.eig(A)
print("\nW = \n",W)
print("\nS = \n",S)
print("\nA2 = \n",(W*S)@np.linalg.inv(W))	#W * diag(S) * W.inv

A = np.array([[-1, 1, 1, 5],
              [-4, 3, 2, -2],
              [1, 0, 3, 1]])

# 奇异值分解
U, S, VT = np.linalg.svd(A)
print("\nU = \n",U)
print("\nS = \n",S)
print("\nV.T = \n",VT)
Sigma = np.zeros(A.shape)
for i,s in enumerate(S):
	Sigma[i,i]=s
print("\nA2 = \n",U@Sigma@VT)	#U * diag(S) * V.T

PCA

主成分分析（Principal Component Analysis，PCA）是非常经典的降维算法。

对于$A\in R^{m\times n}$，它表示$m$维特征空间中的$n$个数据点，但特征的维度$m$过大

方法一：

1. 去中心化，将各特征（如$A_{ij}$）减去它们的均值（如$\frac{1}{n-1}{\sum }_{j=1}^n{A}_{ij}$），得到矩阵$X\in R^{m\times n}$
2. 计算协方差矩阵$C=\frac{1}{n-1}{X}^{T}X\in R^{n\times n}$
3. 做特征值分解，$C=W\Sigma W^T$
4. 选取最大的$k$个特征值，将$W$截取最左边的$k$个特征向量，按行组合成矩阵$P$
5. 做线性变换，$Y=PA\in R^{k\times n}$，这是$k\ll n$维空间中的$n$个数据点

方法二：

1. 去中心化，将各特征（如$A_{ij}$）减去它们的均值（如$\frac{1}{n-1}{\sum }_{j=1}^n{A}_{ij}$），得到矩阵$X\in R^{m\times n}$
2. 做奇异值分解，$X=U\Sigma V^T$
3. 选取最大的$k$个奇异值，将$U$截取最左边的$k$个特征向量，按行组合成矩阵$P$
4. 做线性变换，$Y=PA\in R^{k\times n}$，这是$k\ll n$维空间中的$n$个数据点

PCA的含义：由于截取的$k$个特征值较大，这意味着，在对应的特征向量的方向上的方差较大。由于$P$是由特征向量按行组合的，且这些特征向量彼此单位正交，所以$Y=PA$其实就是将$m$维空间中的$n$个数据点，正交投影到这些特征向量张成的$k$维子空间中（内积就是投影系数），矩阵$P\in R^{k\times m}$是从$m$维空间到$k$维子空间的线性变换。我们认为这$k$维子空间上的投影是消息本身，而另外$m-k$维补空间内的投影则是噪音。

选取合适的主成分个数$k$：
$$1-\frac{{\sum }_{i=1}^k{\Sigma }_{ii}}{{\sum }_{i=1}^n{\Sigma }_{ii}}\le t$$
这里的$t$是误差大小，选取$t=0.01$表示主成分保留了至少$99\%$的原始信息。

代码实现，

def PCA(X,k):
	'X是m*n长矩阵，m维空间中的n个数据点'
	X_mean = np.mean(X,1)			#每一行的均值
	X2 = np.array(X,dtype=float)
	for i,x in enumerate(X2):		#去中心化
		X2[i] = x-X_mean[i]
	U, S, VT = np.linalg.svd(X2)	#SVD
	tmp = 0
	for i in range(k):
		tmp += S[i]
	t = 1 - tmp/np.sum(S)			#误差大小
	P = U[:,:k].T
	return P@X, t					#（投影，误差）