利用SHAP解释Xgboost模型

Xgboost相对于线性模型在进行预测时往往有更好的精度，但是同时也失去了线性模型的可解释性。所以Xgboost通常被认为是黑箱模型。

2017年，Lundberg和Lee的论文提出了SHAP值这一广泛适用的方法用来解释各种模型（分类以及回归），其中最大的受益者莫过于之前难以被理解的黑箱模型，如xgboost和神经网络模型。

本教程中，我们在真实数据集上进行实操，利用SHAP来解释Xgboost模型。

预计学习用时：30分钟。

本教程基于Python 3.6版本、Xgboost 0.82版本以及shap 0.28.5版本。

原创者：东布东 | 修改校对：SofaSofa TeamM |

1. Feature importance

在SHAP被广泛使用之前，我们通常用feature importance或者partial dependence plot来解释xgboost。 feature importance是用来衡量数据集中每个特征的重要性。

简单来说，每个特征对于提升整个模型的预测能力的贡献程度就是特征的重要性。（拓展阅读：随机森林、xgboost中feature importance，Partial Dependence Plot是什么意思？，怎么利用permutation importance来解释xgboost模型）

Feature importance可以直观地反映出特征的重要性，看出哪些特征对最终的模型影响较大。但是无法判断特征与最终预测结果的关系是如何的。

下面这个例子中，我们用2018年足球球员身价数据（请在SofaSofa数据竞赛页面进行数据下载，下载解压后只需要train.csv这个文件）来具体阐述。

# 加载模块
import xgboost as xgb
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt; plt.style.use('seaborn')

# 读取数据，目标变量y是球员的身价（万欧元）
data = pd.read_csv('train.csv')

# 获得当时球员年龄
today = pd.to_datetime('2018-01-01')
data['birth_date'] = pd.to_datetime(data['birth_date'])
data['age'] = np.round((today - data['birth_date']).apply(lambda x: x.days) / 365., 1)

# 选择特征，这里只是举例，未必是最佳组合
# 特征依次为身高（厘米）、潜力、速度、射门、传球、带球、防守、体格、国际知名度、年龄
cols = ['height_cm', 'potential', 'pac', 'sho', 'pas', 'dri', 'def', 'phy', 'international_reputation', 'age']

# 训练xgboost回归模型
model = xgb.XGBRegressor(max_depth=4, learning_rate=0.05, n_estimators=150)
model.fit(data[cols], data['y'].values)

# 获取feature importance
plt.figure(figsize=(15, 5))
plt.bar(range(len(cols)), model.feature_importances_)
plt.xticks(range(len(cols)), cols, rotation=-45, fontsize=14)
plt.title('Feature importance', fontsize=14)
plt.show()

上图中，我们可以看出国际知名度、潜力和年龄是影响球员身价最重要的三个因素。但是这些因素和身价是正相关、负相关还是其他更复杂的相关性，我们无法从上图得知。我们也无法解读每个特征对每个个体的预测值的影响。

2. SHAP value

SHAP的名称来源于SHapley Additive exPlanation。

Shapley value起源于合作博弈论。比如说甲乙丙丁四个工人一起打工，甲和乙完成了价值100元的工件，甲、乙、丙完成了价值120元的工件，乙、丙、丁完成了价值150元的工件，甲、丁完成了价值90元的工件，那么该如何公平、合理地分配这四个人的工钱呢？Shapley提出了一个合理的计算方法（有兴趣地可以查看原论文），我们称每个参与者分配到的数额为Shapley value。

SHAP是由Shapley value启发的可加性解释模型。对于每个预测样本，模型都产生一个预测值，SHAP value就是该样本中每个特征所分配到的数值。 假设第ii个样本为xixi，第ii个样本的第jj个特征为xi,jxi,j，模型对第ii个样本的预测值为yiyi，整个模型的基线（通常是所有样本的目标变量的均值）为ybaseybase，那么SHAP value服从以下等式。

其中f(xi,1)f(xi,1)为xi,jxi,j的SHAP值。直观上看，f(xi,1)f(xi,1)就是第ii个样本中第1个特征对最终预测值yiyi的贡献值，当f(xi,1)>0f(xi,1)>0，说明该特征提升了预测值，也正向作用；反之，说明该特征使得预测值降低，有反作用。

很明显可以看出，与上一节中feature importance相比，SHAP value最大的优势是SHAP能对于反映出每一个样本中的特征的影响力，而且还表现出影响的正负性。

3. SHAP的Python实现

Python中SHAP值的计算由shap这个package实现，可以通过pip install shap安装。

下面我们针对第1节中训练出的模型model，计算其SHAP值。

引用package并且获得解释器explainer。

import shap
# model是在第1节中训练的模型
explainer = shap.TreeExplainer(model)

获取训练集data各个样本各个特征的SHAP值。

因为data中有10441个样本以及10个特征，我们得到的shap_values的维度是10441×1010441×10。

shap_values = explainer.shap_values(data[cols])
print(shap_values.shape)

(10441, 10)

我们也可以获得在第2节中提到的模型的基线ybaseybase。

通过对比发现，我们可以确认基线值就是训练集的目标变量的拟合值的均值。在这里例子中，目标变量是球员的身价（万欧元），也就是球员的平均身价为229万欧元。

y_base = explainer.expected_value
print(y_base)

data['pred'] = model.predict(data[cols])
print(data['pred'].mean())

229.168
229.168

3.1 单个样本的SHAP值

我们可以随机检查其中一位球员身价的预测值以及其特征对预测值的影响。

下面的数据框中第一列是特征名称，第二列是特征的数值，第三列是各个特征在该样本中对应的SHAP值。

# 比如我们挑选数据集中的第30位
j = 30
player_explainer = pd.DataFrame()
player_explainer['feature'] = cols
player_explainer['feature_value'] = data[cols].iloc[j].values
player_explainer['shap_value'] = shap_values[j]
player_explainer

我们知道一个样本中各特征SHAP值的和加上基线值应该等于该样本的预测值。

我们可以做如下的验证。

print('y_base + sum_of_shap_values: %.2f'%(y_base + player_explainer['shap_value'].sum()))
print('y_pred: %.2f'%(data['pred'].iloc[j]))

y_base + sum_of_shap_values: 1332.65
y_pred: 1332.65

shap还提供极其强大的数据可视化功能。下图是对上面数据框的可视化。

蓝色表示该特征的贡献是负数，红色则表示该特征的贡献是正数。最长的红色条是潜力值，球员的潜力值很高，而他的身价也因此增加了1092万；最长的蓝色条是年龄，这个球员年龄较小才20岁出头，尚未到职业巅峰，未来也有诸多不确定性，身价也因此降低了180万元。

shap.initjs()
shap.force_plot(explainer.expected_value, shap_values[j], data[cols].iloc[j])

3.2 对特征的总体分析

除了能对单个样本的SHAP值进行可视化之外，还能对特征进行整体的可视化。

下图中每一行代表一个特征，横坐标为SHAP值。一个点代表一个样本，颜色越红说明特征本身数值越大，颜色越蓝说明特征本身数值越小。

我们可以直观地看出潜力potential是一个很重要的特征，而且基本上是与身价成正相关的。年龄age也会明显影响身价，蓝色点主要集中在SHAP小于0的区域，可见年纪小会降低身价估值，另一方面如果年纪很大，也会降低估值，甚至降低得更明显，因为age这一行最左端的点基本上都是红色的。

shap.summary_plot(shap_values, data[cols])

我们也可以把一个特征对目标变量影响程度的绝对值的均值作为这个特征的重要性。

因为SHAP和feature_importance的计算方法不同，所以我们这里也得到了与第1节不同的重要性排序。

shap.summary_plot(shap_values, data[cols], plot_type="bar")

3.3 部分依赖图Partial Dependence Plot

SHAP也提供了部分依赖图的功能，与传统的部分依赖图不同的是，这里纵坐标不是目标变量y的数值而是SHAP值。

比如下图中，年纪大概呈现出金字塔分布，也就是24到31岁这个年纪对球员的身价是拉抬作用，小于24以及大于31岁的球员身价则会被年纪所累。

shap.dependence_plot('age', shap_values, data[cols], interaction_index=None, show=False)

3.4 对多个变量的交互进行分析

我们也可以多个变量的交互作用进行分析。一种方式是采用summary_plot描绘出散点图，如下：

shap_interaction_values = shap.TreeExplainer(model).shap_interaction_values(data[cols])
shap.summary_plot(shap_interaction_values, data[cols], max_display=4)

我们也可以用dependence_plot描绘两个变量交互下变量对目标值的影响。

shap.dependence_plot('potential', shap_values, data[cols], interaction_index='international_reputation', show=False)