PCA主成分分析python实现
Github源碼:https://github.com/csuldw/MachineLearning/tree/master/PCA
?
PCA(principle component analysis) ,主成分分析,主要是用來降低數據集的維度,然后挑選出主要的特征。原理簡單,實現也簡單。關于原理公式的推導,本文不會涉及,你可以參考下面的參考文獻,也可以去Wikipedia,這里主要關注實現,算是鍛煉一下自己。
本來是在復習LDA的,然后就看到了PCA,就跟著下面這篇文章的步驟,把PCA用Python實現了一遍,具體的思想可以參考這篇文章,講的通俗易懂,主要是有個實例參考,值得擁有!
下面引用的是http://www.cnblogs.com/jerrylead/archive/2011/04/18/2020209.html中關于主成分分析的講解
?
主成分分析(Principal components analysis)-最大方差解釋
???? 在這一篇之前的內容是《Factor Analysis》,由于非常理論,打算學完整個課程后再寫。在寫這篇之前,我閱讀了PCA、SVD和LDA。這幾個模型相近,卻都有自己的特點。本篇打算先介紹PCA,至于他們之間的關系,只能是邊學邊體會了。PCA以前也叫做Principal factor analysis。
1. 問題
???? 真實的訓練數據總是存在各種各樣的問題:
1、 比如拿到一個汽車的樣本,里面既有以“千米/每小時”度量的最大速度特征,也有“英里/小時”的最大速度特征,顯然這兩個特征有一個多余。
2、 拿到一個數學系的本科生期末考試成績單,里面有三列,一列是對數學的興趣程度,一列是復習時間,還有一列是考試成績。我們知道要學好數學,需要有濃厚的興趣,所以第二項與第一項強相關,第三項和第二項也是強相關。那是不是可以合并第一項和第二項呢?
3、 拿到一個樣本,特征非常多,而樣例特別少,這樣用回歸去直接擬合非常困難,容易過度擬合。比如北京的房價:假設房子的特征是(大小、位置、朝向、是否學區房、建造年代、是否二手、層數、所在層數),搞了這么多特征,結果只有不到十個房子的樣例。要擬合房子特征->房價的這么多特征,就會造成過度擬合。
4、 這個與第二個有點類似,假設在IR中我們建立的文檔-詞項矩陣中,有兩個詞項為“learn”和“study”,在傳統的向量空間模型中,認為兩者獨立。然而從語義的角度來講,兩者是相似的,而且兩者出現頻率也類似,是不是可以合成為一個特征呢?
5、 在信號傳輸過程中,由于信道不是理想的,信道另一端收到的信號會有噪音擾動,那么怎么濾去這些噪音呢?
???? 回顧我們之前介紹的《模型選擇和規則化》,里面談到的特征選擇的問題。但在那篇中要剔除的特征主要是和類標簽無關的特征。比如“學生的名字”就和他的“成績”無關,使用的是互信息的方法。
???? 而這里的特征很多是和類標簽有關的,但里面存在噪聲或者冗余。在這種情況下,需要一種特征降維的方法來減少特征數,減少噪音和冗余,減少過度擬合的可能性。
???? 下面探討一種稱作主成分分析(PCA)的方法來解決部分上述問題。PCA的思想是將n維特征映射到k維上(k<n),這k維是全新的正交特征。這k維特征稱為主元,是重新構造出來的k維特征,而不是簡單地從n維特征中去除其余n-k維特征。
2. PCA計算過程
???? 首先介紹PCA的計算過程:
???? 假設我們得到的2維數據如下:
?????
???? 行代表了樣例,列代表特征,這里有10個樣例,每個樣例兩個特征。可以這樣認為,有10篇文檔,x是10篇文檔中“learn”出現的TF-IDF,y是10篇文檔中“study”出現的TF-IDF。也可以認為有10輛汽車,x是千米/小時的速度,y是英里/小時的速度,等等。
?????第一步分別求x和y的平均值,然后對于所有的樣例,都減去對應的均值。這里x的均值是1.81,y的均值是1.91,那么一個樣例減去均值后即為(0.69,0.49),得到
?????
?????第二步,求特征協方差矩陣,如果數據是3維,那么協方差矩陣是
?????
???? 這里只有x和y,求解得
?????
???? 對角線上分別是x和y的方差,非對角線上是協方差。協方差大于0表示x和y若有一個增,另一個也增;小于0表示一個增,一個減;協方差為0時,兩者獨立。協方差絕對值越大,兩者對彼此的影響越大,反之越小。
?????第三步,求協方差的特征值和特征向量,得到
?????
???? 上面是兩個特征值,下面是對應的特征向量,特征值0.0490833989對應特征向量為,這里的特征向量都歸一化為單位向量。
????第四步,將特征值按照從大到小的順序排序,選擇其中最大的k個,然后將其對應的k個特征向量分別作為列向量組成特征向量矩陣。
???? 這里特征值只有兩個,我們選擇其中最大的那個,這里是1.28402771,對應的特征向量是。
?????第五步,將樣本點投影到選取的特征向量上。假設樣例數為m,特征數為n,減去均值后的樣本矩陣為DataAdjust(m*n),協方差矩陣是n*n,選取的k個特征向量組成的矩陣為EigenVectors(n*k)。那么投影后的數據FinalData為
?????
???? 這里是
???? FinalData(10*1) = DataAdjust(10*2矩陣)×特征向量
???? 得到結果是
?????
???? 這樣,就將原始樣例的n維特征變成了k維,這k維就是原始特征在k維上的投影。
???? 上面的數據可以認為是learn和study特征融合為一個新的特征叫做LS特征,該特征基本上代表了這兩個特征。
???? 上述過程有個圖描述:
?????
???? 正號表示預處理后的樣本點,斜著的兩條線就分別是正交的特征向量(由于協方差矩陣是對稱的,因此其特征向量正交),最后一步的矩陣乘法就是將原始樣本點分別往特征向量對應的軸上做投影。
???? 如果取的k=2,那么結果是
?????
???? 這就是經過PCA處理后的樣本數據,水平軸(上面舉例為LS特征)基本上可以代表全部樣本點。整個過程看起來就像將坐標系做了旋轉,當然二維可以圖形化表示,高維就不行了。上面的如果k=1,那么只會留下這里的水平軸,軸上是所有點在該軸的投影。
???? 這樣PCA的過程基本結束。在第一步減均值之后,其實應該還有一步對特征做方差歸一化。比如一個特征是汽車速度(0到100),一個是汽車的座位數(2到6),顯然第二個的方差比第一個小。因此,如果樣本特征中存在這種情況,那么在第一步之后,求每個特征的標準差,然后對每個樣例在該特征下的數據除以。
???? 歸納一下,使用我們之前熟悉的表示方法,在求協方差之前的步驟是:
?????
???? 其中是樣例,共m個,每個樣例n個特征,也就是說是n維向量。是第i個樣例的第j個特征。是樣例均值。是第j個特征的標準差。
???? 整個PCA過程貌似及其簡單,就是求協方差的特征值和特征向量,然后做數據轉換。但是有沒有覺得很神奇,為什么求協方差的特征向量就是最理想的k維向量?其背后隱藏的意義是什么?整個PCA的意義是什么?
3. PCA理論基礎
???? 要解釋為什么協方差矩陣的特征向量就是k維理想特征,我看到的有三個理論:分別是最大方差理論、最小錯誤理論和坐標軸相關度理論。這里簡單探討前兩種,最后一種在討論PCA意義時簡單概述。
3.1 最大方差理論
???? 在信號處理中認為信號具有較大的方差,噪聲有較小的方差,信噪比就是信號與噪聲的方差比,越大越好。如前面的圖,樣本在橫軸上的投影方差較大,在縱軸上的投影方差較小,那么認為縱軸上的投影是由噪聲引起的。
因此我們認為,最好的k維特征是將n維樣本點轉換為k維后,每一維上的樣本方差都很大。
???? 比如下圖有5個樣本點:(已經做過預處理,均值為0,特征方差歸一)
?????
???? 下面將樣本投影到某一維上,這里用一條過原點的直線表示(前處理的過程實質是將原點移到樣本點的中心點)。
?????
???? 假設我們選擇兩條不同的直線做投影,那么左右兩條中哪個好呢?根據我們之前的方差最大化理論,左邊的好,因為投影后的樣本點之間方差最大。
???? 這里先解釋一下投影的概念:
?????
???? 紅色點表示樣例,藍色點表示在u上的投影,u是直線的斜率也是直線的方向向量,而且是單位向量。藍色點是在u上的投影點,離原點的距離是(即或者)由于這些樣本點(樣例)的每一維特征均值都為0,因此投影到u上的樣本點(只有一個到原點的距離值)的均值仍然是0。
???? 回到上面左右圖中的左圖,我們要求的是最佳的u,使得投影后的樣本點方差最大。
???? 由于投影后均值為0,因此方差為:
?????
???? 中間那部分很熟悉啊,不就是樣本特征的協方差矩陣么(的均值為0,一般協方差矩陣都除以m-1,這里用m)。
???? 用來表示,表示,那么上式寫作
??????
???? 由于u是單位向量,即,上式兩邊都左乘u得,
???? 即
???? We got it!就是的特征值,u是特征向量。最佳的投影直線是特征值最大時對應的特征向量,其次是第二大對應的特征向量,依次類推。
???? 因此,我們只需要對協方差矩陣進行特征值分解,得到的前k大特征值對應的特征向量就是最佳的k維新特征,而且這k維新特征是正交的。得到前k個u以后,樣例通過以下變換可以得到新的樣本。
?????
???? 其中的第j維就是在上的投影。
???? 通過選取最大的k個u,使得方差較小的特征(如噪聲)被丟棄。
???? 這是其中一種對PCA的解釋,第二種是錯誤最小化,放在下一篇介紹。
PCA思想
?
主要思想:移動坐標軸,將n維特征映射到k維上(kn),這k維是全新的正交特征。這k維特征稱為主元,是重新構造出來的k維特征,而不是簡單地從n維特征中去除其余n-k維特征。
說到PCA難免會提到LDA(linear discriminate analysis,線性判別分析),以及FA(factor analysis,因子分析)。關于LDA,打算有時間也用代碼實現一遍,下面給出它的主要思想。
LDA思想:最大類間距離,最小類內距離。簡而言之,第一,為了實現投影后的兩個類別的距離較遠,用映射后兩個類別的均值差的絕對值來度量。第二,為了實現投影后,每個類內部數據點比較聚集,用投影后每個類別的方差來度量。
三者的描述如下
以下內容引自?Wikipedia- Linear discriminant analysis
LDA is also closely related to principal component analysis (PCA) and factor analysis in that they both look for linear combinations of variables which best explain the data.[4] LDA explicitly attempts to model the difference between the classes of data. PCA on the other hand does not take into account any difference in class, and factor analysis builds the feature combinations based on differences rather than similarities. Discriminant analysis is also different from factor analysis in that it is not an interdependence technique: a distinction between independent variables and dependent variables (also called criterion variables) must be made.
區別:PCA選擇樣本點投影具有最大方差的方向,LDA選擇分類性能最好的方向。
好了,下面來看下實現源碼!
基本步驟
基本步驟:
- 對數據進行歸一化處理(代碼中并非這么做的,而是直接減去均值)
- 計算歸一化后的數據集的協方差矩陣
- 計算協方差矩陣的特征值和特征向量
- 保留最重要的k個特征(通常k要小于n),也可以自己制定,也可以選擇一個閾值,然后通過前k個特征值之和減去后面n-k個特征值之和大于這個閾值,則選擇這個k
- 找出k個特征值對應的特征向量
- 將m * n的數據集乘以k個n維的特征向量的特征向量(n * k),得到最后降維的數據。
其實PCA的本質就是對角化協方差矩陣。有必要解釋下為什么將特征值按從大到小排序后再選。首先,要明白特征值表示的是什么?在線性代數里面我們求過無數次了,那么它具體有什么意義呢?對一個n*n的對稱矩陣進行分解,我們可以求出它的特征值和特征向量,就會產生n個n維的正交基,每個正交基會對應一個特征值。然后把矩陣投影到這N個基上,此時特征值的模就表示矩陣在該基的投影長度。特征值越大,說明矩陣在對應的特征向量上的方差越大,樣本點越離散,越容易區分,信息量也就越多。因此,特征值最大的對應的特征向量方向上所包含的信息量就越多,如果某幾個特征值很小,那么就說明在該方向的信息量非常少,我們就可以刪除小特征值對應方向的數據,只保留大特征值方向對應的數據,這樣做以后數據量減小,但有用的信息量都保留下來了。PCA就是這個原理。
源碼實現
1.首先引入numpy,由于測試中用到了pandas和matplotlib,所以這里一并加載
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
2.定義一個均值函數
#計算均值,要求輸入數據為numpy的矩陣格式,行表示樣本數,列表示特征
def meanX(dataX):
return np.mean(dataX,axis=0)#axis=0表示按照列來求均值,如果輸入list,則axis=1
3.編寫pca方法,具體解釋參考注釋
"""
參數:
- XMat:傳入的是一個numpy的矩陣格式,行表示樣本數,列表示特征
- k:表示取前k個特征值對應的特征向量
返回值:
- finalData:參數一指的是返回的低維矩陣,對應于輸入參數二
- reconData:參數二對應的是移動坐標軸后的矩陣
"""
def pca(XMat, k):
average = meanX(XMat)
m, n = np.shape(XMat)
data_adjust = []
avgs = np.tile(average, (m, 1))
data_adjust = XMat - avgs
covX = np.cov(data_adjust.T) #計算協方差矩陣
featValue, featVec= np.linalg.eig(covX) #求解協方差矩陣的特征值和特征向量
index = np.argsort(-featValue) #按照featValue進行從大到小排序
finalData = []
if k > n:
print "k must lower than feature number"
return
else:
#注意特征向量時列向量,而numpy的二維矩陣(數組)a[m][n]中,a[1]表示第1行值
selectVec = np.matrix(featVec.T[index[:k]]) #所以這里需要進行轉置
finalData = data_adjust * selectVec.T
reconData = (finalData * selectVec) + average
return finalData, reconData
4.編寫一個加載數據集的函數
#輸入文件的每行數據都以\t隔開
def loaddata(datafile):
return np.array(pd.read_csv(datafile,sep="\t",header=-1)).astype(np.float)
5.可視化結果
因為我將維數k指定為2,所以可以使用下面的函數將其繪制出來:
def plotBestFit(data1, data2):
dataArr1 = np.array(data1)
dataArr2 = np.array(data2)
m = np.shape(dataArr1)[0]
axis_x1 = []
axis_y1 = []
axis_x2 = []
axis_y2 = []
for i in range(m):
axis_x1.append(dataArr1[i,0])
axis_y1.append(dataArr1[i,1])
axis_x2.append(dataArr2[i,0])
axis_y2.append(dataArr2[i,1])
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(axis_x1, axis_y1, s=50, c='red', marker='s')
ax.scatter(axis_x2, axis_y2, s=50, c='blue')
plt.xlabel('x1'); plt.ylabel('x2');
plt.savefig("outfile.png")
plt.show()
6.測試方法
測試方法寫入main函數中,然后直接執行main方法即可:
data.txt可到github中下載:data.txt
#根據數據集data.txt
def main():
datafile = "data.txt"
XMat = loaddata(datafile)
k = 2
return pca(XMat, k)
if __name__ == "__main__":
finalData, reconMat = main()
plotBestFit(finalData, reconMat)
結果展示
最后的結果圖如下:
?
?
?
藍色部分為重構后的原始數據,紅色則是提取后的二維特征!
參考文獻
[1]?http://www.cnblogs.com/jerrylead/archive/2011/04/18/2020209.html?
[2]?Wikipedia- Linear discriminant analysis?
[3]?Wikipedia- Principal_component_analysis
總結
以上是生活随笔為你收集整理的PCA主成分分析python实现的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 播放Wav声音
- 下一篇: python监控桌面捕捉,用Python