1. 程式人生 > >PCA SVD原理詳解及應用

PCA SVD原理詳解及應用

本文分為兩大部分即PCA和SVD,每一部分下又分為原理和應用兩小部分

說明:本文程式碼參考Peter Harrington編寫的Machine Learning in Action,感興趣的小夥伴可以去看一下,筆者認為這本書還不錯

注意:本篇重在說明公式推導,關於具體使用的話python有專門的機器學習庫已經整合,直接用就可以啦,可以在讀完本文的理論部分後再去看筆者另一篇應用了PCA的關於人臉識別的一個簡單例子https://blog.csdn.net/weixin_42001089/article/details/79989788

廢話不多說開始吧

PCA:

(1)原理

(1)XV=\lambda V

(2)XV=V\Sigma \Rightarrow X=V\Sigma V^{-1}=V\Sigma V^{T}

公式(1)是我們熟悉的求特徵值的式子,假設X是m*m,我們根據其求出了\lambda _{1},\lambda _{2},\cdots \lambda _{m}
特徵值以及每個特徵值對應的特徵向量。

---------------------------------------------------------------------------------------------------------------------------------------------------------------

關於怎麼求特徵值這裡簡單舉一個小例子:會的話直接跳過就可以啦

X=\begin{pmatrix} 2 & 1\\ 1& 2 \end{pmatrix}

\begin{vmatrix} 2-\lambda &1 \\1 & 2-\lambda \end{vmatrix}=(2-\lambda )^{2}-1=0\Rightarrow \lambda _{1}=1,\lambda _{2}=3

我們通過上面求得該矩陣的2個特徵值分別為1和3

\bigl(\begin{smallmatrix} 2-1 &1 \\ 1 & 2-1 \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 1 &1 \\ 1 & 1 \end{smallmatrix}\bigr)\Rightarrow \bigl(\begin{smallmatrix} 1 &1 \\ 0 & 0 \end{smallmatrix}\bigr)

上面通過矩陣化簡可以得到解為 k\left ( 1,-1 \right )^{T}

同理3這個特徵值對應\bigl(\begin{smallmatrix} 2-3 & 1\\ 1 & 2-3 \end{smallmatrix}\bigr)= \bigl(\begin{smallmatrix} -1 & 1\\ 1 & -1 \end{smallmatrix}\bigr)\Rightarrow \bigl(\begin{smallmatrix} -1 & 1\\ 0 & 0 \end{smallmatrix}\bigr)的解為:k\left ( 1,1 \right )^{T}

所以特徵向量分別是\left ( 1,-1 \right )^{T} , \left ( 1,1 \right )^{T}

---------------------------------------------------------------------------------------------------------------------------------------------------------------

 

這裡的\Sigma是一個對角矩陣,對角線上值就是所有的特徵值,那麼可以看到X最後可以根據(2)進行分解,一般的習慣VV^{T}=E

所以V^{-1}=V^{T}

------------------------------------------------------------------------------------------------------------------------------------------------------------------

對應到上面的例子中特徵向量應該是:

\left ( \frac{1}{\sqrt{2}}, \frac{-1}{\sqrt{2}} \right )^{T} ,\left ( \frac{1}{\sqrt{2}}, \frac{1}{\sqrt{2}} \right )^{T}

注意當求矩陣單個特徵值時是XV=\lambda V,但當所有特徵值放一起後是XV=V\Sigma而不是XV=\Sigma V這種形式。

即分解結果為:    \bigl(\begin{smallmatrix} 2 &1 \\ 1 & 2 \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} \frac{1}{\sqrt{2}} &\frac{1}{\sqrt{2}} \\ \frac{1}{\sqrt{2}} & \frac{-1}{\sqrt{2}} \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} 3 &0 \\ 0 & 1 \end{smallmatrix}\bigr)\bigl(\begin{smallmatrix} \frac{1}{\sqrt{2}} &\frac{1}{\sqrt{2}} \\ \frac{1}{\sqrt{2}} & \frac{-1}{\sqrt{2}} \end{smallmatrix}\bigr)

---------------------------------------------------------------------------------------------------------------------------------------------------------------------

(2)應用

在實際降維時,往往先要轉化座標軸,使座標軸儘可能的覆蓋多的資料,即找資料最大方差的位置,其往往給出了資料的最重要資訊,比如找到後我們記該方向為x軸,那麼其y軸就是與x軸垂直的方向。

所以思路就是找先將原始矩陣每一個元素減去其自身均值進行標準歸一化,然後求得其協方差,最後再求協方差的特徵值和特徵向量,注意不是求原矩陣的特徵值和特徵向量。

---------------------------------------------------------------------------------------------------------------------------------------------------------------

好了,程式碼也應該很簡單啦:

def pca(datamat,N=3):
    meanVals = mean(datamat,axis=0)
    meanRemoved = datamat - meanVals
    covmat = cov(meanRemoved,rowvar=0)
    eigVals,eigVects = linalg.eig(mat(covmat))
    eigValInd = argsort(eigVals)
    eigValInd = eigValInd[:-(N+1):-1]
    redEigVects = eigVects[:,eigValInd]
    lowDDataMat = meanRemoved*redEigVects
    reconMat = (lowDDataMat *redEigVects.T)+meanVals
    return lowDDataMat ,reconMat

就是在求出協方差的特徵值後進行排序,然後選取前N名(eigValInd)特徵值及其對應的特徵向量(redEigVects)

然後lowDDataMat就是其將原始矩陣降維後,變成的新維度,該維度與原始矩陣相比,行數不變,列數降維到N,但卻保留了原始矩陣所蘊含的大部分資訊

reconMat就是我們用lowDDataMat來重構回原始矩陣,為什麼會是這樣呢?

Cov = low\times red^{T}+meanVals=(meanRemoved\times red)\times red^{T}+meanVals=meanRemoved\times (red\times red^{T})+meanVals=meanRemoved+meanVals=datamat

這裡將 redEigVects 簡寫為red,證明過程使用到了上面所說的VV^{T}=E

重構只是為了驗證一下而已,我們真正想要的正是lowDDataMat即返回的第一個矩陣。即假設datamat是m*m,那麼降維後便是m*N

-------------------------------------------------------------------------------------------------------------------------------------------------------------

但上面分解的前提是X必須為方陣那如果X不是方陣呢即是m*n,這時候SVD隆重登場

---------------------------------------------------------------------------------------------------------------------------------------------------------------------

注意:上面的PCA的原始矩陣雖說可以不是方陣,但我們也不是求原始矩陣的特徵值和特徵向量呀,是求的其協方差,其協方差即cov = E((X-E(X)(X-E(X)^{T})這裡的協方差cov一定是方陣對吧,所以其本質還是對方陣進行的分解,而下面介紹的svd正是直接對不是方陣的情況進行分解

--------------------------------------------------------------------------------------------------------------------------------------------------------------------

SVD

(1)原理

首先看一下SVD矩陣分解公式:

(3)X_{m\times n}= U_{m\times m} \Sigma_{m\times n} V_{n\times n}^{T}

這裡的\Sigma_{m\times n}其實是一個對角矩陣,對角線上面的值就稱為奇異值,U_{m\times m}V_{n\times n}分別稱為左奇異矩陣、右奇異矩陣

-------------------------------------------------------------------------------------------------------------------------------------------------------------------

對比一下方陣求特徵值的公式:

(1)方陣分解那裡對角線是特徵值,這裡是奇異值

(2)方陣分解那裡左右維數相同,互為逆矩陣(VV^{-1}),而這裡是左右奇異矩陣

可以看到這兩種分解在形式上大體相似,只是維數不同

-------------------------------------------------------------------------------------------------------------------------------------------------------------

這裡的U,\Sigma ,V該怎麼求呢?

先來看一下U:

(XX^{T})v=\lambda v 我們求XX^{T}的特徵值和特徵向量,然後將所有的特徵向量展開放到一個矩陣,該矩陣就是U

(X^{T}X)v=\lambda v 我們求X^{T}X的特徵值和特徵向量,然後將所有的特徵向量展開放到一個矩陣,該矩陣就是V

注意這裡我們將U和V中所有特徵向量也進行了歸一化,即 UU^{T}=E , VV^{T}=E

----------------------------------------------------------------------------------------------------------------------------------------------------------------

這裡簡單證明一下:

(4)XX^{T}=U\Sigma V^{T}\times V\Sigma ^{T}U^{T}=U\Sigma ^{2}U^{T}\Rightarrow \left ( XX^{T} \right )U=U\Sigma ^{2}

(5)X^{T}X=V\Sigma ^{T}U^{T}\times U\Sigma V^{T}=V\Sigma ^{2}V^{T}\Rightarrow (X^{T}X)V=V\Sigma ^{2}

通過(4)和(5)可以看到U和V其實分別就是XX^{T}X^{T}X各自所有特徵向量展開合併的矩陣

---------------------------------------------------------------------------------------------------------------------------------------------------------------

那麼\Sigma怎麼求呢,從上面的簡單證明可以明顯看出,其就是特徵值的平方根,

好了最後總結一下,當我們拿到一個X_{m\times n}矩陣需要分解時,首先是求XX^{T}的特徵值和對應的特徵向量,將所有特徵向量合併為一個矩陣,該矩陣就是U,給所有特徵值開平方根就是奇異值,即得到\Sigma,然後求X^{T}X的特徵值和特徵向量,同樣將所有特徵向量合併為一個矩陣,該矩陣就是V(不難看到,其實兩次的特徵最後求得結果應該是一樣的)

 

到這裡我們就將SVD的分解原理敘述完畢啦,回頭看看其實還是挺簡單的對吧,下面來看看SVD的具體用途:

(2)應用

想PCA取前幾個比較大的特徵值一樣,這裡SVD是取前幾個比較大奇異值

用的比較多的應該就是推薦系統裡面的吧

假設現在有一個矩陣,每一行代表一個使用者,每列代表一個物品

  辣條 薯片 可樂 冰紅茶
小花 0 0 2 4
小明 0 0 2 5
小李 3 6 0 0

其中的數字是評價等級1~5 , 0代表沒有評價

我們將該矩陣分解,X_{3\times 4}=U_{3\times 3}\Sigma _{3\times 4}V_{4\times 4}\approx U_{3\times 2}\Sigma _{2\times 2}V_{2\times 4}

這裡假設我們取2個奇異值,那麼便可以近似分解

 

可以看到原本是使用黃色來重構X,現在只使用紅色就可以近似重構X,即達到了降維目的,那麼這裡的U_{3*2}具體是什麼含義呢?

它其實就是將每個使用者對映到二維(原先是4維,分別代表辣條,薯片,可樂,冰紅茶),在這兩維中,小花和小明相似度相近,而小李

這位另一類,這和原始的評分矩陣反應出來的一樣,即小花和小明相似度相近,而小李這位另一類,原先4個維度,現在我們只用2個維度就可以反應出同樣的資訊,這不就是降維了嗎(至於這兩個維度是什麼新的含義呢?我們不必管,或是零食和飲料又或是其他的),總之這裡對於每個使用者來說,將原來的4維降為了2維,但其中還是包含原始原始矩陣的資訊。同理V_{2\times 4}^{T}是對每個物品來說,將原來的3維降為2維,即原始矩陣是這樣的:每個物品下面有3個使用者,而這裡是,每個物品下面有2類人。

其實上面兩種解釋對應著兩種不同的推薦策略:基於使用者(客戶數不變,即行數不變,壓縮列數,歸類物品)和基於物品(物品數不變,即列數不變,壓縮行數,歸類客戶)。

在python中有現成的庫可以用.來進行SVD矩陣分解,其返回的是三個矩陣,分別對應U,\Sigma ,V

--------------------------------------------------------------------------------------------------------------------------------------------------------------

例如:

U,Z,V = linalg.svd(datamat)

這裡的datamat就是我們要分解的矩陣。

現在比如一共求出有10個奇異值,我們取3個來進行降維,即我們要將原始矩陣的列數壓縮到3維(也可以理解為對於每一個使用者來說,原來是對很多商品進行了評價,現在我們要將其總結一下,使其結果是對總結後的3類商品進行的評價),即這裡是求V:

 

X=U\Sigma V^{T}\Rightarrow XV=U\Sigma \Rightarrow U=XV\Sigma ^{-1}

假設原來datamat是15*10,那麼降維後這裡的X_tran就是15*3啦

同理如果是給行降維,即將人進行分類,就是說求V

X=U\Sigma V^{T}\Rightarrow V^{T}=\Sigma ^{-1}U^{T}X\Rightarrow V=X^{T}U(\Sigma ^{-1})^{T}\Rightarrow V=X^{T}U(\Sigma ^{T})^{-1}\Rightarrow V=X^{T}U\Sigma ^{-1}

X_tran = datamat.T*U[:,:3]*mat(eye(4)*U[:3]).I

假設原來datamat是15*10,那麼降維後這裡的X_tran就是3*10啦

--------------------------------------------------------------------------------------------------------------------------------------------------------------

現實中一般是行數對於列數(客戶多,商品少),所以多采用基於商品,即要得到V

假設現在我們要給使用者a推薦物品,那怎麼做呢?

首先我們使用V=X^{T}U\Sigma ^{-1}降維得到V,然後我們從datamat中先篩選出來客戶a還沒有買過的物品標號集合為R1,因為我們要從這裡面挑出一些物品推薦給a,以及其評價過物品的集合R2,最後就是遍歷R1,看其每一個物品和R2中每一個物品的相似度再乘以權重(即評價)

-----------------------------------------------------------------------------------------------------------------------------------------------------------------------

可能上面說的有點不好理解,看一下程式碼吧,就好多啦:

程式碼參考;《機器學習實戰》Peter Harrington著,李銳等譯

def Similarity(datamat,user,simways,item):
    n = shape(datamat)[1]
    simTotal = 0.0
    ratSimTotal = 0.0
    U,Z,V = linalg.svd(datamat)
    X_tran = datamat.T*U[:,:3]*mat(eye(4)*U[:3]).I 
    for j in range(n):
        userRating = datamat[user,j]
        if userRating==0 or j==item:
            continue
        similarity = simways(X_tran[item,:].T,(X_tran[j,:].T)
        simTotal += similarity
        ratSimTotal +=similarity * userRating
     if simTotal=0:
         return 0
     else:
         return ratSimTotal/simTotal

上面的函式有四個輸入,原始矩陣,要給哪個使用者推薦,計算相似度的方法,以及當前推薦物品的id

從中可以看到我們取了3個奇異值

 X_tran = datamat.T*U[:,:3]*mat(eye(4)*U[:3]).I 

對應我們求V的過程,即X_tran就是我們取3個奇異值後對應的V

if userRating==0 or j==item:
            continue

對應部分就是說要計算item和已經評過分的物品的相似度,加入j之前也沒評價過,那比較就沒有意義啦,這裡就相當於上面我們說的處理出R2的過程

similarity = simways(X_tran[item,:].T,(X_tran[j,:].T)

部分就是求兩個向量的相似度

ratSimTotal +=similarity * userRating

部分就是相似度乘以權重,不難想象,如果j的評分高,則userRating高,權重大,此時如果相似度也高,那麼最後結果也高,如果j的評分底,相似度高,相當於帶來的懲罰也越多(和評分差的相似度高,當然要給與高懲罰)

return ratSimTotal/simTotal

部分就是將結果該歸一化到1~5評分等級中

 

好啦!!!!!!!!!!!!!!!!!

有了上面的評分函式,接下來就是在外面遍歷R1就可以啦

def recommend(datamat,user,N=3,simways):
    unratedItems = nonzero(datamat[user,:].A==0)[1]
    if len(unratedItems) == 0:
        return 'you rated everything'
    itemScores = []
    for item in unratedItems:
        estimatedScore = Similarity(datamat,user,simways,item)
        itemScores.append((item, estimatedScore))
    return sorted(itemScores,key=lambda x:x[1],reverse=True)[:N]

這裡的 unratedItems就是我們上面討論的R1集合,N引數就是我們要給該使用者推薦物品的個數

--------------------------------------------------------------------------------------------------------------------------------------------------------

關於上面的simways即相似度的量化方法有很多,比較容易相當的就是歐式距離和餘弦相似度,除此之外的就是皮爾遜相關係數,本文主要講解SVD,這裡就不展開說明相似度的事情啦,網上搜索一下,都很簡單。

 

最後說一下兩者的聯絡,其實svd可以用來做pca的,因為在pca階段我們主要求得是XX^{T}特徵值和特徵向量,對吧(暫不考慮減去均值),而在求svd的左奇異矩陣(U)過程中,也是求得是XX^{T}特徵值和特徵向量對吧,也就是說這裡的SVD奇異值的平方就是pca中的特徵值對吧,其實在SVD的分解過程中,不用像我們上面討論的那樣先求XX^{T}X^{T}X特徵值和特徵向量方法進行分解,而是有很多快速演算法,但我們可以根據SVD的結果進行PCA對吧,其實scikit-learn內部的PCA也是用SVD做的

結束啦!!!!!哪裡有不對還望大佬指正