我們知道當資料維度太大時，進行分類任務時會花費大量時間，因此需要進行資料降維，其中一種非常流行的降維方法叫主成分分析。

Exploratory Data Analysis

鳶尾花資料集：

import numpy as np
from sklearn.datasets import load_iris

iris = load_iris()

type(iris.data) # numpy.ndarray
type(iris.target) # numpy.ndarray

print iris.feature_names
print iris.target_names

X = iris.data
y = iris.target
 

Out:
['sepal length (cm)', 'sepal width (cm)', 'petal length (cm)', 'petal width (cm)']
['setosa' 'versicolor' 'virginica']

4個特徵：花萼長度、花萼寬度、花瓣長度、花瓣寬度
3個分類：山鳶尾、變色鳶尾、維吉尼亞鳶尾

前十個樣本特徵值如下：

[[5.1, 3.5, 1.4, 0.2],
 [4.9, 3. , 1.4, 0.2],
 [4.7, 3.2, 1.3, 0.2],
 [4.6, 3.1, 1.5, 0.2],
 [5. , 3.6, 1.4, 0.2],
 [5.4
 
, 3.9, 1.7, 0.4],
 [4.6, 3.4, 1.4, 0.3],
 [5. , 3.4, 1.5, 0.2],
 [4.4, 2.9, 1.4, 0.2],
 [4.9, 3.1, 1.5, 0.1]]

樣本類別分佈：

print dict(zip(*np.unique(y,return_counts=True)))

Out:
{0: 50, 1: 50, 2: 50}

描述性統計：

print stats.describe(X)

Out:
DescribeResult(
nobs=150L, 
minmax=(array([4.3, 2. , 1. , 0.1]), array
 
([7.9, 4.4, 6.9, 2.5])), 
mean=array([5.84333333, 3.054     , 3.75866667, 1.19866667]), 
variance=array([0.68569351, 0.18800403, 3.11317942, 0.58241432]), 
skewness=array([ 0.31175306,  0.33070281, -0.27171195, -0.10394367]), 
kurtosis=array([-0.57356795,  0.2414433 , -1.3953593 , -1.33524564]))

PCA降維

第一步，計算協方差矩陣，注意計算前需要對資料進行轉置：

cov = np.cov(X.T)
print np.round(cov, decimals=2) # 列印保留2位小數的結果

Out:
[[ 0.69, -0.04,  1.27,  0.52],
 [-0.04,  0.19, -0.32, -0.12],
 [ 1.27, -0.32,  3.11,  1.3 ],
 [ 0.52, -0.12,  1.3 ,  0.58]]

可以看到協方差矩陣是一個對稱矩陣（看到對稱矩陣，立刻就可以想到可對角化了吧，哈哈），對角線表示方差，其它元素表示兩變數間的協方差。

（PS：當資料是中心化（樣本值減去均值）了的，協方差矩陣還可以通過1n−1XTX$\frac{1}{n-1}{X}^{T}X$得到，也就是它的轉置乘以自身再除以樣本數量。）

第二步，特徵值分解（也叫對稱矩陣對角化）：

eig_val,eig_vec = np.linalg.eig(cov)

Out:
array([4.22484077, 0.24224357, 0.07852391, 0.02368303])

array([[ 0.36158968, -0.65653988, -0.58099728,  0.31725455],
       [-0.08226889, -0.72971237,  0.59641809, -0.32409435],
       [ 0.85657211,  0.1757674 ,  0.07252408, -0.47971899],
       [ 0.35884393,  0.07470647,  0.54906091,  0.75112056]])

假設我們把特徵維度從4維降為2維，那我們就選擇特徵值絕對值最大的兩個對應的特徵向量：

indices = np.argsort(np.abs(eig_val))[::-1][:2]

Out:
array([0, 1], dtype=int64)

也就是前兩列特徵向量：

transform_matrix = eig_vec[:, indices]

Out:
[[ 0.36158968, -0.65653988],
 [-0.08226889, -0.72971237],
 [ 0.85657211,  0.1757674 ],
 [ 0.35884393,  0.07470647]]

第三步，將原始資料降維，得到新的資料集：

new_X = np.dot(X, transform_matrix)

Out:
[[ 2.82713597, -5.64133105],
 [ 2.79595248, -5.14516688],
 [ 2.62152356, -5.17737812],
 [ 2.7649059 , -5.00359942],
 [ 2.78275012, -5.64864829],
 [ 3.23144574, -6.06250644],
 [ 2.69045242, -5.23261922],
 [ 2.8848611 , -5.48512908],
 [ 2.62338453, -4.7439257 ],
 [ 2.83749841, -5.20803203]
 ...

這裡需要特別注意的是，這裡的特徵已經不再對應原始資料的特徵{花萼長度、花萼寬度、花瓣長度、花瓣寬度}，經過降維後的資料是在新的座標系下描述的資料。

總結

主成分分析降維方法計算過程分為三步：

1. 計算協方差矩陣
2. 對協方差矩陣進行特徵分解（對角化）
3. 選擇特徵值絕對值最大的特徵值對應的特徵向量作為轉換矩陣，將原始資料降維。