【Python】主成分分析（PCA)をnumpyベースで理解する

主成分分析（Principal Component Analysis, PCA）は、多次元データを圧縮して次元削減することを目的とする統計的手法です。

 PCAは、高次元データの情報を、少数の主成分と呼ばれる新しい特徴量に圧縮することができます。これにより、多次元データをより単純な形式で表すことができます。

 scikit-learnのライブラリーを使えば1行で算出することもできますが、ここでは、numpyを使った実装を通じて、一つずつ計算の過程をみていきます。

 全体感

以下が、numpyを使って主成分分析を行う関数になります。

 ひとつずつ何を計算しているのかをみていきます。

def pca(X, n_components=2):

    ①　X = X- X.mean(axis=0)

    ②   cov = np.cov(X, rowvar=False)

    ③   l, v = np.linalg.eig(cov)

    ④    l_index = np.argsort(l)[::-1]

    ⑤   v_ = v[:, l_index]

    ⑥   components = v_[:,:n_components]

    ⑦  T = np.dot(X, components)

    return T

Xの準備

100×5次元の行列をXとして生成しておきます。

 主成分分析は次元圧縮です。5次元を2次元に圧縮していくのが目的になりま
す。

X = np.random.rand(100,5)

X[:3]

[out]

array([[0.95927047, 0.32503164, 0.16065997, 0.36034609, 0.46927506],

       [0.86403093, 0.02959833, 0.13703161, 0.3636126 , 0.9330399 ],

       [0.65145658, 0.1835818 , 0.91643428, 0.24846823, 0.61529736]])

❶各要素と平均との差を算出

Xの各要素と各次元の平均との差を算出します。

X = X-X.mean(axis=0)

X[:3]

[out]

array([[ 0.40344141, -0.1789226 , -0.35687984, -0.15906571, -0.01172766],

       [ 0.30820188, -0.47435591, -0.38050819, -0.1557992 ,  0.45203718],

       [ 0.09562753, -0.32037244,  0.39889447, -0.27094357,  0.13429464]])

❷共分散行列を生成

共分散行列は、多変量データのk個の変数に対して、各変数の組み合わせの共分散を記述する行列です。

 共分散は、2つの変数がどの程度に同時に変動するかを示します。

 共分散は正の数であれば、2つの変数は正の相関を持ちます。共分散は負の数であれば、2つの変数は負の相関を持ちます。共分散は0であれば、2つの変数は独立であることを示します。

次元数×次元数となります

cov =np.cov(X,rowvar=False)

cov

[out]

array([[ 0.08999058, -0.00547484, -0.0088472 ,  0.00446651,  0.01157695],

       [-0.00547484,  0.09707911, -0.00421609,  0.0080685 , -0.00265034],

       [-0.0088472 , -0.00421609,  0.09334372, -0.00715198,  0.00616764],

       [ 0.00446651,  0.0080685 , -0.00715198,  0.07963141,  0.00369422],

       [ 0.01157695, -0.00265034,  0.00616764,  0.00369422,  0.07095071]])

❸固有値・固有ベクトルを計算

固有値は、行列Aに対する線形変換の「大きさ」を表します。

 固有ベクトルは、行列Aに対する線形変換の「方向」を表します。

 固有値と固有ベクトルは、行列の「特徴」を表すものであり、行列Aに対して、固有値と固有ベクトルを求めることにより、行列Aを解析することができます。

#lに固有値、vに固有ベクトルが入ります

l, v = np.linalg.eig(cov)

l

[out]

array([0.06189794, 0.07396456, 0.08742091, 0.10242342, 0.1052887 ])

v

[out]

array([[-0.40975186, -0.36982006, -0.36271054,  0.69353739,  0.2877212 ],

       [ 0.0064176 , -0.36978536, -0.39478444, -0.63069096,  0.55640971],

       [-0.32303312,  0.04932353, -0.6348379 , -0.22094408, -0.66436537],

       [-0.20261986,  0.84961928, -0.26582834,  0.02539724,  0.40716403],

       [ 0.82864771,  0.04696996, -0.48877649,  0.26790575, -0.02146753]])

`

 ❹❺固有値の大きい順に固有ベクトルを並び替え

argsortを使って、固有値のインデックスを降順にインデックス順に並び替えます・

 [::-1]で降順になります

l_index = np.argsort(l)[::-1]

l_index

[out]

array([4, 3, 2, 1, 0])

v[:,l_index]


[out]

array([[ 0.2877212 ,  0.69353739, -0.36271054, -0.36982006, -0.40975186],

       [ 0.55640971, -0.63069096, -0.39478444, -0.36978536,  0.0064176 ],

       [-0.66436537, -0.22094408, -0.6348379 ,  0.04932353, -0.32303312],

       [ 0.40716403,  0.02539724, -0.26582834,  0.84961928, -0.20261986],

       [-0.02146753,  0.26790575, -0.48877649,  0.04696996,  0.82864771]])

❻n_components個の固有ベクトルを取得

列の次元で最初から、n_components個を取り出す

n_components =2

components = v_[:,:n_components]

components

#5×2の固有ベクトルが取得できる

[out]

array([[ 0.2877212 ,  0.69353739],

       [ 0.55640971, -0.63069096],

       [-0.66436537, -0.22094408],

       [ 0.40716403,  0.02539724],

       [-0.02146753,  0.26790575]])

❼Xに固有ベクトルを掛け合わせることで、Xを2次元に圧縮する

もともとXは100×5の行列でしたが、この操作により、100×2の行列になりま
す。

T = np.dot(X,components)

T[:3]

#3行だけ抽出します。実際は100行あります。

[out]

array([[ 0.18910891,  0.46431532],

       [ 0.0043965 ,  0.71413904],

       [-0.52895739,  0.20934099]])