機器學習之 PCA(主成分分析)

時間 2019-11-08

標籤機器學習 pca 成分分析简体版

原文原文鏈接

以前相關文章目錄：html

下面介紹一種降維算法，即PCA(主成分分析)。python

在機器學習中，有一種問題被稱爲維數災難，在實際機器學習項目中，咱們要處理的的樣本數據的維數多是成千上萬，甚至幾十萬或者更多的，這種狀況下，直接對原始樣本數據進行訓練建模會耗費大量時間，對應的資源消耗是不可接受，這個時候，咱們就須要對數據進行降維，降維固然意味着信息的丟失，不過鑑於實際數據自己經常存在的相關性，咱們能夠想辦法在降維的同時將信息的損失儘可能下降.git

PCA是一種具備嚴格數學基礎而且已被普遍採用的降維方法算法

在具體給出PCA算法以前，先回顧下線性代數中的相關知識，首先，須要強調一點，下面所說的向量，在沒有特殊說明下，都指的是列向量數組

1.向量的內積和投影

兩個維數相同的向量的內積被定義爲：dom

(a_1,a_2,\cdots,a_n)\cdot (b_1,b_2,\cdots,b_n)^\mathsf{T}=a_1b_1+a_2b_2+\cdots+a_nb_n

接下來，看下內積的幾何意義，假設A，B是兩個n維向量，爲了簡單起見咱們假設A和B均爲二維向量，則。則在二維平面上A和B能夠用兩條發自原點的有向線段表示,以下圖：
如今咱們從A點向B所在直線引一條垂線。咱們知道垂線與B的交點叫作A在B上的投影，再設A與B的夾角是a，則投影的矢量長度爲，其中 $|A|=\sqrt{x_1^2+y_1^2}$ 是向量A的模，也就是A線段的標量長度機器學習

咱們知道，內積還可表示爲 $A\cdot B=|A||B|cos(a)$ ，也就是說，A與B的內積等於A到B的投影長度乘以B的模
因此，設向量B的模爲1，則A與B的內積值等於A向B所在直線投影的矢量長度！

基與基變換

咱們知道，一個二維向量能夠對應二維笛卡爾直角座標系中從原點出發的一個有向線段。例以下面這個向量：函數

圖中的那個向量，咱們能夠表示爲（3，2），其中，3，表示的是向量在X軸的投影值爲3，2，表示的是在Y軸的投影值爲2。post

也就是說咱們其實隱式引入了一個定義：以x軸和y軸上正方向長度爲1的向量爲標準。那麼一個向量(3,2)實際是說在x軸投影爲3而y軸的投影爲2。注意投影是一個矢量，因此能夠爲負。學習

因此，向量(x,y)實際上表示線性組合： $$x(1,0)^\mathsf{T}+y(0,1)^\mathsf{T}$$ 而此處的(1,0),(0,1)叫作二維空間中的一組基。也是一組正交基。

咱們之因此默認選擇(1,0)和(0,1)爲基，固然是比較方便，由於它們分別是x和y軸正方向上的單位向量，所以就使得二維平面上點座標和向量一一對應，很是方便。但實際上任何兩個線性無關的二維向量均可以成爲一組基，所謂線性無關在二維平面內能夠直觀認爲是兩個不在一條直線上的向量。

例如，(1,1)和(-1,1)也能夠成爲一組基。通常來講，咱們但願基的模是1，由於從內積的意義能夠看到，若是基的模是1，那麼就能夠方便的用向量點乘基而直接得到其在新基上的座標了！實際上，對應任何一個向量咱們總能夠找到其同方向上模爲1的向量，只要讓兩個份量分別除以模就行了。例如，上面的基能夠變爲$$(\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}})和(-\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}})。$$

如今，咱們想得到(3,2)在新基上的座標，即在兩個方向上的投影矢量值，那麼根據內積的幾何意義，咱們只要分別計算(3,2)和兩個基的內積，不可貴到新的座標爲 $(\frac{5}{\sqrt{2}},-\frac{1}{\sqrt{2}})$ 。下圖給出了新的基以及(3,2)在新基上座標值的示意圖：

另外這裏要注意的是，咱們列舉的例子中基是正交的（即內積爲0，或直觀說相互垂直），但能夠成爲一組基的惟一要求就是線性無關，非正交的基也是能夠的。不過由於正交基有較好的性質，因此通常使用的基都是正交的。

基變換的矩陣表示

經過前面的描述，咱們知道，將(3,2)變換爲新基上的座標，就是用(3,2)去和新的一組基中的每個去作內積運算。咱們能夠用矩陣相乘的形式簡潔的表示這個變換：

\begin{pmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ -1/\sqrt{2} & 1/\sqrt{2} \end{pmatrix} \begin{pmatrix} 3 \\ 2 \end{pmatrix} = \begin{pmatrix} 5/\sqrt{2} \\ -1/\sqrt{2} \end{pmatrix}

推廣到多個向量，假設有m個向量，只要將二維向量按列排成一個兩行m列矩陣，而後用「基矩陣」乘以這個矩陣，就獲得了全部這些向量在新基下的值。例如(1,1)，(2,2)，(3,3)，想變換到剛纔那組基上，則能夠這樣表示：

\begin{pmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ -1/\sqrt{2} & 1/\sqrt{2} \end{pmatrix} \begin{pmatrix} 1 & 2 & 3 \\ 1 & 2 & 3 \end{pmatrix} = \begin{pmatrix} 2/\sqrt{2} & 4/\sqrt{2} & 6/\sqrt{2} \\ 0 & 0 & 0 \end{pmatrix}

由此咱們能夠能夠看到，基變換能夠表示爲矩陣相乘

通常的，若是咱們有M個N維向量，想將其變換爲由R個N維向量表示的新空間中，那麼首先將R個基按行組成矩陣A，而後將向量按列組成矩陣B，那麼兩矩陣的乘積AB就是變換結果，其中AB的第m列爲A中第m列變換後的結果。數學表示爲：

\begin{pmatrix} p_1 & p_2 & \cdots & p_M \end{pmatrix}^\mathsf{T} \begin{pmatrix} a_1 & a_2 & \cdots & a_M \end{pmatrix} = \begin{pmatrix} p_1a_1 & p_1a_2 & \cdots & p_1a_M \\ p_2a_1 & p_2a_2 & \cdots & p_2a_M \\ \vdots & \vdots & \ddots & \vdots \\ p_Ra_1 & p_Ra_2 & \cdots & p_Ra_M \end{pmatrix}

若是R小於N時，這就是將一N維數據變換到更低維度的空間去，所以這種矩陣相乘的表示也能夠表示降維變換。

兩個矩陣相乘的意義是將右邊矩陣中的每一列列向量變換到左邊矩陣中每一行行向量爲基所表示的空間中去

優化目標

經過上面咱們知道，經過選擇不一樣的一組基，能夠將同一組數據給出不一樣的表示，並且若是一組基中基的數量少於向量自己的維數，則能夠達到降維的效果。

那麼，接下來咱們的問題就是如何找到一組最優的基，也就是說，若是咱們有一組N維向量，如今要將它降到K維（K小於N），那麼，咱們該如何選擇K個基才能最大程度的保留原始N維向量的信息。

下面，以一個具體例子來展開，假設咱們的數據有5條記錄，以下：

\begin{pmatrix} 1 & 1 & 2 & 4 & 2 \\ 1 & 3 & 3 & 4 & 4 \end{pmatrix}

每一列爲一條數據，一行爲一個字段，也就是一個特徵。

首先，爲了後續處理方便，咱們將數據進行均值歸零化，也就是將每一個字段內（也就是每一行）全部值都減去字段的均值，處理之後，新的數據每一個字段的均值都爲0。變換後的數據以下：

\begin{pmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{pmatrix}

在座標系中的分佈以下：

如今，咱們要將這些數據降到一維，可是有要儘量的保留原始的信息。如何去處理呢？

經過上面的討論咱們知道，這個問題其實是要在二維平面中選擇一個方向，將全部的數據否投影到這個方向所在的直線上，那麼，如何選擇這個方向（或者說基）才能更多的保留的原始信息呢？一種直觀的見解是：咱們但願使投影后的數值儘量分散

以上圖爲例，能夠看出若是向x軸投影，那麼最左邊的兩個點會重疊在一塊兒，中間的兩個點也會重疊在一塊兒，因而自己四個各不相同的二維點投影后只剩下兩個不一樣的值了，這是一種嚴重的信息丟失，同理，若是向y軸投影最上面的兩個點和分佈在x軸上的兩個點也會重疊。因此看來x和y軸都不是最好的投影選擇。咱們直觀目測，若是向經過第一象限和第三象限的斜線投影，則五個點在投影后仍是能夠區分的。

下面，咱們經過具體的數學方法來表述下這個問題

方差

上面咱們討論出，但願投影后的數值儘量分散，這種分散程度，在數學上能夠用方差來表述。字段a的方差表示爲：

Var = \frac{1}{m} \sum_{i=1}^m{(a_i-\mu)^2}

因爲上面咱們已經對每一個字段進行均值歸零化，因此：

因此，上述咱們的問題被表述爲：尋找一個一維基，使得全部數據變換爲這個基上的座標表示後，方差值最大

協方差

對於上面二維降成一維的問題來講，找到那個使得方差最大的方向就能夠了。不過對於更高維，還有一個問題須要解決。考慮三維降到二維問題。與以前相同，首先咱們但願找到一個方向使得投影后方差最大，這樣就完成了第一個方向的選擇，繼而咱們選擇第二個投影方向。

若是咱們仍是單純只選擇方差最大的方向，很明顯，這個方向與第一個方向應該是「幾乎重合在一塊兒」，顯然這樣的維度是沒有用的，所以，應該有其餘約束條件。從直觀上說，讓兩個字段儘量表示更多的原始信息，咱們是不但願它們之間存在（線性）相關性的，由於相關性意味着兩個字段不是徹底獨立，必然存在重複表示的信息。

數學上能夠用兩個字段的協方差表示其相關性：

因爲已經讓每一個字段均值爲0，則：

能夠看到，在字段均值爲0的狀況下，兩個字段的協方差簡潔的表示爲其內積除以元素數m。

當協方差爲0時，表示兩個字段徹底獨立。爲了讓協方差爲0，咱們選擇第二個基時只能在與第一個基正交的方向上選擇(正交的化，內積爲0，正好和協方差爲0相對應)。所以最終選擇的兩個方向必定是正交的。

至此，咱們獲得了降維問題的優化目標：將一組N維向量降爲K維（K大於0，小於N），其目標是選擇K個單位（模爲1）正交基，使得原始數據變換到這組基上後，各字段兩兩間協方差爲0，而字段的方差則儘量大（在正交的約束下，取最大的K個方差）。

協方差矩陣

咱們看到，最終要達到的目的與字段內方差及字段間協方差有密切關係。所以咱們但願能將二者統一表示，仔細觀察發現，二者都可以表示爲內積的形式，而內積又與矩陣相乘密切相關。因而咱們來了靈感：

假設咱們只有a和b兩個字段，那麼咱們將它們按行組成矩陣X：

X=\begin{pmatrix} a_1 & a_2 & \cdots & a_m \\ b_1 & b_2 & \cdots & b_m \end{pmatrix}

而後咱們用X乘以X的轉置，並乘上係數1/m：

這個矩陣對角線上的兩個元素分別是兩個字段的方差，而其它元素是a和b的協方差。二者被統一到了一個矩陣的。這就是兩個字段間的協方差矩陣

根據矩陣相乘的運算法則，這個結論很容易被推廣到通常狀況：設咱們有m個n維數據記錄，將其按列排成n乘m的矩陣X，設 $C=\frac{1}{m}XX^\mathsf{T}$ ，則C是一個對稱矩陣，其對角線分別個各個字段的方差，而第i行j列和j行i列元素相同，表示i和j兩個字段的協方差

協方差矩陣對角化

根據上述推導，咱們發現要達到優化目前，等價於將協方差矩陣對角化：即除對角線外的其它元素化爲0，而且在對角線上將元素按大小從上到下排列，這樣咱們就達到了優化目的，這樣說可能還不是很明晰，咱們進一步看下原矩陣與基變換後矩陣協方差矩陣的關係：

設原始數據矩陣X對應的協方差矩陣爲C，P是一組基（每一個基都是一個列向量）組成的矩陣，設,則Y爲X對P所包含的那組基作基變換後的數據。設Y的協方差矩陣爲D，則有：

\begin{array}{l l l} D & = & \frac{1}{m}YY^\mathsf{T} \\ & = & \frac{1}{m}(P^TX)(P^TX)^\mathsf{T} \\ & = & \frac{1}{m}P^TXX^\mathsf{T}P \\ & = & P^T(\frac{1}{m}XX^\mathsf{T})P \\ & = & P^TCP \end{array}

如今事情很明白了！咱們要找的P不是別的，而是能讓原始協方差矩陣對角化的P。換句話說，優化目標變成了尋找一個矩陣P，知足 $P^\mathsf{T}CP$ 是一個對角矩陣，而且對角元素按從大到小依次排列，那麼P的前K列就是要尋找的基，用P的前K列組成的矩陣的轉置乘以X就使得X從N維降到了K維並知足上述優化條件。

上文知道，協方差矩陣C是一個是實對稱矩陣，在線性代數上，實對稱矩陣有一系列很是好的性質：
1）實對稱矩陣的特徵值都爲實數
2）實對稱矩陣的全部特徵向量正交。
3）設特徵值 $\lambda$ 重數爲r，則必然存在r個線性無關的特徵向量對應於 $\lambda$ ，所以能夠將這r個特徵向量單位正交化。

由上面兩條可知，一個n行n列的實對稱矩陣必定能夠找到n個單位正交特徵向量，設這n個特徵向量爲 $e_1,e_2,\cdots,e_n$ ，咱們將其按列組成矩陣：

E=\begin{pmatrix} e_1 & e_2 & \cdots & e_n \end{pmatrix}

則對協方差矩陣C有以下結論：

E^\mathsf{T}CE=\Lambda=\begin{pmatrix} \lambda_1 & & & \\ & \lambda_2 & & \\ & & \ddots & \\ & & & \lambda_n \end{pmatrix}

其中 $\Lambda$ 爲對角矩陣，其對角元素爲各特徵向量對應的特徵值（可能有重複）。

到此，咱們就已經找到了咱們須要的矩陣P = E

P是協方差矩陣的特徵向量單位化後按列排列出的矩陣，其中每一列都是C的一個特徵向量。若是設P按照 $\Lambda$ 中特徵值的從大到小，將特徵向量從左到右排列，則用P的前K行組成的矩陣的裝置乘以原始數據矩陣X，就獲得了咱們須要的降維後的數據矩陣Y。

上述，就是整個PCA的數學原理討論

PCA算法

總結一下PCA的算法步驟：

設有m條n維數據。

1）將原始數據按列組成n行m列矩陣X

2）將X的每一行（表明一個屬性字段）進行零均值化，即減去這一行的均值

3）求出協方差矩陣C=\frac{1}{m}XX^\mathsf{T}

4）求出協方差矩陣的特徵值及對應的特徵向量

5）將特徵向量按對應特徵值大小從左到右按列排列成矩陣，取前k列組成矩陣P

6）Y=P^TX即爲降維到k維後的數據

實例

這裏以上文提到的

爲例，咱們用PCA方法將這組二維數據其降到一維。

由於這個矩陣的每行已是零均值，這裏咱們直接求協方差矩陣：

C=\frac{1}{5}\begin{pmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{pmatrix}\begin{pmatrix} -1 & -2 \\ -1 & 0 \\ 0 & 0 \\ 2 & 1 \\ 0 & 1 \end{pmatrix}=\begin{pmatrix} \frac{6}{5} & \frac{4}{5} \\ \frac{4}{5} & \frac{6}{5} \end{pmatrix}

而後求其特徵值和特徵向量，具體求解方法再也不詳述，能夠參考相關資料。求解後特徵值爲：

其對應的特徵向量分別是：

c_1\begin{pmatrix} 1 \\ 1 \end{pmatrix},c_2\begin{pmatrix} -1 \\ 1 \end{pmatrix}

其中對應的特徵向量分別是一個通解，c_1和c_2可取任意實數。那麼標準化後的特徵向量爲：

\begin{pmatrix} 1/\sqrt{2} \\ 1/\sqrt{2} \end{pmatrix},\begin{pmatrix} -1/\sqrt{2} \\ 1/\sqrt{2} \end{pmatrix}

所以咱們的矩陣P是：

P=\begin{pmatrix} 1/\sqrt{2} & -1/\sqrt{2} \\ 1/\sqrt{2} & 1/\sqrt{2} \end{pmatrix}

能夠驗證協方差矩陣C的對角化：

P^\mathsf{T}CP=\begin{pmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ -1/\sqrt{2} & 1/\sqrt{2} \end{pmatrix}\begin{pmatrix} 6/5 & 4/5 \\ 4/5 & 6/5 \end{pmatrix}\begin{pmatrix} 1/\sqrt{2} & -1/\sqrt{2} \\ 1/\sqrt{2} & 1/\sqrt{2} \end{pmatrix}=\begin{pmatrix} 2 & 0 \\ 0 & 2/5 \end{pmatrix}

最後咱們用P的第一行乘以數據矩陣，就獲得了降維後的表示：

Y=\begin{pmatrix} 1/\sqrt{2} \\ 1/\sqrt{2} \end{pmatrix}^T\begin{pmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{pmatrix}=\begin{pmatrix} -3/\sqrt{2} & -1/\sqrt{2} & 0 & 3/\sqrt{2} & -1/\sqrt{2} \end{pmatrix}

代碼實現

import numpy as np
import pandas as pd
複製代碼

## X: 須要降維的原始數據 topNfeat：須要降到的維數
def PCA(X, topNfeat=9999999):
    # 1.原始數據默認都是每一行爲一個樣本數據，每一列爲一個特徵，
    # 因此進行轉置，讓每一列表明一個樣本數據
    X = X.T
    # 2.將數據的每一行（表明一個屬性字段）進行零均值化
    meanValues = np.mean(X, axis=1) #計算每一行的均值
    meanValues = np.mat(meanValues).T #將一個向量轉換成n*1的矩陣
    meanRemoved = X - meanValues   #均值歸零
    
    #3.求出協方差矩陣
    covMat = np.cov(meanRemoved) #cov計算協方差時除的是 （樣本個數-1）,也就是自由度
# covMat = meanRemoved @ meanRemoved.T / (meanRemoved.shape[1])
    print("協方差矩陣：\n",covMat)
    
    #4求出協方差矩陣的特徵值及對應的特徵向量
    eigVals, eigVects = np.linalg.eig(covMat)  #eigVals:特徵值 eigVects：特徵向量
    print("特徵值\n",eigVals)
    print("特徵向量\n",eigVects)
    
    #5將特徵向量按對應特徵值大小從左到右按列排列成矩陣，取前k列組成矩陣
    # argsort函數返回的是數組值從小到大的索引值,參數中加個-號，變爲從大到小
    eigValInd = np.argsort(-eigVals)
    eigValInd = eigValInd[0:topNfeat]  #取出前topNfeat個最大的特徵值所對應的索引
    redEigVects = eigVects[:,eigValInd]  #redEigVects 即爲須要的變換矩陣，即P
    print("變換矩陣：\n",redEigVects)
    
    #6 Y=P^TX即爲降維到k維後的數據
    X_PCA = redEigVects.T @ X
    return X_PCA
    
    
複製代碼

PCA(X, topNfeat = 1)
複製代碼

協方差矩陣：
 [[1.5 1. ]
 [1.  1.5]]
特徵值
 [2.5 0.5]
特徵向量
 [[ 0.70710678 -0.70710678]
 [ 0.70710678  0.70710678]]
變換矩陣：
 [[0.70710678]
 [0.70710678]]





array([[-2.12132034, -0.70710678,  0.        ,  2.12132034,  0.70710678]])
複製代碼

sklearn

sklearn中爲咱們已經封裝好了對應的PCA接口，下面咱們使用PCA對sklearn中自帶的一個手寫數字數據集進行降維。

from sklearn import datasets
複製代碼

digits = datasets.load_digits()
X = digits.data
y = digits.target
X.shape,y.shape
複製代碼

((1797, 64), (1797,))
複製代碼

首先劃分數據集爲訓練集和測試集

from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=666)
X_train.shape
複製代碼

(1347, 64)
複製代碼

首先，在不進行pca降維的狀況下，使用knn算法進行模型的創建和預測

%%time

from sklearn.neighbors import KNeighborsClassifier

knn_clf = KNeighborsClassifier()
knn_clf.fit(X_train, y_train)
knn_clf.score(X_test,y_test)  #準確率爲
複製代碼

Wall time: 69.8 ms
複製代碼

knn_clf.score(X_test,y_test)  #準確率爲
複製代碼

0.9866666666666667
複製代碼

接下來看下 PCA降維

#decomposition(分解)
from sklearn.decomposition import PCA
複製代碼

# n_components 要降到的維數，第一次咱們直接降到2維
pca = PCA(n_components =2)
pca.fit(X_train)
X_train_reducation = pca.transform(X_train)
X_test_reducation = pca.transform(X_test) #注意測試集也須要作降維處理
X_train_reducation.shape
複製代碼

(1347, 2)
複製代碼

pca.explained_variance_ratio_ 
#表明降維後的各主成分的方差值佔總方差值的比例，這個比例越大，則越是重要的主成分。
複製代碼

array([0.14566817, 0.13735469])
複製代碼

%%time 

knn_clf = KNeighborsClassifier()
knn_clf.fit(X_train_reducation, y_train) # #降維後，所用時間明顯減小
複製代碼

Wall time: 2 ms
複製代碼

knn_clf.score(X_test_reducation,y_test) 
#因爲咱們將64維的數據直接降到了兩維，因此信息損失嚴重，準確率只有0.6
複製代碼

0.6066666666666667
複製代碼

#若是咱們直接傳入一個小數值，則表示要保留的信息佔原始信息的比例
#好比下面，保留95%的主要信息
pca = PCA(0.95)
pca.fit(X_train)
X_train_reducation = pca.transform(X_train)
X_test_reducation = pca.transform(X_test)
複製代碼

pca.n_components_  #n_components_能夠查看pca後降到的具體維數
複製代碼

28
複製代碼

knn_clf = KNeighborsClassifier()
knn_clf.fit(X_train_reducation, y_train)
knn_clf.score(X_test_reducation,y_test)   
複製代碼

0.98
複製代碼

能夠看到，準確度基本沒變化，可是維度卻從以前的64維降到28維，從而節約了大量資源。這就是PCA主成分分析法的原理和代碼實現

參考文章： www.cnblogs.com/mikewolf200…

歡迎關注個人我的公衆號 AI計算機視覺工坊，本公衆號不按期推送機器學習，深度學習，計算機視覺等相關文章，歡迎你們和我一塊兒學習，交流。

相關標籤/搜索

每日一句

每一个你不满意的现在，都有一个你没有努力的曾经。