SHAPの論文を読んでみた

# はじめに
SHAPについて調べてみました。

## 参考文献

1. SHAP を提案した論文 <https://arxiv.org/abs/1705.07874>
- shapley 値の解説 <https://dropout009.hatenablog.com/entry/2019/11/20/091450>
- Tree Ensembles に対する shap の近似的な求め方 <https://arxiv.org/abs/1802.03888>
- shap ライブラリ <https://github.com/slundberg/shap>
- shap ライブラリの使い方 <https://orizuru.io/blog/machine-learning/shap/>
- 教科書 <https://christophm.github.io/interpretable-ml-book/shap.html>
- LIMEの論文 <https://www.kdd.org/kdd2016/papers/files/rfp0573-ribeiroA.pdf>
- LIMEの論文要約 <https://www.dendoron.com/boards/40>

参考文献に示した1の論文を中心とし、6と7と8を参考にしてSHAPが何かを説明します。1の論文は用語の説明が抽象化されているので、7の論文の具体例を見た方が分かりやすいです。

## 概要
1. どんなもの？
予測を解釈するための統合フレームワークであるSHAP（SHapley Additive exPlanations）値を提案。

2. 先行研究と比べてどこがすごいの？
特徴量重要度を定義する手法は多数存在するが、手法間の関連性はよく分かっていなかった。
本論文で提案した手法はある理論に基づいて唯一に定まる解を求めることができ、また、新しいクラスは既存の6つの手法を統合した。

3. 技術や手法の"キモ"はどこにある？
（1）加法的特徴重要度尺度の新しいクラスを同定した。
（2）このクラスはユニークな解をもち、望ましいいくつかの特性を持っていることを示した。

## 本論

説明したいモデルを$f(\boldsymbol{x})$とします。

```math
y = f(\boldsymbol{x})
```

説明可能モデルを$g(\boldsymbol{x'})$とします。今から示す手法では、データ全体に対する特徴量の重要度ではなく、各入力$\boldsymbol{x}$に対し、出力である予測$f(\boldsymbol{x})$を説明するための特徴量の重要度を求めます。このような手法をLocal methodsと呼びます。例えばLIMEはLocal methodsです。

説明可能モデルでは、しばしばマップ関数を用いて元のモデルの入力$\boldsymbol{x}$を、説明可能モデルへの入力$\boldsymbol{x'}$へと簡素化します。

```math
\boldsymbol{x} = h_{\boldsymbol{x}}(\boldsymbol{x}')
```
ここでマップ関数は添字$\boldsymbol{x}$があるように、入力値に依存する関数です。この簡素化によって、元のモデルの入力$\boldsymbol{x}$のどの部分が予測に寄与しているかを解釈できるような入力形式$\boldsymbol{x'}$に変換しています。例えば入力がBag of Words (BOW)の場合で、どの単語が予測に寄与しているかを知りたい場合は入力$\boldsymbol{x}$における単語の有無だけを入力とすれば良いです。したがって簡素化は、非ゼロの頻度の値を1に変換する操作になります。また、簡素化された後の入力$\boldsymbol{x'}$は、$h_{\boldsymbol{x}}$で元の入力に戻すことができます。これは$h_{\boldsymbol{x}}$が局所的な入力$\boldsymbol{x}$の情報を持っているからです。

```math
簡素化: [1, 2, 1, 1, 2, 1, 1, 0, 0, 0] \rightarrow [1, 1, 1, 1, 1, 1, 1, 0, 0, 0] \\
h_{\boldsymbol{x}}: [1, 1, 1, 1, 1, 1, 1, 0, 0, 0] \rightarrow [1, 2, 1, 1, 2, 1, 1, 0, 0, 0] = \boldsymbol{x} \odot \boldsymbol{x'}
```
入力が画像で、画像のどの部分が予測に寄与しているか知りたい場合は、隣接している連続部分をパッチ化して、そのパッチの有無を1/0で表現します。$h_{\boldsymbol{x}}$で元の入力空間に戻すときはパッチ部分(1の部分)は元の画像に変換すればよいです。パッチ以外の部分(0の部分)は適当な背景色で戻すのですが、これについては後に説明します。
この辺りの説明は参考文献8のブログが非常に参考になります。

Local methodsは、$\boldsymbol{x}' \simeq \boldsymbol{z}'$となるときに$g(\boldsymbol{z'})\simeq f(h_{\boldsymbol{x}}(\boldsymbol{x}'))$となるように$g$の重みを決定します。

1￥論文では、Additive feature attribution methods と呼ばれる方法で説明可能モデルを構築します。簡素化はバイナリ変数への変換に限定します。したがって説明可能モデルはバイナリ変数の線形関数となります。

```math
g(z') = \phi_0 + \sum_{i=1}^M\phi_i z_i'
```

ここで$z'\in \\{0,1\\}^M$は簡素化された特徴で、$\phi_i\in\mathbb{R}$は特徴の効果を表す係数です。このモデルは、各特徴量の効果が$\phi_i$で表され、すべての特徴量の効果を足し合わせることで元のモデルの出力を近似します。

次の条件を満たすように係数$\phi_i$をユニークに決められます。

- Property 1 (Local accuracy)

  ```math
  \begin{aligned}
  f(\boldsymbol{x}) = g(\boldsymbol{x}') &= \phi_0 + \sum_{i=1}^M\phi_i x_i' \\
  \boldsymbol{x} &= h_\boldsymbol{x}(\boldsymbol{x}') \\
  \end{aligned}
  ```

- Property 2 (Missingness)

  ```math
  x' = 0 \Rightarrow\phi_i=0
  ```

- Property 3 (Consistency)

  $f_x(\boldsymbol{z}') = f(h_\boldsymbol{x}(\boldsymbol{z}'))$と略記し、$z_i’=0$であることを$\boldsymbol{z}'\backslash i$と書くことにする。任意の二つのモデル$f, f'$に対して、

  ```math
  f_\boldsymbol{x}'(\boldsymbol{z}') -   f_\boldsymbol{x}'(\boldsymbol{z}'\backslash i)
  \geq f_\boldsymbol{x}(\boldsymbol{z}') -   f_\boldsymbol{x}(\boldsymbol{z}'\backslash i) \ \ \mathrm{for \ all \ inputs } \ \boldsymbol{z}' \in \{0,1\}^M
  ```

  が成り立つとき、$\phi_i(f', \boldsymbol{x}) \geq \phi_i(f, \boldsymbol{x})$が成り立つ。
これは、特徴量$i$の限界寄与が（他の特徴に関係なく）増加または同じままになるようにモデルが変更された場合、$\phi_i$も増加または同じままになることを要請しています。

以上の三つの条件を課すと、説明可能モデルが唯一に定まることが知られています。

```math
\phi_i(f, \boldsymbol{x}) = \sum_{\boldsymbol{z}'\subseteq \boldsymbol{x}'}\frac{|\boldsymbol{z}'|!(M-|\boldsymbol{z}'|-1)!}{M!}[f_\boldsymbol{x}(\boldsymbol{z}')-f_\boldsymbol{x}(\boldsymbol{z}'\backslash i)]
```

ここで$|\boldsymbol{z}'|$は$1$である成分の数、$\boldsymbol{z}' \subseteq \boldsymbol{x}$は$\boldsymbol{x}'$の非ゼロ成分の部分集合であるすべての$\boldsymbol{z'}$ベクトルを表しています。これは combined cooperative game theory の結果から得られるもので、&\phi_i&は下記のShapley Values と呼ばれるものに一致しています。

```math
\phi_i = \sum_{\mathcal{S}\subseteq \mathcal{N}\backslash i}\frac{1}{\frac{N!}{S!(N-S-1)!}}\left( v\left(\mathcal{S}\cup\{i\}\right) - v\left(\mathcal{S}\right) \right)
```

ここで$S=|\mathcal {S}|$、$N=|\mathcal {N}|$ 、$v$は効用関数で、規格化定数は$\mathcal{N} \backslash i$の並び替えのパターン数です。

以上の結果より、所与の$h_\boldsymbol{x}$に対してAdditive feature attribution methodsは唯一に定まることが分かりました。逆に、Shapley values に基づかない指標は、特徴量重要度が持つべき3つの条件を持たないため、望ましくないと考えられます。

## SHAP (SHapley Additive exPlanation) Values

特徴の重要性の統一的な尺度としてSHAP値を考えます。これはShapley valuesの式において、元のモデル$f$を条件付き期待関数に近似したものです。

```math
f_\boldsymbol{x}(\boldsymbol{z}') \equiv f(h_\boldsymbol{x}(\boldsymbol{z}')) \rightarrow E[f(\boldsymbol{z})|z_S] \\
```
$E[f(\boldsymbol{z})|z_S]$の意味ですが、まず、添字の$S$は$\boldsymbol{z}'$の非ゼロ要素のインデックスを表します。$z_S$は、確率変数$\boldsymbol{z}$に対し、インデックス$S$に対応する成分を入力値$x$の値で固定することを意味します。例えば$\boldsymbol{z}'=[1,1,0,0,...]$であれば、$S=\\{1,2\\}$となり、$E[f(\boldsymbol{z})|z_S]=E[f(\boldsymbol{z})|z_{1,2}]=E[f(\boldsymbol{z})|z_1=x_1, z_2=x_2]$となります。

期待関数への置き換えが必要な理由を説明します。マッピング関数$h_\boldsymbol{x}(\boldsymbol{z}') = z_S$は厳密にはインデックス$S$に含まれない特徴量を欠損させる操作です。というのも、説明可能モデルの入力空間から元の特徴量空間に写像する際に、インデックス$S$に含まれない部分を勝手に他の値に置き換えてしまうと意味を持ってしまう可能性があるからです。例えば画像で考えると、パッチに対応しない部分を全て０で置き換えるべきでしょうか？それとも255で置き換えるべきでしょうか？。一方で重回帰分析の場合は特徴量の欠損は、モデルからその特徴量の項を除くことが特徴量を0とすることと等価なので0に置き換えれば良いですが、ほとんどのモデルはこのように欠損値に自然に対応できません。そこで、マッピング関数はインデックス$S$に含まれない部分を欠損値にするのではなく、データの分布に従ってランダムに値を割り振ることにします。これによってモデルにデータを入力させることが出来ます。これによってモデルの予測値$f$は$z_S$を条件に持つ期待関数になります。

SHapley 値を定義通り計算するとものすごい計算コストがかかるため、計算は困難です。（$\sum_{\boldsymbol{z}'\subseteq \boldsymbol{x}'}$が指数オーダー）。しかし、Additive feature attribution methodsの洞察を組み合わせることで、SHAP値を近似することができます。近似する際には特徴量の独立性およびモデルの線形性を利用します。

```math
\begin{aligned}
E[f(\boldsymbol{z})|z_S] =& E[f(\boldsymbol{z})|z_S]\\
=& E_{\bar z_S}[f(z_S, \bar z_S)] \because \mathrm{特徴量が独立} \\
=& f(z_S, E[\bar z_S]) \because \mathrm{モデルが線形}
\end{aligned}
```
論文中ではSHAPの近似結果として六つの手法が示され、すなわち、これらの手法が統合されたことを意味しています。

- モデルに依存しない近似手法
    - 「Shapley sampling values」
    - 「Kernel SHAP」（新規手法）
- モデル固有の近似手法
    - 「Linear SHAP」
    - 「Low-Order SHAP」
    - 「Max SHAP」（新規手法）
    - 「Deep SHAP」（新規手法）

### Shapley sampling values
論文は[こちら](https://link.springer.com/article/10.1007/s10115-013-0679-x)。有料のため断念。

### Kernel SHAP

Kernel SHAP は 5 つのステップで求められます。

- $z_k'\in \\{0,1\\}^M$, $k\in\{1,...,K\}$(1 = feature present in
  coalition, 0 = feature absent).
- $h_x(z_k')=x_k\cdot z_k'$で特徴を元の空間に戻して、モデル$f$で予測値を得る.
- それぞれの$z_k'$に対して SHAP kernel で weight を求める.
- Fit weighted linear model.
- Return Shapley values $\phi_k$, the coefficients from the linear model.

SHAP kernel は次の式で定義されます。

```math
\pi_\boldsymbol{x}(\boldsymbol{z}') = \frac{M-1}{_MC_{|\boldsymbol{z}'|}|\boldsymbol{z}'|(M-|\boldsymbol{z}'|)}
```

重み付けした損失関数で最適化問題を解きます。

```math
L(f, g, \pi_\boldsymbol{x}) = \sum_{\boldsymbol{z}'\in Z}\left[ f(h_\boldsymbol{x}(\boldsymbol{z}')) - g(\boldsymbol{z}') \right]^2 \pi_\boldsymbol{x}(\boldsymbol{z}')
```

### Linear SHAP

線形モデル$f(\boldsymbol{x})=\sum_{j=1}w_j x_j+b$において、特徴量が独立であれば、SHAP値は次のようになります。導出は[この論文](https://arxiv.org/pdf/1903.10464.pdf)のAppendixに詳細があります。

```math
\begin{aligned}
\phi_0(f, \boldsymbol{x}) &= b \\
\phi_i(f, \boldsymbol{x}) &= w_i(x_i - E[x_i])
\end{aligned}
```

### Low-Order SHAP
（何言ってるのかよくわからん）
線形回帰の場合でもKernel SHAPを使用すると$\mathcal O(2^M+M^3)$だけ計算が必要なので、特徴量を独立にする近似をした方が効率的です（この近似がLow-Order SHAP？）。

### Max SHAP
付属資料があると書いてあるのですが、無い？のでわかりません。査読者には間違っているとか書かれてますね。

### Deep SHAP
Dense層やMaxpool層に対するSHAPの計算ができると主張しており、バックプロパゲーションでSHAPを計算できるということです（非線形の部分で出来るのでしょうか？）。

# おわりに
[決定木ではSHAPの効率的な計算ができる](https://arxiv.org/abs/1802.03888)ようです。LightGBM、Xgboostでは参考文献4のライブラリで簡単にSHAPを出力できます。