标签：... python dfrac 矩阵 降维 pmatrix PCA 数据
背景与原理：
PCA（主成分分析）是将一个数据的特征数量减少的同时尽可能保留最多信息的方法。所谓降维，就是在说对于一个$n$维数据集，其可以看做一个$n$维空间中的点集（或者向量集），而我们要把这个向量集投影到一个$k<n$维空间中，这样当然会导致信息损失，但是如果这个$k$维空间的基底选取的足够好，那么我们可以在投影过程中尽可能多地保留原数据集的信息。
数据降维的目的在于使得数据更直观、更易读、降低算法的计算开销、去除噪声。
接下来我们讨论下如何选取$k$维子空间：
假设原数据集有$m$条数据，每条数据有$n$维，那么可以将其拼成一个$n*m$的矩阵M，而我们想投影到的$k$维空间的一个单位正交基底为$(p_{1},...,p_{k})$，那么我们想把这$m$维向量投影到这个空间中实际上就是进行一次矩阵乘法$\begin{pmatrix} p_{1}\\p_{2}\\...\\p_{k} \end{pmatrix} M$
这个道理是简单易懂的。对于一个向量$\alpha$，其在另一个向量$\beta$方向上的投影是$\dfrac{\alpha \cdot \beta }{|\beta|}$（高中数学）
如果$|\beta|$是一个单位向量，那么这个投影即为$\alpha \cdot \beta=\beta^{T} \alpha$，于是投影到$k$个单位向量为基底的空间中的情况即如上述所示。
因此我们要找到这$k$个单位向量作为基底，然后拼出$P=\begin{pmatrix} p_{1}\\p_{2}\\...\\p_{k} \end{pmatrix}$即可。
那么怎么找呢？我们考虑降维之后我们需要什么，由于我们在降维之后要尽可能多地保留原始信息，因此降维之后的数据要提供最大的信息量，那么这个信息量在这里可以用数据的方差来反映，方差越大，数据的离散程度越高，那么数据的自身特征保留的就越好（个人理解：PCA降维的目的在于突出数据的个体特征，减少信息损失，而如果降维之后数据离散程度低，意味着这些数据全都堆在一起，那数据的特征体现的就不明显了——所有数据全都差不多，这样个体信息保留的就不好了）。
对于一个特征，在$m$组数据中的方差为$\sigma^{2}=\dfrac{1}{m}\sum_{i=1}^{m}(x_{i}-\overline{x})^{2}$，为了便于讨论，我们对所有特征零均值化（即把每个$x_{i}$预先减去$\overline{x}$），这样一个特征的方差即为$\sigma^{2}=\dfrac{1}{m}\sum_{i=1}^{m}x_{i}^{2}$
但是降维过程中只考虑方差是不够的——如果我们发现两个特征之间有很强的线性相关性，那么这两个特征其实差别就不大了，我们当然不需要同时保留这两个特征，因此我们还希望降维之后任意两个特征的协方差（$cov(a,b)=\dfrac{1}{m}\sum_{i=1}^{m}a_{i}b_{i}$，因为已经进行过零均值化了）为零，也就是在说我们选取的子空间的基底一定是正交的。
那么现在的问题就转化成了：对于一个$n$维$m$组数据的$n*m$数据矩阵$X$，我们希望将其投影到$n$维空间的一个$k$维子空间中，因此我们要找到$k$个单位正交基$(p_{1},...,p_{k})$，而如果这$k$个单位正交基构成的矩阵$P=\begin{pmatrix} p_{1}\\p_{2}\\...\\p_{k} \end{pmatrix}$，那么投影过程即为$Y=PX$，$Y$即为降维后所得的$k$维数据集
而结合上述讨论，我们希望$Y$各个特征的方差最大，同时$Y$的两特征的协方差为零，这怎么操作呢？
对于一个$n$维有$m$组数据的$n*m$矩阵$X$，我们考察$C=\dfrac{1}{m}XX^{T}$，那么我们看到如果$X=\begin{pmatrix} x_{11} & x_{12} & ... & x_{1m}\\...\\x_{n1} & x_{n2} &... & x_{nm}\end{pmatrix}$，我们有：
$XX^{T}=\begin{pmatrix} \sum_{i=1}^{m}x_{1i}^{2} & \sum_{i=1}^{m} x_{1i}x_{2i} &...& \sum_{i=1}^{m}x_{1i}x_{ni}\\...\\ \sum_{i=1}^{m}x_{ni}x_{1i} & \sum_{i=1}^{m} x_{ni}x_{2i} &...& \sum_{i=1}^{m}x_{ni}^{2}\end{pmatrix}$
我们称$\dfrac{1}{m}XX^{T}$为协方差矩阵，因为我们看到按照我们上面的解释，这个矩阵是一个实对称矩阵，其主对角线上的元素是一个特征维度的方差，而其余位置上的元素是两个对应特征的协方差！
那么我们的目的是要最大化主对角线上的元素，同时让其余位置上的元素为$0$，那么我们进行的不就是实对称矩阵的正交相似对角化嘛！
形式化地解释一下：我们设$Y$的协方差矩阵为$D$，那么我们希望$D$是一个对角矩阵，同时$D$的主对角线上的元素要尽可能大，那么我们有：
$D=\dfrac{1}{m}YY^{T}=\dfrac{1}{m}(PX)(X^{T}P^{T})=P(\dfrac{1}{m}XX^{T})P^{T}$
那么我们实际进行的不就是把$C=\dfrac{1}{m}XX^{T}$这个协方差矩阵正交相似对角化嘛！
至于我们希望主对角线元素尽可能大，那我们就选取$C$的前$k$大的特征值组成$D$就好了嘛，而此时的$P$就对应于前$k$大的特征值对应的$k$个正交的特征向量构成的矩阵。
那么我们的算法步骤如下：
对于$n$行$m$列的矩阵$X$，我们解释成其有$m$组数据，每组数据有$n$个特征，现在我们欲将其变成$k*m$的矩阵$Y$，表示降维后每组数据只有$k$个特征。
（1）零均值化：对$X$的每个元素，减去自己所在行的均值（即我们是逐特征操作，一行对应于同一个特征，不要搞错这一点）
（2）计算协方差矩阵$C=\dfrac{1}{\textbf{m}}XX^{T}$
（3）对协方差矩阵对角化$C=P\Sigma P^{T}$，找到其单位正交的特征向量$e_{1},...,e_{n}$
（4）选取最大的$k$个特征值对应的特征向量$e_{i_{1}},...,e_{i_{k}}$，拼成一个变换矩阵$P_{k}=\begin{pmatrix} e_{i_{1}}\\e_{i_{2}}\\...\\e_{i_{k}} \end{pmatrix}$
（5）降维后的数据即为$Y=P_{k}X$
（6）如果希望根据降维后的数据集$Y$近似还原原数据集$\hat{X}$，我们有$\hat{X}=P_{k}^{T}Y$（这里的逻辑是如果我们不降维，那么$P_{k}=P$就是一个正交矩阵，那么$P^{T}P=I$，相当于此时数据集没有损失，那么类比这个过程就能导出近似还原方法$\hat{X}=P_{k}^{T}Y$）
代码实现：
import numpy as np
from sympy.matrices import Matrix,GramSchmidt
np.random.seed(1)
x = 7*np.random.rand(100)
y = 0.5*x + 1 + 3*np.random.rand(100)
X = np.hstack([x.reshape(100, 1), x.reshape(100, 1), y.reshape(100, 1), x.reshape(100, 1)])
def centerData(X):
X = X.copy()
X -= np.mean(X, axis=0)
return X
X = centerData(X)print(X[7][2])
C= (np.transpose(X)@X)/100
val,fea=np.linalg.eig(C)
dic=dict()
for i in range(0,4):
dic[val[i]]=fea[:,i]
val=abs(np.sort(-val))
P=np.vstack([dic[val[0]],dic[val[1]]])
Y=X@P.T
reconstruct_X=Y@P
print(reconstruct_X[7][2])
值得注意的问题是这段代码中生成的数据每组数据是一个行向量，每列对应于一个特征，因此所有的计算和上面的推导都构成一个转置。
此外这里使用了numpy里面的linalg.eig方法用来求一个实对称矩阵的特征值和特征向量，返回的val是特征值，fea是特征向量，要特别说明的是不出意外的情况下这里的val都是按从大到小排序的，而fea实际上是一个矩阵，这个矩阵每个列向量对应于一个特征值，因此一定要注意选取的方法。上面代码使用前两个特征向量作为主成分恢复原矩阵，可以看到恢复效果还是不错的。
import numpy as np from sklearn.decomposition import PCA import matplotlib.pyplot as plt from sklearn.datasets import load_iris data=load_iris() X=data.data Y=data.target pca=PCA(n_components=2) X_2d=pca.fit_transform(X) plt.scatter(X_2d[:,0],X_2d[:,1],c=Y) plt.show()
当然，PCA也可以直接使用sklearn里面的包，上述代码加载了经典的鸢尾花数据集，然后进行PCA降维（降成二维，这个n_components参数给出了要降到的维度），然后能清楚看到三个鸢尾花的类别，效果很好。
标签：...,python,dfrac,矩阵,降维,pmatrix,PCA,数据 来源： https://www.cnblogs.com/zhangleo/p/16076052.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。