クラスタリングではクラスタ数を事前に与えることが多いのですが、
クラスタ数も自動的に推定してしまうのがディリクレ過程混合モデルです。

ディリクレ過程混合モデルを動かそうと思うと、
モデル式を書いて、MCMCで数値計算するのがシンプルな方法ですが、
結構面倒なのでちょっと試してみるには少々向いていない。。。

そんな時、簡単に動かすことのできるパッケージを見つけたので、
ここで紹介します。

ディリクレ過程混合モデルとは

過去の記事でディリクレ過程混合モデルについて紹介しているので、
そちらも参照して下さい。

簡単に言ってしまうと、クラスタ数を予め決めずに
データから推定するクラスタリングのためのモデルです。

パッケージのインストール

今回使うパッケージはdirichletprocessというパッケージです。

install.packages("dirichletprocess")

でインストールすることができます。

まずは、パッケージのマニュアルにも乗っている簡単な例で動かして見ます。

library(dirichletprocess)

dp = DirichletProcessMvnormal(as.matrix(faithful %>% scale))
fdp = Fit(dp,500)
plot(fdp)

ここでは、Rの標準データセットであるfailthfulのデータを使って、
クラスタリングを行なっています。

DirichletProcessMvnormalで、
多次元正規分布を混合したディリクレ過程混合モデルを生成して、
Fitでモデルのフィッティングが行われます。

データは行列形式で与えて、各列が一つの変数、
各行が一つのデータサンプルを表します。
ただし、データはscaleで正規化しておきます。

結果は、下の図のようにプロットされます。

ほんの数行で実行できてしまうのでかなりお手軽ですね。

せっかくなので、自分でデータを作って、
もう少しクラス数の多い状況で試してみましょう。

d次元でC個の異なる平均値を持つ正規分布から、
N個のデータを生成します。

多変量正規分布からデータを生成するために、
mvnfastパッケージを利用しているのでインストールしておきましょう。

library(mvnfast)
C = 5
N = C*25
d = 2

y = NULL
for(i in 1:C){
  mu = rnorm(d,sd = 10)
  n = round(N/C)
  y = rbind(y,rmvn(n,mu,diag(1,d)))
}
plot(y)

今の場合、5個の正規分布からデータを生成していて、
例えば、下のようなデータが得られました。

2クラスの場合と同様に、クラスタリングを実行してみます。

dp = DirichletProcessMvnormal(as.matrix(y %>% scale))
fdp = Fit(dp,1000)
plot(fdp)

実行結果は下のようになりました。

真ん中の方と端の方の2つに分類されてます。
確かにそれっぽいといえばそれっぽいのですが、
もう少しクラスを分けて欲しいですよね。

ディリクレ過程混合モデルでは、
クラスタ数はデータから自動的に推定されますが、
その数をある程度調整するパラメータ $\alpha$ があります。
$\alpha$ が大きいほど、多くのクラスタが生まれやすいモデルとなります。

dirichletprocessでは、この $\alpha$ 自体も自動的に調整されるのですが、
$\alpha$ の事前分布を与えることである程度コントロールできます。

$\alpha$ の事前分布はガンマ分布で与えられていて、
デフォルトでは $G(2,4)$ となっています。
この状況では $\alpha$ の値は小さくなりやすいので、
大きな値を取りやすいように設定して実行してみましょう。

dp = DirichletProcessMvnormal(as.matrix(y %>% scale),alphaPriors = c(100,10))
fdp = Fit(dp,1000)
plot(fdp)

DirichletProcessMvnormalの引数alphaPriorsで、
$\alpha$ の事前分布のガンマ分布の形を設定しています。
ここでは $G(100,10)$ の分布を設定してみました。
下の図のように、デフォルトの分布よりも大きな値を取りやすい分布です。

下の図はクラスタリングの結果です。

真ん中の二つはくっ付いてしまっていますが、
それ以外は、かなり綺麗にクラスタリングできていることがわかります。

ここではRで簡単にディリクレ過程混合モデルによる
クラスタリングができるパッケージを紹介しました。

クラスタリングの実行は数行でできてしまうので、
とりあえず試したいという時にぴったりですね。