琥珀色呑んだくれ備忘録

メモとか備忘録とか

random forestを使った解釈性にもとづく推薦システム

[1706.06691] Interpretable Predictions of Tree-based Ensembles via Actionable Feature Tweakingという論文が提案されている。KDD2017でYahoo Researchの研究者によって発表されたもの、とのこと。Ensemble treesによって、あるサンプルがラベルAと予測されたとき、別のラベルBと予測されたかったら元のサンプルのどの特徴量にどのような操作を加えればよいか「改善案」を算出する。手元にあると便利そうだったのでRで実装したので紹介したい。
どんなことができるかの概要は、先日のTokyoRで喋ってきたので、そちらのスライドも参照して欲しい。

www.slideshare.net
大雑把な手順は、下記のとおり:

  1. データのスケーリング
  2. 予測モデル(random forest)の構築
  3. ルール抽出
  4. 各事例への推薦の計算

なお、日本語の解説は、接点QBさんによる丁寧な記事がすでにあるので、そちらを参照されたい。
setten-qb.hatenablog.com

前準備

ライブラリとソースのロード。この記事で使用するコードはgithubに置いてあるので、適当にローカルにおいて読み込む。そのうち時間を作ってパッケージ化したい。
github.com

require(tidyverse)
require(magrittr)
require(randomForest)
devtools::install_github("hoxo-m/pforeach")
require(pforeach)

source("./R/tweak_feature.R")

pforeachに隅々まで依存しているので、別途インストールが必須。

訓練セットとテストセットの準備

例によってirisを使いたいところだが、今回はspamデータセットを使ってspamと分類されるサンプルをnonspamに分類されるように推薦を行う。irisは多クラスだしルールも最小限なんで良いんだけど、アヤメの花びらを切って貼ってもカキツバタにはならないので今回はスキップ。実用上、変更できない説明変数が予測モデルに含まれる場合には、禁止リストをかませるとか適切な層別化するといった工夫が別途必要になる。

# data preparation -------------------------------------------------------------
set.seed(777)

data(spam, package = "kernlab")
dataset <- sample_frac(spam)
n.test <- floor(NROW(dataset) *0.1)

dataset.train <- chop(dataset, n.test)
dataset.test  <- tail(dataset, n.test)

dim(dataset);dim(dataset.train);dim(dataset.test)

chop()はデータを後ろから指定したサイズだけ切り落とすための自作ユーティリティ。

random forestの学習とモデル選択

訓練セットで学習して重要度のプロットと弱学習器の数の評価を行う。この辺りはルーチンなので変わったことはしていないが、ある程度予測精度が担保されていることが前提なので、ここで予測精度がダメな場合は再検討。

# bilding randomForest -----------------------------------------
X <- dataset.train[, 1:(ncol(dataset.train)-1)]
true.y <- dataset.train[, ncol(dataset.train)]

forest.all <- randomForest(X, true.y, ntree=500)
forest.all
par(mfrow=c(1,2))
varImpPlot(forest.all) # to view varImp, x3 & x4 should be removed.
plot(forest.all)
par(mfrow=c(1,1))

f:id:kato-satoshi-0:20180122185241p:plain
この方法は総当たり探索的なアプローチなので、変数(ルール)の数とツリーの数が増えると計算量が爆発的に増える。したがって変数選択とツリーの絞り込みのほか、max_depthを減らす工夫が重要。
今回は上位12変数と100ツリーくらいあれば良さそうなので、適当に切り落とす。

# model shrinkage based on importance -------------------------
top.importance <- forest.all$importance %>% data.frame %>%
  tibble::rownames_to_column(var = "var") %>% 
  arrange(desc(MeanDecreaseGini)) %>% 
  head(12)

dataset.train.fs <- dataset.train %>% select(top.importance$var)
dataset.test.fs  <- dataset.test %>% select(top.importance$var)

データのスケーリングと縮約済みモデルの学習

この手法はハイパーパラメータを一つだけにするため、各変数を標準化する必要がある。実用上、テストセットは訓練セットのスケーリングパラメータで標準化される必要があるのと、推薦後の改善は標準化前の尺度に復元する必要があり、毎回attributeを取ってくるのは面倒なので、それぞれrescale()とdescale()として、ユーティリティ関数を用意した。

# scaling feature-selected data  ---------------------------
X.train <- scale( dataset.train.fs )
X.test  <- rescale( dataset.test.fs, scaled = X.train )

dataset.test.fs[1:6, 1:6]
descale(X.test, scaled = X.train)[1:6, 1:6]
descale(X.test, scaled = X.test)[1:6, 1:6]

ちゃんと復元されていることを確認して欲しい。スケーリングした訓練データでモデルの学習を行う。

forest.rf <- randomForest(X.train, true.y, ntree=100)

forest.all
forest.rf
plot(forest.rf)

f:id:kato-satoshi-0:20180122191903p:plain
変数を減らしたので予測精度が下がっているが、「予測ミス(擬陽性)を起こすサンプルは本質的には改善されて良いのでは」という立場をとって許容範囲と思うことにして目をつぶることにする。

feature tweaking

ここからが提案手法の処理。手順は、random forest予測モデルを学習したのち、(1)ルールを抽出する+変更ルールを算出する(2)個別のサンプルの予測と推薦案を算出する(3)標準化の復元、の1+3ステップとなる。

# (0)random forest予測モデルを学習する
forest.rf <- randomForest(X.train, true.y, ntree=100)

# (1)ルールを抽出する+変更ルールを算出する
es.rf <- set.eSatisfactory(forest.rf, ntree = 30, epsiron = 0.3, resample = TRUE)

# (2)個別のサンプルの予測と推薦案を算出する
tweaked <- tweak(es.rf, newdata= X.test, label.from = "spam", label.to = "nonspam", .dopar = TRUE)

# (3)標準化の復元
dt <- descale.tweakedFeature(tweaked, X.test)
(1)ルールを抽出する+変更ルールを算出する

random forestではgetTree()で弱学習器を取ってこれる。これをつかって各ラベルに分類されるルールセットは、ラベルごとにdata.frameのlistとして整理しなおして抽出する。

# test sampling (for demo)
ep <- getRules.randomForest(forest.rf, k=2, label.to = NULL)
ep %>% str(2)
ep[["spam"]]$path[[1]]

f:id:kato-satoshi-0:20180122194147p:plain
xgboostにも同様にルールセットを取り出すxgb.model.dt.tree()があるので、そのうち対応したい。

抽出されたルールセットと変更ルール(候補案)のセットは対応しているため、抽出と同時に計算してしまう。元論文のアルゴリズムは各サンプルの評価時にこれを毎回計算する構造になっており、素直に実装してしまうと計算量が大変になるので、今から実装しようと思う人は注意。
epsiron は変更の強さを指定するハイパーパラメタ。原理上、どうやっても予測を覆せないサンプルがあり、この量は変更可能な提案が可能なサンプル数に影響する。

es.rf <- set.eSatisfactory(forest.rf, ntree = 30, epsiron = 0.3)

f:id:kato-satoshi-0:20180122200255p:plain
ntreeはルールセットを抽出する弱学習器の数で、省略するとすべてのツリーを使用する。デフォルトでは先頭からntreeで指定したツリーを解析対象にするが、データ不均衡に対処するために訓練セットをアンダーサンプリングして学習した幾つかのモデルをrandomForest::combine(...)で結合したい場合がある。この場合は、先頭から順番に解析してしまうと困るので、resampleオプションをTRUEに指定することでランダムに弱学習器を評価することができる。

es.rf <- set.eSatisfactory(forest.rf, ntree = 30, epsiron = 0.3, resample = TRUE)

f:id:kato-satoshi-0:20180122200705p:plain

(2)個別のサンプルの予測と推薦案を算出する

予測対象のデータセットが label.from と予測されたとき、label.to と予測されるための改善案を計算する。学習データと同じスケールで標準されていないとマトモな推薦が出来ないので注意。

tweaked_ <- tweak(es.rf, newdata= X.test, label.from = "spam", label.to = "nonspam",
                 .dopar = FALSE)

f:id:kato-satoshi-0:20180122202613p:plain
各サンプルについてルールを全部評価するので結構な計算時間がかかるため、並列化オプション.doparはデフォルトではTRUEになっている。

tweaked_ <- tweak(es.rf, newdata= X.test, label.from = "spam", label.to = "nonspam",
                 .dopar = TRUE)
(3)標準化の復元

元のサンプルの操作方向とその量が知りたいのでスケーリングを解除する。元のサンプルの特徴量(original)、予測を変えるための変更後の特徴量(suggest)、予測の変更に必要な努力量=両者の差分(diff)のlistを返す。
>|r|
dt <- descale.tweakedFeature(tweaked, X.test)
dt %>% str(1)
|

f:id:kato-satoshi-0:20180122204026p:plain

個別のサンプルへの推薦パターンの可視化

個別の推薦結果を可視化する。操作量を相対的に比較したいため、標準化されたままのものを可視化している。サンプルによって違う推薦パターンが見えるので面白い。

# plot suggestion for each instance
which(tweaked$predict == "spam")
plot.suggest(tweaked, 4)
plot.suggest(tweaked, 11)

f:id:kato-satoshi-0:20180122203751p:plain
f:id:kato-satoshi-0:20180122203745p:plain
操作量の大きい順に、0でない変数だけ可視化したければ、それぞれ .ordered = TRUEと、.nonzero.only = TRUEを指定する。

plot.suggest(tweaked, 15)
plot.suggest(tweaked, 15, .ordered = TRUE, .nonzero.only = TRUE)

f:id:kato-satoshi-0:20180122203928p:plain
f:id:kato-satoshi-0:20180122203934p:plain

変数重要度の可視化

論文中ではユーザーによる納得感(アンケート調査?)をもとに評価しているが、それだと機械的には評価できないので、「予測サンプル集団において操作回数が多かった変数は重要なんじゃないか」という考えと、「操作された変数の操作の方向性が揃っている変数は、その意味合いを解釈できそう」と考えて、これらを各変数の重要度と思うことにして、それぞれ可視化する。

# Plot population importances
pp <- plot.tweakedPopulation(tweaked, type = "absoluteSum")
pp <- plot.tweakedPopulation(tweaked, type = "direction")
pp <- plot.tweakedPopulation(tweaked, type = "frequency")

type = "frequency"は、サンプル集団において各変数が操作をされた集計値を可視化する。
f:id:kato-satoshi-0:20180122205329p:plain
type = "direction"は、各サンプルの推薦で各変数がどのような操作をされたかを集計して可視化する。相対的な重要性を比較したいので、標準化された状態の推薦量を使用する。
f:id:kato-satoshi-0:20180122210838p:plain
変更量の総和を可視化する"absoluteSum"オプションも作ってみたが、これはあまり有意義と感じなかった。

おわりに

カテゴリ変数を含むデータセットの取り扱い、変更できない変数に対する操作の提案、変更の実現コスト(お金とか労力とかの実リソース)に関するドメイン知識で重みづけ補正する、など、いろいろ工夫の余地があって実際の分析の際にはアイデアの見せどころと思う。
また、この種の感度分析的なアプローチ全般に言えることだが、因果関係の示唆となる関連性を補足しているかもしれないが、かならずしも因果関係そのものを見ているとは限らない点も気を付けたい。