主元分析和奇异值分解

主元分析和奇异值分解

主元分析或者主成分分析，即PCA（Principal Component Analysis）和奇异值分解，即SVD（Singular Value Decomposition）都是常用的对数据进行预处理进行降维的方法，通常我们认为样本的一个特征为一个维度，即特征向量中的维度，若样本的特征非常非常多，那么对于各种算法来说算起来都是很耗时耗内存的，这时就应该降维，同时要保持大部分原有的信息不丢失。

PCA

主成分分析是一种很常见的降维工具，其目标是在样本数据空间中找到新的坐标轴，即新的正交基，使得原有样本能够映射到新的坐标轴上，我们通常认为方差大的数据才是有用的数据，方差小的数据为噪声数据，则在新的空间选取方差较大的坐标轴，将方差小的坐标轴舍弃掉，这就实现了降维

例如数据有A和B两个特征，若A，B的变化是一样的，则完全可以舍弃掉其中一个。

步骤：输入样本集$D=${$x_1,x_2,\cdots ,x_m$}，每个样本有$d$个特征，输出为维数为$\hat{d},\hat{d}<d$的样本

1. 对所有样本进行中心化处理，即使得样本的均值为0。$x_i=x_i-\frac{1}{m}{\sum }_{i=1}^m{x}_i$
2. 计算样本的协方差矩阵$R$，其中$r_{ij}=Cov(x_i,x_j)$
3. 求出协方差矩阵$R$的所有特征值和特征向量
4. 对特征值进行降序排序，取前$\hat{d}$个特征值${\lambda }_1,{\lambda }_2,\cdots ,{\lambda }_{\hat{d}}$及其对应的特征向量$w_1,w_2,\cdots ,w_{\hat{d}}$，这些特征向量便构成映射后的空间的坐标轴
5. 对原有样本进行映射$\hat{D}=D\ast W$，$D$为$m\ast d$矩阵，$W$为$d\ast \hat{d}$矩阵
6. 特征向量$w_i$即映射后的坐标轴$\vec{i}$，其对于原有数据的方差的贡献率为$\frac{{\lambda }_i}{{\sum }_{j=1}^m{\lambda }_j}$

最后的输出维数$\hat{d}$可以根据数据的不同而不同，还可以设置一个阈值，例如$t=90$%，取使下式成立的最小的$\hat{d}$值
$$\frac{{\sum }_{i=1}^{\hat{d}}{\lambda }_i}{{\sum }_{i=1}^d{\lambda }_i}\ge t$$
原理推导：为什么是协方差的特征向量就是新空间的坐标轴呢？特征值越大代表什么？

首先我们的目标是找到新的正交基或者新的坐标轴来映射，并且在该正交基映射后方差最大，即映射后各样本之间的距离足够远。

例如：这个二维数据点映射在蓝色直线上，以该蓝色直线为新的坐标轴$U_1$，这样就可以由2维降到1维，且信息能够最大程度保留。若映射到红色这条线，虽然也是降维，但是数据之间的关系丢失了，信息也丢失了。

所以对于上面的二维数据，我们的目标就是找到一个方向基，使数据映射后的方差有最大值。当然这个只是二维降到一维，若三维降到二维，则还要找到新的方向基，此时就不能仅仅是再找映射后方差最大了，因为方差最大只有一条，那我想方差第二大呢？这也不行，因为这可能只与第一个方向基差一点点，几乎重合，所以我们要想第二个方向基于第一个方向基没有任何相关性，即这两个方向基的协方差为0，形象点就是要正交，加上这个约束条件我们再找方差最大的坐标轴。

• 协方差的定义是$E((x-E(x)(y-E(y)))$

补：这里要补充一下线性代数的理解，矩阵$A$，向量$x$，那么$Ax=y$，相当于$x$在$A$的作用下进行一次变换，更深的理解是以$A$的行向量为新的基，即为新的坐标轴，在这个新空间里，$x$就表示成$y$，这就好理解主成分分析里的映射关系了。

这样我们要找的新的空间来映射，这个新的空间就可以表示为一个矩阵$P$，我们不妨将$P^T$作为新空间的映射矩阵，为什么要转置呢？因为$P^T$的行向量为新的正交基，则$P$的列向量为新的正交基，而我们最后求出来的$P$正是由特征向量组成，特征向量就是新的基。

假设数据集有m个样本，d个特征，一般我们喜欢的样本矩阵是$X$是m*d的，但在这里我们需要$X^T$，即d*m的，这样矩阵$X^T$的每一列是一个样本，这样就将样本能通过映射矩阵进行映射，且有下列变化：
$$\begin{gathered} P^T:\hat{d}\ast d（每一行为一个坐标轴(基),有\hat{d}行，则映射后维度为\hat{d}） \\ {X}^T:d\ast m（每一列为一个样本，每个样本有d个特征） \\ 则P^T{X}^T将数据集进行映射，大小为\hat{d}\ast m \\ 为了符合习惯，我们转置一下大小为m\ast \hat{d}，则样本的特征减少d-\hat{d}个 \\ \hat{X}=(P^T{X}^T{)}^T=XP,X大小为m\ast \hat{d} \end{gathered}$$
那么$\hat{X}$就是我们所要求的降维后的数据集，它应该符合上述的要求，即每个数据的方差最大，且数据间没有相关性，即协方差为0，这怎么做到呢？我们要根据此求出$P$出来。

既涉及数据内方差又涉及数据间协方差，协方差矩阵很完美的能够表示出来了。我们的新数据集$\hat{X}$每一列就是在某一个坐标轴映射上的值，不妨令$\hat{X}=[{\hat{x}}_1,{\hat{x}}_2,\cdots ,{\hat{x}}_{\hat{d}}]$，${\hat{x}}_i$为一个列向量，表示所有样本在坐标轴$U_i$上的映射值，则$\hat{X}$的协方差矩阵为$D,d_{ij}$为$Cov({\hat{x}}_i,{\hat{x}}_j)$，可以看出$D$是一个对称矩阵，大小为$\hat{d}\ast \hat{d}$，对称轴是方差。

• 上述要求等价于这个协方差矩阵$D$对角元素最大，则非对角元素值为0。（条件一）

我们回到开始，PCA步骤一做的是数据的中心化，且这一步是必须的，这有什么用呢？就是为了方便求出协方差矩阵！中心化后的数据均值为0，那么无论是求方差还是求协方差，数据都要减去均值，那么正好减去0就相当于不用动了！

假设矩阵为
$$\left[\begin{array}{cc}{a}_1& b_1\\ a_2& b_2\\ \cdots & \cdots \\ a_d& b_d\end{array}\right]$$
对应我们的样本集$X$，为m*d，即每一行代表一个样本，不同字母代表不同特征，对上面矩阵第一列记为$A$，第二列记为$B$，为什么要记列呢，因为我们协方差矩阵算的是列与列之间的协方差，即不同坐标轴之间的协方差。并且该矩阵已经中心化处理，即每一列的均值为0.

则上述矩阵的协方差矩阵为
$$\left[\begin{array}{cc}{\sum }_{i=1}^d{a}_i^2& {\sum }_{i=1}^d{a}_i\ast b_i\\ {\sum }_{i=1}^d{a}_i\ast b_i& {\sum }_{i=1}^d{b}_i^2\end{array}\right]$$
这个矩阵可以看到$A$的方差为${\sum }_{i=1}^d{a}_i^2$，$B$的方差为${\sum }_{i=1}^d{b}_i^2$，$A$与$B$的协方差为${\sum }_{i=1}^d{a}_i\ast b_i$.

因为中心化之后均值为0那么方差和协方差就是元素值相乘.

并且进一步有公式若矩阵$X$的大小为m*d，则$\frac{1}{d}{X}^{T}X$就是$X$的协方差矩阵（要满足数据中心化）

我们回到条件一，记原有数据集为$X:m\ast d$，降维后数据集$\hat{X}:m\ast \hat{d}$，原有数据集协方差矩阵为$R:d\ast d$，降维后数据集协方差矩阵为$D:\hat{d}\ast \hat{d}$。降维变换矩阵为$P:d\ast \hat{d}$有以下变换：
$$\begin{gathered} D=\frac{1}{\hat{d}}{\hat{X}}^T\hat{X} \\ =\frac{1}{\hat{d}}(XP{)}^T(XP) \\ =\frac{1}{\hat{d}}{P}^T{X}^{T}XP \\ =P^T\frac{1}{\hat{d}}{X}^{T}XP \\ =P^{T}RP \end{gathered}$$
这里发现要想满足条件一，即协方差矩阵$D$对角元素最大，则非对角元素值为0，不就是要矩阵$R$做特征值分解吗？且协方差矩阵$D$对角线元素就是矩阵$R$的特征值，我们只要去最大的$\hat{d}$个特征值，就可以构成$D$，且降维变换矩阵$P$，就是那$\hat{d}$个最大的特征值对于的特征向量呀！

这里还要补充的是$R$是个对称矩阵，大小为d*d，那么肯定有d个不相关的特征向量，正好上述能够成立！

自此PCA就推导出来了，我们有原有数据集的协方差矩阵的特征向量作为变换矩阵。

并且特征值越大，则代表映射后的数据方差越大，效果越好，我们以根据特征值大小来比较贡献，当贡献率到达一定阈值时，新数据能够几乎代表原有数据且维度减少。下图则表示20个特征的原有样本的协方差矩阵可以有20个特征向量和特征值，我们只需取前6大的特征值，就能大幅度将原有样本的方差表示出来，将20降维到6。

代码实现：

def loadDataSet(filename):
    dataSet = []
    fr = open(filename)
    for line in fr.readlines():
        curLine = line.strip().split('\t')
        fltLine = map(float,curLine)
        dataSet.append(list(fltLine))
    return mat(dataSet)
def pca(dataMat,topNfeat=9999999):
    meanVals = mean(dataMat,axis=0) #对列取均值
    meanRemoved = dataMat-meanVals  #数据中心化
    print(shape(meanRemoved))
    covMat = cov(meanRemoved,rowvar=0)  #协方差矩阵，rowvar=0指以一行为一个样本
    print(shape(covMat))
    eigVals,eigVects = linalg.eig(mat(covMat)) #特征值和特征向量
    eigValInd = argsort(eigVals) #argsort是排序后返回index，升序
    eigValInd = eigValInd[:-(topNfeat+1):-1] #从最后一列到topNfeat+1列返回
    redEigVects = eigVects[:,eigValInd] #取出相应的特征向量
    lowDDataMat = meanRemoved*redEigVects #将原数据进行映射
    reconMat = (lowDDataMat*redEigVects.T) + meanVals 
    return lowDDataMat,reconMat

奇异值分解下一篇再讲