PCAとUMAPの組み合わせ(次元削減)

高次元データを扱う場合、UMAPのみで次元削減するのではなくPCAの結果をさらにUMAPで次元削減することでよりよい結果になることがあります。

PCA実行

前回と同様にMNISTのデータセットを使って、PCAで次元削減してみます。

[Google Colaboratory]

1
2
3
4
5
6
pca = PCA(n_components=0.99, random_state=0)
X_pc = pca.fit_transform(digits.data)
df_pca = pd.DataFrame(X_pc, columns=["PC{}".format(i + 1) for i in range(len(X_pc[0]))])
print("主成分の数: ", pca.n_components_)
print("保たれている情報: ", np.sum(pca.explained_variance_ratio_))
display(df_pca.head())

[実行結果(一部略)]

n_componentsに0.99を指定している(1行目)ので、累積寄与率が99%になるPC41までが結果として表示されました。

もともとは64次元でしたので23次元が削減されたことになります。

UMAPとPCA+UMAPの比較

PCAの結果に対してさらにUMAPを実行します。(2行目)

また比較のためにUMAPのみを実行した場合の結果も合わせて表示します。(1行目)

n_neighborsには5, 10, 15を指定します。

[Google Colaboratory]

1
2
create_2d_umap(digits.data, digits.target, digits.target_names, [5,10,15])
create_2d_umap(df_pca, digits.target, digits.target_names, [5,10,15])

[実行結果]

上図がUMAPのみの結果、下図がPCA+UMAPの結果となります。

明確にどちらの分類がよいか判断が難しいところですが、分類結果が変わっていることは確認できます。

PCAとUMAPを組み合わせた手法があるということも覚えておくと検証の幅が広がるかと思います。

Python × AI - UMAP(最適n_neighborsを探索)

n_neighborsは、UMAPの重要なパラメータです。

n_neighborsを大きくするとマクロな構造を反映し、小さくするとミクロな構造を結果に反映することができます。

デフォルト値は15で、2~100の間の値を選択することが推奨されています。

最適n_neighborsを探索(2次元)

最適なn_neighborsを探索する関数を定義します。

n_neighborsを2, 15, 30, 50, 100と設定してUMAPを実行し、結果を2次元に可視化します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def create_2d_umap(target_X, y, y_labels, n_neighbors_list= [2, 15, 30, 50, 100]):
fig, axes = plt.subplots(nrows=1, ncols=len(n_neighbors_list),figsize=(5*len(n_neighbors_list), 4))
for i, (ax, n_neighbors) in enumerate(zip(axes.flatten(), n_neighbors_list)):
start_time = time.time()
mapper = umap.UMAP(n_components=2, random_state=0, n_neighbors=n_neighbors)
Y = mapper.fit_transform(target_X)
for each_label in y_labels:
c_plot_bool = y == each_label
ax.scatter(Y[c_plot_bool, 0], Y[c_plot_bool, 1], label="{}".format(each_label))
end_time = time.time()
ax.legend(loc="upper right")
ax.set_title("n_neighbors: {}".format(n_neighbors))
print("n_neighbors {} is {:.2f} seconds.".format(n_neighbors, end_time - start_time))
plt.show()

create_2d_umap(digits.data, digits.target, digits.target_names)

[実行結果]

n_neighborsの値によって、結果がかなり変化することが分かります。

上記の2次元グラフからは15, 30あたりが良い結果になっているようです。

最適n_neighborsを探索(3次元)

今度は、3次元で最適なn_neighborsを探索する関数を定義します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
def create_3d_umap(target_X, y, y_labels, n_neighbors_list= [2, 15, 30, 50, 100]):
fig = plt.figure(figsize=(5*len(n_neighbors_list),4))
for i, n_neighbors in enumerate(n_neighbors_list):
ax = fig.add_subplot(1, len(n_neighbors_list), i+1, projection="3d")
start_time = time.time()
mapper = umap.UMAP(n_components=3, random_state=0, n_neighbors=n_neighbors)
Y = mapper.fit_transform(target_X)
for each_label in y_labels:
c_plot_bool = y == each_label
ax.scatter(Y[c_plot_bool, 0], Y[c_plot_bool, 1], label="{}".format(each_label))
end_time = time.time()
ax.legend(loc="upper right")
ax.set_title("n_neighbors_list: {}".format(n_neighbors))
print("n_neighbors_list {} is {:.2f} seconds.".format(n_neighbors, end_time - start_time))
plt.show()

create_3d_umap(digits.data, digits.target, digits.target_names)

projectionパラメータ“3d”を設定(4行目)し、3次元のグラフを表示します。

また、n_componentsパラメータを2から3に変更しています。(6行目)

[実行結果]

上記の3次元グラフからは15, 30あたりでの分類結果がよさそうです。


もう少しn_neighborsを変化させて、再度最適値を調べます。

n_neighborsに[10 , 15, 20, 25, 30]を設定して実行します。

[Google Colaboratory]

1
create_3d_umap(digits.data, digits.target, digits.target_names, [10 , 15, 20, 25, 30])

[実行結果]

n_neighborsが10のときに最もうまく分類されています。

UMAPにおいても、n_neighborsの最適値がいくつかというのはデータセットによって違いますので、このように関数化して設定値を変更させながら比較すると最適な設定値を確認しやすくなります。

Python × AI - UMAP(次元削減)

UMAPは、t-SNEと同様の精度がありながら処理速度も速く、4次元以上の圧縮に対応しているという次元削減のトレンドになりつつある手法です。

ただし、t-SNEと同じようにパラメータの調整は必要です。

非常に高次元で大量のデータについても現実的な時間で実行でき、非線形の高次元データを低次元して可視化する手法として主流になっています。

UMAPライブラリのインストール

まずはUMAPライブラリをインストールします。

[Google Colaboratory]

1
!pip3 install umap-learn

[実行結果]

UMAPを実行

UMAPを実行します。(8行目)

また、比較のためにt-SNEも合わせて実行します。(4行目)

前前回読み込んだMNISTデータを使用しています。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
import umap

start_time_tsne = time.time()
X_reduced = TSNE(n_components=2, random_state=0).fit_transform(digits.data)
interval_tsne = time.time() - start_time_tsne

start_time_umap = time.time()
embedding = umap.UMAP(n_components=2, random_state=0).fit_transform(digits.data)
interval_umap = time.time() - start_time_umap

print("tsne : {}s".format(np.round(interval_tsne,2)))
print("umap : {}s".format(np.round(interval_umap,2)))

[実行結果]

t-SNEよりも、UMAPの方が速く終了していることが確認できます。

可視化

t-SNEの結果とUMAPの結果を可視化します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
plt.figure(figsize=(10,8))
plt.subplot(2, 1, 1)
for each_label in digits.target_names:
c_plot_bool = digits.target == each_label
plt.scatter(X_reduced[c_plot_bool, 0], X_reduced[c_plot_bool, 1], label="{}".format(each_label))
plt.legend(loc="upper right")
plt.xlabel("tsne-1")
plt.ylabel("tsne-2")

plt.subplot(2, 1, 2)
for each_label in digits.target_names:
c_plot_bool = digits.target == each_label
plt.scatter(embedding[c_plot_bool, 0], embedding[c_plot_bool, 1], label="{}".format(each_label))
plt.legend(loc="upper right")
plt.xlabel("umap-1")
plt.ylabel("umap-2")
plt.show()

[実行結果]

上図のt-SNEと比較して、下図のUMAPの方がそれぞれのグループがきれいに小さくまとまっていて、明確に分類されています。

つまり、より適切に情報を落とさず次元削減できているということになります。

UMAPは、どんなデータにも適用でき複雑なデータの関連性が視覚的に分かりやすくなるので、非常に有効な手法です。

次元削減アルゴリズムの選び方

次元削減をする際のアルゴリズムの選び方としては、まずPCAUMAPの2つを試してみてより分類結果が確認しやすい方を選択するのがよいでしょう。

その2つの分類結果がいまいちであれば、t-SNEなど他のアルゴリズムを検討してみましょう。

Python × AI - t-SNE(最適なPerplexity探索)

t-SNEにとって重要なパラメータであるPerplexityの最適値を調べます。

Perplexityとは、どれだけ近傍の点を考慮するかを決めるためのパラメータであり、データの局所的な特性と全体的な特性のどちらをより考慮するかというバランスを表します。

デフォルトは30であり、5から50の間の値を選択することが推奨されています。

複数のPerplexityを設定して結果を確認することが、基本的なアプローチになります。

最適なPerplexityを探索 (2次元)

最適なPerplexityを調べるための関数を定義し、結果を2次元で表示します。

(前回読み込んだMNISTデータを使用しています)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import time
def create_2d_tsne(target_X, y, y_labels, perplexity_list= [2, 5, 30, 50, 100]):
fig, axes = plt.subplots(nrows=1, ncols=len(perplexity_list),figsize=(5*len(perplexity_list), 4))
for i, (ax, perplexity) in enumerate(zip(axes.flatten(), perplexity_list)):
start_time = time.time()
tsne = TSNE(n_components=2, random_state=0, perplexity=perplexity)
Y = tsne.fit_transform(target_X)
for each_label in y_labels:
c_plot_bool = y == each_label
ax.scatter(Y[c_plot_bool, 0], Y[c_plot_bool, 1], label="{}".format(each_label))
end_time = time.time()
ax.legend()
ax.set_title("perplexity: {}".format(perplexity))
print("perplexity {} is {:.2f} seconds.".format(perplexity, end_time - start_time))
plt.show()

create_2d_tsne(digits.data, digits.target, digits.target_names)

定義した関数の引数には次の3パラメータを設定します。

  • 元データ
  • ラベル名
  • ラベル名のユニークリスト

Perplexityにそれぞれ2,5,30,50,100を設定(2行目)し、t-SNEを実行(6,7行目)します。

それぞれの結果を2次元で可視化したものが下記になります。

[実行結果]

2次元ではPerplexityが30,50の時に、うまく分類されていることが分かりました。

t-SNEの結果を3次元表示

今度は3次元で、最適なPerplexityを調べるための関数を定義します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
def create_3d_tsne(target_X, y, y_labels, perplexity_list= [2, 5, 30, 50, 100]):
fig = plt.figure(figsize=(5*len(perplexity_list),4))
for i, perplexity in enumerate(perplexity_list):
ax = fig.add_subplot(1, len(perplexity_list), i+1, projection="3d")
start_time = time.time()
tsne = TSNE(n_components=3, random_state=0, perplexity=perplexity)
Y = tsne.fit_transform(target_X)
for each_label in y_labels:
c_plot_bool = y == each_label
ax.scatter(Y[c_plot_bool, 0], Y[c_plot_bool, 1], label="{}".format(each_label))
end_time = time.time()
ax.legend()
ax.set_title("Perplexity: {}".format(perplexity))
print("perplexity {} is {:.2f} seconds.".format(perplexity, end_time - start_time))
plt.show()

create_3d_tsne(digits.data, digits.target, digits.target_names)

projectionパラメータ“3d”を設定(4行目)し、3次元のグラフを表示します。

また、n_componentsパラメータを2から3に変更しています。(6行目)

[実行結果]

3次元の表示からは、Perplexityが30の方がより特徴ごとに分かれているようです。

今回はPerplexityのデフォルト値である30でうまく分類できることが分かりましたが、データセットによって最適なPerplexityの設定値が異なりますので、今回実施したように複数の結果を作成し比較することをお勧めします。

Python × AI - t-SNE(次元削減)

t-SNE(ティースニー)は、2次元または3次元への圧縮に特化しているアルゴリズムです。

(次元数4以上の場合の結果は保証されていません)

MNISTデータの読み込み

0~9の手書き数字の画像データセットMNIST(エムニスト)を読み込みます。

画像データはベクトル化すると高次元になるため、次元削減アルゴリズムがとても有効です。

[Google Colaboratory]

1
2
3
4
from sklearn.datasets import load_digits
digits = load_digits()
print(digits.data.shape)
print(digits.data)

[実行結果]

1797個の画像データがあり、各データは 8 × 8 = 64 個の数値配列、つまり64次元のデータセットになります。

この64次元のデータセットをt-SNEで2次元に削減し可視化します。

データの先頭から100文字を画像として表示するコードは下記の通りです。

[Google Colaboratory]

1
2
3
4
5
import matplotlib.pyplot as plt 
fig, axes = plt.subplots(10, 10, figsize=(8, 8),subplot_kw={"xticks":[], "yticks":[]})
for i, ax in enumerate(axes.flat):
ax.imshow(digits.images[i], cmap="binary", interpolation="nearest")
ax.text(0, 0, str(digits.target[i]))

[実行結果]

手書き数字の画像データが確認できました。

PCAで次元削減

まずはPCAで次元削減を実行し(2行目)、可視化を行います。

[Google Colaboratory]

1
2
3
4
5
6
7
from sklearn.decomposition import PCA
X_reduced = PCA(n_components=2).fit_transform(digits.data)
for each_label in digits.target_names:
c_plot_bool = digits.target == each_label
plt.scatter(X_reduced[c_plot_bool, 0], X_reduced[c_plot_bool, 1], label="{}".format(each_label))
plt.legend()
plt.show()

[実行結果]

数字ごとに分類されて欲しかったのですが、うまく分類できていないようです。

これはPCAが非線形データに対応できていないためです。

t-SNEで次元削減

次にt-SNEで次元削減を実行し(2行目)、可視化を行います。

[Google Colaboratory]

1
2
3
4
5
6
7
from sklearn.manifold import TSNE
X_reduced = TSNE(n_components=2, random_state=0).fit_transform(digits.data)
for each_label in digits.target_names:
c_plot_bool = digits.target == each_label
plt.scatter(X_reduced[c_plot_bool, 0], X_reduced[c_plot_bool, 1], label="{}".format(each_label))
plt.legend()
plt.show()

[実行結果]

t-SNEは、非線形データに対応できるためうまく分類されています。

1,3,8,9がやや混じっているように見えますが、これは人の目から見ても形が似ている数字があるためしかたがないように思えます。

入力データが複雑になればなるほどPCAはうまく次元削減できなくなることが多いのですが、t-SNEは複雑なデータでもかなり高精度な次元削減を行うことができます。

ただし、次のようなデメリットもあるので注意が必要です。

  • 処理に時間がかかる。
  • 4次元以上には不向き。
  • パラメータ調整が必要になる。

Python × AI - Isomap(非線形データの次元削減)

前回まで実行していたPCAは、データが多次元正規分布に従うことを仮定しているため、非線形データに対してはうまく動作しないという問題があります。

この問題に対処するため、非線形データの変換を行う手法がいくつかあります。

今回はその中でIsomap(Isometric mapping)を試してみます。

Isomapは多様体上の距離を測定し、多次元尺度構成法で表現する手法です。

多次元尺度構成法は近いもの同士は近くに配置し、遠いものは遠くに配置する手法で、Isomapは近いもの同士をより考慮した手法になります。

ムーンデータの取得

まずは非線形データとして、ムーンデータを取得します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn import preprocessing, decomposition, manifold
from sklearn import datasets
from sklearn.decomposition import PCA
X,Y = datasets.make_moons(n_samples=200, noise=0.05, random_state=0)
sc=preprocessing.StandardScaler()
sc.fit(X)
X_norm=sc.transform(X)
plt.figure(figsize=(10,3))
plt.scatter(X[:,0],X[:,1], c=Y)
plt.xlabel("x")
plt.ylabel("y")

[実行結果]

Isomap実行

取得したムーンデータに対してIsomapを実行します。(4~8行目)

比較のためPCAも実行しています。(1~2行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
pca = PCA(n_components=2)
X_reduced = pca.fit_transform(X_norm)

isomap_5 = manifold.Isomap(n_neighbors=5, n_components=2)
X_isomap_5 = isomap_5.fit_transform(X_norm)

isomap_10 = manifold.Isomap(n_neighbors=10, n_components=2)
X_isomap_10 = isomap_10.fit_transform(X_norm)

Isomapは全てのサンプル間の距離を計算しますが、細かく計算するのはサンプルごとに最も近いN個のサンプル間の距離だけです。

可視化した結果が良くない場合は、Nを変更して再び試行することになります。

このNに対応するのがn_neighborsパラメータになります。

4行目ではn_neighborsを5に、7行目ではn_neighborsを10に設定しています。(上記のソースコード)

Isomap実行結果の可視化

順番に、PCAの結果、Isomapでn_neighborsが5の結果、Isomapでn_neighborsが10の結果を可視化します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
plt.figure(figsize=(10,6))
plt.subplot(3, 1, 1)
plt.scatter(X_reduced[:, 0], X_reduced[:, 1], c=Y)
plt.xlabel("pca-1")
plt.ylabel("pca-2")

plt.subplot(3, 1, 2)
plt.scatter(X_isomap_5[:,0],X_isomap_5[:,1], c=Y)
plt.xlabel("isomap_n5-1")
plt.ylabel("isomap_n5-2")

plt.subplot(3, 1, 3)
plt.scatter(X_isomap_10[:,0],X_isomap_10[:,1], c=Y)
plt.xlabel("isomap_n10-1")
plt.ylabel("isomap_n10-2")
plt.show

[実行結果]

一番下のIsomapでn_neighborsが10の結果が、一番きれいに分かれていて新たな軸が作成されていることが確認できました。

Python × AI - 寄与率(次元削減数を探索)

寄与率(explained_variance_ratio_)をもとに有効なPC数を探索します。

(前回記事で取得したワインデータと主成分分析結果を使います。)

寄与率を確認

寄与率は、主成分がどの程度、元データの情報を保持しているかを表します。

各固有値を固有値で割ったものが寄与率になります。

固有値と同じくPC1がもっとも大きく、次第に小さくなります。

寄与率を表示するソースコードは下記のようになります。

[Google Colaboratory]

1
pd.DataFrame(np.round(pca.explained_variance_ratio_,2), index=["PC{}".format(x + 1) for x in range(len(df_pca.columns))], columns=["寄与率"])

[実行結果]

PC1が36%、PC2が19%となっておりこの2つを合わせた寄与率は55%で、PC1とPC2だけで半分以上説明できるということになります。

累積寄与率を可視化

縦軸が累積寄与率で、横軸がPCのグラフを表示します。(4行目)

また累積寄与率が90%になるまでの主成分を判断するために、90%の基準線も表示します。(2,8行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
import matplotlib.ticker as ticker
line = np.full(14, 0.9)
plt.gca().get_xaxis().set_major_locator(ticker.MaxNLocator(integer=True))
plt.plot([0] + list( np.cumsum(pca.explained_variance_ratio_)), "-o")
plt.xlabel("PC")
plt.ylabel("cumulative contribution rate")
plt.yticks( np.arange(0, 1.1, 0.1))
plt.plot(line, "s-")
plt.grid()
plt.show()

[実行結果]

PC8のところで90%の基準線を超えていますので、PC8までを有効な次元数だと判断することができます。

より次元を削減したい場合は70%を基準することもあり、その場合はPC4までが有効な次元数となります。

累積寄与率の指定

PCAのパラメータであるn_components0~1を指定すると累積寄与率が設定でき、指定した累積寄与率を超過するまでの主成分を返してくれます。

n_components0.9(90%)を指定して実行します。(5行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
sc = preprocessing.StandardScaler()
X = df_wine.iloc[:, 1:]
X_norm=sc.fit_transform(X)

pca = PCA(n_components=0.9, random_state=0)
X_pc = pca.fit_transform(X_norm)
df_pca = pd.DataFrame(X_pc, columns=["PC{}".format(i + 1) for i in range(len(X_pc[0]))])
print("主成分の数: ", pca.n_components_)
print("保たれている情報: ", round(np.sum(pca.explained_variance_ratio_),2))
display(df_pca.head())

[実行結果]

累積寄与率が90%を超過するPC8までが結果として表示されました。

累積寄与率は92%となっています。

寄与率をあらかじめ決めている場合や、大まかに有効な次元数を確認したい場合は便利な手法となります。

PCAのデータ標準化

PCAにおいてデータの標準化は基本的に実施することが推奨されていますが、ノイズデータが多い場合正しく軸をとれないことがあります。

そのため、標準化するパターン標準化しないパターンの両方を実施して結果の良い方を使用するこをお勧めします。

Python × AI - スクリープロット(次元削減数を探索)

PCAは、教師あり学習の説明変数としての用途もあります。

その際には、次元を減らしつつも元の情報をなるべく保持しているPCまで使う必要があります。

今回は、有効なPC数を探索していきます。

ワインデータの取得

アイリスデータより次元数の多いワインデータを使用します。

ワインデータは列名が入っていないので、列名を定義しています。(2行目)

[Google Colaboratory]

1
2
3
4
df_wine=pd.read_csv("https://archive.ics.uci.edu/ml/machine-learning-databases/wine/wine.data", header=None)
df_wine.columns = ["class", "Alcohol", "Malic acid", "Ash", "Alcalinity of ash","Magnesium", "Total phenols", "Flavanoids", "Nonflavanoid phenols","Proanthocyanins", "Color intensity", "Hue","OD280/OD315 of diluted wines", "Proline"]
display(df_wine.shape)
display(df_wine.head())

[実行結果]

正解データが入っているため、14次元サンプル数が178件のデータになります。

ワインデータの可視化

取得したワインデータを散布図行列で可視化します。

[Google Colaboratory]

1
2
3
from pandas import plotting
plotting.scatter_matrix(df_wine.iloc[:, 1:], figsize=(8, 8), c=list(df_wine.iloc[:, 0]), alpha=0.5)
plt.show()

[実行結果]

次元数が多いため、とても判断が難しくなっています。

だいたい3~4種類に分かれているような気がしなくもありません。

PCA

PCAを実行します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
from sklearn.decomposition import PCA
from sklearn import preprocessing
import numpy as np
sc=preprocessing.StandardScaler()
X = df_wine.iloc[:, 1:]
X_norm=sc.fit_transform(X)

pca = PCA(random_state=0)
X_pc = pca.fit_transform(X_norm)
df_pca = pd.DataFrame(X_pc, columns=["PC{}".format(i + 1) for i in range(len(X_pc[0]))])
print("主成分の数: ", pca.n_components_)
print("保たれている情報: ", round(np.sum(pca.explained_variance_ratio_),2))
display(df_pca.head())

[実行結果]

主成分の数がPC1からPC13の13個あります。

パラメータのn_componentsを指定しない場合、主成分は元データの次元数と一致します。

全次元を対象にしたので保たれている情報は1.0(100%)になります。

固有値

固有値を確認して、PC1からどこまでが有効な主成分かを探索してみましょう。

[Google Colaboratory]

1
pd.DataFrame(np.round(pca.explained_variance_, 2), index=["PC{}".format(x + 1) for x in range(len(df_pca.columns))], columns=["固有値"])

[実行結果]

固有値(explained_variance)は、主成分の分散のことで主成分の情報量の大きさを表します。

PC1がもっとも大きく、次第に小さくなっていきます。

標準化している場合には、固有値が1.0以上のものを使うというのがもっともシンプルな判断基準となります。

今回の結果ですと、PC3までを使うことになります。

スクリープロットで可視化

次にスクリープロットを確認します。

スクリープロットは、固有値を最大から最小まで降順でプロットしたグラフのことで、崖(スクリー:scree)のような形になります。

固有値からスクリープロットを作成(2行目)し、グラフ表示するソースは以下の通りです。

固有値1.0の基準線も合わせて表示しています。(1,3行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
line = np.ones(14)
plt.plot(np.append(np.nan, pca.explained_variance_), "s-")
plt.plot(line, "s-")
plt.xlabel("PC")
plt.ylabel("explained_variance")
plt.xticks( np.arange(1, 14, 1))
plt.grid()
plt.show()

[実行結果]

固有値が、ある段階から急に小さな値となって以降は安定する箇所までを利用することが望ましいとされています。

上記のグラフだとPC7あたりかと思いますが、使用用途に応じて±2程度で調整する必要があります。

Python × AI - 主成分と元データの相関(PCA)

前回作成した主成分を解析します。

PCAの結果

PCAの結果から取得できる値は次の通りです。

  • 固有ベクトル : components_
    PCAでデータのばらつきが大きい方向に軸を取り直した結果のベクトル。
    各主成分と元データとの相関関係(-1~1)を意味し、元のデータと主成分の影響度合いを表す。
  • 主成分得点 : explained_variance
    固有ベクトルと元データをかけ合わせた値。
     主成分得点 = 元データ × 固有ベクトル
  • 固有値 : explained_variance_
    固有ベクトルの方向に沿ったデータの分散の大きさ。
    固有値が大きい固有ベクトルほど、データの分散をよく説明しており、データの重要な特徴を捉えている。
    データを標準化している場合、各PCは1以上あれば元データより情報をもっていることになり、4次元データであれば全ての合計は4になる。
  • 寄与率 : explained_variance_ratio_
    固有値から算出した、データ特徴の捉え度合い。
     寄与率 = 固有値 ÷ 固有値の合計

相関図

元データとの関係を固有ベクトル(components)から確認します。

固有ベクトルのヒートマップを表示します。(4~12行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
import seaborn as sns
fig = plt.figure(figsize=(8,6))
ax = fig.add_subplot(111)
sns.heatmap(pca.components_,
cmap="Blues",
annot=True,
annot_kws={"size": 14},
fmt=".2f",
xticklabels=["SepalLength", "SepalWidth", "PetalLength", "PetalLength"],
yticklabels=["PC1", "PC2", "PC3", "PC4"],
ax=ax)
plt.show()

[実行結果]

上図は、各主成分(縦軸)と元データ(横軸)との相関関係を表しています。

PC1をみると、SepalLength、PetalLength、PetalWidthの3つが平均として大きな値となっています。

PC2については、SepalWidthが強く影響しています。

このように主成分がどんな内容の軸になったのかは、PCAの結果を見て判断する必要があります。

Python × AI - 主成分分析(PCA)

主成分分析(PCA:Principal Component Analysis)は、多次元データのもつ情報をできるだけ損なわずに低次元とする方法です。

次元削減で最も簡単な方法であり、広い分野で使われています。

アイリスデータの読み込み

まずアイリスデータを読み込みます。

結果の確認用にtarget_nameに正解の花の名称を追加しています。(5~8行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
import pandas as pd
from sklearn.datasets import load_iris
iris = load_iris()
df=pd.DataFrame(iris.data, columns=iris.feature_names)
df["target"] = iris.target
df.loc[df["target"]==0, "target_name"] = "setosa"
df.loc[df["target"]==1, "target_name"] = "versicolor"
df.loc[df["target"]==2, "target_name"] = "virginica"
df.head()

[実行結果]

アイリスデータの散布図行列

読み込んだアイリスデータを散布図行列で可視化します。

hueパラメータに名称”target_name”を設定することで色をつけています。(2行目)

[Google Colaboratory]

1
2
import seaborn as sns
sns.pairplot(df, vars=df.columns[:4], hue="target_name")

[実行結果]

品種ごとに色分けされていて分類されていることは把握できますが、次元が多いため解釈するのが大変です。

アイリスの3次元立体図

花の品種ごとにpetal_width以外の情報を取得します。(5~6行目)

その後、3次元での可視化を行います。(7~12行目)

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
fig = plt.figure(figsize=(8,6))
ax = fig.add_subplot(1, 1, 1, projection="3d")
for c in df["target_name"].unique():
ax.scatter(df.iloc[:, 0][df["target_name"]==c], df.iloc[:, 1][df["target_name"]==c] , df.iloc[:, 2][df["target_name"]==c], label=c)
ax.set_title("iris 3D")
ax.set_xlabel("sepal_length")
ax.set_ylabel("sepal_width")
ax.set_zlabel("petal_length")
ax.legend(loc=2, title="legend", shadow=True)
plt.show()

[実行結果]

3品種に分かれているように見えますが、petal_widthの情報は全く無視してしまっています。

この情報が特徴量のある重要なデータだったという可能性もあります。

PCAを使用すれば、この4次元データを2次元データで表現することが可能になります。

PCA実行

PCAを実行し(3~4行目)、結果である主成分得点を表示します(5~8行目)。

[Google Colaboratory]

1
2
3
4
5
6
7
8
from sklearn.decomposition import PCA
import numpy as np
pca = PCA(random_state=0)
X_pc = pca.fit_transform(df.iloc[:, 0:4])
df_pca = pd.DataFrame(X_pc, columns=["PC{}".format(i + 1) for i in range(len(X_pc[0]))])
print("主成分の数: ", pca.n_components_)
print("保たれている情報: ", np.sum(pca.explained_variance_ratio_))
display(df_pca.head())

[実行結果]

横軸が主成分(PC:Principal Component)で、縦軸が各サンプルを表します。

主成分(PC)とは、データを要約(縮小)したあとの新しい合成変数で、第1主成分(PC1)に最も多くの情報が集まっていて第2主成分(PC2)以降にだんだんと情報が小さくなります。

PC1とPC2を可視化

多くの情報が集まっているPC1PC2を可視化します。

[Google Colaboratory]

1
sns.scatterplot(x="PC1", y="PC2", data=df_pca, hue=df["target_name"])

[実行結果]

4次元あったデータを2次元で可視化することができました。

うまく3種類に分類されています。

PCAの用途はいくつかあるのですが、多次元データの特徴を低次元に次元削減し可視化する手段としてとても有効です。