? 本文主要基于同名的兩篇外文參考文獻 A Tutorial on Principal Component Analysis 。

? PCA，亦即主成分分析，主要用于對特征進行降維。如果數(shù)據(jù)的特征數(shù)非常多，我們可以認(rèn)為其中只有一部分特征是真正我們感興趣和有意義的，而其他特征或者是噪音，或者和別的特征有冗余。從所有的特征中找出有意義的特征的過程就是降維，而PCA是降維的兩個主要方法之一(另一個是LDA).

? Jonathon Shlens的論文中舉了一個物理學(xué)中測試?yán)硐肭闆r下彈簧振動的例子，非常生動，詳見[1](中文翻譯見[5])。

? 我們首先看一下給定一個代表數(shù)據(jù)記錄的矩陣A，如果計算其主成分P，并如何利用P得到降維后的數(shù)據(jù)矩陣B，然后介紹一下這個計算過程背后的原理，最后會有在Python中實現(xiàn)PCA和在Weka中調(diào)用PCA算法的實例。

? 1. 計算過程：

? 假設(shè)我們有n條數(shù)據(jù)記錄，每條記錄都是m維，我們可以把這些數(shù)據(jù)記錄表示成一個n*m矩陣A。

? 對矩陣A的每一列，求出其平均值，對于A中的每一個元素，減去該元素所在列的平均值，得到一個新的矩陣B。

? 計算矩陣Z=B ^T B/(n-1)。其實m*m維矩陣Z就是A矩陣的協(xié)方差矩陣。

? 計算矩陣Z的特征值D和特征向量V,其中D是1*m矩陣，V是一個m*m矩陣，D中的每個元素都是Z的特征值，V中的第i列是第i個特征值對應(yīng)的特征向量。

? 下面，就可以進行降維了，假設(shè)我們需要把數(shù)據(jù)由m維降到k維，則我們只需要從D中挑選出k個最大的特征向量，然后從V中挑選出k個相應(yīng)的特征向量，組成一個新的m*k矩陣N。

? N中的每一列就是A的主成分(Principal Component). 計算A*N得到n*k維矩陣C，就是對源數(shù)據(jù)進行降維后的結(jié)果，沒條數(shù)據(jù)記錄的維數(shù)從m降到了k。

? 2. 原理

? 要對數(shù)據(jù)進行降維的主要原因是數(shù)據(jù)有噪音，數(shù)據(jù)的軸(基)需要旋轉(zhuǎn)，數(shù)據(jù)有冗余。

? (1) 噪音

? 上圖是一個記錄彈簧振動的二維圖。我們發(fā)現(xiàn)沿正45度方向數(shù)據(jù)的方差比較大，而沿負(fù)45度方向數(shù)據(jù)的方差比較小。通常情況下，我們都認(rèn)為方差最大的方向記錄著我們感興趣的信息，所以我們可以保留正45度方向的數(shù)據(jù)信息，而負(fù)45度方向的數(shù)據(jù)信息可以認(rèn)為是噪音。

? (2) 旋轉(zhuǎn)

? 在線性代數(shù)中我們知道，同一組數(shù)據(jù)，在不同的基下其坐標(biāo)是不一樣的，而我們一般認(rèn)為基的方向應(yīng)該與數(shù)據(jù)方差最大的方向一致，亦即上圖中的基不應(yīng)該是X,Y軸，而該逆時針旋轉(zhuǎn)45度。

? (3) 冗余

? 上圖中的a，c分別代表沒有冗余和高度冗余的數(shù)據(jù)。在a中，某個數(shù)據(jù)點的X,Y軸坐標(biāo)值基本上是完全獨立的，我們不可能通過其中一個來推測另外一個，而在c中，數(shù)據(jù)基本上都分布在一條直線上，我們很容易從一個坐標(biāo)值來推測出另外一個坐標(biāo)值，所以我們完全沒有必要記錄數(shù)據(jù)在X,Y兩個坐標(biāo)軸上的值，而只需要記錄一個即可。數(shù)據(jù)冗余的情況跟噪音比較相似，我們只需要記錄方差比較大的方向上的坐標(biāo)值，方差比較小的方向上的坐標(biāo)值可以看做是冗余(或噪音).

? 上面三種情況，歸結(jié)到最后都是要求出方差比較大的方向(基)，然后在求數(shù)據(jù)在這個基下的坐標(biāo)，這個過程可以表示為：

? PX=Y。

? 其中k*m矩陣P是一個正交矩陣，P的每一行都對應(yīng)一個方差比較大的基。m*n矩陣X和k*n矩陣Y的每一列都是一個數(shù)據(jù)(這一點和1中不同,因為這是兩篇不同的論文,表示方式不一樣，本質(zhì)上是一樣的)。

? X是原數(shù)據(jù)，P是一個新的基，Y是X在P這個新基下的坐標(biāo)，注意Y中數(shù)據(jù)記錄的維數(shù)從m降到了k，亦即完成了降維。

? 但是我們希望得到的是一個什么樣的Y矩陣呢？我們希望Y中每個基下的坐標(biāo)的方差盡量大，而不同基下坐標(biāo)的方差盡量小，亦即我們希望C _Y =YY ^T /(n-1)是一個對角線矩陣。

? C _Y ? =YY ^T /(n-1)=P(XX ^T )P ^T /(n-1)

? 令A(yù)=XX ^T ，我們對A進行分解：A=EDE ^T

? 我們?nèi)=E ^T ，則C _Y =E ^T AE/(n-1)=E ^T EDE ^T E/(n-1),因為E ^T =E ^-1 ，所以C _Y =D/(n-1)是一個對角矩陣。

? 所以我們應(yīng)該取P的每一行為A的特征向量，得到的Y才會有以上性質(zhì)。

? 3. 實現(xiàn)

? 1. PCA的Python實現(xiàn)

? 需要使用Python的科學(xué)計算模塊numpy

      
         1
      
      
        import
      
      
         numpy as np

      
      
         2
      
      
         3
      
         mat=[(2.5,0.5,2.2,1.9,3.1,2.3,2,1,1.5,1.1),    (2.4,0.7,2.9,2.2,3.0,2.7,1.6,1.1,1.6,0.9
      
        )]

      
      
         4
      
      
        #
      
      
        轉(zhuǎn)置,每一行是一條數(shù)據(jù)
      
      
         5
      
         data=
      
        np.matrix(np.transpose(mat))

      
      
         6
      
         data_adjust=data-
      
        mean

      
      
         7
      
      
        #
      
      
        求協(xié)方差矩陣
      
      
         8
      
         covariance=np.transpose(data_adjust)*data_adjust/9

      
         9
      
      
        #
      
      
        求得協(xié)方差矩陣的特征值和特征向量
      
      
        10
      
         eigenvalues,eigenvectors=
      
        np.linalg.eig(covariance)         

      
      
        11
      
         feature_vectors=
      
        np.transpose(eigenvectors)

      
      
        12
      
      
        #
      
      
        轉(zhuǎn)換后的數(shù)據(jù)
      
      
        13
      
         final_data=feature_vectors*np.transpose(data_adjust)
    

? 2. 在Weka中調(diào)用PCA：

      
        import
      
      
         java.io.FileReader as FileReader

      
      
        import
      
      
         java.io.File as File

      
      
        import
      
      
         weka.filters.unsupervised.attribute.PrincipalComponents as PCA

      
      
        import
      
      
         weka.core.Instances as Instances

      
      
        import
      
      
         weka.core.converters.CSVLoader as CSVLoader

      
      
        import
      
      
         weka.filters.Filter as Filter



      
      
        def
      
      
         main():
    
      
      
        #
      
      
        使用Weka自帶的數(shù)據(jù)集cpu.arff
      
      
    reader=FileReader(
      
        '
      
      
        DATA/cpu.arff
      
      
        '
      
      
        )
    data
      
      =
      
        Instances(reader)
    pca
      
      =
      
        PCA()
    pca.setInputFormat(data)
    pca.setMaximumAttributes(
      
      5
      
        )
    newData
      
      =
      
        Filter.useFilter(data,pca)
    
      
      
        for
      
       n 
      
        in
      
      
         range(newData.numInstances()):
        
      
      
        print
      
      
         newData.instance(n)

      
      
        if
      
      
        __name__
      
      ==
      
        '
      
      
        __main__
      
      
        '
      
      
        :
    main()
      
    

? 參考文獻：

? [1]. Jonathon Shlens. A Tutorial on Principal Component Analysis.

? [2]. Lindsay I Smith.? A Tutorial on Principal Component Analysis.

? [3]. 關(guān)于PCA算法的一點學(xué)習(xí)總結(jié)

? [4]. 主成分分析PCA算法 ?原理解析

? [5]. 主元分析(PCA)理論分析及應(yīng)用

PCA算法 原理與實現(xiàn)