【基础方法】主成分分析（Principal Component Analysis）

我们在做图像处理的时候，为了避免提取的特征维数太高使得计算速度变慢，可以用下面几种方法降低特征的维数：

在一开始的时候对图片预处理，降低图片本身的维数。
在特征提取的过程中对特征进行处理，降低特征的的维数。（比如对LBP特征进行直方图统计，就有均衡模式的LBP这个方法来降低特征的维数）
采用一些降维德方法。（比如我们要介绍的PCA）

数学理论

PCA降维的思想方法是找到一个投影矩阵，这个矩阵对应了一个维数更低的空间，样本经过和这个投影矩阵运算，可以变成维数更小的矩阵，从而达到投影的目的。至于怎么这过程中的理论推导，用到了很多线性代数的知识，我深觉不能讲的更好，推荐一个现有的教程：PCA数学原理，我从中受益颇多。
然后我们来划一下重点。

问：首先我们的目的是什么呢？

答：是为了降低矩阵的维度。

问：那么怎么降低维度呢？

答：我们通过将现有的特征投影到维数更低的空间的基上，将维数高的向量空间用维数低的向量空间的基表示。

问：那么怎么求的我们想要的基呢？

答：不是随随便便一组基就可以的，我们需要找到在这组基上所有向量的离散度都能最大，即方差最大化，而且我们希望这些基是线性无关的，即它们的协方差为0。我们在计算协方差矩阵（$C=\frac{1}{m}XX^T$）的时候，会发现它的对角线上是对应的是不同维度的方差，而其他位置是一些维度的协方差，而且它是一个对称的矩阵。那么我们通过实矩阵对角化的方式，可以得到这个协方差矩阵对应的特征向量矩阵中的一部分就是我们所要的基。根据特征值的大小，我们将对应的特征矩阵进行排序，得到前k个特征向量组成的矩阵（每个特征向量为一个行向量），即我们所要求得的投影矩阵。

在这里补充一个链接，其中有个回答我觉得很有意思，大概就是关于和奶奶讲PCA的，讲的很浅显，很有意思，网上也有中文版。

算法步骤

假设有$m$条$n$维数据。

将原始数据按列组成$m$行$n$列矩阵$X$。
将$X$的每一行（代表一个属性字段）进行零均值化，即求得每一行的均值，并将每个元素减去对应行的均值。
求出均值化后的$X$的协方差矩阵$C$，其中$C=\frac{1}{m}XX^T$。
求出协方差矩阵$C$的特征值及对应的特征向量。
将特征向量按对应特征值大小从上到下按行排列成矩阵，取前$k$行组成投影矩阵$P$。
将矩阵$X$和$P$进行运算，得到新矩阵$Y$，其中$Y=PX$。

具体实现

matlab实现

我简单实现了一下PCA，学习到了很多新的matlab函数。

    function  [X, P]  = PCA( XMat, k, percent )

    if nargin <1
        help PCA;
    end
    if nargin < 2
        k=0;
        percent=0.9;
    end
    if nargin<3
        percent=0.9;
    end
    
    % 对所有特征进行归一化
    avg=mean(XMat);
    X_norm = bsxfun(@minus, XMat, avg);
    
    % 求协方差矩阵
    C= cov(X_norm);
    
    % 求协方差矩阵的特征值和特征向量
    [fvec,fvalue]=eig(C);
    
    % 将特征向量按对应特征值大小从上到下按行排列成矩阵
    fvalue=diag(rot90(rot90(fvalue)));
    fvec=rot90(fvec);
    
    % 如果未给出k值，则计算k值
    sumper=0;
    t=sum(fvalue);
    if k==0
        while sumper<percent
            k=k+1;
            sumper=sumper+fvalue(k)/t;
        end
    end
    
    P=fvec(1:k,:);
    X=(P*XMat')';
    end

matlab自带函数princomp()

我在matlab命令行中执行help命令时，得到了如下反馈：

    princomp Principal Components Analysis (PCA). 
    princomp will be removed in a future release. Use PCA instead. 
 
    Calls to princomp are routed to PCA.
 
    [coeff, score, latent, tsquare] = princomp(x,econFlag)
  
    See also: PCA

    
    
    princomp 的参考页

我更改目录用help pca试了一下，发现matlab本身就有pca函数，而对princomp()函数的调用实际上是调用pca()函数，所以这里关于函数的调用也只有一句说明：

    [coeff, score, latent, tsquare] = princomp(x,econFlag)

我在网上看了一些介绍，大致了解了下：

x：待处理矩阵，每列是一个维度，每行是一个样本。
coeff：原矩阵所对应的协方差阵V的所有特征向量组成的矩阵，即变换矩阵或称投影矩阵。
score：对主分的打分，也就是说原矩阵在主成分空间的表示。
latent：协方差矩阵的特征值。
tsquare：对每个样本点Hotelling的T方统计量。

跑了一下我写的和matlab自带的pca的效率对比，默默在一旁流泪。

python实现

    import numpy as np

    def pca(XMat, k):
        avg = np.mean(XMat, axis=0)#col average
        m, n = np.shape(XMat)#m: the number of samples  n:the size of vector
        avgs = np.tile(avg, (m, 1))#create a matrix of m*1, each row is a column average
        Xadjust=XMat-avgs
        Xcov=np.cov(Xadjust.T)#cov matrix
        feaVec, feaVal = np.linalg.eig(Xcov) # feature vector and feature value
        index=np.argsort(-feaVal)
        if k>n:
            print("k must lower than feature number")
            return
        selectVec = np.matrix(feaVec.T[index[:k]])
        finalData = Xadjust * selectVec.T
        return selectVec, finalData