PCA原理與實踐

標籤：通過   stc   float   current   ogr   pyw   gray   cell   slax   

　　在對資料進行預先處理時，我們經常會遇到資料的維數非常之大，如果不進行相應的特徵處理，那麼演算法的資源開銷會很大，這在很多情境下是我們不能接受的。而對於資料的若干維度之間往往會存在較大的相關性，如果能將資料的維度之間進行相應的處理，使它們在保留最大資料資訊的同時降低維度之間的相關性，就可以達到降維的效果。PCA（主成分分析）便是利用這樣的概念將資料對應到新的維度空間中，選擇最重要的幾個成分作為新空間向量的基，這樣在新的座標空間中，資料既可以保留大部分的資料資訊又可以達到降維的效果。在機器學習實戰中對於PCA的描述太過簡單，只是利用代碼可能並不能為我們詳細理解其背後的演算法原理，通過學習一些部落格，其中有一篇部落格http://blog.codinglabs.org/articles/pca-tutorial.html，這篇部落格通俗易懂讓我從中學習到了很多，現在將學習到的知識總結如下：　

總體思路：

　　資料最初所在空間的維度數目是由資料的特徵數量決定的，我們想要把它映射到新的空間中，這個空間的維度數目是由資料最重要的特徵所決定並且當資料對應後能夠儘可能的保留資料資訊，而且資料的分離程度較明顯。這裡就用到了數學上方差和共變數的概念，PCA演算法的主要理論依據便是基於這些概念進行一步步求解最佳化的。

演算法原理：

　　因為在http://blog.codinglabs.org/articles/pca-tutorial.html這篇部落格中原理解釋的非常棒，所以以下是基於這篇這部落格的原理講解，感謝博主！

向量內積以及矩陣相乘的幾何意義：

　　兩個向量內積從公式上可以看出最終得到了一個實數，也就是說內積運算是將兩個向量映射為一個實數。

　　那麼內積的幾何意義是什麼呢？其實它是將向量A投影到了B上的向量長度。假設A和B是兩個n維向量，我們知道n維向量可以等價表示為n維空間中的一條從原點發射的有向線段，為了簡單起見我們假設A和B均為二維向量，則A=(x1,y1)，B=(x2,y2)。則在二維平面上A和B可以用兩條發自原點的有向線段表示，見：

現在我們從A點向B所在直線引一條垂線。我們知道垂線與B的交點叫做A在B上的投影，再設A與B的夾角是a，則投影的向量長度為|A|cos(a)，其中

是向量A的模，也就是A線段的標量長度。注意這裡我們專門區分了向量長度和標量長度，標量長度總是大於等於0，值就是線段的長度；而向量長度可能為負，其絕對值是線段長度，而符號取決於其方向與標準方向相同或相反。同時向量內積還有一個通用的公式：

A?B=|A||B|cos(a)

如果我們假設B的模為1，即讓|B|=1，那麼就變成了：A?B=|A|cos(a)，也就是說設向量B的模為1，則A與B的內積值等於A向B所在直線投影的向量長度

矩陣運算的幾何意義：

如果我們確定了N維空間中的一組基，就可以利用矩陣相乘將向量在該空間中進行線性變換。現在我們有一個向量，要確定它在由N個基向量確定的空間中的具體位置，就可以利用該向量在各個基向量上的投影來確定。通常基向量都是單位正交基，這樣我們可以抽取出一個方向上的基向量與該向量做內積運算，利用上一節的內積幾何意義不就得到了向量在基向量上的投影，因此可以利用矩陣形式來做向量的基變換。用一個例子來說明：

在二維空間中，一般情況下我們是以（0,1）和（1,0）作為一組基的，但同時(1,1)和(-1,1)也可以成為一組基。一般來說，我們希望基的模是1，因為從內積的意義可以看到，如果基的模是1，那麼就可以方便的用向量點乘基而直接獲得其在新基上的座標了！因此將其單位化上面的基可以變為和。現在，我們想獲得(3,2)在新基上的座標，即在兩個方向上的投影向量值，那麼根據內積的幾何意義，我們只要分別計算(3,2)和兩個基的內積，不難得到新的座標為。給出了新的基以及(3,2)在新基上座標值的：

矩陣形式：

即將右邊矩陣中的每一列列向量變換到左邊矩陣中每一行行向量為基所表示的空間中去。

PCA演算法原理：

　　我們在做PCA時希望能夠得到資料在哪個方向上的分離程度最大並且能夠保留資料的資訊儘可能多。在數學上方差能夠很好的表徵資料的離散程度，方差越大資料越分離，我們越容易辨認。而共變數能夠很好的表徵資料在不同維度上的相關性。因此我們在處理資料時先對資料去均值化，將均值變為0，然後將資料矩陣乘以矩陣的轉置得到了我們需要的共變數矩陣。在該共變數矩陣中，主對角線上的元素為方差，其餘元素為共變數。

　　現在所要做的事情就是要選擇出方差最大的向量方向作為基向量，以此為基礎在選擇與第一個基向量正交並且方差第二大的向量作為第二維度基向量，以此類推。通常前r個方向上的基向量就能夠最大化的近似未經處理資料，這裡r遠遠小於n，從而達到降維的效果。因此整個邏輯就是最佳化共變數矩陣，將其對角化操作。由於共變數矩陣是個實對稱矩陣故而其必然能夠對角化，而且它的不同特徵值所對應的特徵向量都正交。由線性代數的知識我們可以很容易的得到共變數矩陣的特徵值和特徵向量，特徵值表示特徵向量的重要程度，特徵向量則表示了資料方差變化的方向。對特徵值按降序排列後，將其所對應的特徵向量構成了新的座標空間中的基向量。這樣將未經處理資料重構於新的座標空間，根據矩陣相乘的意義就可以做出映射後的資料。

　　這裡證明一下為什麼我們需要的空間向量基就是共變數矩陣的特徵向量。設未經處理資料矩陣X對應的共變數矩陣為C，而P是一組基按行組成的矩陣，設Y=PX，則Y為X對P做基變換後的資料。設Y的共變數矩陣為D，我們推導一下D與C的關係：

　　我們要找的P不是別的，而是能讓原始共變數矩陣對角化的P，即就是未經處理資料共變數的特徵向量。換句話說，最佳化目標變成了尋找一個矩陣P，滿足PCPT是一個對角矩陣，並且對角元素按從大到小依次排列，那麼P的前r行就是要尋找的基，用P的前r行組成的矩陣乘以X就使得X從N維降到了r維並滿足上述最佳化條件。

機器學習實戰應用執行個體：

1. 將資料轉換成前N個主成分

   1 def loadData(path):2     fr = open(path)3     dataSet= []4     for line in fr.readlines():5         currentData = line.strip().split(‘ ‘)6         dataSet.append(currentData)7     dataArr = [map(float,lineData) for lineData in dataSet]8     return mat(dataArr

    

 #PCA
 1 def PCA(dataMat,featureNum): 2     dataMean = np.mean(dataMat , axis = 0)            #求均值 3     dataMat = dataMat -dataMean                       #資料去均值化 4     covDataMat = np.cov(dataMat,rowvar=0)             #求共變數 5     featVal,featVec = linalg.eig(covDataMat)          #求特徵值和特徵向量 6     featIndex = argsort(featVal)                      #對特徵值排序後取前k個特徵向量 7     featIndex=list(reversed(list(featIndex))) 8     featIndex = featIndex[:featureNum:1] 9     featVec = mat(featVec[:,featIndex])    10     newDataMat = dataMat * featVec                    #降維11     originMat = (newDataMat*featVec.T)+dataMean       #重構未經處理資料12     return newDataMat,originMat       

#畫圖
1 def plotData(dataMat,originMat):2     fig =plt.figure()3     ax = fig.add_subplot(111)4     ax.scatter(dataMat[:,0].flatten().A[0],dataMat[:,1].flatten().A[0],marker=‘^‘,s=90)5     ax.scatter(originMat[:,0].flatten().A[0],dataMat[:,1].flatten().A[0],s=50,c=‘red‘)6     plt.show()

　　2.針對機器學習實戰提供的有關半導體資料，含有590個特徵。我們想要做的是，這些特徵中有多少主特徵，把這些主特徵提取出來。由於資料中含有很多缺失值是以NaN標識的，所以第一步要用均值來代替這些資料。

1 def replaceNan(dataMat):2     m = shape(dataMat)[1]3     for i in range(m):4         dataOfMean = mean(dataMat[nonzero(~isnan(dataMat[:,i]))[0],i])5         dataMat[nonzero(isnan(dataMat[:,i]))[0],i] = dataOfMean6     return dataMat

稍微改動下PCA程式在裡面加入判別，如果前k個主成分的累積方差佔總方差的百分比大於98%就停止，選擇這K個特徵向量：

def PCA(dataMat,featureNum):    dataMean = np.mean(dataMat , axis = 0)    dataMat = dataMat -dataMean    covDataMat = np.cov(dataMat,rowvar=0)    featVal,featVec = linalg.eig(covDataMat)    featIndex = argsort(featVal)    featIndex=list(reversed(list(featIndex)))    featIndex = featIndex[:featureNum:1]    sumOfFeatVal = sum(featVal)    addtotal = 0.0    for k in range(featureNum):        addtotal +=featVal[featIndex[k]]        percent= addtotal / sumOfFeatVal        if percent>0.98:            print ‘the number of %d has occupied %f‘ % (k, percent)            print ‘the best number of feature is ‘,k            break        print ‘the number of %d has occupied %f‘%(k,percent)    fig = plt.figure()    ax = fig.add_subplot(111)    percentFeat = featVal[featIndex][0:20] / sumOfFeatVal    x = arange(20)    ax.plot(x, percentFeat, ‘o-‘, c=‘r‘)    plt.xlabel(‘the number of Princple feature‘)    plt.ylabel(‘var percent(%)‘)    plt.grid()    plt.show()    featVec = mat(featVec[:,featIndex])    newDataMat = dataMat * featVec    originMat = (newDataMat*featVec.T)+dataMean    return newDataMat,originMat,featVal

 

最後的統計結果：

the number of 0 has occupied 0.592541the number of 1 has occupied 0.833779the number of 2 has occupied 0.925279the number of 3 has occupied 0.948285the number of 4 has occupied 0.962877the number of 5 has occupied 0.968065the number of 6 has occupied 0.971291the number of 7 has occupied 0.974438the number of 8 has occupied 0.977069the number of 9 has occupied 0.979382the number of 10 has occupied 0.981557the best number of feature is  10

參考：

1.部落格http://blog.codinglabs.org/articles/pca-tutorial.html

2.《機器學習實戰》

PCA原理與實踐