PCA——主成分分析

　　PCA（Principal Components Analysis）主成分分析是一个简单的机器学习算法，利用正交变换把由线性相关变量表示的观测数据转换为由少量线性无关比变量表示的数据，实现降维的同时尽量减少精度的损失，线性无关的变量称为主成分。大致流程如下：

　　首先对给定数据集（数据是向量）进行规范化，使得数据集的平均值为0，方差为1（规范化是为了使数据散布在原点附近，而不是远离原点的某块区域，便于后面的计算）。之后对每个数据进行正交变换，把数据投影到几个少量的相互正交的方向（这些方向构成了数据空间的一个子空间）上。数据在每个方向上都有对应的坐标，而用这些方向和对应的坐标（坐标×方向的累加）就能近似表示原来高维的数据。因此这些方向的求解就是PCA的关键，这些方向有如下性质：

　　如果再由这些坐标通过对应的方向映射回原来的数据，精度损失是同等方向数量的方向集合（或者叫同维度的子空间吧）中最小的，而数据集在各个方向上的坐标的方差之和是同等方向数量的方向集合中最大的，也正对应着方差越大表示这个方向上保存数据的信息量越大（方差越小所有数据越集中在一个点上信息量当然越小）。数据集在这些方向的上的坐标的方差，从大到小排序，就把这每一个方向称为第一主成分、第二主成分……

证明

　　接下来证明什么样的方向是主成分，即什么样的方向集合能保存原数据集更多的信息，并进一步说明数据集在主成分上的坐标的方差之和最大。下面把数据映射为各个方向上的坐标称为编码（即降维），再反映射回来称为解码。

　　定义矩阵$Xin R^{m imes n}$为数据集，由$m$个$n$维行向量数据${x_1,x_2,...,x_m}$组成。求得数据均值与方差：

$displaystyleoverline{x}=sumlimits_{i=1}^{m}frac{x_i}{m}$

$displaystyle s_i=frac{1}{m-1}sumlimits_{j=1}^{m}left(x_{ji}-{overline{x}}_i ight)^2,i=1,2,ldots,n$

　　规范化样本：

$displaystyle x_{ij}^ast=frac{x_{ij}-overline{x_j}}{sqrt{s_j}},i=1,2,ldots,m,j=1,2,ldots,n$

　　获得规范化后的数据集$X^*$。随后PCA将把$X^*$编码成矩阵$Y^*in R^{m imes k},k<n$，通过右乘列向量单位正交的方向矩阵$Din R^{n imes k}$得到编码后的矩阵$Y^*$：

$Y^* = X^*cdot D$

　　而再解码回来是再乘上$D^T$：

$hat{X^*} = Y^*cdot D^{T}$

　　为了最小化精度损失，容易想到，优化$D$来最小化$hat{X^*}$与$X$之差的Frobenius范数：

$ egin{aligned} &D = displaystyle ext{arg}minlimits_D{||hat{X^*} - X^*||^2_F}\ &D = ext{arg}minlimits_D{|| X^*DD^T - X^*||^2_F}\ &D = ext{arg}minlimits_D{ ext{Tr}[(X^*DD^T - X^*)^T(X^*DD^T - X^*)]}\ &D = ext{arg}minlimits_D{ ext{Tr}[DD^TX^{*T}X^*DD^T-DD^TX^{*T}X^*-X^{*T}X^*DD^T+X^{*T}X^*]}\ end{aligned} $

　　去掉与$D$不相关的$X^{*T}X^*$项，再由迹的性质，矩阵乘法循环右移迹的值不变：

 $ egin{aligned} &D = ext{arg}minlimits_D{ ext{Tr}[X^*DD^TDD^TX^{*T}-X^*DD^TX^{*T}-D^TX^{*T}X^*D]}\ &D = ext{arg}minlimits_D{ ext{Tr}[X^*DD^TX^{*T}-X^*DD^TX^{*T}-D^TX^{*T}X^*D]}\ &D = ext{arg}maxlimits_D{ ext{Tr}[D^TX^{*T}X^*D]}\ end{aligned} $

　　容易发现当$D$的列向量取$X^{*T}X^*$的前$k$大特征值对应的特征向量时，有这个迹的值最大，等于前$k$大特征值之和。此时压缩的精度损失最少。我们又容易发现$X^{*T}X^*$实际上就是$X$行数据每个元素之间的协方差矩阵（没有除以$m-1$），再通过计算可以发现，映射之后的$Y^*$行向量的每个元素（也就是数据在每个方向上的坐标）的方差就是协方差矩阵的最大特征值，正印证了坐标方差越大信息量越大保存的信息越多。

　　另外，这个映射矩阵$D$的列向量在人脸识别中就是所谓的特征脸。

手动代码实现

　　实验代码：

import matplotlib.pyplot as plt
import pylab                #原本在jupyter里才能显示的图片，可以用窗口显示
import numpy  as np 

def PCA_img(img,n):
    t = np.mean(img,0)                                          #以行向量为压缩单元，求行平均
    s = np.std(img,0)                                           #求行标准差
    img=(img-t)/s                                               #每行减去平均值进行规范化
    v,vecters = np.linalg.eig(np.dot(img.T,img)/(len(img)-1))   #生成规范化矩阵后，转置乘自身再除以（行数-1），即为行向量各个元素的协方差矩阵，求协方差矩阵的特征值与特征向量
    indv = v.argsort()                                          #将特征值排序，从小到大
    vecters = vecters[:,indv]                                   #将特征向量按特征值大小排序
    vecters = vecters[:,-n:]                                    #选择前n大特征值对应的特征向量，作为映射矩阵。在人脸识别中就是特征脸
    img = np.dot(np.dot(img,vecters),vecters.T)                 #规范化矩阵先乘映射矩阵再乘映射矩阵的转置，也就是先压缩再解压缩
    img = np.multiply(img,s)+t                                  #把减掉的均值加回去变为原矩阵
    return img
        
img=plt.imread("aaa.jpg") #读取图片  

img.flags["WRITEABLE"] = True
img = img.astype(float) 
n = 50
img[:,:,0] = PCA_img(img[:,:,0],n)
img[:,:,1] = PCA_img(img[:,:,1],n)
img[:,:,2] = PCA_img(img[:,:,2],n)
fig = plt.figure()           
ax = fig.add_subplot(111)    
ax.imshow(img/255)  
print(img)           
pylab.show()              

　　原图是1200×1920的图：

　　压缩为1200×50后再解压回去：

sklearn包PCA

　　sklearn中计算PCA比直接使用numpy矩阵计算快得多，这是因为它里面实现了比较高效的随机算法，虽然每次算出来的值不一样，但是相差不大。代码如下：

import sklearn.decomposition as dc 

pca = dc.PCA(n_components = 10)       #规定PCA要取前k大特征值对应的特征向量的个数，即要映射到几维，PCA默认压缩行，比如100*100压缩到100*10
pca.fit(A)                            #将要用于压缩的矩阵传入，用来“学习”压缩矩阵
D = pca.components_                   #获取“学习”到的编码矩阵（压缩矩阵）D。原矩阵A经过AD^T会得到压缩后的矩阵（这里的D是上面推算的D^T）
C = pca.transform(B)                  #传入新的矩阵B，用“学习”到的编码矩阵进行压缩，即BD^T。在人脸特征提取中，同一个人的不同人脸照片，用同一个压缩矩阵压缩后的坐标向量通常相似，因此可以用于人脸识别的降维