PCA原理簡述

原貼： http://www.fengchang.cc/post/89

以前讀書看過，沒有總結，忘得差很少了，最近讀書又看了一遍，寫總結以下。

 

仍是老方法，先創建直覺，提供一個比較好的創建直覺的資料： A TUTORIAL ON PRINCIPAL COMPONENT ANALYSIS 個人圈評​

 

我就引用這裏面的例子來講明。

上圖是一個彈簧實驗的示意圖。假設咱們經過三個攝影機 A B C來記錄連續若干時刻的球的位置，目標是求出球的運動方程。首先彈簧自己會沿着x方向伸縮，可是因爲現實世界的偏差的緣由，在y方向上依然可能會有輕微的擾動。在已經肯定好x, y, z的座標軸的狀況下，這個問題倒也不難，難就難在，在現實世界的不少問題中，並不存在這樣的設定好的座標軸給你參考，因此基本上你收集到的信息多是比較雜亂的，好比上圖中，咱們假設每一個攝像機記錄兩個維度的信息A記錄(Ax, Ay)， B記錄(Bx, By), C記錄(Cx, Cy)， 注意這裏x, y 只表示兩個維度，不表明在哪一個座標軸，那麼結果能夠認爲就是一個攝像機記錄一個觀察角度下球的向量，最終的一條記錄數據長這樣（Ax,Ay,Bx,By,Cx,Cy）,這種數據咱們能夠想象的一個缺陷是，它有6個維度，可是對實際彈簧球的運動來講，實在太多了，且沒有重點，由於彈簧只在一個方向上運動爲主，那麼若是咱們拿這樣的數據集去作機器學習中的訓練，這個訓練數據其實質量是至關差的。

 

又比如我收集學生的以下6個維度的信息：年齡，性別，智商，膚色，是否直髮，眼睛顏色，來預測其人物理考試的成績會落在那個分數段，可想而知，智商對結果的影響必定是最大的，而是否直髮這個維度的信息可能根本就可有可無。那麼PCA要作的一件事就是經過「降維」來處理訓練數據。

 

問題來了，怎麼降呢？若是是上述物理考試成績預測的問題，直接把是否直髮這樣的維度刪掉，就是一種處理方案，這樣就把6維訓練數據，降到5維了。可是這麼粗暴的方法，並不適合全部的問題，好比上述彈簧球的問題就不適合，由於你並不能明確知道到底哪一個維度的數據是可有可無的，甚至實際上可能某一個維度的信息，其實是由另外兩個維度的信息疊加合成的（能夠參考牛頓物理學中基本的運動的合成），這個時候問題就比較複雜了。

 

那麼PCA是怎麼作的呢？簡單來講就是線性代數當中的基變換。

 

姿式前提：基：假設線性子空間的基B={v1,v2,...,vk}, 向量 a = v1c1+v2c2+...+vkck，那麼(c1,c2,...,ck)即爲a基於B的座標。

 

這裏不得不先說一下，理解PCA的前提是基本的線性代數知識，涉及到方陣，矩陣對角化，特徵值，特徵向量，這幾個概念理解得越好，對PCA的理解就越水到渠成。我本人只是在大學數學學過線性代數，沒有多精通，因此也是先複習了一下這幾個概念，因此自認爲對PCA的理解尚未達到多深的水平。

 

假設咱們的原始訓練數據是一個m乘n矩陣，其中m是觀察的原始維度，n是訓練樣本的個數, 因此至關於每一列是一個觀察樣本，每一行是一個維度的全部觀察集合，那麼針對上面彈簧實驗的問題，假如咱們在記錄了200個數據時刻的數據，那麼訓練數據就是6*200的矩陣， 我就記錄爲X吧，這樣的矩陣，若是左乘一個矩陣P,例如P的尺寸爲5*6, 那麼結果就變成一個5*200的矩陣，先無論目的爲什麼，這樣左乘，實際上就是對原始數據作了一個基變換，且從6維降到了5維。

 

而後咱們但願達到一個目標，若是在基變換後的訓練數據，各個維度之間的相關度都很低，那就至關的好，爲何呢？由於在你胡亂觀察數據的時候，收集來的信息多是有很大冗餘的，好比你既收集了一我的的出生年月，又收集了他的年齡，又收集了他的屬相，而後存爲3個維度的信息，這很明顯就有冗餘，由於這三者是相關的，從出生年月這個維度，徹底能夠推導出另外兩個維度。這種維度之間的冗餘性，就是PCA要解決的問題，經過線性變換，達到一個效果，讓維度之間儘量獨立，減小冗餘性。至於爲何線性變換能夠達到這個效果，不在本文的討論範圍內了， 我也理解得不深。但這就是咱們的優化目標，那麼下面就須要用一種正式的數學公式來衡量這個冗餘性，PCA是這樣作的：

假設咱們的訓練數據叫X（就以咱們上面說的6*200的訓練數據爲例），那麼經過上圖中的公式定義了一個叫covariance matrix的矩陣，這個矩陣裏面的每個元素，都表明着兩個維度之間的協方差，協方差這個東西，就是定義兩個變量相關性的的一個量化指標，這個指標的定義跟什麼均值，方差之類的指標是一個level的，我這裏就再也不贅述了。總之上圖就能夠定義維度與維度之間的協方差（若是你們對比協方差的標準公式定義可能會有一點困惑，就是均值哪去了？實際上是這樣的，這裏的數據都作了預處理，每一個原始數據都已經減去了均值，這樣處理後的均值都變成0了，因此看起來均值就不見了）。

 

接下來這個協方差矩陣就是一個有表明性的東西，它的每一個元素都表明了兩個維度之間的相關性，假設第i行第j列的元素a<i,j>,就表明了第i維數據跟第j維數據的相關性。那麼還記得咱們的優化目標嗎？就是要經過線性變換，讓不一樣維度之間的相關性儘量小，對應到這個協方差矩陣，就是但願該矩陣對角線之外的相關性都是0，只有對角線上的元素非0。這裏就要用到線性代數了，首先根據定義，這個協方差矩陣必定是個方陣，方陣的話，若是能對角化就說明，咱們老是有辦法找到這樣一個線性變換，讓它變成對角矩陣，這不就達到咱們的優化目標了嗎？

插入一個矩陣可對角化的定義：若是存在一個 可逆矩陣 P 使得 P  −1AP 是對角矩陣，則它就被稱爲可對角化的

 

再梳理一下，能夠這樣描述PCA作的事，對原始訓練數據矩陣（記爲X），要尋找一個線性變換矩陣P,使得PX對應的協方差矩陣爲對角陣。找到這樣一個P以後，PX就是咱們要變換後的結果。

 

這應該是能一句話說明的最簡單描述了，可是這裏P的維度實際上是要注意的，由於理論上可能有多種，好比對於上面彈簧的訓練數據（6*200的訓練數據），P多是5*6,也多是4*6或者3*6的尺寸，總之就是k*6,這個k就是咱們降維後想要達到的維度，一般是人爲指定的，至於如何選這個k，就是另外一個話題了，之後有空再說。

 

下面說一下在指定了k以後，如何解這個P，基本上就是所謂矩陣對角化那一套了，我這裏不說底層細節，由於其實展開了說挺難的，涉及到的都是線性代數的運算知識，包括特徵值，大致就是把矩陣的特徵值全找出來，而後挑出k個最大的特徵值，找出對應的k個特徵向量，這k個特徵向量拼起來就是矩陣P。

 

下面給一個python sklearn的代碼實現，用的是Iris數據集，長這樣，也很少，就150條記錄。

比較簡單，並且沒有底層的實現細節，由於sklearn已經封裝好了，不過經過debug的過程，能夠看看矩陣的維度是怎麼變化的，加深對PCA該過程的理解。

 

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62

#!/usr/bin/python # -*- coding: utf-8 -*-   """ ========================================================= PCA example with Iris Data-set =========================================================   Principal Component Analysis applied to the Iris dataset.   See `here <https://en.wikipedia.org/wiki/Iris_flower_data_set>`_ for more information on this dataset.   """ print (__doc__)     # Code source: Gaël Varoquaux # License: BSD 3 clause   import  numpy as np import  matplotlib.pyplot as plt from  mpl_toolkits.mplot3d  import  Axes3D     from  sklearn  import  decomposition from  sklearn  import  datasets   np.random.seed( 5 )   centers  =  [[ 1 ,  1 ], [ - 1 ,  - 1 ], [ 1 ,  - 1 ]] iris  =  datasets.load_iris() X  =  iris.data y  =  iris.target   fig  =  plt.figure( 1 , figsize = ( 4 ,  3 )) #clear figure plt.clf() ax  =  Axes3D(fig, rect = [ 0 ,  0 , . 95 ,  1 ], elev = 48 , azim = 134 )   # clear axis plt.cla() pca  =  decomposition.PCA(n_components = 3 ) pca.fit(X) X  =  pca.transform(X)   for  name, label  in  [( 'Setosa' ,  0 ), ( 'Versicolour' ,  1 ), ( 'Virginica' ,  2 )]:      ax.text3D(X[y  = =  label,  0 ].mean(),                X[y  = =  label,  1 ].mean()  +  1.5 ,                X[y  = =  label,  2 ].mean(), name,                horizontalalignment = 'center' ,                bbox = dict (alpha = . 5 , edgecolor = 'w' , facecolor = 'w' )) # Reorder the labels to have colors matching the cluster results y  =  np.choose(y, [ 1 ,  2 ,  0 ]).astype(np. float ) ax.scatter(X[:,  0 ], X[:,  1 ], X[:,  2 ], c = y, cmap = plt.cm.nipy_spectral,             edgecolor = 'k' )   ax.w_xaxis.set_ticklabels([]) ax.w_yaxis.set_ticklabels([]) ax.w_zaxis.set_ticklabels([])   plt.show()

剛纔說了，對線性代數以上若干概念的理解直接關乎對PCA理解的深度，附知乎帖子一篇：  如何理解矩陣特徵值​

 

代碼原貼： PCA example​

 

 

再附一個言簡意賅的解釋，這時看當更明白其真意了。