Python
scikit-learn
clustering
textmining
tf-idf

text miningや方向情報の特徴抽出に便利な,超球上に分布する点をクラスタリングする方法

More than 1 year has passed since last update.

初めに

研究の関係で,3次元のベクトルデータを機械学習の入力とする必要がありました.
データの形式としては,数秒ごとに長さがばらばらの3次元ベクトルを取得します.

ベクトルの方向のみに着目するとし,合成ベクトルの長さを1へ正規化します.
このように大量の時系列で得られた方向データをどのように特徴量抽出すれば良いか,というのが本内容です.
結論から言って,von Mises-Fisher Distribution (vMFD)からパラメータを最尤推定し,Mixture of vMFD (movMFD)で超球上の方向データをクラスタリングします.
このモデル化された値が,方向データにおけるある種の特徴量と言えます.
自然言語処理をされている方はわかると思いますが,これは文書中の単語頻出指標tf-idfに関係が深い手法です.
また物理空間で「ベクトルの向き」が時間変化するようなシステムに活用が可能です.

今回参考にした本は「異常検知と変化検知」の7章,「方向データの異常検知」です.

また,使用したプログラムはclara-labsのsphereclusterです.
さらにこちらの論文( "Clustering on the Unit Hypersphere using von Mises-Fisher Distributions")の詳細説明を一部使用しています.
あくまでD次元のベクトルデータを特徴量抽出することが目的ですが,sphereclusterの例題に沿って,特徴量抽出した後にクラスタリングするところまで全部やってしまいます.
このパッケージがあったおかげで,作業効率が一気に上がりました.本当にありがとうございます.

von Mises-Fisher Distributionや方向情報を活用するのは今回が初めてなので,発展的な意見や認識間違いがあればコメントをいただけるとありがたいです.

sphereclusterを使ってみる

パッケージの概要

sphereclusterはこちらの論文で紹介されているアルゴリズムを元に作成されたpython packageです.
主にSpherical K-means (spkmeans)とmovMFの2つのアルゴリズムから構成されており,movMFにはsoft-, hard-の2種類が存在します.
まずはこのsphereclusterをインストールします.

pip install spherecluster

簡単です.ちなみにsklearnやnumpy等のパッケージに依存しています.

Algorithms

Spherical Clustering

KMeansの概念を球面上へ拡張した手法.
入力空間が球面として定義されているだけで,計算方法はKMeansと同等.

vMFDについて

vMFDは超球上に分布する点を,正規分布のように確率分布として表現する一般的な表現方法です.wikipedia参照

Screen Shot 2017-01-27 at 13.26.54.png

数式を見てみると,入力xに対して[$\mu$, $\kappa$]という2つのパラメータで表現されていることがわかります.
$\mu$は分布している点の平均方向,$\kappa$は$\mu$周りの集中度です.
つまり下図のように,集中度$\kappa$が大きくなると平均方向$\mu$(実線)から離れた位置まで点が存在するというイメージです.

wikiのイメージ図 $\kappa$と分散の関係($\theta$は平均方向からの角度)
Screen Shot 2017-01-26 at 8.09.47.png

また,数式を見てわかる通り,正則化定数の計算には第1種変形ベッセル関数という複雑な関数を含んでおり,理解するのに少々時間がかかりそうです.

今回使用するmovMFDはこのvMFDの混合分布(Mixture Distribution)です.
本記事では数式の理解ではなく,実際の使用例に注目して解説していきます.

movMFD

混合分布モデル(Mixture Model)とvMFDで構成された混合分布の尤度をEMアルゴリズムで解く.

Screen Shot 2017-01-27 at 7.43.35.png

よくある混合ガウシアンモデルと同様のアプローチ.
つまりこのように球面上に独立したvMFDとしてクラスタリングされることになる.

sphere_w_clusters.png

Example1: small-mix

2次元のダミーデータをKMeans, Spherical KMeans, soft-movMF, hard-movMFの4通りでクラスタリングします.
実際のプログラムはこちらにあります.

初めにライブラリのインストールをします.

import numpy as np
from matplotlib import pyplot as plt
from sklearn.cluster import KMeans
from sklearn import metrics

import sys
sys.path.append('../spherecluster')

from spherecluster import SphericalKMeans
from spherecluster import VonMisesFisherMixture
from spherecluster import sample_vMF

そして次の2つのvMFに基づいたデータを生成します.
$\mu_0$=[-0.251, -0.968], $\kappa_0$=8
$\mu_1$=[0.399, 0.917], $\kappa_1$=2
データ数num_points_per_classは100です.

# Generate small-mix dataset
mu_0 = np.array([-0.251, -0.968])
mu_1 = np.array([0.399, 0.917])
mus = [mu_0, mu_1]
kappa_0 = 8 # concentration parameter
kappa_1 = 2 # concentration parameter
kappas = [kappa_0, kappa_1]
num_points_per_class = 100

X_0 = sample_vMF(mu_0, kappa_0, num_points_per_class)
X_1 = sample_vMF(mu_1, kappa_1, num_points_per_class)
X = np.zeros((2 * num_points_per_class, 2))
X[:num_points_per_class, :] = X_0
X[num_points_per_class:, :] = X_1
labels = np.zeros((2 * num_points_per_class, ))
labels[num_points_per_class:] = 1

そしてKMean,Spherical K-Means,soft-movMF,hard-movMFへフィッティングします.

この時に,次のようなエラーが発生する場合があります.

ValueError: Data l2-norm must be 1

movMFの関数内では入力データを各セットずつ固有値計算しているのですが,元々のsample_vMFが生成するデータには5%前後の確率で固有値の絶対値が1e-4より大きくなることがあります.この時エラー処理でプログラムが停止するため,ここでは無理やり固有値計算をしてエラーの原因をはじくことにします.

次のプログラムをフィッティング前に実行しておくと,エラーが無くなります.
(ただしエラー原因のデータも消えるので,サンプル総数も減少します)

import scipy.sparse as sp
dummy = []
n_samples = X.shape[0]
for ee in range(n_samples):
    if sp.issparse(X):
        n = sp.linalg.norm(X[ee, :])
    else:
        n = np.linalg.norm(X[ee, :])

    if np.abs(n - 1.) > 1e-4:
        dummy.append([ee, n])

X = np.delete(X, np.array(dummy)[:, 0], 0)

実行結果は次の通り.

small_mix_2d.png

元々のデータが綺麗に2グループへ分かれているため,手法ごとの目立った違いは見当たりません.

Example2: small-mix_3d

Example1と同様の内容で,3次元データに対してクラスタリングをしてみます.
元のコードはこちら.

スクリプトを見ると,Examaple1とほぼ同じことがわかります.
実行してみた結果がこちら.

figure_2.png

ちなみに,推定結果がこのようになっています.

Screen Shot 2017-01-27 at 10.58.06.png
少なくとも手元で生成したデータに関しては,kmeans, spkmeansと比べて良いパフォーマンスを発揮しているようです.

データ紹介

今回取り扱うデータを紹介します.

Screen Shot 2017-01-26 at 6.56.25.png

idは各ノードに紐づいた値で,今回は固定値とします.[x, y, z]は3次元ベクトル,azimuth, angleは置いといて,normは[x, y, z]を合成した「長さ」です.
こちらのデータが3000サンプル程度あるとします.

EDA

ひとまず簡単な可視化だけ済ませます.

各変数のヒストグラムと相関図は次の通り.
なんとなくパラメータごとに特徴があるように見えますが,具体的にどのようにしてその特徴を抽出すれば良いかまではわかりません.

download (4).png

3Dで可視化してみた結果がこちら.
ますますよくわかりません.

download (5).png

sphereclusterのexample(make a sphere graph)を参考にし,上の散布図から方向データだけ抽出し3Dの球面上へプロットしてみます.

全データの球面上プロット 前半と後半の球面上プロット
download.png download.png

なんとなく特徴が出てきました.
一応時間的な影響を考慮し,データのタイムスタンプを前半(青色)と後半(赤色)へ分けてプロットもしてみましが,特に前半と後半で時間的変化はなさそうです.

movMFへフィッティングする

いよいよmovMFをテストデータへフィッティングします.

###############################################################################
# Mixture of von Mises Fisher clustering (soft)
vmf_soft = VonMisesFisherMixture(n_clusters=1, posterior_type='soft', n_init=20)
vmf_soft.fit(X)

print "center:{}".format(vmf_soft.cluster_centers_)
print "kappa:{}".format(vmf_soft.concentrations_[0])

結果は次の通り.

center:[[-0.82651971  0.17847126  0.53386626]]
kappa:3.91851479297

無事3次元ベクトルから方向とその分散度合いを,確率密度分布でモデル化することができました.
例えば,Xの次元数を文書中に出現する単語の数と一致させると,D次元の超球上に分布する点をmovMFで分布モデル化,そしてクラスタリングをすることが可能となります.