🐍
KernelSHAPをフルスクラッチで実装しながら理解する

2024/12/03に公開
 はじめに一般的に、機械学習モデルにおいてはモデルの複雑さと解釈性はトレードオフの関係にあるとされています。これは、モデルが複雑化するにつれて、「なぜそのような予測結果を出力したのか」を人間が解釈することが難しくなるためです。モデルの解釈が難しい状況は、人間にとって望ましくない影響を引き起こす場合があります。
たとえば、2018年に Amazon の AI 採用ツールが女性に対する偏見を示したために利用停止となったこと[1]が報じられました。また、再犯率を予測する機械学習モデルが黒人に対する人種的偏見を示した[2]ケースも存在します。
このように、機械学習モデルは時として人間にとって好ましくないバイアスを学習し、誤った判断を下すことがあります。この問題に対応するためには、開発者がモデルを十分に解釈し、その動作を理解することが不可欠です。
このような背景から、近年では機械学習モデルの解釈に関する研究が活発に行われています。
本記事では、機械学習モデルの解釈手法の中でも代表的な SHAP (SHapley Additive exPlanations)について、関連論文を読み進めながら簡易的な実装をフルスクラッチで行った記録をまとめていきます。
本記事はあくまで学習の記録を目的としているので、元の論文や実装等を参照しながら読み進めていただければ幸いです。
また、内容に誤りがあれば、後学のためにコメントをいただけますと幸いです。

 Shapley ValueSHAPは、協力ゲーム理論におけるShapley Valueを機械学習に援用することで、各予測結果に対する各特徴量の寄与度を示す手法です。
ここではまず、SHAPのベースとなるShapley Valueについて説明します。

 Shapley Valueとはあるゲームがあり、複数人で参加でき、なんらかの報酬を獲得できるとします。
このゲームに複数人で参加した場合、報酬をプレイヤー全員に均等に分配すると、貢献度の大きいプレイヤーが不利益を被り、貢献度の小さいプレイヤーが利益を得る状況が発生してしまう可能性があります。
そのため、各プレイヤーの貢献度に応じて、公正に報酬を分配する仕組みが求められます。
しかし、プレイヤーそれぞれの貢献度をどのように算出すればよいでしょうか？
ここで、プレイヤーの貢献度を定量的に算出する方法の一つとして Shapley Value があります。Shapley Valueでは、それぞれのプレイヤーがゲームに参加した場合、参加しなかった場合と比較して、どれだけ報酬が増加したかを考慮して貢献度を算出します。
Shapley Value は以下の式で定義されます。


\phi_i = \sum_{S \subseteq N \backslash \{i\}}P_r(S)(\upsilon (S\cup\{i\}) - \upsilon (S))
S: 集合Nの部分集合(プレイヤーiは含まれない)

\upsilon (S): チームで協力してゲームを行ったときの報酬

\upsilon (S\cup\{i\}): にプレイヤーが加わった場合の報酬

P_r(S): チームSが組まれる確率

具体的に下の例について、実際にShapley Valueを計算してみます。
プレイヤーはA、B、Cの3人。

A、B、Cが一人で得られる利益は、それぞれ2、3、5であるとする。

また、

AとBが提携したときに得られる利益は15であり、

BとCが提携したときに得られる利益は10であり、

AとCが提携したときに得られる利益は20であり、

A、B、Cの3人が提携したときに得られる利益は30であるとする。

このとき、A、B、CそれぞれのShapley Valueを求めたい。
このとき、プレイヤーAのShapley Valueはこのように計算できます。

\begin{aligned}
\phi_A &= \frac{1}{3}(2-0) + \frac{1}{6}(15-3) + \frac{1}{6}(20-5) + \frac{1}{3}(30-10) \\
       &= \frac{71}{6} \\
       &\simeq 11.83
\end{aligned}
同様にプレイヤーB, CについてもShapley Valueを求めると下のようになり、貢献度はA→C→Bの順に高いといえます。


プレイヤー
Shapley Value


A
11.83

B
7.33

C
10.83


 Shapley Valueの機械学習への応用Shapley Valueにおけるゲームおよび報酬、プレイヤーを下のように読み替えることで機械学習に援用することができます。
任意の入力に対する機械学習モデルの予測をゲーム、予測値の期待値と任意の予測値の差を報酬、任意の入力における各特徴量の値をプレイヤーとして扱い、考えうる全ての特徴量の組み合わせに対する予測についてその平均値を取ることで、各特徴量の貢献度を算出します。
このとき、2つの特徴量X_1とX_2を持つモデルがあり、予測モデルfに対してそれぞれの特徴量のShapley Value、\phi_{X_1}, \phi_{X_2}を計算すると、特徴量X_1が予測に強く貢献している場合、\phi_{X_1}は大きくなり、特徴量X_2の貢献が小さい場合、\phi_{X_2}は小さくなります。
このように、Shapley Valueを機械学習に援用することで、モデルの特定の予測結果に対して、どの特徴量がどのような影響を与えているかを直感的に理解することが可能になります。

 SHAP
 Shapley Valueの課題Shapley Valueを機械学習に援用する上での課題として計算量の多さがあります。

Shapley Valueを計算するためには、すべての特徴量の組み合わせに対して評価を行う必要があり、特徴量の数がnである場合、Shrapley Valueを計算するために必要な特徴量の組み合わせは、2^n通り存在します。したがって、Shapley Valueの計算量は特徴量の数に対して指数的に増加します。
これを避けるために、考案された手法がSHAPです。

 SHAPSHAPの基本的なアイデアは、モデルの予測値f(x)を以下のように各特徴量の貢献度の和として分解することです。Lundbergらは、これを望ましい3つの性質とともにAdditive Feature Attribution Methodとしてまとめています。(Additive Feature Attribution Methodsについてはこちらの論文を参照してください。)


f(x) = \phi_0 + \sum_{i=1}^n \phi_i

f(x): モデルの予測値

\phi_0: すべての特徴量が存在しない場合のモデルの出力

\phi_i: 特徴量iのSHAP Value

SHAPでは、各特徴量の貢献度をShapley Valueの考え方を用いることで算出しており、特徴量 i のSHAP Value \phi_i は、以下の式で定義されます。


\phi_i = \sum_{S \subseteq N \setminus \{i\}} \frac{|S|! \, (|N| - |S| - 1)!}{|N|!} \left[ f(S \cup \{i\}) - f(S) \right]

N = \{1, 2, \dots, n\}: すべての特徴量の集合

S \subseteq N \setminus \{i\}: 特徴量 i を除いた部分集合

f(S): 部分集合 S の特徴量を用いたモデルの予測値

|S|: 部分集合 S の要素数

この式では、すべての可能な順序で特徴量がモデルに追加される際の期待される寄与を計算していますが、SHAPでは、なんらかの近似アルゴリズムを使用してSHAP Valueを効率的に計算します。
近似方法には様々ありますが、すべてのモデルに対して適用可能であり、広く利用されている手法としてKernel SHAPがあります。

Kernel SHAPでは、特徴量の部分集合をランダムにサンプリングすることで、Shapley Valueを近似的に計算し、計算コストを大幅に削減します。これにより、すべての特徴量の組み合わせを列挙する必要がなくなり、現実的な計算時間で特徴量の寄与を評価することが可能になります。

 Kernel SHAPKernel SHAPでは、Shapley Valueの算出を重み付き線形回帰の問題として定式化します。

g(Z^{\prime}) = \phi_0 + \sum_{i=1}^{M}\phi_iz_i^{\prime}
ここで、重み\pi_xは以下のように定義されます。

\pi_x(z^{\prime}) = \frac{(M - 1)}{\binom{M}{|z^{\prime}|} |z^{\prime}| (M - |z^{\prime}|)}
重み\pi_xを用いて、線形モデルgを損失関数Lで最適化することによって学習します。

L(f, g, \pi_x) = \sum_{z^{\prime} \in Z} \left[ f(h_x(z^{\prime})) - g(z^{\prime}) \right]^2 \pi_x(z^{\prime})
このとき、線形モデルgの係数\phi_iがSHAP Valuesとなります。

 Kernel SHAPの実装ここまでで、Kernel SHAPの理論的な部分について触れてきました。
ここからは、実際にKernel SHAPをフルスクラッチで実装していきます。

 モデルの作成まずは、SHAPを使って解釈するモデルを訓練します。
データはscikit-learnのmake_regression()を用いて生成し、モデルはSVMを用います。
from sklearn.datasets import make_regression
from sklearn.svm import SVR

# 乱数シードの固定
np.random.seed(42)

# データの生成
X, y = make_regression(n_samples=1000, n_features=5, noise=0.1, random_state=42)

# モデルの訓練
svr_model = SVR(kernel='rbf')
svr_model.fit(X, y)
また、訓練したモデルにデータを入力して、予測値を返す関数を定義しておきます。
def model(x: np.ndarray) -> np.ndarray:
    """
    予測関数 f(x)

    Parameters:
    - x (numpy.ndarray): 入力データ行列

    Returns:
    - y_pred (numpy.ndarray): モデルの予測値
    """

    return svr_model.predict(x)

 部分集合のサンプリング特徴量の組み合わせをサンプリングします。
ここでは、特徴量の数を入力としてサンプリングされた特徴量の組み合わせの行列Zを生成する関数を定義します。
def sample_coalitions(self, M: int) -> np.ndarray:
    """
    特徴量の組み合わせをサンプリングする

    Parameters:
    - M (int): 特徴量の総数

    Returns:
    - Z (np.ndarray): サンプリングされた特徴量の組み合わせ行列
    """

    num_samples = self.num_samples

    # 各特徴量を含む(1)か欠損(0)かをランダムに決定
    Z = np.random.randint(2, size=(num_samples, M))

    return Z

 データの生成次に、サンプリングした特徴量の組み合わせをもとに説明モデルgに入力するデータx^{\prime}を生成します。
x_baselineはベースライン値と呼ばれる欠損した特徴量を補完するために使用される値です。
通常特徴量の平均値や0が用いられます。
def generate_data(self, Z: np.ndarray, x: np.ndarray) -> np.ndarray:
    """
    サンプリングした特徴量の組み合わせに基づいてデータを生成する

    Parameters:
    - Z (np.ndarray): 特徴量の組み合わせ行列
    - x (np.ndarray): 解析対象のデータ点

    Returns:
    - X (np.ndarray): 生成されたデータ行列
    """

    x_baseline = self.x_baseline
    X = Z * x + (1 - Z) * x_baseline

    return X

 重みの計算Zを用いて重み\pi_xを算出します。
def kernel_shap_weight(self, Z: np.ndarray, M: int) -> np.ndarray:
    """
    各特徴量の組み合わせに対して重みを計算する

    Parameters:
    - Z (numpy.ndarray): 特徴量の組み合わせ行列
    - M (int): 特徴量の総数

    Returns:
    - weights (numpy.ndarray): 各組み合わせに対応する重み
    """

    Nz = np.sum(Z, axis=1)
    M = float(M)
    weights = (M - 1) / (comb(M, Nz) * Nz * (M - Nz))

    return weights

 SHAP Valueを計算SHAP Valueを算出します。
generate_data()によって生成されたデータをモデル関数に入力することでサンプリングされた組み合わせに対する予測値を得ます。
得られた予測値とkernel_shap_weight()によって算出した重みを用いて重みつき線形回帰を解いた上で、その係数をSHAP Valueとして出力します。
def explain(self, x: np.ndarray) -> np.ndarray:
    """
    入力データ x に対するSHAP Valueを計算する

    Parameters:
    - x (numpy.ndarray): 解析対象のデータ点

    Returns:
    - phi (numpy.ndarray): 各特徴量のSHAP Value
    """

    M = x.shape[0]
    Z = self.sample_coalitions(M)
    X = self.generate_data(Z, x)
    y = self.model(X)
    weights = self.kernel_shap_weight(Z, M)
    X_reg = np.concatenate([np.ones((Z.shape[0], 1)), Z], axis=1)

    # 重み付き線形回帰を解く
    W_sqrt = np.sqrt(weights)
    X_weighted = X_reg * W_sqrt[:, np.newaxis]
    y_weighted = y * W_sqrt
    theta = np.linalg.lstsq(X_weighted, y_weighted, rcond=None)[0]

    # Shapley値は theta の 1 行目以降
    phi = theta[1:]

    return phi

 実装の全体像ここまでの実装を全て組み合わせると以下のようになります。
class KernelSHAP:
    def __init__(self, model, x_baseline: np.ndarray, num_samples: int):
        """
        Kernel SHAPクラスの初期化

        Parameters:
        - model (callable): 入力としてデータの配列を受け取り予測値の配列を返す関数
        - x_baseline (numpy.ndarray): ベースライン値(Ex. 特徴量の平均値やゼロ), 欠損した特徴量を補完するために使用
        - num_samples (int, optional): サンプリングする特徴量の組み合わせ数
        """

        self.model = model
        self.x_baseline = x_baseline
        self.num_samples = num_samples

    def sample_coalitions(self, M: int) -> np.ndarray:
        """
        特徴量の組み合わせ z' をサンプリングする

        Parameters:
        - M (int): 特徴量の総数

        Returns:
        - Z (numpy.ndarray): サンプリングされた特徴量の組み合わせ行列
        """

        num_samples = self.num_samples
        # 各特徴量を含む(1)か欠損(0)かをランダムに決定
        Z = np.random.randint(2, size=(num_samples, M))
        # 全ての特徴量が欠損または全て存在する場合を除外
        Z = Z[(Z.sum(axis=1) != 0) & (Z.sum(axis=1) != M)]

        return Z

    def generate_data(self, Z: np.ndarray, x: np.ndarray) -> np.ndarray:
        """
        サンプリングした特徴量の組み合わせに基づいてデータ x' を生成する

        Parameters:
        - Z (numpy.ndarray): 特徴量の組み合わせ行列
        - x (numpy.ndarray): 解析対象のデータ点

        Returns:
        - X (numpy.ndarray): 生成されたデータ行列
        """
        x_baseline = self.x_baseline
        X = Z * x + (1 - Z) * x_baseline
        return X

    def kernel_shap_weight(self, Z: np.ndarray, M: int) -> np.ndarray:
        """
        各特徴量の組み合わせに対して重み π_x(z') を計算する

        Parameters:
        - Z (numpy.ndarray): 特徴量の組み合わせ行列
        - M (int): 特徴量の総数

        Returns:
        - weights (numpy.ndarray): 各組み合わせに対応する重み
        """
        Nz = np.sum(Z, axis=1)
        M = float(M)
        weights = (M - 1) / (comb(M, Nz) * Nz * (M - Nz))
        return weights

    def explain(self, x: np.ndarray) -> np.ndarray:
        """
        入力データ x に対するSHAP Valueを計算する

        Parameters:
        - x (numpy.ndarray): 解析対象のデータ点

        Returns:
        - phi (numpy.ndarray): 各特徴量のSHAP Value
        """
        M = x.shape[0]
        Z = self.sample_coalitions(M)
        X = self.generate_data(Z, x)
        y = self.model(X)
        weights = self.kernel_shap_weight(Z, M)
        # 線形回帰のためにデザイン行列を作成
        X_reg = np.concatenate([np.ones((Z.shape[0], 1)), Z], axis=1)
        # 重み付き線形回帰を解く
        W_sqrt = np.sqrt(weights)
        X_weighted = X_reg * W_sqrt[:, np.newaxis]
        y_weighted = y * W_sqrt
        theta = np.linalg.lstsq(X_weighted, y_weighted, rcond=None)[0]
        # Shapley値は theta の 1 行目以降
        phi = theta[1:]
        return phi

 モデルの解釈とSHAPパッケージとの比較ここでは、上で作成したクラスとSHAPパッケージの出力を比較してみます.
まず、作成したクラスをモデルに対して適用します。
# ベースライン値の定義
x_baseline = np.mean(X, axis=0)

# 解析対象のデータ点の選択
# ここではデータセットの最初のデータ点を使用
x = X[0]

# Kernel SHAPのインスタンスを作成
kernel_shap = KernelSHAP(model, x_baseline, num_samples=100)

# SHAP Value値を計算
shap_values_custom = kernel_shap.explain(x)
# [ 32.710083    22.85739521   5.64838539 -29.21839368   6.46725837]
次にSHAPパッケージをモデルに対して適用します。
# ベースライン値を2次元配列に変換
x_baseline_shap = x_baseline.reshape(1, -1)
x_shap = x.reshape(1, -1)

# KernelExplainerのインスタンスを作成
explainer = shap.KernelExplainer(model, x_baseline_shap)

# SHAP Valueを計算
shap_values_shap = explainer.shap_values(x_shap, nsamples=1000)
# [ 32.84550755  23.16058231   6.20604045 -28.56450713   6.82698873]
出力されたSHAP Valueを比較してみます。
青色のバーが今回実装したクラスの出力、オレンジ色のバーがSHAPパッケージによる出力を示しています。

SHAPパッケージと近しい値が得られていることが分かります。
わずかな差分は、サンプリング方法の違いや線形回帰における正則化の有無などによるものであるかと思われます。

 おわりに本記事のコードはこちらから確認できます。
実装を通じて、計算量を削減するための工夫が図られているものの、実際の分析にSHAPを適用する場合は、次元数の少ないデータセット対して適用するなど限定的な状況に限られるという印象を受けました。
SHAPに関しては、KernelSHAP以外にもTreeSHAP[3]やDeepSHAPなどの特定のモデルに対してのみ適用可能な方法や、Shapley Flow[4]などのSHAPを拡張した手法なども存在するため、それらについても理解を深めていければと思います。
ここまで読んでいただきありがとうございました。

 参考文献https://arxiv.org/abs/1705.07874
https://github.com/shap/shap
https://hacarus.github.io/interpretable-ml-book-ja/shap.html
https://www.rand.org/content/dam/rand/pubs/papers/2021/P295.pdf

脚注
https://www.reuters.com/article/us-amazon-com-jobs-automation-insight/amazon-scraps-secret-ai-recruiting-tool-that-showed-bias-against-women-idUSKCN1MK08G/ ↩︎
https://www.propublica.org/article/machine-bias-risk-assessments-in-criminal-sentencing ↩︎
https://www.nature.com/articles/s42256-019-0138-9 ↩︎
https://proceedings.mlr.press/v130/wang21b/wang21b.pdf ↩︎
はじめに

Shapley Value

Shapley Valueとは

Shapley Valueの機械学習への応用

SHAP

Shapley Valueの課題

SHAP

Kernel SHAP

Kernel SHAPの実装

モデルの作成

部分集合のサンプリング

データの生成

重みの計算

SHAP Valueを計算

実装の全体像

モデルの解釈とSHAPパッケージとの比較

おわりに

参考文献

Discussion
