主成分分析(PCA)による教師なし次元削減

  • 分散が最大となる方向を見つけ出し、元の次元と同じかそれよりも低い次元の新しい部分空間へ射影する。
  • 部分空間の直交軸(主成分)を分散が最大となる方向とみなすことができる。
  • PCAの方向はデータのスケーリングに対して敏感なため、特徴量の標準化が必要な場合がある。

PCAのステップ

  1. $ d $次元のデータセットを標準化
  2. 標準化したデータセットの共分散行列を作成する
  3. 共分散行列を固有ベクトルと固有値に分解する
  4. 最も大きい$ k $個の固有値に対応する$ k $個の固有ベクトルを選択する
  5. 上位$ k $個の固有ベクトルから射影行列$ {\boldsymbol W} $を作成する
  6. 射影行列$ {\boldsymbol W} $を使って$ d $次元の入力データセット$ {\boldsymbol X} $を作成し、新しい$ k $次元の特徴部分空間を取得する

以下ではまずscikit-learnの機能を使わずに主成分分析を行う。


In [2]:
# データの取得と前処理(標準化)
import pandas as pd
df_wine = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/wine/wine.data', header=None)

from sklearn.cross_validation import train_test_split
from sklearn.preprocessing import StandardScaler

X, y = df_wine.iloc[:, 1:].values, df_wine.iloc[:, 0].values
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)
sc = StandardScaler()
X_train_std = sc.fit_transform(X_train)
X_test_std = sc.transform(X_test)


C:\Users\tono\Anaconda3\lib\site-packages\sklearn\cross_validation.py:44: DeprecationWarning: This module was deprecated in version 0.18 in favor of the model_selection module into which all the refactored classes and functions are moved. Also note that the interface of the new CV iterators are different from that of this module. This module will be removed in 0.20.
  "This module will be removed in 0.20.", DeprecationWarning)

共分散行列の固有対を求める

特徴量のペアごとの共分散を共分散行列として保持する。

共分散行列を$ \Sigma $で表す。固有ベクトル$ {\boldsymbol v} $、固有値$ \lambda $との関係は以下の通り

$$ \Sigma{\boldsymbol v} = \lambda{\boldsymbol v} $$

下記np.linalg.eigでWineデータセットの共分散行列の固有対を取得する。


In [3]:
import numpy as np
cov_mat = np.cov(X_train_std.T)
eigen_vals, eigen_vecs = np.linalg.eig(cov_mat)
print(eigen_vals)
print(eigen_vecs)


[ 4.8923083   2.46635032  1.42809973  1.01233462  0.84906459  0.60181514
  0.52251546  0.08414846  0.33051429  0.29595018  0.16831254  0.21432212
  0.2399553 ]
[[  1.46698114e-01  -5.04170789e-01   1.17235150e-01  -2.06254611e-01
    1.87815947e-01   1.48851318e-01  -1.79263662e-01  -5.54687162e-02
   -4.03054922e-01  -4.17197583e-01   2.75660860e-01   4.03567189e-01
    4.13320786e-04]
 [ -2.42245536e-01  -2.42168894e-01  -1.49946576e-01  -1.30489298e-01
   -5.68639776e-01   2.69052764e-01  -5.92636731e-01   3.32731614e-02
   -1.01833706e-01   2.17101488e-01  -8.13845005e-02  -1.52474999e-01
   -8.78560762e-02]
 [ -2.99344215e-02  -2.86984836e-01  -6.56394387e-01  -1.51536318e-02
    2.99209426e-01   9.33386061e-02   6.07334578e-02  -1.00618575e-01
    3.51841423e-01   1.28549846e-01  -1.29751275e-02   1.68376064e-01
   -4.52518598e-01]
 [ -2.55190023e-01   6.46871827e-02  -5.84282337e-01   9.04220851e-02
    4.12499478e-02   1.01342392e-01   2.50323869e-01   5.61658566e-02
   -5.00457282e-01   4.73344124e-02   9.89088030e-02  -6.70902927e-02
    4.86169765e-01]
 [  1.20797723e-01  -2.29953850e-01  -8.22627466e-02   8.39128346e-01
    2.71971315e-02  -1.12567350e-01  -2.85240559e-01   9.58423947e-02
    8.37391743e-02  -2.78918776e-01  -9.59297663e-02  -1.02396856e-01
    1.14764951e-01]
 [  3.89344551e-01  -9.36399132e-02  -1.80804417e-01  -1.93179478e-01
   -1.40645426e-01  -1.22248798e-02   5.31455344e-02  -4.21265116e-01
    1.35111456e-01  -2.80985650e-01   2.83897644e-01  -6.18600153e-01
    9.45645138e-02]
 [  4.23264856e-01  -1.08862204e-02  -1.42959330e-01  -1.40459548e-01
   -9.26866486e-02   5.50345182e-02   7.98994076e-02   8.47224703e-01
    3.36016514e-03  -3.91442963e-02   1.16729207e-01  -1.39680277e-01
   -1.00444099e-01]
 [ -3.06349555e-01  -1.87021637e-02  -1.72234753e-01  -3.37332618e-01
    8.58416771e-02  -6.95340883e-01  -2.97371718e-01   1.66256803e-01
    1.90120758e-01  -2.78622194e-01  -3.96566280e-02   1.63323514e-03
    2.00128778e-01]
 [  3.05722194e-01  -3.04035180e-02  -1.58362102e-01   1.14752900e-01
   -5.65105241e-01  -4.98354410e-01   2.02519133e-01  -1.66197468e-01
   -1.76029939e-01   1.48539457e-01   8.60602743e-02   3.88568490e-01
   -1.39942067e-01]
 [ -9.86919131e-02  -5.45270809e-01   1.42421708e-01  -7.87857057e-02
   -1.32346052e-02  -1.59452160e-01   3.97364107e-01   3.96173606e-02
   -2.14930670e-01  -4.10240865e-03  -5.71651893e-01  -3.08345904e-01
   -1.15349466e-01]
 [  3.00325353e-01   2.79243218e-01  -9.32387182e-02  -2.41740256e-02
    3.72610811e-01  -2.16515349e-01  -3.84654748e-01  -1.05383688e-01
   -5.17259438e-01   1.97814118e-01  -1.98844532e-01  -2.00456386e-01
   -3.02254353e-01]
 [  3.68211538e-01   1.74365000e-01  -1.96077407e-01  -1.84028641e-01
   -8.93796748e-02   2.35172361e-01  -8.62903341e-02  -9.95055559e-02
    1.36456039e-01  -2.38138151e-01  -6.50869713e-01   2.84100327e-01
    3.18414303e-01]
 [  2.92597130e-01  -3.63154608e-01   9.73171134e-02  -5.67677845e-02
    2.17529485e-01  -1.05621383e-01  -1.30298291e-01  -1.60661776e-02
    1.67758429e-01   6.37350206e-01   7.12377082e-02   3.75546771e-02
    5.03247839e-01]]
  • 上記固有ベクトルから、データに含まれる大半の情報を含んでいるもの(主成分)を選択する。
  • 固有値が固有ベクトルの大きさを表しているため、固有値の大きい順に並べ替えて判定する。

In [4]:
tot = sum(eigen_vals)
var_exp = [(i/tot) for i in sorted(eigen_vals, reverse=True)]
cum_var_exp = np.cumsum(var_exp)
import matplotlib.pyplot as plt
plt.bar(range(1,14), var_exp, alpha=0.5, align='center', label='individual exp var')
plt.step(range(1,14), cum_var_exp, where='mid', label='cum exp var')
plt.ylabel('ex')
plt.xlabel('prin')
plt.legend(loc='best')
plt.show()


特徴変換

  • Wineデータセットを新しい主成分軸に変換する

In [9]:
# 固有値の大きいものから固有対を並び替える
eigen_pairs = [(np.abs(eigen_vals[i]), eigen_vecs[:,i]) for i in range(len(eigen_vals))]
eigen_pairs.sort(key=lambda k: k[0], reverse=True)
# 最も大きい2つの固有値に対応する固有ベクトルを集める
w = np.hstack((eigen_pairs[0][1][:, np.newaxis], eigen_pairs[1][1][:, np.newaxis]))
print('Matrix W:\n',w)


Matrix W:
 [[ 0.14669811 -0.50417079]
 [-0.24224554 -0.24216889]
 [-0.02993442 -0.28698484]
 [-0.25519002  0.06468718]
 [ 0.12079772 -0.22995385]
 [ 0.38934455 -0.09363991]
 [ 0.42326486 -0.01088622]
 [-0.30634956 -0.01870216]
 [ 0.30572219 -0.03040352]
 [-0.09869191 -0.54527081]
 [ 0.30032535  0.27924322]
 [ 0.36821154  0.174365  ]
 [ 0.29259713 -0.36315461]]

In [13]:
# 124 x 13 の特徴量を124 x 2 に変換
X_train_pca = X_train_std.dot(w)
print(X_train_std.shape)
print(w.shape)
print(X_train_pca.shape)


(124, 13)
(13, 2)
(124, 2)

In [16]:
colors = ['r', 'b', 'g']
markers = ['s', 'x', 'o']

for l, c, m in zip(np.unique(y_train), colors, markers):
    plt.scatter(X_train_pca[y_train==l, 0], X_train_pca[y_train==l, 1], c=c, label=l, marker=m)
plt.xlabel('PC 1')
plt.xlabel('PC 2')
plt.legend(loc='lower left')
plt.show()



In [17]:
# scikit-learnのPCAクラスを使用した場合
from matplotlib.colors import ListedColormap

def plot_decision_regions(X, y, classifier, resolution=0.02):

    # setup marker generator and color map
    markers = ('s', 'x', 'o', '^', 'v')
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')
    cmap = ListedColormap(colors[:len(np.unique(y))])

    # plot the decision surface
    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution),
                           np.arange(x2_min, x2_max, resolution))
    Z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)
    Z = Z.reshape(xx1.shape)
    plt.contourf(xx1, xx2, Z, alpha=0.4, cmap=cmap)
    plt.xlim(xx1.min(), xx1.max())
    plt.ylim(xx2.min(), xx2.max())

    # plot class samples
    for idx, cl in enumerate(np.unique(y)):
        plt.scatter(x=X[y == cl, 0], y=X[y == cl, 1],
                    alpha=0.8, c=cmap(idx),
                    edgecolor='black',
                    marker=markers[idx], 
                    label=cl)

In [19]:
from sklearn.linear_model import LogisticRegression
from sklearn.decomposition import PCA

# 主成分数2で初期化
pca = PCA(n_components=2)
lr=LogisticRegression()
X_train_pca = pca.fit_transform(X_train_std)
X_test_pca = pca.transform(X_test_std)
lr.fit(X_train_pca, y_train)

plot_decision_regions(X_train_pca, y_train, classifier=lr)
plt.xlabel('PC1')
plt.ylabel('PC2')
plt.legend(loc='lower left')
plt.show()



In [ ]: