クラスタリングのクラスタ数を決めるのにはNbClust::NbClust()が便利

第88回R勉強会@東京にてさけのわデータプロジェクト | さけのわに関するLTをさせて頂きました。

好みの日本酒を呑みたい！〜さけのわデータで探す自分好みの酒🍶〜 from Takashi Kitano

好みの銘柄に近い銘柄を探す主旨から銘柄間の距離行列を求め、そこから階層型クラスタリングを行っていますが、クラスタ数をいくつにするか悩んだので試行錯誤の経緯が参考になれば幸い。

デンドログラムの描画

階層型クラスタリングなのでまずはデンドログラムを目視で見てみる。

## フレーバーチャート
charts <-
  jsonlite::read_json("https://muro.sakenowa.com/sakenowa-data/api/flavor-charts",
                      simplifyVector = TRUE) %>% 
  purrr::pluck("flavorCharts")

## 銘柄の距離行列
brands.dist.mat <-
  charts %>% 
  # 行名を削除
  tibble::remove_rownames() %>% 
  # 行名に銘柄IDを設定
  tibble::column_to_rownames(var = "brandId")  %>% 
  # コサイン類似度を計算
  proxy::dist(method = "cosine") 

brands.dist.mat %>% 
  # ウォード法で階層型クラスタリング
  hclust(method = "ward.D2") %>% 
  # デンドログラムを描画
  factoextra::fviz_dend()

少なくともクラスタ数4でも良さそう。6〜9のどれにするかは悩ましい・・・

{NbClust}で最適なクラスタ数を求める

こういう時に便利なのが以下の書籍で紹介されていた{NbClust}パッケージ。Best Number of Clustersということで最適なクラスタ数を提案してくれる。

アクセンチュアのプロフェッショナルが教えるデータ・アナリティクス実践講座

作者:アクセンチュアアナリティクス
発売日: 2016/05/31
メディア: 単行本（ソフトカバー）

以下でも使い方が紹介されている。

impsbl.hatenablog.jp

さけのわデータのフレーバーチャードで試してみる。indexにはkl, ch, hartigan, ccc, scott, marriot, trcovw, tracew, friedman, rubin, cindex, db, silhouette, duda, pseudot2, beale, ratkowsky, ball, ptbiserial, gap, frey, mcclain, gamma, gplus, tau, dunn, hubert, sdindex, dindex, sdbwといった指標を指定できるがallを指定すると計算量が大きいGAP, GAP, Gplus, Tauを除いて計算してくれる。（GAP, GAP, Gplus, Tauを含める場合はalllongを指定する。（手元の環境では終わらんかった・・・））

cls.hc.nb <- 
  NbClust::NbClust(
    data = charts[, -1], 
    method = "ward.D2", 
    index = "all", 
    min.nc = 2, 
    max.nc = 9
  )

f:id:tak95:20200920002609p:plain

クラスタ数3もしくは6が良さそう。樹形図をみると3だと少なすぎるように感じるので6で分割する。

brands.dist.mat %>% 
  # ウォード法で階層型クラスタリング
  factoextra::hcut(method = "ward.D2", k = 6) %>% 
  # デンドログラムを描画
  factoextra::fviz_dend()

f:id:tak95:20200920003647p:plain

良さそう！

非階層クラスタリングで試す

今回は距離行列を求めた経緯もあり階層型クラスタリングを行ったが{NbClust}はkmeansのような非階層クラスタリングでも活用できる。

cls.km.nb <- 
  NbClust::NbClust(
    data = scale(charts[, -1]), 
    method = "kmeans", 
    index = "all"
  )

f:id:tak95:20200920011516p:plain

クラスタ数は3もしくは6, 8が良さそう。

factoextra::fviz_cluster()でそれぞれの配置を見てみる。

cls.km <-
  c(3, 6, 8) %>% 
  purrr::map(function(x) {
    charts[, -1] %>% 
      scale() %>% 
      kmeans(centers = x) %>% 
      factoextra::fviz_cluster(data = charts[, -1])
  })