説明可能なAI④（SHAP値の散布図）

SHAP値の散布図

dependence_plotは、特定の説明変数とSHAP値の散布図で、相関関係を確認する場合に便利です。

前回記事で上位に表示された、LSTATをdependence_plotで表示してみます。

[Google Colaboratory]

1
2
3
4
shap.dependence_plot(ind="LSTAT",
                     interaction_index=None,
                     shap_values=shap_values,
                     features=X_test)

indオプションで確認したい説明変数(LSTAT)を指定しています。（1行目）

[実行結果]

横軸に説明変数、縦軸に同じ説明変数のSHAP値をプロットしています。

説明変数のSHAP値と、相関関係がみられるほど、予測への影響も高くなります。

上のグラフからはLSTATが低いほどSHAP値が高く、予測の結果に大きく影響を与えることが確認できます。


またinteraction_indexを指定することで、色に対して別の説明変数を指定することができます。

interaction_indexにRMを指定して実行してみます。（2行目）

[Google Colaboratory]

1
2
3
4
shap.dependence_plot(ind="LSTAT",
                     interaction_index="RM",
                     shap_values=shap_values,
                     features=X_test)

[実行結果]

LSTATが高くなると、RMが低くなる傾向が確認できます。

各説明変数の関係性や、どのように予測に影響しているかを確認する場合にdependence_plotはとても有効です。

説明可能なAI③（SHAPの可視化）

summary_plot(bar)

summary_plotを使うと、どの説明変数が大きく影響していたかを可視化することができます。

大局的に結果を確認する場合に便利です。

[Google Colaboratory]

1
2
3
4
shap.summary_plot(shap_values=shap_values,
                  features=X_test, 
                  plot_type="bar", 
                  max_display=5)

plot_typeに“bar”を指定することで、各説明変数を貢献度順に確認することができます。（3行目）

max_displayは上位項目の表示数で、今回は上位5項目まで表示しています。（4行目）

[実行結果]

横軸は平均SHAP値、縦軸は説明変数の項目になります。

縦軸の上位項目ほどモデルへの貢献度が高いことを表しており、今回のモデルではRM、LSTATの貢献度が高いことが確認できます。

summary_plot(dot)

次に、plot_typeに“dot”を指定して実行してみます。（3行目）

[Google Colaboratory]

1
2
3
4
shap.summary_plot(shap_values=shap_values,
                  features=X_test,
                  plot_type="dot",
                  max_display=5)

[実行結果]

各点がデータで、横軸にSHAP値、縦軸に説明変数の項目、色は説明変数の大小を表しています。

RMが高くなるほど予測値も高くなり、RMが低くなるほど予測値も低くなる傾向があります。

一方、LSTATが高ければ予測値は低くなり、LSTATが低ければ予測値は高くなるという傾向を視覚的に確認することができます。

説明可能なAI②（SHAP）

SHAPのインストール

SHAPライブラリのインストールを行います。

[Google Colaboratory]

1
!pip install shap

[実行結果]

正常にインストールができました。

SHAPモデルの作成

回帰モデルを引数にして、SHAPモデルを作成します。（2行目）

[Google Colaboratory]

1
2
3
import shap
explainer = shap.TreeExplainer(tree_reg)
explainer

前回記事で用意したのは決定木モデルなので、決定木モデル用のshap.TreeExplainerクラスを使っています。

このクラスは決定木以外にも、ランダムフォレストやXgBoostなどで利用できます。

[実行結果]

SHAPモデルが作成できました。

SHAP値の確認

SHAP値を確認します。

[Google Colaboratory]

1
2
shap_values = explainer.shap_values(X_test)
shap_values

[実行結果]

SHAP値は、入力したデータセットと同じ次元の要素数になり、値が大きいほど予測への影響が大きくなります。

行方向に見れば特定の予測に各説明変数がどれくらい貢献したかを確認できます。

列方向に見れば予測全体でその説明変数がどれくらい貢献したかを確認できます。


SHAP値には、グラフを描画する仕組みが用意されています。

次回は各グラフを表示してみます。

説明可能なAI①（SHAP）

SHAP

SHAPは、学習済みモデルにおいて各説明変数が予測値にどのような影響を与えたかを貢献度として定義して算出するモデルです。

各データごとに結果を出力し、可視化することができます。

前準備

前準備として、ボストンデータセットを用いた回帰モデル（決定木）を作成し、予測結果を確認します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
import pandas as pd
from sklearn.datasets import load_boston
from sklearn.model_selection import train_test_split
from sklearn.tree import DecisionTreeRegressor
 
boston = load_boston()
df = pd.DataFrame(boston.data,columns=boston.feature_names)
df["MEDV"] = boston.target
X= df[boston.feature_names]
y = df[["MEDV"]]
X_train, X_test, y_train, y_test = train_test_split(X, y,test_size=0.3,random_state=0)
 
print(len(X_train))
display(X_train.head(1))
print(len(X_test))
display(X_test.head(1))
 
tree_reg = DecisionTreeRegressor(max_depth=3, random_state=0).fit(X_train,y_train)

[実行結果]

以上で、回帰系の決定木モデルが作成できました。

重要度

作成したモデルの説明変数ごとに重要度を表示します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
import matplotlib.pyplot as plt
import numpy as np
 
features = X_train.columns
importances = tree_reg.feature_importances_
indices = np.argsort(importances)
 
plt.figure(figsize=(6,6))
plt.barh(range(len(indices)), importances[indices], color="b", align="center")
plt.yticks(range(len(indices)), features[indices])
plt.show()

feature_importances_を参照し、説明変数ごとの重要度を取得しています。（5行目）

[実行結果]

このモデルでは、RM、LSTATなどの重要度が高く、予測に強く影響していることが確認できます。

予測値の算出

予測値を算出します。

[Google Colaboratory]

1
2
3
X_test_pred = X_test.copy()
X_test_pred["pred"] = np.round(tree_reg.predict(X_test), 2)
X_test_pred.describe()[["RM","LSTAT","CRIM","DIS","PTRATIO","pred"]]

テストデータで予測値を算出し、結果を説明変数とマージしています。（2行目）

重要度の高い説明変数の上位5項目と、予測値(pred)を表示しています。（3行目）

[実行結果]

RMは8.7付近が最大値となっており、predの最大は50になっています。

予測値の表示

最も重要度の高かった説明変数であるRMでソートして、結果を確認してみます。

[Google Colaboratory]

1
X_test_pred.sort_values("RM")

[実行結果]

predの最大値が50だったので、RMが高いほどpredも高く、RMが低いほどpredも低く出ている傾向が見られます。

2番目に重要度が高かったLSTATはその逆で、LSTATが高いとpredは低く、LSTATが低ければpredが高くなっているようです。

このようにfeature_importances_は、モデル作成時にどのような説明変数が重要であるかを知るために大局的な指標となります。

一方SHAPは、作成したモデルの各説変数がどのように予測に寄与してしるかを知るための局所的な指標となります。

次回はSHAPを実装して予測結果を確認していきます。

分類モデルの評価⑦（まとめ）

今回は、これまで実行してきたいろいろな分類モデルをまとめて評価します。

PR曲線の可視処理を関数化

まず、前回実行したPR曲線の可視化処理を関数化します。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
def plot_pr_curve(y_true,proba):
    precision, recall, thresholds = precision_recall_curve(y_true, proba[:,0], pos_label=0)
    auc_score = auc(recall, precision)

    plt.figure(figsize=(12, 4))
    plt.subplot(1,2,1)

    plt.plot(recall, precision,label=f"PR Curve (AUC = {round(auc_score,3)})")
    plt.plot([0,1], [1,1], linestyle="--", color="red", label="Ideal Line")

    tg_thres = [0.3,0.5,0.8]
    for thres in tg_thres:
    tg_index = np.argmin(np.abs(thresholds - thres))
    plt.plot(recall[tg_index], precision[tg_index], marker = "o",markersize=10, label=f"Threshold = {thres}")

    plt.legend()
    plt.title("PR curve")
    plt.xlabel("Recall")
    plt.ylabel("Precision")
    plt.grid()

    plt.subplot(1,2,2)

    plt.plot(np.append(thresholds, 1), recall, label = "Recall")
    plt.plot(np.append(thresholds, 1), precision, label = "Precision")
    plt.xlabel("Thresholds")
    plt.ylabel("Score")
    plt.grid()
    plt.legend()

    plt.show()

スケーリング

決定木以外のアルゴリズムのために、データのスケーリングを行います。

[Google Colaboratory]

1
2
3
4
5
from sklearn.preprocessing import StandardScaler

scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

各モデルの定義

辞書型で各モデルの定義を行います。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
from sklearn.linear_model import LogisticRegression
from sklearn.svm import SVC
from sklearn.neighbors import KNeighborsClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier

models = {"Logistic Regression":LogisticRegression(), 
          "Linear SVM":SVC(kernel="linear",probability=True,random_state=0),
          "Kernel SVM":SVC(kernel="rbf",probability=True,random_state=0),
          "K Neighbors":KNeighborsClassifier(),
          "Decision Tree":DecisionTreeClassifier(max_depth=3,random_state=0),
          "Random Forest":RandomForestClassifier(max_depth=3,random_state=0)}

SVCクラスでpredict_probaを使用するために、probabilityをTrueにします。

[Google Colaboratory]

1
data_set = {"Train":[X_train_scaled,y_train],"Test":[X_test_scaled,y_test]}

事前準備は以上で終了です。

各モデルの評価

各モデルごとに以下の処理を行います。

  1. モデルの構築・学習（4行目）
  2. 予測（11行目）
  3. データセットごとに分類レポートを出力（13～16行目）
    output_dictをTrueにすることで辞書型で出力し、それをデータフレームにしています。
  4. テストデータに関してのPR曲線の可視化（20行目）

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
for model_name in models.keys():
    print()
    print(f"{model_name} Score Report")
    model = models[model_name].fit(X_train_scaled, y_train)

    for data_set_name in data_set.keys():

        X_data = data_set[data_set_name][0]
        y_true = data_set[data_set_name][1]

        y_pred = model.predict(X_data)

        score_df = pd.DataFrame(classification_report(y_true, y_pred, output_dict=True))
        score_df["model"] = model_name
        score_df["type"] = data_set_name
        display(score_df)

        if data_set_name == "Test":
            proba = model.predict_proba(X_data)
            plot_pr_curve(y_true, proba)

各モデルの結果は次のようになりました。

[ロジスティック回帰の評価結果]

[線形SVMの評価結果]

[カーネルSVMの評価結果]

[K近傍法の評価結果]

[決定木の評価結果]

[ランダムフォレストの評価結果]

スコアを見ると、各アルゴリズムの評価に大きな差はないようです。

決定木やロジスティック回帰のF1値が高く出ているので、今回のデーセットに対しては、シンプルなアルゴリズムでもうまく分類できていると言えます。

分類レポートのスコアは全て閾値が50%として出力されているので、PR曲線などを見て閾値を調整すると評価結果が変わります。

いろいろな閾値を試してみると良いかもしれません。

分類モデルの評価⑥（AUC）

AUC (Area Under the Curve)

AUCとは、PR曲線からモデルの性能を測るための指標です。

PR曲線より下の部分の面積を表す指標で、AUCが1に近いほどモデルの性能が良いことになります。

AUCの算出

AUCを算出するためには、aucメソッドに再現率と適合率を渡します。（2行目）

[Google Colaboratory]

1
2
from sklearn.metrics import auc
auc(recall, precision)

[実行結果]

算出されたAUCは0.95となりました。

閾値を切り口とした可視化

最適な閾値を確認しやすくするため、横軸に閾値を表示したグラフを作成します。

[Google Colaboratory]

1
2
3
4
5
6
7
plt.plot(np.append(thresholds, 1), recall, label = "Recall")
plt.plot(np.append(thresholds, 1), precision, label = "Precision")
plt.xlabel("Thresholds")
plt.ylabel("Score")
plt.grid()
plt.legend()
plt.show()

[実行結果]

閾値ごとに再現率と適合率の変化を確認しやすくなりました。

どの閾値を最適とするかは、どの程度の偽陽性や陽性未検出を許容できるかに依存しますので、状況に応じて閾値を設定することになります。

分類モデルの評価⑤（PR曲線）

PR曲線

PR曲線(Precision-Recall Curve)は、モデルの評価精度に使用されるとともに、最適な閾値を調べる時にも用いられる手法です。

PR曲線は縦軸に適合率、横軸に再現率の値をとり、閾値の変化による適合率と再現率のトレードオフ関係を表現します。

[トレードオフ関係]

  • 閾値を上げる
    陽性の判定をより厳しく行う。
    予測の正確性が上がる。
  • 閾値を下げる
    偽陽性が増える。
    予測の網羅性は上がる。

PR曲線の算出

PR曲線を引くためには適合率、再現率、閾値をそれぞれ取得する必要があります。

precision_recall_curveを使うとそれらの算出値をまとめて取得することができます。（4行目）

[Google Colaboratory]

1
2
3
4
5
6
7
8
from sklearn.metrics import precision_recall_curve
from sklearn.metrics import auc

precision, recall, thresholds = precision_recall_curve(y_test,pred_proba_test[:,0], pos_label=0)

print(precision[:3])
print(recall[:3])
print(thresholds[:3])

[実行結果]

各閾値(thresholds)における適合率(precision)と再現率(recall)を算出することができました。

PR曲線の可視化

matplotlibを使って、PR曲線の可視化を行います。（1行目）

参考として、３点閾値（30%、50%、80%）をプロットしています。（3～6行目）

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
plt.plot(recall, precision,label="PR Curve")

tg_thres = [0.3,0.5,0.8]
for thres in tg_thres:
    tg_index = np.argmin(np.abs(thresholds - thres))
    plt.plot(recall[tg_index], precision[tg_index], marker = "o",markersize=10, label=f"Threshold = {thres}")

plt.plot([0,1], [1,1], linestyle="--", color="red", label="Ideal Line")

plt.legend()
plt.title("PR curve")
plt.xlabel("Recall")
plt.ylabel("Precision")
plt.grid()
plt.show()

[実行結果]

閾値が上がるほど適合率が上がり、再現率が下がっていることが確認できます。

分類モデルの評価④（確信度）

前回記事までの分類予測では、陽性／陰性という最終的な予測結果のみを見てきましたが、「その予測がどれだけの確信度をもって行われたのか」という情報も重要になります。

確信度

predict_probaメソッドに予測対象データを渡すと、確率の形式で予測結果を出力することができます。

[Google Colaboratory]

1
2
3
4
5
pred_proba_train = rf_cls.predict_proba(X_train)
pred_proba_test = rf_cls.predict_proba(X_test)

print(pred_proba_train[:5])
print(pred_proba_test[:5])

[実行結果]

出力された配列には、陽性に分類される確率と陰性に分類される確率がそれぞれ格納されています。

デフォルトでは、50%を上回っているカテゴリを最終的な予測結果とするようになっていますが、このように確率を算出することで60%以上であれば陽性とするなど、新たな閾値を設定することも可能です。

確信度の可視化

領域ごとの確信度を可視化してみます。

[Google Colaboratory]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
import numpy as np

x_min, x_max = X["mean radius"].min() - 0.5, X["mean radius"].max() + 0.5
y_min, y_max = X["mean texture"].min() - 0.5, X["mean texture"].max() + 0.5

step = 0.5
x_range = np.arange(x_min, x_max, step)
y_range = np.arange(y_min, y_max, step)
xx, yy = np.meshgrid(x_range, y_range)

Z = rf_cls.predict_proba(np.c_[xx.ravel(), yy.ravel()])[:, 0]
Z = Z.reshape(xx.shape)

plt.figure(figsize=(10,7))
plt.contourf(xx, yy, Z, alpha=0.8,cmap=plt.cm.coolwarm)
plt.colorbar()
plt.scatter(X_train["mean radius"], X_train["mean texture"], c=y_train,marker="o", edgecolors="k", cmap=plt.cm.coolwarm_r, label="Train")
plt.scatter(X_test["mean radius"], X_test["mean texture"], c=y_test, marker="^",edgecolors="k", cmap=plt.cm.coolwarm_r, label="Test")
plt.xlabel("mean radius")
plt.ylabel("mean texture")
plt.legend()
plt.show()

[実行結果]

背景のメッシュ部分が陽性に分類される確信度を表しています。

各変数（mean／radius／mean texture）の最小値から最大値の範囲内で、それぞれ0.5刻みでデータを新たに生成し、それぞれのデータを説明変数としてpredict_probaで算出した確信度を可視化しています。

確信度を可視化することで、どの領域で誤分類が起こりやすいかを把握することができます。

分類モデルの評価③（分類レポート）

前回までは、正解率や適合率、再現率、F1値などの評価指標を個別に確認してきましたが、scikit-learnにはそれらの指標をまとめて出力してくれる関数が用意されています。

分類レポート

classification_reportを使って、分類レポートを出力します。

[Google Colaboratory]

1
2
3
4
5
6
from sklearn.metrics import classification_report

print("Train Score Report")
print(classification_report(y_train,y_train_pred))
print("Test Score Report")
print(classification_report(y_test,y_test_pred))

[実行結果]

悪性を陽性としたケースだけでなく、良性を陽性としたケースでの各評価指標も算出しています。

また、カテゴリごとの指標のマクロ平均(macro avg)や、サンプル数(support)に応じて重み付けを行った平均(weighted avg)も出力されています。

分類レポートは、モデル精度の全体像を確認するのに非常に有効ですので、分類の評価フェーズで積極的に活用していきましょう。

分類モデルの評価②（適合率／再現率／F1値）

適合率

適合モデルは、モデルがどれだけ正確に予測できいるかを示す指標です。

適合率の算出には、precision_scoreを使用します。

[Google Colaboratory]

1
2
3
4
from sklearn.metrics import precision_score

print(f"訓練データ適合率：{precision_score(y_train,y_train_pred,pos_label=0)}")
print(f"テストデータ適合率：{precision_score(y_test,y_test_pred,pos_label=0)}")

適合率の算出に使用する分母は陽性と予測したサンプル数であり、偽陽性件数を減らすことが適合率の向上につながります。

また混同行列の上部合計に対して、正解がどのくらいだったを示しています。

[実行結果]

テストデータの適合率は87.93%となっています。

前回記事の混同行列を見ると真陽性（左上）の５１件と、偽陽性（右上）の７件を足すと５８件となりこれが分母となります。

５８件のうち正しく予測できた５１件を分子として計算すると、今回のテストデータ適合率0.879と一致します。

再現率

再現率は、モデルがどれだけ網羅的に予測できているかを示す指標です。

再現率の算出には、recall_scoreを使用します。

[Google Colaboratory]

1
2
3
4
from sklearn.metrics import recall_score

print(f"訓練データ再現率：{recall_score(y_train,y_train_pred,pos_label=0)}")
print(f"テストデータ再現率：{recall_score(y_test,y_test_pred,pos_label=0)}")

再現率の算出に使用する分母は実際に陽性であるサンプル数であり、偽陽性件数は考慮されません。

偽陽性が増えても多くの悪性を検出したい場合は、再現性を重視した方がよいことになります。

また混同行列の左部合計に対して、正解がどのくらいだったを示しています。

[実行結果]

混同行列を見ると真陽性（左上）の５１件と、偽陰性（左下）の１２件を足すと６３件となりこれが分母となります。

６３件のうち正しく予測できた５１件を分子として計算すると、今回のテストデータ再現率0.809と一致します。

再現率(80.9%)は適合率(87.9%)と比較して低い結果となっています。

F1値

F1値は適合率と再現率 両方のバランスを考慮した指標です。

適合率と再現率の調和平均をとったものがF1値となります。

F1値の算出には、f1_scoreを使用します。

[Google Colaboratory]

1
2
3
4
from sklearn.metrics import f1_score

print(f"訓練データF1値：{f1_score(y_train,y_train_pred,pos_label=0)}")
print(f"テストデータF1値：{f1_score(y_test,y_test_pred,pos_label=0)}")

[実行結果]

適合率と再現率の間をとったような結果(84.29%)になりました。

F1値は、分類モデルの総合的な評価指標として使用されることが多いです。

ただし運用で使えるモデルかどうかを測るためには、F1値だけではなく、適合率や再現率などの個別指標にも注意する必要があります。