PCA演算法

1. 去平均值，即每一位特徵減去各自的平均值
2. 計算新矩陣的協方差矩陣 設$X=(X_1, X_2…X_N)^T $，在鳶尾花例子裡N=4,會生成一個4*4的協方差矩陣 稱矩陣 $C=(c_{ij})_{n \times n}=\begin{pmatrix} c_{11} & c_{12} & \cdots & c_{1n} \\ c_{21} & c_{22} & \cdots & c_{2n} \\ \cdots & \cdots & \cdots & \cdots \\ c_{n1} & c_{n2} & \cdots & c_{nn} \end{pmatrix}$
   其中$c_{ij}=Cov(X_i, X_j)$

備註:協方差衡量兩個隨機變數的相關性，當協方差為0是認為兩個隨機變數是獨立的。 $Cov(X, Y) = E[(X-E[X])(Y-E(Y))]$ 3. 計算協方差矩陣的特徵值與特徵向量 4. 對特徵值從大到小排序，保留大的特徵值對應的特徵向量 5. 將原始資料矩陣乘特徵向量矩陣（減去忽略的特徵值）得到新的特徵空間 詳見程式例子

PCA程式實現（基於sklearn）

備註: 在本例中我們將鳶尾花資料的四個特徵通過PCA降低成兩個，即原有的150 * 4個特徵降維成150 * 2，之後再進行貝葉斯分析。

# -*- coding: utf-8 -*-
import numpy  as np
from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from sklearn.naive_bayes import GaussianNB 
from sklearn.pipeline import Pipeline


def iris_type(s):

    it = {b'Iris-setosa': 0, b'Iris-versicolor': 1, b'Iris-virginica': 2}
    return it[s]

def _test_pca():

    """
    採用鳶尾花卉資料集進行特徵降維
    """

    data = np.loadtxt("./data/iris.data", dtype=float, delimiter=',', converters={4: iris_type})
    x, y = np.split(data, (4,), axis=1)

    pca = PCA(n_components=4)
    pca.fit(x)

    """
    協方差矩陣 特徵值: [4.22484077 0.24224357 0.07852391 0.02368303]
    協方差矩陣 特徵所佔比重: [0.92461621 0.05301557 0.01718514 0.00518309]
    協方差矩陣 特徵向量: 
    [[ 0.36158968 -0.08226889  0.85657211  0.35884393]
    [ 0.65653988  0.72971237 -0.1757674  -0.07470647]
    [-0.58099728  0.59641809  0.07252408  0.54906091]
    [ 0.31725455 -0.32409435 -0.47971899  0.75112056]]
    """

    print('協方差矩陣 特徵值: {0}'.format(pca.explained_variance_))
    print('協方差矩陣 特徵所佔比重: {0}'.format(pca.explained_variance_ratio_))
    print('協方差矩陣 特徵向量: {0}'.format(pca.components_))

    #保留兩個特徵
    v = pca.components_[[0,1],:]


    """
    [[ 0.36158968 -0.08226889  0.85657211  0.35884393]
    [ 0.65653988  0.72971237 -0.1757674  -0.07470647]]
    """
    print(v)

    #重新構造特徵資料 x * v 稱為一個150 * 2的矩陣

    x_new = np.dot(x, np.transpose(v))

    #利用新的特徵資料來做分類
    x_train, x_test, y_train, y_test = train_test_split(x_new, y, test_size=0.3, random_state=1)


    #如果這裡採用MultinomialNB引數，則預設認為概率分佈為多項式分佈，在鳶尾花例子中會報錯：
    #ValueError: Input X must be non-negative
    model = Pipeline([
        ('sc', StandardScaler()),
        ('clf', GaussianNB())])

    model.fit(x_train, y_train)

    y_test_pre = model.predict(x_test)      

    # 訓練集上的預測結果
    y_test = y_test.reshape(-1)
    result = (y_test_pre == y_test)

    acc = np.mean(result)
    
    print('準確度: %.2f%%' % (100 * acc)) #93.33%

    pass


"""
說明：

SVD程式碼實現，對應的筆記《12.主成分分析(PCA)》

作者：fredric

日期：2018-9-12

"""
if __name__ == "__main__":

    _test_pca()