主成分分析のお勉強

やりたいこと

三次元点群に最もフィットする主軸を導出したい!

教科書

• 主成分分析の典型的な例題   点群の主軸を求める - 理数アラカルト -

• http://www.sist.ac.jp/~kanakubo/research/statistic/shuseibun_bunseki.html

• 主成分分析を使ってバウンティボックスを作る - Pashango’s Blog <- わかりやすい

• ２直線の交点を求める公式 - 具体例で学ぶ数学

実装コード

#coding:utf-8
import numpy
from sklearn.decomposition import PCA

N = 20
t = numpy.random.rand(N)
x = t+numpy.random.rand(N)*0.3
y = t+numpy.random.rand(N)*0.3
z = t+numpy.random.rand(N)*0.3

pca_o = PCA(n_components=3)
pca_o.fit(zip(x,y,z))
pca = pca_o.components_

v  = []
mx = numpy.mean(x)
my = numpy.mean(y)
mz = numpy.mean(z)

for tx,ty,tz in zip(x,y,z):
    ax,ay,az = numpy.dot(pca,[(tx-mx),(ty-my),(tz-mz)])
    v.append( dict(ax=ax,ay=ay,az=az,x=tx,y=ty,z=tz) )

maxx = max(v,key=lambda x:x['ax'])
minx = min(v,key=lambda x:x['ax'])
maxy = max(v,key=lambda x:x['ay'])
miny = min(v,key=lambda x:x['ay'])
maxz = max(v,key=lambda x:x['az'])
minz = min(v,key=lambda x:x['az'])

#グラフ描画
from matplotlib.pyplot import *
from mpl_toolkits.mplot3d import Axes3D

fig = figure()
ax = Axes3D(fig)
ax.scatter3D(x,y,z)

# 表示範囲の設定
ax.set_xlim(0, 1.5)
ax.set_ylim(0, 1.5)
ax.set_zlim(0, 1.5)

#基底ベクトルの描画
v1x,v1y,v1z = pca[0]
ax.plot([mx+v1x,mx-v1x],[my+v1y,my-v1y],[mz+v1z,mz-v1z])
v2x,v2y,v2z = pca[1]
ax.plot([mx+v2x,mx-v2x],[my+v2y,my-v2y],[mz+v2z,mz-v2z])
v3x,v3y,v3z = pca[2]
ax.plot([mx+v3x,mx-v3x],[my+v3y,my-v3y],[mz+v3z,mz-v3z])

ax.plot([maxx['x']], [maxx['y']], [maxx['z']], "*", color="b", ms=10, mew=0.5)
ax.plot([minx['x']], [minx['y']], [minx['z']], "*", color="b", ms=10, mew=0.5)
ax.plot([maxy['x']], [maxy['y']], [maxy['z']], "*", color="r", ms=10, mew=0.5)
ax.plot([miny['x']], [miny['y']], [miny['z']], "*", color="r", ms=10, mew=0.5)
ax.plot([maxz['x']], [maxz['y']], [maxz['z']], "*", color="g", ms=10, mew=0.5)
ax.plot([minz['x']], [minz['y']], [minz['z']], "*", color="g", ms=10, mew=0.5)

show()

ブルーの線が主軸です！

☆がそれぞれの軸での最大・最小の点!!

これを使って何をするかはお楽しみ♪