【模式识别】PCA原理、推导及实现

主成分分析

主成分分析(Principal Component Analysis, PCA) 是一种线性的特征提取方法,同时也是一种线性的降维(dimensionality reduction)方法.

PCA的原理推导可以从两方面进行(egin{cases} 最大化方差法 \最小平方误差end{cases}),其中最小平方误差法是通过求解一个线性函数实现对样本中每个点的拟合,所以设定的优化目标函数是到样本中所有元素的距离的总和最小;最大化方差法就是找到一组基向量,使得当原数据变换到以该组基向量组成的坐标中时可以保留最多的原始信息,而原始信息的多少是通过方差来衡量。

下文通过最大化方差法来推导步骤和原理。

1.前提知识:矩阵的线性变换

以坐标(1,1)(2,2)(3,3)为例,若将这三个数据点变换到以((frac 1{sqrt{2}},frac 1{sqrt{2}}))((-frac {1}{sqrt{2}},frac {1}{sqrt{2}}))为基的坐标中,则相应的坐标该是多少?(这种也叫在新的基空间中的投影值是多少?)

在原坐标空间中,基底是(1,0)(0,1),所以这几个数据点的坐标表示如下:

[egin{pmatrix} 1&0\ 0&1end{pmatrix}egin{pmatrix} 1&2&3\1&2&3end{pmatrix}=egin{pmatrix} 1&2&3\1&2&3end{pmatrix} ]

当空间的基底变为((frac 1{sqrt{2}},frac 1{sqrt{2}}))((-frac {1}{sqrt{2}},frac {1}{sqrt{2}}))时,则原数据点的坐标在此空间中可以表示为:

[egin{pmatrix} frac 1{sqrt{2}}&-frac 1{sqrt{2}}\ frac 1{sqrt{2}}&frac 1{sqrt{2}}end{pmatrix}egin{pmatrix} 1&2&3\1&2&3end{pmatrix}=egin{pmatrix} 0&0&0\frac 2{sqrt{2}}&frac 4{sqrt{2}}&frac 6{sqrt{2}}end{pmatrix} ]

从上述的计算可以看出,某原始数据点(用列向量表示)左乘一个基矩阵(用行向量表示),所得结果就是原始数据点在新的基坐标空间中的坐标表示(用列向量表示)或者叫投影值。

其实用向量内积来看的话,某数据点向量(坐标)与空间中的基向量分别作内积,则与每一个基向量都会得到一个内积的标量,而把这几个标量组合起来表示成向量形式,就是该数据点在这个基向量空间中的坐标表示(投影)

将以上形式推广,若有m个数据向量(列向量),有n个基向量,则这m个数据向量在这n个基向量的空间中的坐标表示(投影值)可以这样计算:

[egin{pmatrix} K_1\ K_2\K_3\vdots\K_nend{pmatrix}egin{pmatrix}P_1& P_2&P_3&cdots&P_m end{pmatrix}=egin{pmatrix} K_1P_1&K_1P_2&cdots&K_1P_m\ K_2P_1&K_2P_2&cdots&K_2P_m\K_3P_1&K_3P_2&cdots&K_3P_m\vdots&vdots&ddots&vdots\K_nP_1&K_nP_2&cdots&K_nP_mend{pmatrix}=A ]

若将上述计算所得矩阵记为A,则m个数据在新的n维空间中的坐标(投影)为矩阵A的列向量。

由上述也可以得到矩阵相乘的物理解释:两个矩阵相乘的意义是就是将右边矩阵中的每一个列向量变换到左边矩阵中以每一行行向量为基所表示的空间中去,也即矩阵乘法就是一种线性变换。

上述的矩阵相乘是一种线性变换,如果将左乘的基矩阵的基的个数减少,那么将原数据变换到该空间时,相应的 维数也会减少,这就是降维。注意:这里的减少基的个数是将Ki减少,而不是将Ki所代表的行向量的维数减少,要是行向量的维数减少了那就不能和原数据作内积了。

2.PCA算法推导

  • 问题的提出:减少哪几个基可以使数据尽可能多的保留原始信息?或者用什么来衡量原始信息保留的多少?

  • 在最大可分性中,用原数据左乘基矩阵后所得新的坐标数据(即投影值)的分散程度来表示。我们希望投影(线性变换)后的投影值可以尽可能的分散,因为要是有投影值(坐标值)重合了,那就代表着该数据点的信息被重复表示了,浪费掉了。从信息熵角度来看,就是熵越大,数据所含的信息也就越多。

  • 那么用什么来表示数据的分散程度呢?

    • 单个向量可以用方差来表示,多个向量可以用协方差来表示。有关方差和协方差细节见【方差协方差协方差矩阵】小节。

      当协方差为0的时候,也即表示两个变量之间不相关。在PCA降维时,就是选择几个基,当原始数据变换到该组基上时,使得各变量之间的协方差为零而变量自身的方差却很大。为了让协方差为0,选择第二个基的时候,与第一个基正交,第三个与第二个正交,这样两两正交的基,最后的相关性就是0。

      综上,如果要使某n维向量数据降维为一维,也即要找到一个基向量,那么就应该找到那个使得投影值的方差最大的基向量;如果要使某n维向量数据降维为k维,也即要找到k个基向量,那么就应该找到使得协方差矩阵中对角线上元素值最大,并且其他位置元素为零的k个基向量,因为协方差对角线上的元素就是投影后数据自己的方差,要保证最大,这样才可以尽可能多的表示原始数据信息,而其他位置的元素都是与其他向量之间的协方差要为零,因为不为零的话,两个向量之间相关,有重复表示的信息。

下面在计算协方差等时,是已经零均值化了的,零均值化问题见【方差协方差协方差矩阵】

记原始数据矩阵为P,记我们要寻找的一组基向量组成的矩阵为K,记将原始数据P经过矩阵K降维后所得投影矩阵为X,即X=KP。而我们所希望的使协方差矩阵对角线上的元素最大化,其余位置元素为0的协方差矩阵,就是投影后的矩阵X的协方差矩阵cov(X),而(cov(X)=frac 1{m-1}XX^T(零均值化后))所以我们的优化目标就是使得矩阵(frac 1{m-1}XX^T)的非对角元素为0,对角线上元素从大到小的排列

  • 那么什么样的矩阵K与原始数据矩阵P相乘得到的X的协方差可以有这样的效果呢?
    • 记原始数据矩阵P的协方差为cov(P),则有(cov(P)=frac 1{m-1}PP^T),看下面的推导:

[ egin{aligned} cov(X)&=frac 1{m-1}XX^T\ & =frac 1{m-1}KP(KP)^T \ & =frac 1{m-1}KPP^TK^T\ & =Kfrac 1{m-1}PP^TK^T\ & = Kcov(P)K^T \ end{aligned} ]

如上所见,我们要找的由基向量组成的矩阵K,就是使原数据的协方差矩阵变为投影后数据的协方差矩阵的变换矩阵。而投影后矩阵的协方差我们所期望的是,它的非对角元素全为零,对角线元素的值从大到小的排列(因为我们选出的基向量要使得投影后数据的方差尽可能的大嘛,而对角线元素就是方差)。故我们要找的矩阵K就是要使协方差矩阵cov(P)对角化,并且使对角线元素按从大到小的顺序排列的矩阵。即我们可直接求出协方差矩阵cov(P)的特征值和特征向量,由它的特征值组成的矩阵就是一个对角阵,而对应的特征向量即可作为变换矩阵K。所以有下面的分解因为协方差矩阵cov(P)是实对称矩阵,所以根据【特征值分解】就有正交矩阵Q和对角矩阵(Lambda)使得式

[Lambda=Q^Tcov(P)Q ]

成立,其中的Q就是由协方差矩阵cov(P)求出的特征向量组成的矩阵,(Lambda)为由特征值组成的对角矩阵。故选择k个基向量,就是选择矩阵cov(P)的最大的几个特征值对应的特征向量。 所以有如下步骤。

3.PCA算法步骤

设原始数据有m条,每条数据有n个特征,则可以得到PCA降维的算法步骤:

  1. 输入要降维的原始数据矩阵P,行为n(特征个数),列为m(样本个数);
  2. 计算每一行的均值,组合成向量(mu);[注意对于矩阵P来说,它的每一行都是不同的样本对应的相同位置的特征,我们求均值时,是求相同特征类型的均值,而不是求某一个样本的所有特征的均值,这一点在[方差协方差协方差矩阵]小节有说明.总的来说我们就是要求协方差矩阵,而求协方差矩阵就要这样求特征的均值,而不是样本的均值]
  3. 求取协方差矩阵(cov(P)=frac 1{m-1}(P-hatmu)(P-hatmu)^T);(若已经零均值化了,则不用减去均值,看下文补充)
  4. 求协方差矩阵cov(P)的特征值和特征向量;
  5. 将所求出的特征向量单位化,并将特征向量根据其所对应的特征值的大小按从大到小排序;
  6. 选取前k个特征值对应的特征向量组成矩阵K;
  7. (X=K^T(P-hatmu))即为降为k维后的数据。(选出的特征向量的个数即为降维后的维数,由矩阵的线性变化那部分知识知道,数据矩阵P要和特征向量相乘,所以要转置K)

补充:未零均值化的协方差计算:具体的也可看【方差协方差协方差矩阵】

[设P=egin{pmatrix}P_{11}& P_{12}&cdots&P_{1m}\ P_{21}&P_{22}&cdots&P_{2m}\ vdots&vdots&&vdots\ P_{n1}&P_{n2}&cdots&P_{nm}end{pmatrix} =egin{pmatrix}P_{1}\ P_{2}\ vdots\ P_{n}end{pmatrix}  设hatmu=egin{pmatrix}mu_{1}&mu_{1}&cdots&mu_{1}\ mu_{2}&mu_{2}&cdots&mu_{2}\ vdots&vdots&&vdots\ mu_{n}&mu_{n}&cdots&mu_{n}end{pmatrix}\ 则计算协方差矩阵,cov(P)=frac 1{m-1}(P-hatmu)(P-hatmu)^T\ =frac1{m-1}egin{pmatrix}P_{11}-mu_1& P_{12}-mu_1&cdots&P_{1m}-mu_1\ P_{21}-mu_2&P_{22}-mu_2&cdots&P_{2m}-mu_2\ vdots&vdots&&vdots\ P_{n1}-mu_n&P_{n2}-mu_n&cdots&P_{nm}-mu_nend{pmatrix}egin{pmatrix}P_{11}-mu_1& P_{21}-mu_2&cdots&P_{n1}-mu_n\ P_{12}-mu_1&P_{22}-mu_2&cdots&P_{n2}-mu_n\ vdots&vdots&&vdots\ P_{1m}-mu_1&P_{2m}-mu_2&cdots&P_{nm}-mu_nend{pmatrix}\ =frac1{m-1}egin{pmatrix}sum_{i=1}^m(P_{1i}-mu_1)(P_{1i}-mu_1)&sum_{i=1}^m(P_{1i}-mu_1)(P_{2i}-mu_2)&cdots&sum_{i=1}^m(P_{1i}-mu_1)(P_{ni}-mu_n)\sum_{i=1}^m(P_{2i}-mu_2)(P_{1i}-mu_1)&sum_{i=1}^m(P_{2i}-mu_2)(P_{2i}-mu_2)&cdots&sum_{i=1}^m(P_{2i}-mu_2)(P_{ni}-mu_n)\ vdots&vdots&&vdots\sum_{i=1}^m(P_{ni}-mu_n)(P_{1i}-mu_1)&sum_{i=1}^m(P_{ni}-mu_n)(P_{2i}-mu_2)&cdots&sum_{i=1}^m(P_{ni}-mu_n)(P_{ni}-mu_n) end{pmatrix}\ =egin{pmatrix}cov(P_1,P_1)&cov(P_1,P_2)&cdots&cov(P_1,P_n)\ cov(P_2,P_1)&cov(P_2,P_2)&cdots&cov(P_2,P_n)\vdots&vdots&ddots&vdots\ cov(P_n,P_1)&cov(P_n,P_2)&cdots&cov(P_n,P_n)end{pmatrix} ]

4.PCA的一些补充

总的来说,PCA就只需要计算训练样本的协方差矩阵的特征值和特征向量,然后由这些特征向量得到一组基向量构成新的坐标系,然后将特征向量在这些基向量上投影即可达到降维的目的。

  • 协方差矩阵计算得到的特征值和特征向量都是实数吗?要是有复数怎么比大小?
    • 因为协方差矩阵是一个实对称矩阵,所以特征值都是实数,并且由于协方差矩阵是半正定矩阵,所以其所有的特征值都是大于等于0的实数,其特征向量也是实数。此外,上文也提到了,在选择基向量(也即后文中的特征向量)的时候,为了让协方差非对角元素为零,所以第一个基向量与第二个正交,以此类推,选的所有的基向量是两两正交的,也即最后的K矩阵是个正交矩阵。(或者说,相同特征值对应的特征向量可以正交化处理,所以是正交的),也即由主成分分析法得到的新的坐标系是一个直角坐标系。
  • 如果将所有的特征向量都选上了,那么经过PCA变换就没了降维,但是还仍有去除相关性的方法,因为PCA变换后的数据的协方差是一个对角阵。所以说经过PCA变换就会消除样本集在各个特征之间的相关性。

  • 用降维后的数据还原原始数据:

    • 由上文也可以知道,最后我们选好的基矩阵K是一个正交矩阵,所以从变换(X=K^T(P-hatmu))求出P为:

[ K^TP-K^Thatmu=X\ K^TP=X+K^Thatmu\ P=K(X+K^Thatmu) =KX+hatmu ]

  • PCA降维,降到多少维度才合适?
    • 因为投影值的方差等于特征值,所以当一个特征值为0,那么它所对应的特征向量在保持原始数据的分布信息方面完全不起作用。同理,当特征值也非常小的时候,我么也有认为该特征向量起的作用也很小,所以,为了保留一个合理的比例的方差,常用的经验法则是:若累计的特征值之和超过所有的特征值之和的90%,那么就在此停止。有时也用0.85和0.95.

5.PCA的MATLAB实现

function [new_data,check] = PCA(data,dimension)
%% PCA降维和去相关性的实现
%输入:data为要降维的数据矩阵,行应该为特征个数,列为样本个数;
       %dimension为要降到多少维度,降得维度一般为保持求出的特征值之和的90%;
%输出:输出new_data为降维后的数据;
       %输出check为降维后的数据的协方差矩阵,计算所得应该是个对角阵,用来检查计算过程的准确性;
%%
%计算均值,注意是计算行向量的均值
[row,column]=size(data);
mu=zeros(row,1);
for i=1:row
    mu(i)=mean(data(i,:));
end

%减去均值,data_mu代表零均值后的数据
data_mu=zeros(row,column);
for i=1:row
    for j=1:column
        data_mu(i,j)=data(i,j)-mu(i);
    end
end

%计算协方差矩阵,cov_data为零均值后的数据的协方差矩阵
cov_data=(data_mu*data_mu')/(length(data_mu)-1);
%计算协方差矩阵的特征值与特征向量
[eigenvector,eigenvalue]=eig(cov_data);%eigenvector,eigenvalue分别为特征向量(列)和特征值
%按照特征值从大到小的顺序排列,选出对应的特征向量
[~,index]=sort(diag(eigenvalue),'descend');
basevector=eigenvector(:,index(1:dimension));
%则可以得到降维后的数据new_data
new_data=basevector'*data_mu;
%计算一下降维后的数据的协方差是不是对角阵
check=cov(new_data');
end
原文地址:https://www.cnblogs.com/LENMOD/p/13597371.html