ノート/ノート

訪問者数　7678　　　　　　最終更新　2015-07-22 (水) 12:48:41

Pythonでクラスタリング （2015-07-22） †

↑

x-meansという論文 †

クラスタリングのk-means法は広く使われ、パッケージも各言語環境で作られている。Pythonの 場合、scipyにあるkmeansと、多分その導出系のscikit-learnにあるKMeansが利用可能になっている。

k-means法はあらかじめクラスタ数を決めて分類するので、クラスタ数が不明なときにそのままは対応できない。 それに対して、BIC（ベイズ情報量規準）やAIC（赤池情報量規準）などを使って、最適なモデルを 選択しようという考え方があり、その選択のロジックを含めたのがxmeansパッケージらしい。

xmeanについてはhttp://web-salad.hateblo.jp/entry/2014/07/19/200347

↑

Scikit-learnでのk-meansパッケージ †

Scikit-learnのkmeansパッケージの使い方は、numpyのarrayを、次元別に用意したものを 組合せるのが良いらしい。

x_means = xmeans.XMeans(random_state = 1).fit(np.c_[x,y,z])

但し、x, y, zはnumpyのarrayデータ。np.c_は

c=array([1,2,3]), d=array([4,5,6])のとき、np.c_[c,d]=array([[1,4],[2,5],[3,6]])

のような動作をする。

↑

サンプル（matplotlibによる表示を含む） †

２次元の場合

#!/usr/bin/env python
# -*- coding: utf-8 -*-
'''
xmeansを、(1)外部からライブラリとして使う。(2)３次元でのクラスタリング。
xmeanについてはhttp://web-salad.hateblo.jp/entry/2014/07/19/200347
'''

import sys, string, codecs
import numpy as np
from mpl_toolkits.mplot3d import Axes3D
import xmeans
import matplotlib.pyplot as plt
 
# データの準備
x = np.array([np.random.normal(loc, 0.1, 20) for loc in np.repeat([1,2], 2)]).flatten()
y = np.array([np.random.normal(loc, 0.1, 20) for loc in np.tile([1,2], 2)]).flatten()
# 実際の効果は　xdom = np.array([1,1,2,2])
#               ydom = np.array([1,2,1,2])

# クラスタリングの実行
x_means = xmeans.XMeans(random_state = 1).fit(np.c_[x,y]) #それぞれの要素のペアを作る
#     c=array([1,2,3]), d=array([4,5,6])のとき、np.c_[c,d]=array([[1,4],[2,5],[3,6]])

print(x_means.labels_)
print(x_means.cluster_centers_)
print(x_means.cluster_log_likelihoods_)
print(x_means.cluster_sizes_)

# 結果をプロット
#plt.rcParams["font.family"] = "Hiragino Kaku Gothic Pro"
plt.rcParams["font.family"] = "VL Gothic"
plt.scatter(x, y, c = x_means.labels_, s = 30)
plt.scatter(x_means.cluster_centers_[:,0], x_means.cluster_centers_[:,1], c = "r", marker = "+", s = 100)
plt.xlim(0, 3)
plt.ylim(0, 3)
plt.title(u"改良x-means法の実行結果  参考: 石岡(2000)")
plt.show()
# plt.savefig("clustering.png", dpi = 200)

３次元の場合

#!/usr/bin/env python
# -*- coding: utf-8 -*-
'''
xmeansを、(1)外部からライブラリとして使う。(2)３次元でのクラスタリング。
xmeanについてはhttp://web-salad.hateblo.jp/entry/2014/07/19/200347
３次元にしたので、mplot3dについてはhttp://qiita.com/mojaie/items/c993fbb3aa63d0001c1c
'''

import sys, string, codecs
import numpy as np
from mpl_toolkits.mplot3d import Axes3D
import xmeans

import matplotlib.pyplot as plt

# データの準備
xdom = np.array([1,1,1,1,2,2,2,2])
ydom = np.array([1,1,2,2,1,1,2,2])
zdom = np.array([1,2,1,2,1,2,1,2])
x = np.array([np.random.normal(loc, 0.1, 200) for loc in xdom]).flatten()
y = np.array([np.random.normal(loc, 0.1, 200) for loc in ydom]).flatten()
z = np.array([np.random.normal(loc, 0.1, 200) for loc in zdom]).flatten()
# クラスタリングの実行
x_means = xmeans.XMeans(random_state = 1).fit(np.c_[x,y,z]) 
#  c=array([1,2,3]), d=array([4,5,6])のとき、np.c_[c,d]=array([[1,4],[2,5],[3,6]])

print(x_means.labels_)
print(x_means.cluster_centers_)
print(x_means.cluster_log_likelihoods_)
print(x_means.cluster_sizes_)

# 結果をプロット
fig = plt.figure()
ax = Axes3D(fig)
# 軸ラベルの設定
ax.set_xlabel("X-axis")
ax.set_ylabel("Y-axis")
ax.set_zlabel("Z-axis")
# 表示範囲の設定
ax.set_xlim(0, 3)
ax.set_ylim(0, 3)
ax.set_zlim(0, 3)

# グラフ描画色分け編
'''
labelsはクラスタ番号。それによって配列x, y, zから値を抜きたい
 x = np.arange(9);   s = np.array([0,1,0,1,0,1,0,1,0])
 np.select([s==0],[x])
とすると、抜ける。結果はarray([0, 0, 2, 0, 4, 0, 6, 0, 8])
'''

labels = x_means.labels_
number_of_clusters = len(x_means.cluster_sizes_)
colors = ["#ff0000", "#00ff00", "#0000ff", "#cccc00", "#00cccc", "#cc00cc", "#ccccc0", "#c0cccc", "#ccc0cc", "#0c0c0c"]

if number_of_clusters > len(colors):
  print "number of clusters %d exceeds number of colors" % number_of_clusters
  quit()
for i in range(number_of_clusters):    # クラスタごとに色を変えてプロットする
  ux = np.select([labels==i], [x])     # ラベルの値がiの要素だけを[x]から選ぶ
  uy = np.select([labels==i], [y])
  uz = np.select([labels==i], [z])
  ax.plot(ux, uy, uz, "o", color=colors[i], ms=4, mew=0.5)

plt.show()

得られた絵は、

正規分布分散0.1正規分布分散0.27
3dxmeans_0.1.png3dxmeans_0.27.png
添付ファイル: file3dxmeans_0.1.png 637件 [詳細] file3dxmeans_0.27.png 1255件 [詳細]
トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2015-07-22 (水) 12:48:41 (1237d)