?PCA主成分分析原理分析和Matlab實現方法（三）

【尊重原創，轉載請注明出處】http://blog.csdn.net/guyuealian/article/details/68487833 ? ? ??網上關于PCA（主成分分析）原理和分析的博客很多，本博客并不打算長篇大論推論PCA理論，而是用最精簡的語言說明鄙人對PCA的理解，并在最后給出用Matlab計算PCA過程的三種方法，方便大家對PCA的理解。 ? ? ? ?PS：本博客所有源代碼，都可以在附件中找到下載：http://download.csdn.net/detail/guyuealian/9799160 ? ? ? 關于PCA原理的文章，可參考： ? ? ? [1]http://blog.csdn.net/guyuealian/article/details/68483384? ? ? ? [2]http://blog.csdn.net/guyuealian/article/details/68483213 ? ? ? [3]張錚的《精通Matlab數字圖像處理與識別》

一、PCA原理簡要說明

? ? ? PCA算法主要用于降維，就是將樣本數據從高維空間投影到低維空間中，并盡可能的在低維空間中表示原始數據。PCA的幾何意義可簡單解釋為：

? ? ? 0維-PCA：將所有樣本信息都投影到一個點，因此無法反應樣本之間的差異；要想用一個點來盡可能的表示所有樣本數據，則這個點必定是樣本的均值。
? ? ? 1維-PCA：相當于將所有樣本信息向樣本均值的直線投影；
? ? ? 2維-PCA：將樣本的平面分布看作橢圓形分布，求出橢圓形的長短軸方向，然后將樣本信息投影到這兩條長短軸方向上，就是二維PCA。（投影方向就是平面上橢圓的長短軸方向）；
? ? ? 3維-PCA：樣本的平面分布看作橢圓形分布，投影方法分別是橢圓球的赤道半徑a和b，以及是極半徑c（沿著z軸）；

? ? ?PCA簡而言之就是根據輸入數據的分布給輸入數據重新找到更能描述這組數據的正交的坐標軸，比如下面一幅圖，對于那個橢圓狀的分布，最方便表示這個分布的坐標軸肯定是橢圓的長軸短軸而不是原來的x ，y軸。

? ? ?那么如何求出這個長軸和短軸呢？于是線性代數就來了：我們需要先求出這堆樣本數據的協方差矩陣，然后再求出這個協方差矩陣的特征值和特征向量，對應最大特征值的那個特征向量的方向就是長軸(也就是主元)的方向，次大特征值的就是第二主元的方向，以此類推。

?

二、PCA實現方法

? ? ? 實現PCA的方法，?可【1】直接調用Matlab工具箱princomp( )函數實現，也可【2】?自己實現PCA的過程，當然也可以【3】使用快速PCA算法的方法。 （1）方法一：[COEFF SCORE latent]=princomp(X)
? ? 參數說明：
? ? 1）COEFF 是主成分分量，即樣本協方差矩陣的特征向量；
? ? 2）SCORE主成分，是樣本X在低維空間的表示形式，即樣本X在主成份分量COEFF上的投影?，若需要降k維，則只需要取前k列主成分分量即可
? ? 3）latent：一個包含樣本協方差矩陣特征值的向量； 實例：假設有8個樣本，每個樣本有4個特征（屬性），使用PCA方法實現降維（k維，k小于特征個數4），并提取前2個主成份的特征，即將原始數據從4維空間降維到2維空間。 %% 樣本矩陣X，有8個樣本，每個樣本有4個特征，使用PCA降維提取k個主要特征（k<4） k=2; %將樣本降到k維參數設置 X=[1 2 1 1; %樣本矩陣3 3 1 2; 3 5 4 3; 5 4 5 4;5 6 1 5; 6 5 2 6;8 7 1 2;9 8 3 7] %% 使用Matlab工具箱princomp函數實現PCA [COEFF SCORE latent]=princomp(X) pcaData1=SCORE(:,1:k) %取前k個主成分 運行結果： X =1 2 1 13 3 1 23 5 4 35 4 5 45 6 1 56 5 2 68 7 1 29 8 3 7 COEFF =0.7084 -0.2826 -0.2766 -0.58460.5157 -0.2114 -0.1776 0.81110.0894 0.7882 -0.6086 0.01530.4735 0.5041 0.7222 -0.0116 SCORE =-5.7947 -0.6071 0.4140 -0.0823-3.3886 -0.8795 0.4054 -0.4519-1.6155 1.5665 -1.0535 1.2047-0.1513 2.5051 -1.3157 -0.77180.9958 -0.5665 1.4859 0.77751.7515 0.6546 1.5004 -0.61442.2162 -3.1381 -1.6879 -0.13055.9867 0.4650 0.2514 0.0689 latent =13.21512.95501.50690.4660 pcaData1 =-5.7947 -0.6071-3.3886 -0.8795-1.6155 1.5665-0.1513 2.50510.9958 -0.56651.7515 0.65462.2162 -3.13815.9867 0.4650 （2）方法二：自己編程實現
? ??PCA的算法過程，用一句話來說，就是“將所有樣本X減去樣本均值m，再乘以樣本的協方差矩陣C的特征向量V，即為PCA主成分分析”，其計算過程如下：
? ??[1].將原始數據按行組成ｍ行ｎ列樣本矩陣X（每行一個樣本，每列為一維特征）
? ??[2].求出樣本X的協方差矩陣C和樣本均值m；（Matlab可使用cov()函數求樣本的協方差矩陣C，均值用mean函數）
? ??[3].求出協方差矩陣的特征值D及對應的特征向量V；（Matlab可使用eigs()函數求矩陣的特征值D和特征向量V）
? ??[4].將特征向量按對應特征值大小從上到下按行排列成矩陣，取前k行組成矩陣P；（eigs()返回特征值構成的向量本身就是從大到小排序的）
? ??[5].Y=(X-m)×P即為降維到k維后的數據； ? ? PS：關于協方差矩陣，很多人有點郁悶，有些教程用協方差矩陣，而有些資料是用散步矩陣計算，其實協方差矩陣和散步矩陣就是一個倍數關系：協方差矩陣C×（n-1）=散步矩陣S。我們可以使用Matlab工具箱的cov函數求協方差矩陣C。 %% 自己實現PCA的方法 [Row Col]=size(X); covX=cov(X); %求樣本的協方差矩陣（散步矩陣除以(n-1)即為協方差矩陣） [V D]=eigs(covX); %求協方差矩陣的特征值D和特征向量V meanX=mean(X); %樣本均值m %所有樣本X減去樣本均值m，再乘以協方差矩陣（散步矩陣）的特征向量V，即為樣本的主成份SCORE tempX= repmat(meanX,Row,1); SCORE2=(X-tempX)*V %主成份：SCORE pcaData2=SCORE2(:,1:k) 運行結果： SCORE2 =-5.7947 0.6071 -0.4140 0.0823-3.3886 0.8795 -0.4054 0.4519-1.6155 -1.5665 1.0535 -1.2047-0.1513 -2.5051 1.3157 0.77180.9958 0.5665 -1.4859 -0.77751.7515 -0.6546 -1.5004 0.61442.2162 3.1381 1.6879 0.13055.9867 -0.4650 -0.2514 -0.0689pcaData2 =-5.7947 0.6071-3.3886 0.8795-1.6155 -1.5665-0.1513 -2.50510.9958 0.56651.7515 -0.65462.2162 3.13815.9867 -0.4650 ? ? ?對比方法一和方法可知，主成份分量SCORE和SCORE2的絕對值是一樣的（符號只是相反方向投影而已，不影響分析結果），其中pcaData是從SCORE提取前2列的數據，這pcaData就是PCA從4維降到2維空間的數據表示形式，pcaData可以理解為：通過PCA降維，每個樣本可以用2個特征來表示原來4個特征了。

（3）方法三：使用快速PCA方法 ? ? ?PCA的計算中最主要的工作量是計算樣本協方差矩陣的本征值和本征向量。假設樣本矩陣X的大小為n ×d (n個d?維樣本特征向量)，則樣本散布矩陣(協方差矩陣) S 將是一個d×d的方陣，故當維數d較大時計算復雜度會非常高。例如當維數d=10000，S是一個10 000 ×10 000的矩陣，此時如果采用上面的princomp函數計算主成份，Matlab通常會出現內存耗盡的問題（out of memory）， 即使有足夠多的內存，要得到S的全部本征值可能也要花費數小時的時間。 ? ? 快速PCA的方法相關理論，可以參考張錚的《精通Matlab數字圖像處理與識別》第12章（P307），PDF可在附件下載：http://download.csdn.net/detail/guyuealian/9799160

? ? fastPCA函數用來對樣本矩陣A進行快速主成分分析和降維(降至k維)，其輸出pcaA為維后的k維樣本特征向量組成的矩陣，每行一個樣本，列數k為降維后的樣本特征維數，相當于princomp函數中的輸出SCORE， 而輸出V為主成分分量，相當于princomp函數中的輸出COEFF。

?

%% 使用快速PCA算法實現的方法 [pcaData3 COEFF3] = fastPCA(X, k ) ? ? 其中fastPCA函數的代碼實現如下：

?

function [pcaA V] = fastPCA( A, k ) % 快速PCA % 輸入：A --- 樣本矩陣，每行為一個樣本 % k --- 降維至 k 維 % 輸出：pcaA --- 降維后的 k 維樣本特征向量組成的矩陣，每行一個樣本，列數 k 為降維后的樣本特征維數 % V --- 主成分向量 [r c] = size(A); % 樣本均值 meanVec = mean(A); % 計算協方差矩陣的轉置 covMatT Z = (A-repmat(meanVec, r, 1)); covMatT = Z * Z'; % 計算 covMatT 的前 k 個本征值和本征向量 [V D] = eigs(covMatT, k); % 得到協方差矩陣 (covMatT)' 的本征向量 V = Z' * V; % 本征向量歸一化為單位本征向量 for i=1:kV(:,i)=V(:,i)/norm(V(:,i)); end % 線性變換（投影）降維至 k 維 pcaA = Z * V; % 保存變換矩陣 V 和變換原點 meanVec
運行結果為： pcaData3 =-5.7947 -0.6071-3.3886 -0.8795-1.6155 1.5665-0.1513 2.50510.9958 -0.56651.7515 0.65462.2162 -3.13815.9867 0.4650COEFF3 =0.7084 -0.28260.5157 -0.21140.0894 0.78820.4735 0.5041
? ? 源代碼和附件下載地址：http://download.csdn.net/detail/guyuealian/9799160 如果你覺得該帖子幫到你，還望貴人多多支持，鄙人會再接再厲，繼續努力的~

總結

以上是生活随笔為你收集整理的PCA原理分析和Matlab实现方法（三）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。