データサイエンスの基本：Sklearnで主成分分析

主成分分析(Principal component analysis PCA)

主成分分析(Principal component analysis PCA)は最小の一つのデータ次元減少のアルゴリズムです。発見者の KARL PEARSONが1901年、パソコン時代の前に発表した論文でこの手法を発表し、現在それは一般的なデータサイエンスの手法になりました　（Pearson 1901）。そして、無意識に私たちの脳の中でも、同じようなことをしています。

例：「散歩中に猫を見つけると、その猫は屋根の上に飛び上がって寝そべりました。 貴方はそれが可愛いと思います」

脳では反射された光を読み込んで画像に整えて、その画像で猫という動物が認識されました。それから画像の変化を猫の行動として認識して、X、Y、Z軸に移動を反映して物体を追跡しました。そしてこの画像の情報は「可愛い」という概念を表します自分の脳の中で、情報を吸収し、整え、大事なことを取得しました。そのようなことは数学でPCAが行います。

PCAの場合は脳のように無意識でなく、データがどの方向に最大分散しているか、で見つけます。そして、新しい軸に集中します。その結果、分散が小さい方向を縮めて、重要な情報だけが残ります。このポストでは、PythonでSKlearnというモジュールを使用し、PCAを行います

この投稿のコードはすべて、gitHubにアップし、GoogleColabでノートブックがあり、自由に使用してください

Sklearn PCA

SklearnではPCAはクラスの関数に存在し、行うまえに機械学習のモデルのように定義しなければなりません

sklearn.decomposition.PCA(n_components=None, copy=True, whiten=False, svd_solver='auto', tol=0.0, iterated_power='auto', random_state=None)

from sklearn.decomposition import PCA
n_components=2
pca = PCA(n_components=n_components, random_state=42)

定義した後で、データを入り込めます。しかし、このデータは正規化しなければならないです。ない場合は結果が可笑しくなります。今回のデータはfetch_olivetti_facesという顔のデータセットです。

pca.fit(X) , Where X must be normalized in a shape of  (n_samples, n_features) 

from sklearn.datasets import fetch_olivetti_faces
data=data_dic["images"] 
ids=data_dic["target"]
def normalize(X):
 mu=np.mean(X, axis=0) #Mean 
 std=np.std(X,axis=0) #Standard deviataion 
 #Deal with std vales of zero
 std_filled = std.copy()
 std_filled[std==0] = 1. 
 Xbar = ((X-mu)/std_filled)  #Normalized so mean is zero 
 return Xbar,mu,std

dim=data.shape #Data shape 
Xbar,mu,std= normalize(data/255)  #Normailzed data 
X=Xbar.reshape(-1,dim[1]*dim[2]) 
print(X.shape)  #Fit data n_samples (Imanges), n_features(pixels)
pca.fit(X)

PCAのモデルにデータを入り込んでから、PCAの主成分を抜け出れます。下記のコードでは、主成分を抽出して、結果をプロットします。

pca_components=pca.components_.reshape(2,dim[1],dim[2])

fig,axs= plt.subplots(1,len(pca_components))

for i,ax in enumerate(axs):
 ax.imshow(pca_components[i])
 ax.set_title(f"PCA component {i+1}")

二つの主成分だけですが、顔の構造を見えます。しかし、ちょっと怖いですね。幽霊みたいです。次に、データは主成分の軸に映しましょう。それのために、pca.transform(x)という変化する関数を使用します。

pca_features=pca.transform(X)
print(pca_features.shape) #(400,2)
df=pd.DataFrame(pca_features, columns=['PC1', 'PC2'])
df["ids"]=ids
df.plot.scatter(x="PC1",y="PC2",c="ids",cmap="jet",figsize=(10,10))

このプロットは、それぞれの画像の固有値を示して、もっとお洒落になるために、色によって画像のIDを表していました。その結果は、もともとの64x64のピクセルのデータが1x1の固有値のデータに減少しました。しかし、減少度は何ですか？　pca.explained_variance_ratio_ （分散の割合）を使用し、減少度を表示できます

print("pca.explained_variance_ratio_ :",pca.explained_variance_ratio_)
print(f"information retained {round(np.sum(pca.explained_variance_ratio_)*100,2)} %")
#pca.explained_variance_ratio_ : [0.26876438 0.12342374]
#information retained 39.22 %

分散の割合は、主成分ともともとのデータにより分散を上げます。今回の場合は、二つの主成分があり、それを加えば、もっと画像から残っている情報の割合を表せます。そして、SklearnのPCAでは、もし残っている情報の割合の主成分が欲しい場合は、PCAのn_componentsでは、整数でなく、少数に記入すれば、これをできます。

pca2 = PCA(n_components=0.66)
pca2.fit(X)
print('infromation retained: ', np.sum(pca2.explained_variance_ratio_))
print('n_compoments: ', pca2.n_components_) 
#--------------------------------
#infromation retained:  0.6605544
#n_compoments:  10

66％の情報の場合は、十の主成分があり、この主成分を表示しましょう。

pca2_components=pca2.components_.reshape(pca2.n_components_,dim[1],dim[2])

fig,axs= plt.subplots(2,len(pca2_components)//2,figsize=(12,12))
axs=axs.reshape(-1)
for i,ax in enumerate(axs):
 ax.imshow(pca2_components[i])
 ax.set_title(f"PCA component {i+1}")
fig.tight_layout(h_pad=-25)

これはまた怖いですね。しかし、今回は、もっと顔の特徴を見えます。鼻、眼鏡、髭などをちょっと見えますね。それは、情報量が上げたことのおかげです。PCAのよりデータ量を減少す手法は、激しいパソコンの計算にはとても役に立ちます。例えば、機械学習アルゴリズムには、データ量を減れば、あまり精度を落ちらずに、学習のスピードを上げれます。

すべてのコード

この投稿のコードはすべて、gitHubにアップし、GoogleColabでノートブックがあり、自由に使用してください

import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
%matplotlib inline
import random
import pandas as pd
from sklearn.datasets import fetch_olivetti_faces
#=========================#
# Load the faces datasets
data_dic= fetch_olivetti_faces()
data=data_dic["images"] 
ids=data_dic["target"]

#=========================#
from sklearn.decomposition import PCA
n_components=2
pca = PCA(n_components=n_components, random_state=42)


#=========================#

def normalize(X):
 mu=np.mean(X, axis=0) #Mean 
 std=np.std(X,axis=0) #Standard deviataion 
 #Deal with std vales of zero
 std_filled = std.copy()
 std_filled[std==0] = 1. 
 Xbar = ((X-mu)/std_filled)  #Normalized so mean is zero 
 return Xbar,mu,std

dim=data.shape #Data shape 
Xbar,mu,std= normalize(data/255)  #Normailzed data 
X=Xbar.reshape(-1,dim[1]*dim[2]) 
print(X.shape)  #Fit data n_samples (Imanges), n_features(pixels)
pca.fit(X)


#=========================#

pca_components=pca.components_.reshape(2,dim[1],dim[2])
print(pca_components.shape)
fig,axs= plt.subplots(1,len(pca_components))

for i,ax in enumerate(axs):
 ax.imshow(pca_components[i])
 ax.set_title(f"PCA component {i+1}")

#=========================#

pca_features=pca.transform(X)
print(pca_features.shape)
df=pd.DataFrame(pca_features, columns=['PC1', 'PC2'])
df["ids"]=ids
df.plot.scatter(x="PC1",y="PC2",c="ids",cmap="jet",figsize=(10,10))

#=========================#

print("pca.explained_variance_ratio_ :",pca.explained_variance_ratio_)
print(f"information retained {round(np.sum(pca.explained_variance_ratio_)*100,2)} %")


#=========================#

pca2 = PCA(n_components=0.66)
pca2.fit(X)
print('infromation retained: ', np.sum(pca2.explained_variance_ratio_))
print('n_compoments: ', pca2.n_components_) 

#=========================#
pca2_components=pca2.components_.reshape(pca2.n_components_,dim[1],dim[2])

fig,axs= plt.subplots(2,len(pca2_components)//2,figsize=(12,12))
axs=axs.reshape(-1)
for i,ax in enumerate(axs):
 ax.imshow(pca2_components[i])
 ax.set_title(f"PCA component {i+1}")
fig.tight_layout(h_pad=-25)

#=========================#