続・機械学習モデルを解釈する方法 SHAP value

前回の話

kaggleの中に、Machine Learning for Insights Challengeという4日間の講座がある。
後半2日間はSHAP valueが題材だったので、SHAP valueについてまとめる。
Machine Learning for Insights Challengeの内容、前半2日間の内容については前回のエントリを参照。
linus-mk.hatenablog.com

ちなみに今気づいたのだが、この4日間講座はkaggle Learnの講座一覧の中に加わっている。以下のサイトから各日の講義や演習問題に行ける。
Machine Learning Explainability | Kaggle
……って、タイトルに『Explainability（説明可能性）』って書いてあるじゃん! やっぱり解釈可能性の話じゃん！
最初の講座のときはExplainabilityなんて一言も書いてなくて、謎の「insight」って書いてあったのに……

SHAP valueとは何か

SHAP value(SHapley Additive exPlanationsの略) は、それぞれの予想に対して、「それぞれの特徴量がその予想にどのような影響を与えたか」を算出するものである。

1つのインスタンスを指定すると、このような図ができる。（講座ページから引用）

f:id:soratokimitonoaidani:20181027155817p:plain — SHAP value の例

赤色の矢印は予測値の増加を表し、青色の矢印は予測値の減少を表している。
また、SAHP valueでは「ベースラインの値」との比較を実施している。実際には、「ベースライン」というのは（データセット全体の）期待値である。

用いたデータセットはPredict FIFA 2018 Man of the Match | Kaggleであり、これはサッカーの試合のデータである。
「各試合でのチームの成績データを入力として、そのチームの選手がその試合でMan of the Match（サッカーなどスポーツの試合で最も活躍した選手）を獲得したかを予想する」
1つのインスタンス（ある試合での片方のチームの成績データ）を入力すると、モデルが「このチームがMan of the Matchを獲得する確率は0.70だ」と判定したが、その内訳を特徴量ごとに分解して示してくれる。

例えば、赤色の最長の矢印は「Goal Scored = 2」と書いてある。この矢印の長さは約0.11である。
そのチームの得点が、ベースラインの値（データの平均値）の代わりに2であったことで、Man of the Match獲得確率は約0.11ほど増加した、という風に読める。
逆に、青色の最長の矢印を見ると、そのチームのBall Possession % (ボールのポゼッション、支配率)が、ベースラインの値の代わりに38%であったことで、Man of the Match獲得確率は約0.06ほど減少した、という風に読みとれる。

応用（どういう目的で使えるか）

次のような応用に使うことができる。

銀行がある人にお金を貸すべきでないと機械学習モデルは言っており、銀行はそれぞれのローン拒否の根拠を説明する法的な必要がある
医療提供者は、それぞれの患者がある疾病にかかるリスクをどの要因が高めているのかを特定したい。狙いを定めた健康介入によって、それらのリスク要因に直接対処できるようにしたいからである。

数学的な背景については省略

1日目のPermutation Importanceや2日目のPartial Dependence Plot (PDP)では理論的な説明があった。しかしSHAP valueについては、数学的・理論的な仕組みにほとんど触れていない。

感覚的に言えば、SHAP valueは「特定のインスタンスとベースラインとの差分」を、それぞれの特徴量に関して「分割して」くれる。どうやって差分を「うまく分割」しているのかがアルゴリズムの鍵なのだが、その仕組については講座の中では言及していない。

4日目の方では、単純な例として
$y = 4x_1+2x_2$

を考えて、 $x_1$ が0(ベースライン)から2に変化したら、 $x_1$ に対するSHAP_valueは8になるよ、と書いてある。
……いや、そりゃ、 $x_1$ と $x_2$ の影響が分離されているモデルなら各変数の寄与を求めるのは簡単である。
モデルが仮に
$y = 4x_1+2x_2 +x_1x_2$
だったら、 $x_1$ がyにいくら寄与したかを判定するのは難しくなる。
$(x_1, x_2)$ が(0, 0)から(2, 3)に変化したらyは0から20 (4×2+2×3+2×3) に変化するわけだが、ではこの20のうち $x_1$ による変化はいくつだろうか。これはそんなに自明ではない。
ましてや、機械学習モデルのような複雑な式の場合にはどうやってそれぞれの特徴量の寄与を計算するんだろうか?

個人的には背後の理論はとても気になるところだが、理論を追うよりも、手を動かしてコードを書けるようになることに注力したいので
今は数学の理論は一旦おいておこう。

オプション：Summary Plots 特徴量ごとに分布を見る

SHAP valueのライブラリにはいくつかの描画方法があり、その1つには特徴量ごとに分布を見られるsummary_plot関数がある。

以下の図は講座ページからの引用。 f:id:soratokimitonoaidani:20181027173512p:plain

縦にそれぞれの特徴量が並んでいる。
1つの点が1つのインスタンスを表す。
点の色は、その特徴量の値が大きいか小さいかを表す。
横の位置は、そのインスタンスの特徴量のSHAP valueの値を示す。

その他の色々な描画方法については以下を参照。
SHAPライブラリのGitHub：GitHub - slundberg/shap: A unified approach to explain the output of any machine learning model.
SHAPライブラリの公式ドキュメント：Explainers — SHAP latest documentation

例：Irisデータセットの場合

最後に、試しにIrisデータセットに対してRandom Forestを適用して、SHAP valueを表示させてみよう。

from sklearn.datasets import load_iris
from sklearn.ensemble import RandomForestClassifier

iris = load_iris()
rnd_clf = RandomForestClassifier(n_estimators=500, n_jobs=-1, random_state=42)
rnd_clf.fit(iris["data"], iris["target"])

row_to_show = 76
data_for_prediction = iris["data"][row_to_show]  # use 1 row of data here. Could use multiple rows if desired
data_for_prediction = pd.Series(data_for_prediction, index=iris["feature_names"])
data_for_prediction_array = data_for_prediction.values.reshape(1, -1)

rnd_clf.predict_proba(data_for_prediction_array)

# → array([[0.   , 0.966, 0.034]])

76は適当なサンプル。（クラス1で、正しいクラスの確率が高くなりすぎないものを適当に選んでいる）
今回の76番はクラス1のサンプルである。そして学習したモデルの予測結果は、「確率0.966でクラス1」となっている。すなわち、正しく予測できている。
では、このサンプルに対して正しく学習できたのはどの特徴量のおかげだろうか？

Irisデータの特徴量は4次元であるが、最初の2次元の特徴量はあまり重要ではない。
pythonで散布図行列を描くでIrisデータの散布図を描いているのを見ると、
最初の2つの特徴量を見ても、'versicolor'と'virginica' (クラス1と2)は分離できないことが分かる。
大雑把に言えば「最初の2つの特徴量を与えられてもクラス1だとは分かりません。だけど、後半の2つの特徴量を見るとクラス1だと予想できます」という結果になるはずだ。
以上の予想を踏まえて、SHAP valueを表示させてみる。

import shap  # package used to calculate Shap values

# Create object that can calculate shap values
explainer = shap.TreeExplainer(rnd_clf)

# Calculate Shap values
shap_values = explainer.shap_values(data_for_prediction)

target_class = iris["target"][row_to_show]
shap_values[target_class]
# → array([0.02789272, 0.00855236, 0.23464138, 0.36383353])

iris["data"][row_to_show]
# → array([6.8, 2.8, 4.8, 1.4])

shap.initjs()
shap.force_plot(explainer.expected_value[target_class], shap_values[target_class], data_for_prediction)

f:id:soratokimitonoaidani:20181028225048p:plain

確かに前半2つの特徴量のSHAP valueは小さく、後半2つの特徴量のSHAP valueは大きい。
「最初の2つの特徴量を与えられてもクラス1だとは分かりません。だけど、後半の2つの特徴量を見るとクラス1だと予想できます」という事前の予想と合致した出力になっている。