【機械学習入門】「PCA主成分分析；多変量解析」の軸変換を3Dプロットで見える化♬

今回はわかりやすいと思うので、ちょっと理論から入ることとする。
データを$x_i=(x_{i1},…,x_{id})^T$, $(i=1,…,N)$、ここで$d$は独立な測定項目の次元である。
データ行列を $X=(x_i,…,x_{N})^T$ 、原点を平均 $\bar{x} = \sum_{i=1}^{N}x_{i}$としたデータ行列を $\bar{X}=(x_i−\bar{x},…,x_{N}−\bar{x})^{T}$とすると、PCA主成分分析は$S=\frac{1}{N}\bar{X}^{T} \bar{X}$という共分散行列を使って、以下のようないわゆる固有値問題となる。

```math
Su=λu
```
を解いて固有値を求めて$S$を対角化し、新たな主軸（分散の大きい順に第一主成分、第二主成分。。。）を求めることである。
ここでこの固有ベクトルが新たな主軸の軸方向を示しており、この固有ベクトルを利用すると、

```math
P^{-1}SP
```
と対角化できる。そして、元の方程式は以下のように変形でき、変換テンソルPは固有ベクトル$u_1,u_2,u_3...$で以下のように書けることが分かる。

```math
P^{-1}SP P^{-1}u=λP^{-1}u　\\
P=(u_1,u_2,u_3...)
```
つまり、このPを使って、新たな主軸（基底という）における各座標が

```math
P^{-1}u
```
で求められることが分かる。
また、この固有値を$s_k$とすると、各成分の寄与は以下のようになる。

```math
c_k=\frac{s_k}{V_{total}}
```
ここで重要なことは、上記の共分散行列の各行列要素は比例していないことが前提となる。すなわち測定項目間の相関が無いかあるいはあっても小さいことが前提条件になる。

今回はこの対角化の様子を多変量解析についてグラフ化して見える化することを試みる。
手法は前回の３D見える化と同じであるが、上式の新たな座標上の点への座標変換$P^{-1}u$を利用して示そうと思う。
【参考】
・[[Python][Scikit-learn]主成分分析を用いた次元削減、主成分ベクトルを用いた予測と線形回帰による予測の比較](https://dev.classmethod.jp/statistics/pythonscikit-learn-pca1/)
・[主成分分析と固有値問題@Aidemy Blog](https://blog.aidemy.net/entry/2017/10/19/222941)
・[固有値問題・固有空間・スペクトル分解@武内修＠筑波大](https://dora.bk.tsukuba.ac.jp/~takeuchi/?%E7%B7%9A%E5%BD%A2%E4%BB%A3%E6%95%B0II%2F%E5%9B%BA%E6%9C%89%E5%80%A4%E5%95%8F%E9%A1%8C%E3%83%BB%E5%9B%BA%E6%9C%89%E7%A9%BA%E9%96%93%E3%83%BB%E3%82%B9%E3%83%9A%E3%82%AF%E3%83%88%E3%83%AB%E5%88%86%E8%A7%A3)
###やったこと
・最初の一歩２次元
・三次元でGifアニメーション
・六次元を３DGifアニメーション
###今回のコードは以下に置きました
・[MachineLearning/pca_exaple.py](https://github.com/MuAuan/MachineLearning/blob/master/pca_exaple.py)
・[MachineLearning/out_gif.py](https://github.com/MuAuan/MachineLearning/blob/master/out_gif.py)
###・最初の一歩２次元
説明はこの2次元の場合が分かりやすい。

```py
import numpy as np
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
```
データはscikit_learnの例から以下のとおりのものを使う。
【参考】
・[sklearn.decomposition.PCA](https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html)

```py
X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
```
まず、2次元のPCA分析を実施する。
そして、その結果として、主成分方向のベクトルv_pca0と第二主成分方向のベクトルv_pca1が以下のとおり、得られる。

```py
pca = PCA(n_components=2)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
```
このベクトルはそれぞれ以下のとおり求められる。

```py
X_pca0=[0,1*v_pca0[0]]
Y_pca0=[0,1*v_pca0[1]]
X_pca1=[0,1*v_pca1[0]]
Y_pca1=[0,1*v_pca1[1]]
```
元々の空間で入力座標をプロットする。
さらに、上記で求めた第一主成分(赤)と第二主成分(青)のベクトルをプロットする。

```py
plt.scatter(X[:,0],X[:,1],c="g",marker = "*")
plt.plot(X_pca0,Y_pca0,c="r",marker = "o")
plt.plot(X_pca1,Y_pca1,c="b",marker = "o")
plt.savefig('k-means/pca_example/pca_data_scatter.jpg')
plt.pause(1)
plt.close()
```
第一主成分が赤い線の方向、第二主成分が青い線の方向である。
![pca_data_scatter.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/a5a8b0d1-799c-4f23-6fe1-315ceddae7cb.jpeg)
次に固有値ベクトルを利用して、測定点Xと主成分ベクトルYを変換して、変換された座標上にプロットする。この変換はscikit-learnではPCAの関数pca.transform()を利用して以下のとおり変換できる。

```py
Y=np.array([[0,0],[v_pca0[0],v_pca0[1]],[v_pca1[0],v_pca1[1]]])
pcatran_X = pca.transform(X)
pcatran_Y = pca.transform(Y)
print("pcatran_Y={}".format(pcatran_Y))
```
そして、変換された座標上に以下のとおりプロットする。

```py
plt.scatter(pcatran_X[:,0],pcatran_X[:,1],c="g",marker = "*")
plt.scatter(pcatran_Y[:,0],pcatran_Y[:,1],c="r",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[0,1]],[pcatran_Y[2,0],pcatran_Y[2,1]],c="b")
plt.plot([pcatran_Y[2,0],pcatran_Y[2,1]],[pcatran_Y[0,0],pcatran_Y[0,1]],c="r")  #,marker = "o")
plt.savefig('k-means/pca_example/pca_tran_scatter.jpg')
plt.pause(1)
plt.close()
```
まともな方向にプロットされ、新しい座標軸も第一主成分及び第二主成分方向になった。
![pca_tran_scatter.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/9b84ddfc-40cb-7d30-90a5-79876dbaec78.jpeg)
次に、これらの成分の寄与を以下のように求める。

```py
s=2
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca2d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()
```
![pca2d_contribution_PCA12_plotSn2.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/64793aa3-07aa-52dd-696a-bce6e36b4a9e.jpeg)
###・三次元でGifアニメーション
次元が増えるだけで、同じ流れで進みます。
データは、以下のようにmake_blobsで生成しました。

```py
from sklearn.datasets import make_blobs
X, Y = make_blobs(random_state=8,
                  n_samples=100, 
                  n_features=3, 
                  cluster_std=4,
                  centers=3)
```
PCA分析の成分は3とします。そして、3つの成分ベクトルが求まります。
そして、それぞれの方向のベクトルを定義します。

```py
pca = PCA(n_components=3)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
v_pca2=pca.components_[2]
X_pca0=[0,10*v_pca0[0]]
Y_pca0=[0,10*v_pca0[1]]
Z_pca0=[0,10*v_pca0[2]]
X_pca1=[0,10*v_pca1[0]]
Y_pca1=[0,10*v_pca1[1]]
Z_pca1=[0,10*v_pca1[2]]
X_pca2=[0,10*v_pca2[0]]
Y_pca2=[0,10*v_pca2[1]]
Z_pca2=[0,10*v_pca2[2]]
```
座標軸は３次元なので、以下のように立体的に描画します。


```py
def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.plot(X_pca0,Y_pca0,Z_pca0, c='b', marker='o', alpha=0.5, label='1st pricipal')  
    ax.plot(X_pca1,Y_pca1,Z_pca1, c='r', marker='o', alpha=0.5, label='2nd principal') 
    ax.plot(X_pca2,Y_pca2,Z_pca2, c='g', marker='o', alpha=0.5, label='3rd principal') 
    ax.scatter(X[:, 0], X[:, 1], X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/pca3d/keiba3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)
```
![pca3d512_72_n360.gif](https://qiita-image-store.s3.amazonaws.com/0/233744/418a4574-5efe-78a3-ab63-34157377f506.gif)
３次元でも同じようにYを定義して、一括変換して新たな座標上のベクトルを求め、スナップを書いてみます。

```py
Y=np.array([[0,0,0],[10*v_pca0[0],10*v_pca0[1],10*v_pca0[2]],[10*v_pca1[0],10*v_pca1[1],10*v_pca1[2]],[10*v_pca2[0],10*v_pca2[1],10*v_pca2[2]]])
print(Y)

pcatran_X = pca.transform(X)
print("pcatran_X={}".format(pcatran_X))
pcatran_Y = pca.transform(Y)
print("pcatran_Y={}".format(pcatran_Y))

plt.scatter(pcatran_X[:,0],pcatran_X[:,1],c="g",marker = "*")
plt.scatter(pcatran_Y[:,0],pcatran_Y[:,1],c="r",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[1,0]],[pcatran_Y[0,1],pcatran_Y[1,1]],c="b",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[2,0]],[pcatran_Y[0,1],pcatran_Y[2,1]],c="r",marker = "o")
plt.savefig('k-means/pca_example/pca_tran_scatter3.jpg')
plt.pause(1)
plt.close()
```
![pca_tran_scatter3.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/502f08c6-2d1f-3ee1-3a32-7f9d7277e7cc.jpeg)
次に新しい座標軸上に３dプロットします。

```py
def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.scatter(pcatran_X[:, 0], pcatran_X[:, 1], pcatran_X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('tran_x')
    ax.set_ylabel('tran_y')
    ax.set_zlabel('tran_z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/tran_pca3d/pca3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)
```
今回は変換された軸上にプロットしているので、軸の記載はしていません。
![tran_pca3d512_72_n360.gif](https://qiita-image-store.s3.amazonaws.com/0/233744/c2319217-ddc5-0f3c-884c-d2fdd23c8226.gif)
最後に寄与のグラフは以下のとおりになっています。

```py
s=3
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2,3], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca3d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()
```
第一主成分とそん色なく第二主成分や第三主成分の広がりがあるように見えますが、寄与は以下のとおり、７０％，１９％，１１％程度になっています。
![pca3d_contribution_PCA12_plotSn3.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/5a5bd6db-21e0-b9f0-c753-bd8f7553ade9.jpeg)
###・六次元を３DGifアニメーション
次元、広がり、そして３中心（n_features=6,cluster_std=4,centers=3）で生成しました。

```py
s=6
X, Y = make_blobs(random_state=8,
                  n_samples=100, 
                  n_features=6, 
                  cluster_std=4,
                  centers=3)
```
上と同じように、固有ベクトルはとりあえず第三主成分まで求めます。

```py
pca = PCA(n_components=6)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
v_pca2=pca.components_[2]

X_pca0=[0,10*v_pca0[0]]
Y_pca0=[0,10*v_pca0[1]]
Z_pca0=[0,10*v_pca0[2]]
X_pca1=[0,10*v_pca1[0]]
Y_pca1=[0,10*v_pca1[1]]
Z_pca1=[0,10*v_pca1[2]]
X_pca2=[0,10*v_pca2[0]]
Y_pca2=[0,10*v_pca2[1]]
Z_pca2=[0,10*v_pca2[2]]
```
そして、以下のように三次元プロットをします。

```py
def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.plot(X_pca0,Y_pca0,Z_pca0, c='b', marker='o', alpha=0.5, label='1st pricipal')  
    ax.plot(X_pca1,Y_pca1,Z_pca1, c='r', marker='o', alpha=0.5, label='2nd principal') 
    ax.plot(X_pca2,Y_pca2,Z_pca2, c='g', marker='o', alpha=0.5, label='3rd principal') 
    ax.scatter(X[:, 0], X[:, 1], X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/pca6d/keiba3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)
```
![pca6d_pca3d512_72_n360.gif](https://qiita-image-store.s3.amazonaws.com/0/233744/3f2d3bb7-0793-f4cc-0832-6605bd51426a.gif)
今回も同じように固有ベクトルで変換し、３ｄプロットします。

```py
pcatran_X = pca.transform(X)
print("pcatran_X={}".format(pcatran_X))

def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.scatter(pcatran_X[:, 0], pcatran_X[:, 1], pcatran_X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('tran_x')
    ax.set_ylabel('tran_y')
    ax.set_zlabel('tran_z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/tran_pca6d/pca3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)
```
![tran_pca6d_pca3d512_72_n360.gif](https://qiita-image-store.s3.amazonaws.com/0/233744/17a90c78-edc6-5731-c9cb-a3d728cab270.gif)
以下のようにそれぞれの成分の寄与を見ます。

```py
print("v_pca0={},v_pca1={},v_pca2={}".format(v_pca0,v_pca1,v_pca2))
s=3
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2,3,4,5,6], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca6d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()
```
第三成分までで、４５％，２５％，１０％で８０％程度をカバーしています。
![pca6d_contribution_PCA12_plotSn3.jpg](https://qiita-image-store.s3.amazonaws.com/0/233744/91b27ff9-5b8d-b499-ed52-d81c58293b20.jpeg)
###まとめ
・固有値問題としてのPCA分析の理論を見た
・二次元、三次元、そして六次元の場合について座標変換の様子を見た
・こうして次元圧縮の様子が分かった

・これをベースにさらなるクラスタリング手法と利用を見ていく
###おまけ
【参考】
・[機械学習 〜 データセット生成 〜](https://qiita.com/fujin/items/bb82d77b0b08c107f819)


```py
import numpy as np
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
from pandas.tools import plotting # 高度なプロットを行うツールのインポート
import pandas as pd
from mpl_toolkits.mplot3d import Axes3D

X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])

pca = PCA(n_components=2)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
X_pca0=[0,1*v_pca0[0]]
Y_pca0=[0,1*v_pca0[1]]
X_pca1=[0,1*v_pca1[0]]
Y_pca1=[0,1*v_pca1[1]]

plt.scatter(X[:,0],X[:,1],c="g",marker = "*")
plt.plot(X_pca0,Y_pca0,c="r",marker = "o")
plt.plot(X_pca1,Y_pca1,c="b",marker = "o")
plt.savefig('k-means/pca_example/pca_data_scatter.jpg')
plt.pause(1)
plt.close()
print(v_pca0[0],v_pca0[1])
Y=np.array([[0,0],[v_pca0[0],v_pca0[1]],[v_pca1[0],v_pca1[1]]])
print(Y)

pcatran_X = pca.transform(X)
pcatran_Y = pca.transform(Y)
print("pcatran_Y={}".format(pcatran_Y))

plt.scatter(pcatran_X[:,0],pcatran_X[:,1],c="g",marker = "*")
plt.scatter(pcatran_Y[:,0],pcatran_Y[:,1],c="r",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[0,1]],[pcatran_Y[2,0],pcatran_Y[2,1]],c="b")
plt.plot([pcatran_Y[2,0],pcatran_Y[2,1]],[pcatran_Y[0,0],pcatran_Y[0,1]],c="r")  #,marker = "o")
plt.savefig('k-means/pca_example/pca_tran_scatter.jpg')
plt.pause(1)
plt.close()

print("v_pca0={},v_pca1={}".format(v_pca0,v_pca1))
s=2
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca2d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()

print("pca.singular_values_={}".format(pca.singular_values_))  

from sklearn.datasets import make_regression
#X, Y, coef = make_regression(random_state=12, n_samples=100, n_features=3, n_informative=2, noise=10.0, bias=-0.0, coef=True)

from sklearn.datasets import make_blobs
X, Y = make_blobs(random_state=8,
                  n_samples=100, 
                  n_features=3, 
                  cluster_std=4,
                  centers=3)

pca = PCA(n_components=3)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
v_pca2=pca.components_[2]
X_pca0=[0,10*v_pca0[0]]
Y_pca0=[0,10*v_pca0[1]]
Z_pca0=[0,10*v_pca0[2]]
X_pca1=[0,10*v_pca1[0]]
Y_pca1=[0,10*v_pca1[1]]
Z_pca1=[0,10*v_pca1[2]]
X_pca2=[0,10*v_pca2[0]]
Y_pca2=[0,10*v_pca2[1]]
Z_pca2=[0,10*v_pca2[2]]

def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    #for kx, ky,kz, name in zip(kmeans_model.cluster_centers_[:,0], kmeans_model.cluster_centers_[:,1], kmeans_model.cluster_centers_[:,2], labels):
        #ax.text(kx, ky,kz, name, alpha=0.8, size=20)
    ax.plot(X_pca0,Y_pca0,Z_pca0, c='b', marker='o', alpha=0.5, label='1st pricipal')  
    ax.plot(X_pca1,Y_pca1,Z_pca1, c='r', marker='o', alpha=0.5, label='2nd principal') 
    ax.plot(X_pca2,Y_pca2,Z_pca2, c='g', marker='o', alpha=0.5, label='3rd principal') 
    ax.scatter(X[:, 0], X[:, 1], X[:, 2],  c="b", marker='o', alpha=1)
    #ax.scatter(kmeans_model.cluster_centers_[:,0], kmeans_model.cluster_centers_[:,1], kmeans_model.cluster_centers_[:,2], c = "b", marker = "*", s = 100)
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/pca3d/keiba3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)

print(v_pca0[0],v_pca0[1],v_pca0[2])
Y=np.array([[0,0,0],[10*v_pca0[0],10*v_pca0[1],10*v_pca0[2]],[10*v_pca1[0],10*v_pca1[1],10*v_pca1[2]],[10*v_pca2[0],10*v_pca2[1],10*v_pca2[2]]])
print(Y)

pcatran_X = pca.transform(X)
print("pcatran_X={}".format(pcatran_X))
pcatran_Y = pca.transform(Y)
print("pcatran_Y={}".format(pcatran_Y))

plt.scatter(pcatran_X[:,0],pcatran_X[:,1],c="g",marker = "*")
plt.scatter(pcatran_Y[:,0],pcatran_Y[:,1],c="r",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[1,0]],[pcatran_Y[0,1],pcatran_Y[1,1]],c="b",marker = "o")
plt.plot([pcatran_Y[0,0],pcatran_Y[2,0]],[pcatran_Y[0,1],pcatran_Y[2,1]],c="r",marker = "o")
plt.savefig('k-means/pca_example/pca_tran_scatter3.jpg')
plt.pause(1)
plt.close()

def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.scatter(pcatran_X[:, 0], pcatran_X[:, 1], pcatran_X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('tran_x')
    ax.set_ylabel('tran_y')
    ax.set_zlabel('tran_z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/tran_pca3d/pca3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)

print("v_pca0={},v_pca1={},v_pca2={}".format(v_pca0,v_pca1,v_pca2))
s=3
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2,3], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca3d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()

print("pca.singular_values_={}".format(pca.singular_values_))  

s=6

X, Y = make_blobs(random_state=8,
                  n_samples=100, 
                  n_features=6, 
                  cluster_std=4,
                  centers=3)


pca = PCA(n_components=6)
pca.fit(X)
v_pca0=pca.components_[0]
v_pca1=pca.components_[1]
v_pca2=pca.components_[2]

X_pca0=[0,10*v_pca0[0]]
Y_pca0=[0,10*v_pca0[1]]
Z_pca0=[0,10*v_pca0[2]]
X_pca1=[0,10*v_pca1[0]]
Y_pca1=[0,10*v_pca1[1]]
Z_pca1=[0,10*v_pca1[2]]
X_pca2=[0,10*v_pca2[0]]
Y_pca2=[0,10*v_pca2[1]]
Z_pca2=[0,10*v_pca2[2]]

def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    #for kx, ky,kz, name in zip(kmeans_model.cluster_centers_[:,0], kmeans_model.cluster_centers_[:,1], kmeans_model.cluster_centers_[:,2], labels):
        #ax.text(kx, ky,kz, name, alpha=0.8, size=20)
    ax.plot(X_pca0,Y_pca0,Z_pca0, c='b', marker='o', alpha=0.5, label='1st pricipal')  
    ax.plot(X_pca1,Y_pca1,Z_pca1, c='r', marker='o', alpha=0.5, label='2nd principal') 
    ax.plot(X_pca2,Y_pca2,Z_pca2, c='g', marker='o', alpha=0.5, label='3rd principal') 
    ax.scatter(X[:, 0], X[:, 1], X[:, 2],  c="b", marker='o', alpha=1)
    #ax.scatter(kmeans_model.cluster_centers_[:,0], kmeans_model.cluster_centers_[:,1], kmeans_model.cluster_centers_[:,2], c = "b", marker = "*", s = 100)
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/pca6d/keiba3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)

print(v_pca0[0],v_pca0[1],v_pca0[2])

pcatran_X = pca.transform(X)
print("pcatran_X={}".format(pcatran_X))

def show_with_angle(angle):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.view_init(elev = 10., azim = angle)
    ax.scatter(pcatran_X[:, 0], pcatran_X[:, 1], pcatran_X[:, 2],  c="b", marker='o', alpha=1)
    ax.set_xlabel('tran_x')
    ax.set_ylabel('tran_y')
    ax.set_zlabel('tran_z')
    ax.legend()
    plt.pause(0.01)
    plt.savefig('k-means/pca_example/tran_pca6d/pca3_PCA3d_angle_'+str(angle)+'.jpg')
    plt.close()
 
for angle in range(0, 360, 5):
    show_with_angle(angle)

print("v_pca0={},v_pca1={},v_pca2={}".format(v_pca0,v_pca1,v_pca2))
s=3
print("pca.explained_variance_ratio_{}".format(pca.explained_variance_ratio_))
plt.bar([1, 2,3,4,5,6], pca.explained_variance_ratio_, align = "center",label="pca.explained_variance_ratio_={}".format(pca.explained_variance_ratio_))
plt.title("pca.explained_variance_ratio_")
plt.xlabel("components")
plt.ylabel("contribution")
plt.legend()
plt.savefig('k-means/pca_example/pca6d_contribution_PCA12_plotSn'+str(s)+'.jpg')
plt.pause(1)
plt.close()

print("pca.singular_values_={}".format(pca.singular_values_))  
```