matlab pca基礎知識

PCA的一些基本資料

最近因為最人臉表情識別，提取的gabor特征太多了，所以需要用PCA進行對提取的特征進行降維。

 

本來最早的時候我沒有打算對提取的gabor特征進行降維，但是如果一個圖像時64*64,那么使用五個尺度八個方向的gabor濾波器進行濾波，這樣提取的特征足足有64*64*5*8這么多，如果圖像稍微大一點，比如128*128的圖像，那么直接提取的特征就會幾十萬，所以不降維的話直接用SVM訓練分類器是非常困難的。

 

所以在這段時間我就學習了一下PCA降維的基本原理和使用方法，網上給出的資料都比較亂，而且很不清楚，經過這幾天的學習和測試，終於把調理弄清楚了，給大家分享一下，下面只是我對於PCA的個人理解，肯定有不對的地方，還請各位大牛多多指教。

 

下面先給出一下PCA的資料地址，都是我收集的：

http://hi.baidu.com/yicomrdztxbeiwd/item/913f28c05cf7ebc4994aa06f

http://blog.sciencenet.cn/blog-265205-544681.html

http://blog.csdn.net/mpbchina/article/details/7384425

http://blog.sina.com.cn/s/blog_6833a4df0100pvk7.html

http://stackoverflow.com/questions/4991343/matlab-principal-component-analysis-eigenvalues-order

http://stackoverflow.com/questions/10400230/what-is-score-in-princomp

http://www.mathworks.com/matlabcentral/newsreader/view_thread/152608

http://stats.stackexchange.com/questions/27572/matlab-princomp-latent

http://www.nlpca.org/pca-principal-component-analysis-matlab.html

http://www.matlabsky.com/thread-11751-1-1.html

http://stackoverflow.com/questions/10818718/principal-component-analysis

http://www.mathworks.cn/cn/help/stats/princomp.html

http://www.mathworks.cn/cn/help/stats/pca.html#bti6n7k-2

http://lovelittlebean.blog.163.com/blog/static/116582186201181213911729/

http://www.ilovematlab.cn/thread-54493-1-1.html

http://www.ilovematlab.cn/forum.php?mod=viewthread&tid=146626

http://www.ilovematlab.cn/forum.php?mod=viewthread&tid=204069

http://www.ilovematlab.cn/forum.php?mod=viewthread&tid=54600

http://search.discuz.qq.com/s/aa8585553/princomp+%E9%99%8D%E7%BB%B4.html

http://www.ilovematlab.cn/thread-68796-1-1.html

http://www.ilovematlab.cn/thread-209229-1-1.html

http://www.ilovematlab.cn/thread-209229-1-1.html

http://blog.sina.com.cn/s/blog_61c0518f0100f4mi.html

http://blog.csdn.net/haitao111313/article/details/7875392

http://media.cs.tsinghua.edu.cn/~ahz/digitalimageprocess/chapter11/chapt11_ahz.htm

http://hi.baidu.com/845777018/item/7438e555df1138404fff2011

http://en.wikipedia.org/wiki/Principal_component_analysis

http://baike.baidu.com/view/852194.htm

http://wenku.baidu.com/view/bd9284fcfab069dc51220107.html

http://wenku.baidu.com/view/c0bde56da98271fe910ef9b8.html

http://wenku.baidu.com/view/9f69930790c69ec3d5bb75d3.html

http://www.ilovematlab.cn/thread-54600-1-1.html

http://www.cnblogs.com/sunwufan/archive/2011/08/31/2159952.html

http://zhidao.baidu.com/question/416895922.html

上面的網址都是一些pca原理啊，實現什么的介紹。

 

具體的PCA的算法的理論基礎呢，我這里就不詳細說了，因為我也沒有看具體詳細，所以如果想要徹底的弄明白PCA的工作原來，還是請到wiki上看吧，寫的非常清晰，我因為臨時用一下，就寫個大致的原理就可以了。

 

ＰＣＡ原理：

 

PCA的原理就是將原來的樣本數據投影到一個新的空間中，相當於我們在矩陣分析里面學習的將一組矩陣映射到另外的坐標系下。通過一個轉換坐標，也可以理解成把一組坐標轉換到另外一組坐標系下，但是在新的坐標系下，表示原來的原本不需要那么多的變量，只需要原來樣本的最大的一個線性無關組的特征值對應的空間的坐標即可。

比如，原來的樣本是30*1000000的維數，就是說我們有30個樣本，每個樣本有1000000個特征點，這個特征點太多了，我們需要對這些樣本的特征點進行降維。那么在降維的時候會計算一個原來樣本矩陣的協方差矩陣，這里就是1000000*1000000，當然，這個矩陣太大了，計算的時候有其他的方式進行處理，這里只是講解基本的原理，然后通過這個1000000*1000000的協方差矩陣計算它的特征值和特征向量，最后獲得具有最大特征值的特征向量構成轉換矩陣。比如我們的前29個特征值已經能夠占到所有特征值的99%以上，那么我們只需要提取前29個特征值對應的特征向量即可。這樣就構成了一個1000000*29的轉換矩陣，然后用原來的樣本乘以這個轉換矩陣，就可以得到原來的樣本數據在新的特征空間的對應的坐標。30*1000000 * 1000000*29 = 30 *29， 這樣原來的訓練樣本每個樣本的特征值的個數就降到了29個。

 

一般來說，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降維的函數的說明，就是princomp，這個函數的說明也是用的樣本的個數多余特征的維數。后來經過試驗是證實，證實了那個浙大的博士的認識是錯誤的，pca降維肯定不需要樣本的個數大於特征的維數，要不然還降維個什么意思。比如我有30*1000000的特征矩陣，那么降維后肯定是每個樣本在新的空間中的表示的特征維數不超過30.

 

迷惑二

   

另外一個迷惑，在最初剛開始做的時候，就是為什么這么大的數據，比如30*1000000直接就降到了30*29，這不是減少的數據有點太多了么，會不會對性能造成影響。之所以有這個迷惑，是因為最初並不了解pca的工作方式。 pca並不是直接對原來的數據進行刪減，而是把原來的數據映射到新的一個特征空間中繼續表示，所有新的特征空間如果有29維，那么這29維足以能夠表示非常非常多的數據，並沒有對原來的數據進行刪減，只是把原來的數據映射到新的空間中進行表示，所以你的測試樣本也要同樣的映射到這個空間中進行表示，這樣就要求你保存住這個空間坐標轉換矩陣，把測試樣本同樣的轉換到相同的坐標空間中。

 

    有些同學在網上發帖子問對訓練樣本降維以后，怎么對測試樣本降維，是不是還是使用princomp這個函數進行降維，這個是錯誤的。如果你要保證程序運行正常，就要保證訓練樣本和測試樣本被映射到同一個特征空間，這樣才能保證數據的一致性。

 

迷惑三

 

網上有不同的pca降維的代碼，每個代碼也實現的不一樣，那么對於同一個數據是否是pca降維以后都是獲得相同的數據呢，也就是說不管你用哪種方式進行pca降維，不管你是從哪里下載到的或者自己根據算法實現的pca降維，同樣的矩陣降維以后的數據是否一致？這個我個人認為，不同的算法最后導致的pca降維的數據肯定不一致。因為pca降維以后，只是把原來的數據映射到新的特征空間，所以如果你的算法不同，那么選擇的協方差矩陣肯定就不同，最后獲得的轉換矩陣肯定也不一樣。那么訓練樣本和測試樣本和不同的轉換矩陣相乘以后最終肯定會獲得不同的降維坐標。所以使用不同的算法應該最后不會有相同的坐標結果，這個也是我一直實驗的結果，我也使用了matlab自帶的princomp降維，並且使用相同的數據使用網上下載的一些降維方法進行降維，得到的數據都不一致。

 

比如說princomp這個matlab自帶的函數，在降維之前就將每一個樣本減去了一個所有樣本的平均值，也可能有很多樣本沒有減去平均值。princomp這里使用一行表示一個樣本，每行包括這個樣本的所有的特征值。而網上大部分都是每一列表示一個樣本，這樣這一列的所有行都表示這個樣本的特征值。網上的程序使用列表示樣本是有一定好處的，比如我的樣本是1000000*30，總共有30個訓練樣本，每個樣本的特征值個數是1000000，那么這個矩陣獲得的協方差矩陣是30*30，計算起來非常的方便，不想30*1000000這樣的矩陣獲得到的協方差矩陣式1000000*1000000，直接就內存溢出了，不過matlab有自己的實現方式，巧妙的解決了這個問題。

 

 

pca的實現（matlab）

 

我在網上看了很多pca降維的例子，都大同小異，原理差不多，都是活的原來矩陣的協方差矩陣，然后計算協方差矩陣的特征值和特征向量，最后通過特征向量的根據特征值由大到小的排序進行KL變換神馬的獲得一個轉換矩陣。

 

1. matlab自帶的實現方式

 

　PCA在matlab中的實現舉例

 

　　以下資料來自matlab的help，翻譯和注解部分由筆者添加：(重點部分添加了翻譯！)

 

　　princomp-----函數名稱

 

　　Principal component analysis (PCA) on data

 

　　Syntax------函數調用語法

 

　　[COEFF,SCORE] = princomp(X)

 

　　[COEFF,SCORE,latent] = princomp(X)

 

　　[COEFF,SCORE,latent,tsquare] = princomp(X)

 

　　[...] = princomp(X,'econ')

 

　　Description -----函數描述

 

　　COEFF = princomp(X) performs principal components analysis (PCA) on the n-by-p data matrix X, and returns the principal component coefficients, also known as loadings. Rows of X correspond to observations, columns to variables. COEFF is a p-by-p matrix, each column containing coefficients for one principal component. The columns are in order of decreasing component variance.

 

　　在n行p列的數據集X上做主成分分析。返回主成分系數。X的每行表示一個樣本的觀測值，每一列表示特征變量。COEFF是一個p行p列的矩陣，每一列包含一個主成分的系數，列是按主成分變量遞減順序排列。(按照這個翻譯很難理解，其實COEFF是X矩陣所對應的協方差陣V的所有特征向量組成的矩陣，即變換矩陣或稱投影矩陣，COEFF每列對應一個特征值的特征向量，列的排列順序是按特征值的大小遞減排序，后面有具體例子解釋，見說明1)

 

　　princomp centers X by subtracting off column means, but does not rescale the columns of X. To perform principal components analysis with standardized variables, that is, based on correlations, use princomp(zscore(X)). To perform principal components analysis directly on a covariance or correlation matrix, use pcacov.

 

　　計算PCA的時候，MATLAB自動對列進行了去均值的操作，但是並不對數據進行規格化，如果要規格化的話，用princomp(zscore(X))。另外，如果直接有現成的協方差陣，用函數pcacov來計算。

 

　　[COEFF,SCORE] = princomp(X) returns SCORE, the principal component scores; that is, the representation of X in the principal component space. Rows of SCORE correspond to observations, columns to components.

 

　　返回的SCORE是對主分的打分，也就是說原X矩陣在主成分空間的表示。SCORE每行對應樣本觀測值，每列對應一個主成份(變量)，它的行和列的數目和X的行列數目相同。

 

　　[COEFF,SCORE,latent] = princomp(X) returns latent, a vector containing the eigenvalues of the covariance matrix of X.

 

　　返回的latent是一個向量，它是X所對應的協方差矩陣的特征值向量。

 

　　[COEFF,SCORE,latent,tsquare] = princomp(X) returns tsquare, which contains Hotelling's T2 statistic for each data point.

 

　　返回的tsquare，是表示對每個樣本點Hotelling的T方統計量(我也不很清楚是什么東東)。

 

　　The scores are the data formed by transforming the original data into the space of the principal components. The values of the vector latent are the variance of the columns of SCORE. Hotelling's T2 is a measure of the multivariate distance of each observation from the center of the data set.

 

　　所得的分(scores)表示由原數據X轉變到主成分空間所得到的數據。latent向量的值表示SCORE矩陣每列的方差(見說明2)。Hotelling的T方是用來衡量多變量間的距離，這個距離是指樣本觀測值到數據集中心的距離。

 

　　When n <= p, SCORE(:,n:p) and latent(n:p) are necessarily zero, and the columns of COEFF(:,n:p) define directions that are orthogonal to X.

 

　　[...] = princomp(X,'econ') returns only the elements of latent that are not necessarily zero, and the corresponding columns of COEFF and SCORE, that is, when n <= p, only the first n-1. This can be significantly faster when p is much larger than n.

 

　　當維數p超過樣本個數n的時候，用[...] = princomp(X,'econ')來計算，這樣會顯著提高計算速度

 

　　Examples--舉例

 

　　(上面說了那么多廢話，看了還不一定懂，還不如舉例容易理解，下面樣本數據集為ingredients，這個數據集是matlab自帶的)

 

　　Compute principal components for the ingredients data in the Hald data set, and the variance accounted for by each component.

 

　　load hald; %載入matlab內部數據

 

　　[pc,score,latent,tsquare] = princomp(ingredients); %調用pca分析函數

 

　　ingredients,score,pc,latent,tsquare %顯示得到的結果

 

　　ingredients =

 

　　7 26 6 60

 

　　1 29 15 52

 

　　11 56 8 20

 

　　11 31 8 47

 

　　7 52 6 33

 

　　11 55 9 22

 

　　3 71 17 6

 

　　1 31 22 44

 

　　2 54 18 22

 

　　21 47 4 26

 

　　1 40 23 34

 

　　11 66 9 12

 

　　10 68 8 12

 

　　score =

 

　　36.8218 -6.8709 -4.5909 0.3967

 

　　29.6073 4.6109 -2.2476 -0.3958

 

　　-12.9818 -4.2049 0.9022 -1.1261

 

　　23.7147 -6.6341 1.8547 -0.3786

 

　　-0.5532 -4.4617 -6.0874 0.1424

 

　　-10.8125 -3.6466 0.9130 -0.1350

 

　　-32.5882 8.9798 -1.6063 0.0818

 

　　22.6064 10.7259 3.2365 0.3243

 

　　-9.2626 8.9854 -0.0169 -0.5437

 

　　-3.2840 -14.1573 7.0465 0.3405

 

　　9.2200 12.3861 3.4283 0.4352

 

　　-25.5849 -2.7817 -0.3867 0.4468

 

　　-26.9032 -2.9310 -2.4455 0.4116

 

　　pc =

 

　　-0.0678 -0.6460 0.5673 0.5062

 

　　-0.6785 -0.0200 -0.5440 0.4933

 

　　0.0290 0.7553 0.4036 0.5156

 

　　0.7309 -0.1085 -0.4684 0.4844

 

　　latent =

 

　　517.7969

 

　　67.4964

 

　　12.4054

 

　　0.2372

 

　　tsquare =

 

　　5.6803

 

　　3.0758

 

　　6.0002

 

　　2.6198

 

　　3.3681

 

　　0.5668

 

　　3.4818

 

　　3.9794

 

　　2.6086

 

　　7.4818

 

　　4.1830

 

　　2.2327

 

　　2.7216

 

　　%下面我們來做一個驗證

 

　　%下面為計算ingredients協方差矩陣：

 

　　cov_ingredients=cov(ingredients)

 

　　cov_ingredients =

 

　　34.6026 20.9231 -31.0513 -24.1667

 

　　20.9231 242.1410 -13.8782 -253.4167

 

　　-31.0513 -13.8782 41.0256 3.1667

 

　　-24.1667 -253.4167 3.1667 280.1667

 

　　%下面為計算ingredients所對應的協方差矩陣(也就是cov_ingredients矩陣)的特征值和特征

 

　　%向量，下面的矩陣V為特征向量，D為特征值(對比上面的latent)組成的對角線矩陣

 

　　[V,D] = eig(cov_ingredients)

 

　　V =

 

　　0.5062 0.5673 0.6460 -0.0678

 

　　0.4933 -0.5440 0.0200 -0.6785

 

　　0.5156 0.4036 -0.7553 0.0290

 

　　0.4844 -0.4684 0.1085 0.7309

 

　　D =

 

　　0.2372 0 0 0

 

　　0 12.4054 0 0

 

　　0 0 67.4964 0

 

　　0 0 0 517.7969

 

　　%說明1：對比一下矩陣V和矩陣pc，現在很容易明白為什么COEFF是按列遞減順序排列的

 

　　% 了！(V中第三列與pc中倒數第三列差個負號，學過線性代數的人都知道這沒問題)

 

　　%下面再驗證一下說明2

 

　　diag(cov(score))

 

　　ans =

 

　　517.7969

 

　　67.4964

 

　　12.4054

 

　　0.2372

 

　　%說明2：以上結果顯示latent確實表示SCORE矩陣每列的方差，517.7969表示第一列方差

 

　　下面做圖表示結果：

 

　　上面說了半天還沒有達到我們終極想要的，其實我們要的是由函數[pc,score,latent,tsquare] = princomp(ingredients)它所產生的pc和latent。由latent可以算出降維后的空間所能表示原空間的程度，只要這個累積的值大於95%就行了。

 

　　The following command and plot show that two components account for 98% of the variance:

 

　　cumsum(latent)./sum(latent)

 

　　ans =

 

　　0.86597

 

　　0.97886

 

　　0.9996

 

　　1

 

　　%由以上ans值可以看出前兩個主成分就能表示原空間的97.886%,所以取pc中的前兩列可

 

　　%做主成分變換矩陣tranMatrix = pc(:,1:2)。則從原來的4維空間降到2維空間。對任意一個

 

　　%原空間樣本,例如a=(7 ,26 ,6 ,60)變到低維空間的表達式為a1 = a*tranMatrix。(當然你也可

 

　　%以取pc中的前三列，由原來的4維空間變到3維空間)

 

　　biplot(pc(:,1:2),'Scores',score(:,1:2),'VarLabels',...

 

　　{'X1' 'X2' 'X3' 'X4'})

 

上面這個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') returns only the elements of LATENT that are
    not necessarily zero, i.e., when N <= P, only the first N-1, and the
    corresponding columns of COEFF and SCORE.  This can be significantly
    faster when P >> N.

 

latent是返回的按降序排列的特征值，根據這個你可以手動的選擇降維以后的數據要選擇前多少列。

cumsum(latent)./sum(latent)

，通過這樣計算特征值的累計貢獻率，一般來說都選擇前95%的特征值對應的特征向量，還是原來的矩陣30*1000000，如果你計算得到前25個特征值的累計貢獻率已經超過99.9%，那么就完全可以只要降維后的數據的前25列。

 

tsquared是個什么東西我也不知道。。。不過貌似很少有人能用到，網絡上也沒有神馬資料，各位如果需要用的再查閱吧，一般情況下也用不到。

 

 

如果你需要對測試樣本降維，一般情況下，使用matlab自帶的方式，肯定需要對測試樣本減去一個訓練樣本均值，因為你在給訓練樣本降維的時候減去了均值，所以測試樣本也要減去均值，然后乘以coeff這個矩陣，就獲得了測試樣本降維后的數據。比如說你的測試樣本是1*1000000，那么乘上一個1000000*29的降維矩陣，就獲得了1*29的降維后的測試樣本的降維數據。

 

princomp(x)使用的行表示一個樣本，每行的所有的列數據都是這個樣本的特征值。降維以后比如是30*29，那么每一行就是降維以后的數據。每個樣本有29個特征值。

 

 

2. 一個自實現的pca降維方式

 

下面是來自mpb同學的一個自實現的例子，很牛的一個人，我們本科同學。

原文地址：http://blog.csdn.net/mpbchina/article/details/7384425

 

下面引用原文內容：

%訓練
%Lx=X'*X
clear;
clc;
train_path='..\Data\TrainingSet\';
phi=zeros(64*64,20);
for i=1:20
path=strcat(train_path,num2str(i),'.bmp');
Image=imread(path);
Image=imresize(Image,[64,64]);
phi(:,i)=double(reshape(Image,1,[])');
end;
%mean
mean_phi=mean(phi,2);
mean_face=reshape(mean_phi,64,64);
Image_mean=mat2gray(mean_face);
imwrite(Image_mean,'meanface.bmp','bmp');
%demean
for i=1:19
X(:,i)=phi(:,i)-mean_phi;
end
Lx=X'*X;
tic;
[eigenvector,eigenvalue]=eigs(Lx,19);
toc;
%normalization
for i=1:19
%K-L變換
UL(:,i)=X*eigenvector(:,i)/sqrt(eigenvalue(i,i));
end
%display Eigenface
for i=1:19
Eigenface=reshape(UL(:,i),[64,64]);
figure(i);
imshow(mat2gray(Eigenface));
end

得到的均值圖像mean_face：

前19個最大主元對應的“特征臉”：

測試：
測試用樣本：

 

%使用測試樣本進行測試
clc;
test_path='..\Data\TestingSet\';
error=zeros([1,4]);
for i=1:4
path=strcat(test_path,num2str(i),'.bmp');
Image=imread(path);
Image=double(imresize(Image,[64,64]));
phi_test=zeros(64*64,1);
phi_test(:,1)=double(reshape(Image,1,[])');
X_test=phi_test-mean_phi;
Y_test=UL'*X_test;
X_test_re=UL*Y_test;
Face_re=X_test_re+mean_phi;
calculate error rate
e=Face_re-phi_test;


%%display figure
Face_re_2=reshape(Face_re(:,1),[64,64]);
figure(i);

imshow(mat2gray(Image));
title('Original');
figure(10+i);
imshow(mat2gray(Face_re_2));
title('Reconstruct');
error(1,i)=norm(e);

%dispaly error rate
error_rate=error(1,i);
display(error_rate);
end

重建出的測試樣本與原樣本的對比：

 

四副測試樣本的重建誤差分別為：
1.4195e+003
1.9564e+003
4.7337e+003
7.0103e+003

可見測試樣本為人臉的樣本的重建誤差顯然小於非人臉的重建誤差。

 

 上面的降維的例子中，每一列表示一個樣本，這樣就一共有4096*20的待降維矩陣，然后對這個矩陣降維，請注意，如果采用列表示一個樣本，那么獲得的降維矩陣，是一個4096*19的矩陣，然后用這個降維矩陣對測試樣本和訓練樣本降維，我們的測試樣本是4096*1的矩陣，降維的時候這樣：

Y_test=UL'*X_test;  

 

UL是計算獲得降維矩陣，UL' （對UL進行轉至）獲得的19*4096的矩陣，19*4096 * 4096*1，就獲得了19*1的降維后的數據。

 

如果是使用matlab自帶的princomp進行降維，那么得到的coeff就是降維矩陣，使用測試樣本,這里的訓練樣本和測試樣本都要轉換成行向量，每一行表示一個樣本，測試樣本是1*4096，降維矩陣是 4096*29，那么就是  用待降維的樣本 x乘上降維矩陣 ， x * coeff ,注意這兩種不同的樣本表示方法中降維的使用，降維矩陣的不同位置。這樣降維后獲得1*4096 * 4096*29 = 1*29 的降維后的數據。

 

通過 上面的自己實現的pca降維的代碼，還可以對降維后的數據進行重建，獲得重建后的圖像，上面的程序中已經給出了。下面給出一個通過princomp降維后再對降維后的數據進行重建的程序。

 

通過 princomp降維后的數據進行重建

 

clear;
clc;
train_path='E:\TrainingSet\angry\positive\';
images = dir('E:\TrainingSet\angry\positive\*.bmp');
phi=zeros(30,64*64);

% 加載樣本圖像到 30*(64*64)的矩陣中，每一行代表一幅圖像
for i=1:30
    path=strcat(train_path,images(i).name);
    Image=imread(path);
    Image=imresize(Image,[64,64]);
    phi(i,:)=double(reshape(Image,1,[]));
end;

% 計算平均臉，並保存用以查看
mean_phi=mean(phi,1);
mean_face=reshape(mean_phi,64,64);
Image_mean=mat2gray(mean_face);
imwrite(Image_mean,'meanface2.bmp','bmp');

% 使用matlab自帶的pca進行降維
[coeff, score, latent, TSQUARED] = princomp(phi,'econ');

%display Eigenface
for i=1:29
    Eigenface=reshape(coeff(:,i),[64,64]);
    figure(i);
    imshow(mat2gray(Eigenface));
end

% 進行測試
%使用測試樣本進行測試
clc;
test_path='E:\BIT\code\FER\meanface.bmp';
error=zeros([1,4]);

Image=imread(test_path);
Image=double(imresize(Image,[64,64]));
phi_test=zeros(1,64*64);
phi_test(1,:)=double(reshape(Image,1,[])); % 讀入的測試圖像保存為一行，行向量
X_test=phi_test-mean_phi; % 檢測訓練樣本的平均臉
Y_test=X_test*coeff;  % 進行降維<span style="background-color: rgb(248, 248, 248);"></span>

X_test_re=Y_test*coeff'; % 重構
Face_re=X_test_re+mean_phi;
%calculate error rate
e=Face_re-phi_test;

%%display figure
Face_re_2=reshape(Face_re(1,:),[64,64]);
figure(i);

imshow(mat2gray(Image));
title('Original');
figure(10+i);
imshow(mat2gray(Face_re_2));
title('Reconstruct');
error(1,i)=norm(e);

%dispaly error rate
error_rate=error(1,i);
display(error_rate);

 

上面的程序關鍵處都有注釋，應該挺好理解的。

 

 

關於網絡上的一些解釋個人理解（僅供大家參考理解）

 

1. 

 

原文地址：http://www.cnblogs.com/sunwufan/archive/2011/08/31/2159952.html

 

原文：

最近看了些主成分分析，混跡Matlab論壇，翻了n多帖子，對princomp函數有了些了解。

在此只講一些個人理解，並沒有用術語，只求通俗。

貢獻率：每一維數據對於區分整個數據的貢獻，貢獻率最大的顯然是主成分，第二大的是次主成分......

[coef,score,latent,t2] = princomp(x);（個人觀點）：

x：為要輸入的n維原始數據。帶入這個matlab自帶函數，將會生成新的n維加工后的數據（即score）。此數據與之前的n維原始數據一一對應。

score：生成的n維加工后的數據存在score里。它是對原始數據進行的分析，進而在新的坐標系下獲得的數據。他將這n維數據按貢獻率由大到小排列。（即在改變坐標系的情況下，又對n維數據排序）

latent：是一維列向量，每一個數據是對應score里相應維的貢獻率，因為數據有n維所以列向量有n個數據。由大到小排列（因為score也是按貢獻率由大到小排列）。

coef：是系數矩陣。通過cofe可以知道x是怎樣轉換成score的。

 

則模型為從原始數據出發：
score= bsxfun(@minus,x,mean(x,1))*coef;(作用：可以把測試數據通過此方法轉變為新的坐標系)
逆變換：
x= bsxfun(@plus,score*inv(coef),mean(x,1))

例子：

%%
%清屏
clear
%%
%初始化數據
a=[-14.8271317103068,-3.00108550936016,1.52090778549498,3.95534842970601;-16.2288612441648,-2.80187433749996,-0.410815700402130,1.47546694457079;-15.1242838039605,-2.59871263957451,-0.359965674446737,1.34583763509479;-15.7031424565913,-2.53005662064257,0.255003254103276,-0.179334985754377;-17.7892158910100,-3.32842422986555,0.255791146332054,1.65118282449042;-17.8126324036279,-4.09719527953407,-0.879821957489877,-0.196675865428539;-14.9958877514765,-3.90753364293621,-0.418298866141441,-0.278063876667954;-15.5246706309866,-2.08905845264568,-1.16425848541704,-1.16976057326753;];
x=a;
%%
%調用princomp函數
[coef,score,latent,t2] = princomp(x);
score
%測試score是否和score_test一樣
score_test=bsxfun(@minus,x,mean(x,1))*coef;
score_test

latent=100*latent/sum(latent)%將latent總和統一為100，便於觀察貢獻率
pareto(latent);%調用matla畫圖

上圖是通過自帶函數繪制，當貢獻率累加至95%，以后的維數會不在顯示，最多只顯示10維。

下面用自己編寫的表示：

之前的錯誤認識：

1.認為主成分分析中latent顯示的貢獻值是原始數據的，其實是加工后的數據的。解釋:對原始數據既然選擇PCA方法，那么計算機認為原始數據每維之間可能存在關聯，你想去掉關聯、降低維數。所以采用這種方法的。所以計算機並不關心原始數據的貢獻值，因為你不會去用了，用的是加工后的數據（這也是為什么當把輸入數據每一維的順序改變后，score、latent不受影響的原因）。

2.認為PCA分析后自動降維，不對。PCA后會有貢獻值，是輸入者根據自己想要的貢獻值進行維數的改變，進而生成數據。（一般大家會取貢獻值在85%以上，要求高一點95%）。

3.PCA分析，只根據輸入數據的特征進行主成分分析，與輸出有多少類型，每個數據對應哪個類型無關。如果樣本已經分好類型，那PCA后勢必對結果的准確性有一定影響，我認為對於此類數據的PCA，就是在降維與准確性間找一個平衡點的問題，讓數據即不會維數多而使運算復雜，又有較高的分辨率。

我的個人見解：這篇文章中的解釋挺靠譜的，可以用來參考。第二點其實matlab的輸出結果score這個數據已經是降維后的數據，不過大家可以根據自己的需要取前多少列的數據。

 

 

2。

 

原文地址：http://www.ilovematlab.cn/thread-54600-1-1.html

 

部分原文：

回復 8# 5342245 的帖子                                                                                                                                                                                                                                                                                                                        設原始數據為X，先不做任何預處理。
[coef,score,latent,t2] = princomp(X);
則那些參數的底層算法大體過程如下：
x0 = bsxfun(@minus,X,mean(X,1)); %x0為將X去均值后的數據。
[coef,ignore] = eig(x0'*x0); 這就是coef的由來。  【當然最終的還有排序什么亂七八糟的。。】
scroe = x0*coef % 這就是score的由來，就是一個簡單的線性變換，將原來的X的坐標轉換到主成分空間中的坐標。僅此而已

則模型為從原始數據出發：
score = bsxfun(@minus,X,mean(X,1))*coef;

逆變換：
X = bsxfun(@plus,score*inv(coef),mean(X,1))


以上這些你可以自己驗證，看是否正確。
關於你的第三問。對於每一個主成分，就看coef的相應的列就能知道原始的變量那個對該主成分貢獻大了啊。。

上面是沒有預處理的。如果加了可逆的預處理。則原始數據亦可從預處理后的數據表示出。進而 bla bla....
===============這回夠通俗易懂吧。。O(∩_∩)O
PS：pca算法流程，你熟悉嗎？只要知道那個算法過程。這些都不難理解啊。。
建議您看看書把pca算法流程再過一遍。。否則別人再怎么說也沒用。。。

 

我的個人見解：

這里我想說的是，再對測試樣本進行降維的時候，一定要減去訓練樣本的均值，使用訓練樣本得到的轉換矩陣，保證訓練樣本和測試樣本轉換到相同的樣本空間中，這樣才有意思。大家有時間可以去看看英文的資料，說的都比較詳細。再用測試樣本減去均值以后，就可以進行轉換了。

 

很多同學可能在開始的時候和我一樣，都是不知道如果對測試樣本進行降維，很多人就選擇了還是使用princomp這個函數處理測試樣本，那么這樣測試樣本被映射到一個新的空間中，和原來的訓練樣本完全不是在一個空間，一點意義都沒有，還是要使用測試樣本減去均值，然后乘上訓練樣本降維的時候獲得降維矩陣，轉換到相同的空間中。

 

 

 

基本的對pca的認識就都說完了，比較亂，沒有條理，不過如果認真看下來的話，應該還是可以理解的。目前網上沒有關於pca的綜合的介紹個注意事項，說以我就把我的經驗和大家分享一下，還望文明轉載，轉載聲明出處。我也沒有對pca進行詳細的學習，肯定有不正確的地方，還請大家多多指教，共同探討。