法律や金融などの分野では説明性が重視されるケースがあるぞ
どうやったら説明性を付与できますか。
shapというライブラリを使用すればその問題点が解消できるぞ
目次
環境構築
環境構築はGoogle Colabで行います。
Google Colab上でpipで必要なライブラリをインストールします。
https://github.com/slundberg/shap
下記のコードをベースに処理をみていきます
pipで必要なライブラリをインストールします。
!pip install shap
shapに用意されているサンプルデータで線形モデルを学習します。今回はボストンにある物件の価格を予測するモデルの学習になります。
import pandas as pd
import shap
import sklearn
# a classic housing price dataset
X,y = shap.datasets.boston()
X100 = shap.utils.sample(X, 100)
# a simple linear model
model = sklearn.linear_model.LinearRegression()
model.fit(X, y)
モデルの出力と各特徴量の相関係数を導出します。
print("Model coefficients:\n")
for i in range(X.shape[1]):
print(X.columns[i], "=", model.coef_[i].round(4))
下記のような出力になります。
Model coefficients:
CRIM = -0.108
ZN = 0.0464
INDUS = 0.0206
CHAS = 2.6867
NOX = -17.7666
RM = 3.8099
AGE = 0.0007
DIS = -1.4756
RAD = 0.306
TAX = -0.0123
PTRATIO = -0.9527
B = 0.0093
LSTAT = -0.5248
最も正の相関がある特徴量RMに関して図にして見てみます。灰色の横軸は価格の期待値で縦軸はRMの期待値です。青色の線はRMが大きくなればなるほど家の価格が大きくなることを示しています。
shap.plots.partial_dependence("RM", model.predict, X100, ice=False, model_expected_value=True, feature_expected_value=True)

では負の相関がある特徴量NOXを見てみます。RMとは逆の傾向が確認できます。
shap.plots.partial_dependence("NOX", model.predict, X100, ice=False, model_expected_value=True, feature_expected_value=True)

SHAP値
下記リンクの内容が分かりやすかったです。
SHAP値を用いたモデルの説明
各特徴量がモデルにどの程度、影響するかを表す際にSHAP値を用いつのでその値も確認してみます。
# compute the SHAP values for the linear model
explainer = shap.Explainer(model.predict, X100)
shap_values = explainer(X)
shap.plots.scatter(shap_values[:,"RM"])

同様に”NOX”も見てみます。

モデルの出力に各特徴量がどの程度貢献しているかを把握することもできます。
下記の図は各特徴量を入力していき、現在のモデルの出力との差がどの程度あるか把握できます。
差の絶対値が大きいほどモデルの出力に貢献していることを把握できます。
shap.plots.waterfall(shap_values[sample_ind], max_display=14)

GAMモデル
先程までは線形モデルではモデルの性能が十分でないケースがあるので、この説明能力を保ちながらモデルの性能を上げたいケースがあります。
その性質を持つGAMのモデルを使用します。
`shap.Explainer`にモデルとデータを入力するだけで良いので非常に簡単に行えます。interpretライブラリが必要なので導入します。
! pip install interpret
# fit a GAM model to the data
import interpret.glassbox
model_ebm = interpret.glassbox.ExplainableBoostingRegressor()
model_ebm.fit(X, y)
# explain the GAM model with SHAP
explainer_ebm = shap.Explainer(model_ebm.predict, X100)
shap_values_ebm = explainer_ebm(X)
# make a standard partial dependence plot with a single SHAP value overlaid
fig,ax = shap.partial_dependence_plot(
"RM", model_ebm.predict, X, model_expected_value=True,
feature_expected_value=True, show=False, ice=False,
shap_values=shap_values_ebm[sample_ind:sample_ind+1,:]
)

同様に特徴量”RM”に対するSHAP値を確認してみます。
shap.plots.scatter(shap_values_ebm[:,"RM"])

同様にどの特徴量が貢献度が高いかを見てみます。
shap.plots.waterfall(shap_values_ebm[sample_ind], max_display=14)

XGboost
次により複雑なモデルであるXgboostを試してみます。
基本的に今までと同様に行えます。
# train XGBoost model
import xgboost
model_xgb = xgboost.XGBRegressor(n_estimators=100, max_depth=2).fit(X, y)
# explain the GAM model with SHAP
explainer_xgb = shap.Explainer(model_xgb, X100)
shap_values_xgb = explainer_xgb(X)
# make a standard partial dependence plot with a single SHAP value overlaid
fig,ax = shap.partial_dependence_plot(
"RM", model_xgb.predict, X, model_expected_value=True,
feature_expected_value=True, show=False, ice=False,
shap_values=shap_values_ebm[sample_ind:sample_ind+1,:]
)

shap.plots.scatter(shap_values_xgb[:,"RM"])

ロジスティック回帰
データを少し変えて、国勢調査のデータを使用します。
こちらのデータは年間$50K以上稼いでいるかどうかがTrueとFalseで設定されているデータ・セットになります。
https://archive.ics.uci.edu/ml/datasets/adult
まずデータを取得してモデルを学習します。
# a classic adult census dataset price dataset
X_adult,y_adult = shap.datasets.adult()
# a simple linear logistic model
model_adult = sklearn.linear_model.LogisticRegression(max_iter=10000)
model_adult.fit(X_adult, y_adult)
def model_adult_proba(x):
return model_adult.predict_proba(x)[:,1]
def model_adult_log_odds(x):
p = model_adult.predict_log_proba(x)
return p[:,1] - p[:,0]
各特徴量との関係を導出するためSHAP値を導出します。
background_adult = shap.maskers.Independent(X_adult, max_samples=100)
explainer = shap.Explainer(model_adult_proba, background_adult)
shap_values_adult = explainer(X_adult[:1000])
今までとは違い、線形の性質が出てきません。これは出力が確率のため、SHAPの計算の前提にある加算できる性質が考慮されていないためです。
shap.plots.scatter(shap_values_adult[:,"Age"])

そのため、対数に変えて見てみます。対数に変えて計算をしてみます。
# compute the SHAP values for the linear model
explainer_log_odds = shap.Explainer(model_adult_log_odds, background_adult)
shap_values_adult_log_odds = explainer_log_odds(X_adult[:1000])
対数にすると加算性が担保されるので線形な性質を確認できました。
shap.plots.scatter(shap_values_adult_log_odds[:,"Age"])

XGBoost詳細分析
先程の国勢調査のデータのデータを使用してより詳細に分析します。
同様にモデルの学習とSHAP値を導出します。
# train XGBoost model
model = xgboost.XGBClassifier(n_estimators=100, max_depth=2).fit(X_adult, y_adult)
# compute SHAP values
explainer = shap.Explainer(model, background_adult)
shap_values = explainer(X_adult)
# set a display version of the data to use for plotting (has string values)
shap_values.display_data = shap.datasets.adult(display=True)[0].values
平均SHAP値を導出します。年収が$50kを超えるのに関係のある特徴量はRelationshipのようです。
shap.plots.bar(shap_values)

最大SHAP値を導出します。年収が$50kを超えるのに関係のある特徴量はCaptial Gainのようです。これは当たり前の結果で示唆が得にないです。

各特徴量のSHAP値の分布と各特徴量の値を確認します。SHAP値と各特徴量には関連があるように見えます。
shap.plots.beeswarm(shap_values)

SHAP値の絶対値を出して、各特徴量の分布を見てみます。
shap.plots.beeswarm(shap_values.abs, color="shap_red")

SHAP値と”Age”の関係性を見てみます。
shap.plots.scatter(shap_values[:,"Age"])

”Age”と”Capital Gain”の関係も見てみます。
shap.plots.scatter(shap_values[:,"Age"], color=shap_values[:,"Capital Gain"])

各特徴量でどれがもっとも貢献度が高いかを導出できます。
まずSHAP値でクラスタリングします。
clustering = shap.utils.hclust(X_adult, y_adult)
可視化をしてみます。SHAP値に貢献する50%の特徴量は”Relationship”と”Material Status”になり、他の値は削除されます。
shap.plots.bar(shap_values, clustering=clustering)

SHAP値に貢献する80%の特徴量はclustering_cutoffに0.8を設定します。
“Relationship”と”Material Status”に”Education-Num”と”Occupation”に加えられています。
shap.plots.bar(shap_values, clustering=clustering, clustering_cutoff=0.8)

自然言語処理への適応
映画のレビューデータ・セットを使用して試します。
どの単語が予測の貢献度が高いかを可視化することができます。
! pip install transformers nlp
下記のコードでモデルの学習からSHAP値の導出も行っています。
import transformers
import nlp
import torch
import numpy as np
import scipy as sp
import shap
# load a BERT sentiment analysis model
tokenizer = transformers.DistilBertTokenizerFast.from_pretrained("distilbert-base-uncased")
model = transformers.DistilBertForSequenceClassification.from_pretrained(
"distilbert-base-uncased-finetuned-sst-2-english"
).cuda()
# define a prediction function
def f(x):
tv = torch.tensor([tokenizer.encode(v, padding='max_length', max_length=500, truncation=True) for v in x]).cuda()
outputs = model(tv)[0].detach().cpu().numpy()
scores = (np.exp(outputs).T / np.exp(outputs).sum(-1)).T
val = sp.special.logit(scores[:,1]) # use one vs rest logit units
return val
# build an explainer using a token masker
explainer = shap.Explainer(f, tokenizer)
# explain the model's predictions on IMDB reviews
imdb_train = nlp.load_dataset("imdb")["train"]
shap_values = explainer(imdb_train[:10], fixed_context=1)
文章中の各単語のSHAP値を導出しています。絶対が大きい単語が赤くなるようになっています。
shap.plots.text(shap_values[0])

各単語におけるSHAP値の絶対値の平均を導出します。’impressive’がもっとも貢献している単語になっています。
shap.plots.bar(shap_values.abs.mean(0))

各単語におけるSHAP値の合計値を導出します。これは文章によく使用される単語が高くなっているため、あまり意味のない分析になっているようです。
shap.plots.bar(shap_values.abs.sum(0))
