PCA原理:
PCA的原理就是将原来的样本数据投影到一个新的空间中,相当于我们在矩阵分析里面学习的将一组矩阵映射到另外的坐标系下。通过一个转换坐标,也可以理解成把一组坐标转换到另外一组坐标系下,但是在新的坐标系下,表示原来的原本不需要那么多的变量,只需要原来样本的最大的一个线性无关组的特征值对应的空间的坐标即可。
下面是百度百科中对pca降维的一段解释,还是挺清晰的:
“对于一个训练集,100个对象模板,特征是10维,那么它可以建立一个100*10的矩阵,作为样本。求这个样本的协方差矩阵,得到一个10*10的协方差矩阵,然后求出这个协方差矩阵的特征值和特征向量,应该有10个特征值和特征向量,我们根据特征值的大小,取前四个特征值所对应的特征向量,构成一个10*4的矩阵,这个矩阵就是我们要求的特征矩阵,100*10的样本矩阵乘以这个10*4的特征矩阵,就得到了一个100*4的新的降维之后的样本矩阵,每个特征的维数下降了。
当给定一个测试的特征集之后,比如1*10维的特征,乘以上面得到的10*4的特征矩阵,便可以得到一个1*4的特征,用这个特征去分类。”
pca的实现(matlab)
我在网上看了很多pca降维的例子,都大同小异,原理差不多,都是活的原来矩阵的协方差矩阵,然后计算协方差矩阵的特征值和特征向量,最后通过特征向量的根据特征值由大到小的排序进行KL变换神马的获得一个转换矩阵。
1.matlab自带的实现方式
PCA在matlab中的实现举例
以下资料来自matlab的help,翻译和注解部分由笔者添加:(重点部分添加了翻译!)
princomp-----函数名称
Principalcomponentanalysis(PCA)ondata
Syntax------函数调用语法
[COEFF,SCORE]=princomp(X)
[COEFF,SCORE,latent]=princomp(X)
[COEFF,SCORE,latent,tsquare]=princomp(X)
[]=princomp(X,'econ')
Description-----函数描述
上面这个matlab函数的说明呢,只是引用百度百科,也可以看看matlab的函数说明,但是多少还是有点难懂。
我把我的理解简单的说说。
[COEFF,SCORE,LATENT,TSQUARED]=PRINCOMP(X)
上面这个函数,coeff矩阵是返回的转换矩阵,也就是把样本转换到新的空间中的准换矩阵,这个准换矩阵式比较大的,比如你的降维矩阵式30*100000,那么这个准换矩阵一般都是10000*29的维数。
score是原来的样本矩阵在新的坐标系中的表示,也就是原来的样本乘上转换矩阵,但是还不是直接乘,要减去一个样本的均值。将原来的数据转换到新的样本空间中的算法是这样实现的:
x0=bsxfun(@minus,x,mean(x,1));
score=x0*coeff;
然后就会得到和[COEFF,SCORE,LATENT,TSQUARED]=PRINCOMP(X)输出一样的score数据。同时这个也是原来的样本矩阵降维后的结果,如果使用降维后的数据就使用这个数据。一般情况下,如果你的每个样本的特征维数远远大于样本数,比如30*1000000的维数,princomp要加上'econ',就是princomp(x,'econ')这样使用,可以很大程度的加快计算速度,而且不会内存溢出,否则会经常报内存溢出。
[]=PRINCOMP(X,'econ')returnsonlytheelementsofLATENTthatare
notnecessarilyzero,,whenN=P,onlythefirstN-1,andthe
fasterwhenPN.
latent是返回的按降序排列的特征值,根据这个你可以手动的选择降维以后的数据要选择前多少列。
cumsum(latent)./sum(latent),通过这样计算特征值的累计贡献率,一般来说都选择前95%的特征值对应的特征向量,还是原来的矩阵30*1000000,如果你计算得到前25个特征值的累计贡献率已经超过99.9%,那么就完全可以只要降维后的数据的前25列
如果你需要对测试样本降维,一般情况下,使用matlab自带的方式,肯定需要对测试样本减去一个训练样本均值,因为你在给训练样本降维的时候减去了均值,所以测试样本也要减去均值,然后乘以coeff这个矩阵,就获得了测试样本降维后的数据。比如说你的测试样本是1*1000000,那么乘上一个1000000*29的降维矩阵,就获得了1*29的降维后的测试样本的降维数据。
princomp(x)使用的行表示一个样本,每行的所有的列数据都是这个样本的特征值。降维以后比如是30*29,那么每一行就是降维以后的数据。每个样本有29个特征值。
2.一个自实现的pca降维方式
1.%训练
2.%Lx=X'*X
3.clear;
4.clc;
5.train_path='..\Data\TrainingSet\';
6.phi=zeros(64*64,20);
7.fori=1:20
8.path=strcat(train_path,num2str(i),'.bmp');
9.Image=imread(path);
10.Image=imresize(Image,[64,64]);
11.phi(:,i)=double(reshape(Image,1,[])');
12.;
13.%mean
14.mean_phi=mean(phi,2);
15.mean_face=reshape(mean_phi,64,64);
16.Image_mean=mat2gray(mean_face);
17.imwrite(Image_mean,'','bmp');
18.%demean
19.fori=1:19
20.X(:,i)=phi(:,i)-mean_phi;
21.
22.Lx=X'*X;
23.tic;
24.[eigenvector,eigenvalue]=eigs(Lx,19);
25.toc;
26.%normalization
27.fori=1:19
28.%K-L变换
29.UL(:,i)=X*eigenvector(:,i)/sqrt(eigenvalue(i,i));
30.
31.%displayEigenface
32.fori=1:19
33.Eigenface=reshape(UL(:,i),[64,64]);
34.figure(i);
35.imshow(mat2gray(Eigenface));
36.
得到的均值图像mean_face:

前19个最大主元对应的“特征脸”:

测试:
测试用样本:

37.%使用测试样本进行测试
38.clc;
39.test_path='..\Data\TestingSet\';
40.error=zeros([1,4]);
41.fori=1:4
42.path=strcat(test_path,num2str(i),'.bmp');
43.Image=imread(path);
44.Image=double(imresize(Image,[64,64]));
45.phi_test=zeros(64*64,1);
46.phi_test(:,1)=double(reshape(Image,1,[])');
47.X_test=phi_test-mean_phi;
48.Y_test=UL'*X_test;
49.X_test_re=UL*Y_test;
50.Face_re=X_test_re+mean_phi;
51.calculateerrorrate
52.e=Face_re-phi_test;
53.
54.
55.%%displayfigure
56.Face_re_2=reshape(Face_re(:,1),[64,64]);
57.figure(i);
58.
59.imshow(mat2gray(Image));
60.title('Original');
61.figure(10+i);
62.imshow(mat2gray(Face_re_2));
63.title('Reconstruct');
64.error(1,i)=norm(e);
65.
66.%dispalyerrorrate
67.error_rate=error(1,i);
68.display(error_rate);
69.
重建出的测试样本与原样本的对比:

四副测试样本的重建误差分别为:
1.4195e+003
1.9564e+003
4.7337e+003
7.0103e+003
可见测试样本为人脸的样本的重建误差显然小于非人脸的重建误差。
下面先给出一下PCA的资料地址,供大家参考