可解释AI简介：Shapley值
这是介绍如何使用Shapley值解释机器学习模型的入门教程。Shapley值是合作博弈论中广泛使用的方法，具有理想的特性。本教程旨在帮助建立对如何计算和解释基于Shapley的机器学习模型解释的坚实理解。我们将采用实践性的动手方法，使用``shap`` Python包逐步解释更复杂的模型。这是一个不断更新的文档，并作为``shap`` Python包的介绍。因此，如果您有反馈或贡献，请打开一个问题或拉取请求，以使本教程更好！
大纲
解释线性回归模型
解释广义加性回归模型
解释一个非加性提升树模型
解释线性逻辑回归模型
解释一个非加性提升树逻辑回归模型
处理相关输入特征
解释一个变压器NLP模型
解释线性回归模型
在使用Shapley值解释复杂模型之前，了解它们在简单模型中的工作原理是有帮助的。最简单的模型类型之一是标准线性回归，因此我们在下面的`加利福尼亚住房数据集 <https://www.dcc.fc.up.pt/~ltorgo/Regression/cal_housing.html>`__上训练了一个线性回归模型。该数据集包含1990年加利福尼亚州的20,640个房屋区块，我们的目标是根据8个不同的特征预测房屋中位数价格的自然对数：
MedInc - 街区组中的收入中位数
HouseAge - 街区组中的房屋中位年龄
AveRooms - 每个家庭的平均房间数
AveBedrms - 每户家庭的平均卧室数量
人口 - 街区组人口
AveOccup - 家庭成员平均数量
纬度 - 街区组纬度
经度 - 街区组经度
[1]:
import sklearn
import shap
# a classic housing price dataset
X, y = shap.datasets.california(n_points=1000)
X100 = shap.utils.sample(X, 100) # 100 instances for use as the background distribution
# a simple linear model
model = sklearn.linear_model.LinearRegression()
model.fit(X, y)
[1]:
LinearRegression()
检查模型系数
理解线性模型最常见的方式是检查为每个特征学习的系数。这些系数告诉我们当我们改变每个输入特征时，模型输出会变化多少：
[2]:
print("Model coefficients:\n")
for i in range(X.shape[1]):
print(X.columns[i], "=", model.coef_[i].round(5))
Model coefficients:
MedInc = 0.45769
HouseAge = 0.01153
AveRooms = -0.12529
AveBedrms = 1.04053
Population = 5e-05
AveOccup = -0.29795
Latitude = -0.41204
Longitude = -0.40125
虽然系数在告诉我们改变输入特征值时会发生什么方面非常有用，但它们本身并不是衡量特征整体重要性的好方法。这是因为每个系数的值取决于输入特征的尺度。例如，如果我们用分钟而不是年份来衡量房屋的年龄，那么HouseAge特征的系数将变为0.0115 / (365∗24∗60) = 2.18e-8。显然，自房屋建成以来的年数并不比分钟数更重要，但其系数值却大得多。这意味着系数的幅度不一定是在线性模型中衡量特征重要性的好方法。
使用部分依赖图的更完整画面
要理解模型中某个特征的重要性，有必要了解该特征的变化如何影响模型的输出，以及该特征值的分布情况。为了对线性模型进行可视化，我们可以构建一个经典的偏依赖图，并在x轴上显示特征值的分布作为直方图：
[3]:
shap.partial_dependence_plot(
"MedInc",
model.predict,
X100,
ice=False,
model_expected_value=True,
feature_expected_value=True,
)
上图中的灰色水平线表示模型应用于加利福尼亚住房数据集时的预期值。垂直的灰色线表示中位收入特征的平均值。注意，蓝色的部分依赖图线（当我们固定中位收入特征为某个给定值时，模型输出的平均值）总是通过两条灰色预期值线的交点。我们可以将这个交点视为部分依赖图相对于数据分布的“中心”。当我们接下来转向Shapley值时，这种中心化的影响将变得清晰。
从部分依赖图中读取SHAP值
基于Shapley值的机器学习模型解释的核心思想是利用合作博弈论中的公平分配结果，将模型输出 \(f(x)\) 的信用分配给其输入特征。为了将博弈论与机器学习模型联系起来，需要将模型的输入特征与游戏中的玩家匹配，并将模型函数与游戏规则匹配。由于在博弈论中，玩家可以选择加入或不加入游戏，因此我们需要一种方式让特征“加入”或“不加入”模型。定义特征“加入”模型的最常见方式是，当我们知道该特征的值时，就说该特征“加入了模型”，而当我们不知道该特征的值时，就说该特征没有加入模型。为了在只有特征子集 \(S\) 是模型的一部分时评估现有模型 \(f\)，我们使用条件期望值公式对其他特征进行积分。这种公式可以有两种形式：
或
在第一种形式中，我们知道S中特征的值，因为我们*观察*它们。在第二种形式中，我们知道S中特征的值，因为我们*设定*它们。一般来说，第二种形式通常更可取，这不仅因为它告诉我们如果我们干预并改变模型的输入，模型将如何表现，而且因为它更容易计算。在本教程中，我们将完全关注第二种表述。我们还将使用更具体的术语“SHAP值”来指代应用于机器学习模型条件期望函数的Shapley值。
SHAP 值的计算可能非常复杂（一般来说它们是 NP-难的），但线性模型非常简单，我们可以直接从部分依赖图中读取 SHAP 值。当我们解释一个预测 \(f(x)\) 时，特定特征 \(i\) 的 SHAP 值就是模型输出的期望值与特征值 \(x_i\) 处的部分依赖图之间的差异：
[4]:
# compute the SHAP values for the linear model
explainer = shap.Explainer(model.predict, X100)
shap_values = explainer(X)
# make a standard partial dependence plot
sample_ind = 20
shap.partial_dependence_plot(
"MedInc",
model.predict,
X100,
model_expected_value=True,
feature_expected_value=True,
ice=False,
shap_values=shap_values[sample_ind : sample_ind + 1, :],
)
经典部分依赖图与SHAP值之间的紧密对应关系意味着，如果我们绘制整个数据集中特定特征的SHAP值，我们将精确地描绘出该特征的部分依赖图的均值中心化版本：
[5]:
shap.plots.scatter(shap_values[:, "MedInc"])
Shapley值的加性特性
Shapley值的一个基本性质是，它们总是总和为所有玩家在场时的游戏结果与没有玩家在场时的游戏结果之间的差异。对于机器学习模型来说，这意味着所有输入特征的SHAP值将总是总和为基线（预期）模型输出与当前模型输出之间的差异，用于解释预测。最直观的方式是通过瀑布图来观察这一点，瀑布图从房价的背景先验期望 \(E[f(X)]\) 开始，然后逐个添加特征，直到达到当前模型输出 \(f(x)\)：
[6]:
# the waterfall_plot shows how we get from shap_values.base_values to model.predict(X)[sample_ind]
shap.plots.waterfall(shap_values[sample_ind], max_display=14)
解释一个加性回归模型
线性模型的部分依赖图与SHAP值有如此紧密联系的原因在于，模型中的每个特征都是独立于其他特征处理的（效果只是简单相加）。我们可以在放松直线线性要求的同时保持这种加性特性。这就产生了广为人知的广义加性模型（GAMs）类。虽然有许多方法可以训练这些类型的模型（如将XGBoost模型设置为深度1），但我们使用InterpretML的可解释增强机器，这些机器是专门为此设计的。
[7]:
# fit a GAM model to the data
import interpret.glassbox
model_ebm = interpret.glassbox.ExplainableBoostingRegressor(interactions=0)
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(
"MedInc",
model_ebm.predict,
X100,
model_expected_value=True,
feature_expected_value=True,
show=False,
ice=False,
shap_values=shap_values_ebm[sample_ind : sample_ind + 1, :],
)
[8]:
shap.plots.scatter(shap_values_ebm[:, "MedInc"])
[9]:
# the waterfall_plot shows how we get from explainer.expected_value to model.predict(X)[sample_ind]
shap.plots.waterfall(shap_values_ebm[sample_ind])
[10]:
# the waterfall_plot shows how we get from explainer.expected_value to model.predict(X)[sample_ind]
shap.plots.beeswarm(shap_values_ebm)
## 解释一个非加性增强树模型
[11]:
# 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(
"MedInc",
model_xgb.predict,
X100,
model_expected_value=True,
feature_expected_value=True,
show=False,
ice=False,
shap_values=shap_values_xgb[sample_ind : sample_ind + 1, :],
)
[12]:
shap.plots.scatter(shap_values_xgb[:, "MedInc"])
[13]:
shap.plots.scatter(shap_values_xgb[:, "MedInc"], color=shap_values)
## 解释线性逻辑回归模型
[14]:
# 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]
需要注意的是，解释线性逻辑回归模型的概率在输入中并不是线性的。
[15]:
# make a standard partial dependence plot
sample_ind = 18
fig, ax = shap.partial_dependence_plot(
"Capital Gain",
model_adult_proba,
X_adult,
model_expected_value=True,
feature_expected_value=True,
show=False,
ice=False,
)
如果我们使用SHAP来解释线性逻辑回归模型的概率，我们会看到强烈的交互效应。这是因为线性逻辑回归模型在概率空间中不是加性的。
[16]:
# compute the SHAP values for the linear model
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])
Permutation explainer: 1001it [00:58, 14.39it/s]
[17]:
shap.plots.scatter(shap_values_adult[:, "Age"])
如果我们解释模型的对数几率输出，我们会看到模型输入和模型输出之间存在完美的线性关系。记住你所解释的模型的单位是什么，以及解释不同的模型输出可能会导致对模型行为的截然不同的看法，这一点很重要。
[18]:
# 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])
Permutation explainer: 1001it [01:01, 13.61it/s]
[19]:
shap.plots.scatter(shap_values_adult_log_odds[:, "Age"])
[20]:
# make a standard partial dependence plot
sample_ind = 18
fig, ax = shap.partial_dependence_plot(
"Age",
model_adult_log_odds,
X_adult,
model_expected_value=True,
feature_expected_value=True,
show=False,
ice=False,
)
## 解释非加性提升树逻辑回归模型
[21]:
# train XGBoost model
model = xgboost.XGBClassifier(n_estimators=100, max_depth=2).fit(
X_adult, y_adult * 1, eval_metric="logloss"
)
# 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
The use of label encoder in XGBClassifier is deprecated and will be removed in a future release. To remove this warning, do the following: 1) Pass option use_label_encoder=False when constructing XGBClassifier object; and 2) Encode your labels (y) as integers starting with 0, i.e. 0, 1, 2, ..., [num_class - 1].
98%|===================| 31839/32561 [00:12<00:00]
默认情况下，SHAP 条形图将取数据集中所有实例（行）上每个特征的平均绝对值。
[22]:
shap.plots.bar(shap_values)
但是平均绝对值并不是创建特征重要性全局度量的唯一方法，我们可以使用任意数量的变换。这里我们展示了如何使用最大绝对值来突出资本收益和资本损失特征，因为它们具有不频繁但影响幅度大的特点。
[23]:
shap.plots.bar(shap_values.abs.max(0))
如果我们愿意处理稍微复杂一些的情况，我们可以使用蜂群图来总结每个特征的SHAP值的整个分布。
[24]:
shap.plots.beeswarm(shap_values)
通过取绝对值并使用纯色，我们在这两种图表的复杂性之间取得了折中：条形图和全蜂群图。请注意，上面的条形图仅是下面蜂群图中显示的值的汇总统计数据。
[25]:
shap.plots.beeswarm(shap_values.abs, color="shap_red")
[26]:
shap.plots.heatmap(shap_values[:1000])
[27]:
shap.plots.scatter(shap_values[:, "Age"])
[28]:
shap.plots.scatter(shap_values[:, "Age"], color=shap_values)
[29]:
shap.plots.scatter(shap_values[:, "Age"], color=shap_values[:, "Capital Gain"])
[30]:
shap.plots.scatter(shap_values[:, "Relationship"], color=shap_values)
## 处理相关特征
[31]:
clustering = shap.utils.hclust(X_adult, y_adult)
[32]:
shap.plots.bar(shap_values, clustering=clustering)
[33]:
shap.plots.bar(shap_values, clustering=clustering, clustering_cutoff=0.8)
[34]:
shap.plots.bar(shap_values, clustering=clustering, clustering_cutoff=1.8)
## 解释一个变压器NLP模型
这展示了如何将SHAP应用于具有高度结构化输入的复杂模型类型。
[35]:
import datasets
import numpy as np
import scipy as sp
import torch
import transformers
# 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 = datasets.load_dataset("imdb")["train"]
shap_values = explainer(imdb_train[:10], fixed_context=1, batch_size=2)
2022-06-15 14:43:09.022292: I tensorflow/core/platform/cpu_feature_guard.cc:142] This TensorFlow binary is optimized with oneAPI Deep Neural Network Library (oneDNN) to use the following CPU instructions in performance-critical operations: AVX2 AVX512F FMA
To enable them in other operations, rebuild TensorFlow with the appropriate compiler flags.
2022-06-15 14:43:09.731330: I tensorflow/core/common_runtime/gpu/gpu_device.cc:1510] Created device /job:localhost/replica:0/task:0/device:GPU:0 with 8395 MB memory: -> device: 0, name: NVIDIA GeForce RTX 2080 Ti, pci bus id: 0000:15:00.0, compute capability: 7.5
2022-06-15 14:43:09.732184: I tensorflow/core/common_runtime/gpu/gpu_device.cc:1510] Created device /job:localhost/replica:0/task:0/device:GPU:1 with 9631 MB memory: -> device: 1, name: NVIDIA GeForce RTX 2080 Ti, pci bus id: 0000:21:00.0, compute capability: 7.5
Reusing dataset imdb (/home/slundberg/.cache/huggingface/datasets/imdb/plain_text/1.0.0/2fdd8b9bcadd6e7055e742a706876ba43f19faee861df134affd7a3f60fc38a1)
Partition explainer: 9it [00:25, 2.80s/it]Token indices sequence length is longer than the specified maximum sequence length for this model (720 > 512). Running this sequence through the model will result in indexing errors
Partition explainer: 11it [00:33, 4.21s/it]
[36]:
# plot a sentence's explanation
shap.plots.text(shap_values[2])
[37]:
shap.plots.bar(shap_values.abs.mean(0))
Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray
[38]:
shap.plots.bar(shap_values.abs.sum(0))
有更多有用示例的想法吗？我们鼓励提交增加此文档笔记本的拉取请求！