R与机器学习系列|15.可解释的机器学习算法（Interpretable Machine Learning）（上）

苹果醋3

于 2024-09-26 19:04:10 发布

阅读量578

点赞数 8

分类专栏：面试题汇总与解析文章标签： nginx 运维 java spring boot mysql

本文链接：https://blog.csdn.net/pingguocu3/article/details/142569713

版权

面试题汇总与解析专栏收录该内容

30 篇文章 0 订阅

订阅专栏

Machine learning is changing the world!

Machine learning is changing the world!

在之前的章节中（见 个人微信公众号连载），我们学习了如何训练多种不同形式的高级机器学习模型。 由于机器学习算法模型的内部结构的负责性，它们通常被认为是“黑箱”模型。然而，正是由于它们的复杂性，它们通常可以更准确地预测非线性、微弱或一些罕见的现象。不幸的是，以上过程通常以丧失可解释性为代价。幸运的是，近年来在机器学习算法模型的可解释性方面（可解释性机器学习，Interpretable Machine Learning）已经取得了一定进展，以帮助解释机器学习模型。本章将重点介绍常见的可解释性机器学习（IML）算法。

IML的想法来源

在机器学习任务中，仅仅确定一个能够优化预测性能的机器学习模型是远远不够的。具有可解释性及可信性是一个性能较好的模型的标志，也是别人采用我们模型的前提条件。当我们在应用和嵌入越来越复杂的预测模型和机器学习算法时，不论是我们模型开发者还是模型应用者，都需要一些方法来解释和理解模型的结果，以保证模型可以应用于相关决策。

现在，可解释性的引入使我们能够从最先进的机器学习模型中提取关键见解和可操作信息。这些进展使我们能够回答以下问题：

哪些最重要的客户属性影响其行为？
这些属性与行为之间的关系如何？
多个属性是否相互作用，导致不同客户之间的不同行为？
我们为什么预计某个客户会做出特定的决策？
我们基于预测结果做出的决策是否公平可靠？

解释模型的方法可以根据回答上述示例问题的范围，大致分为提供全局或局部解释两种方式。全面了解在全局范围内训练的整个模型相当重要，同时也要放大到数据或预测的局部区域，并得出解释。能够回答这些问题并提供两个层次的解释对于任何机器学习项目的成功和接受、采用、嵌入和正确利用至关重要。

全局解释（Global interpretation）

全局解释性是基于对其特征的整体观点以及它们如何影响底层模型结构的整体观点来理解该模型如何做出预测。它回答了哪些特征相对具有影响力，这些特征如何影响响应变量，以及可能发生的潜在交互作用。全局模型解释性有助于理解响应变量与各个特征（或其子集）之间的关系。可以说，全面的全局模型解释性在实践中很难实现。任何超过少数特征的模型将难以完全理解，因为我们无法一次理解整个模型结构。
虽然全局模型解释性通常难以实现，但我们有机会至少可以在模块级别上理解一些模型。这通常围绕着了解哪些特征最具有影响力（通过特征重要性），然后关注最具影响力的变量如何驱动模型输出（通过特征效应）。尽管我们可能无法完全理解具有数百个特征的模型，但实际上，只有十几个变量对于驱动模型的性能起到了真正的重要作用。而且，我们有可能对这十几个变量如何影响模型的性能有很清晰的理解。

局部解释（Local interpretation）

全局解释性方法帮助我们理解输入与响应变量的整体关系，但在某些情况下可能会产生误导（例如，在出现强烈的交互作用时）。尽管某个特征可能会影响整体模型的预测准确性，但并不意味着该特征对于特定观测值或者一组观测值的预测值具有最大的影响。局部解释帮助我们理解哪些特征影响了特定观测值的预测响应。局部解释不仅可以帮助我们回答客户可能会做什么，还可以解释模型对于特定观测值的具体预测原因。
局部解释有三种主要方法：
局部可解释模型无关解释（Local interpretable model-agnostic explanations， LIME）
Shapley值
局部逐步过程
以上三种方法下文会有系统介绍。当然，这些方法的目标是相同的：解释哪些变量在预测一组观测值上最具影响力。

模型特定 vs. 模型无关（Model-specific vs. model-agnostic）

还有一点很重要，那就是解释模型的方法可以分为模型特定和模型无关两种。在前面的章节中，我们看到的许多特征重要性理解方法都是模型特定的。例如，在线性模型中，我们可以使用t统计量的绝对值作为衡量特征重要性的指标（不过当线性模型涉及交互项和转换时，情况会变得复杂）。另一方面，随机森林可以记录数据的OOB部分的预测准确性，然后在对每个预测变量进行置换后重复同样的操作，并计算两个准确性之间的差异，然后在所有树中进行平均，并通过标准误差进行归一化。这些模型特定的解释工具局限于各自的模型类别。使用模型特定的方法也有一些优点，因为它们与模型性能更紧密相关，可能能够更准确地纳入预测变量之间的相关性结构。然而，这些方法也存在一些缺点。例如，许多机器学习算法（例如堆叠集成算法）没有固有的方法来衡量特征重要性。
此外，在不同模型之间比较特征的模型特定重要性非常困难，因为这种情况下是在比较不同的测量指标（例如，在线性模型中是(t)-统计量的大小，而在随机森林中是预测准确度的降低）。在模型无关的方法中，模型被视为一个“黑箱”。将可解释性与具体模型分离使我们能够轻松地在不同模型之间比较特征的重要性。

置换特征重要性

概念

在之前的章节中，我们展示了一些针对特征重要性的模型特定方法（例如，对于线性模型，我们使用了(t)-统计量的绝对值）。然而，对于支持向量机（SVM）来说，我们必须依赖于一种模型无关的方法，该方法是基于布雷曼（Breiman，2001）提出的用于随机森林的排列特征重要性测量方法，并由费舍尔（Fisher）、鲁丁（Rudin）和多米尼奇（Dominici）（2018）进一步扩展。

置换特征重要性是一种衡量特征对模型预测的重要性的方法。它的基本思想是在训练数据中对某个特征进行随机置换，然后观察模型预测性能的变化。如果某个特征对模型的预测结果有重要影响，那么在置换该特征后，模型的性能会明显下降。
具体步骤如下：
训练模型：首先，使用原始的训练数据来训练机器学习模型。
计算基准性能：在训练完成后，使用模型对测试数据进行预测，并计算模型的性能指标（如准确率、均方根误差等），作为基准性能。
置换特征：接下来，随机选择一个特征，并在训练数据中将该特征的值进行随机置换。
计算置换性能：使用置换后的数据进行模型预测，并计算性能指标。
计算特征重要性：将置换后的性能与基准性能进行比较，通过计算性能的变化来衡量特征的重要性。如果性能变化较大，说明该特征对模型预测有重要影响，其重要性较高；反之，如果性能变化较小，说明该特征对模型预测影响较小，其重要性较低。
通过反复执行以上步骤，可以得到每个特征的重要性排序，从而帮助理解模型对不同特征的依赖程度。这种方法在解释复杂模型的预测结果以及选择最重要特征方面具有重要的应用价值。

实现

置换特征重要性计算可以通过DALEX、iml和vip包来实现，每个包都有其独特的优势。
iml包提供了FeatureImp()函数，使用置换法计算一般预测模型的特征重要性。FeatureImp()函数允许用户指定通用损失函数或从预定义列表中选择（例如，loss = "mse"表示均方误差）。用户还可以指定重要性是作为原始模型误差与置换后的模型误差之差还是比值来衡量。用户还可以指定在对每个特征进行排列时使用的重复次数，以帮助稳定过程中的变异性。
DALEX包也通过variable_importance()函数提供置换特征重要性分数。与iml::FeatureImp()类似，该函数允许用户指定损失函数以及如何计算重要性分数（例如，使用差异或比值）。variable_importance()函数还提供了一个选项，在对数据计算重要性之前对训练数据进行抽样（默认使用n_sample = 1000），这可以加快计算速度。
vip包专门关注变量重要性图（VIPs），并提供多种模型特定和模型无关的的方法来计算变量重要性，包括置换法。使用vip包，我们可以使用自定义的损失函数（或从预定义列表中选择），执行蒙特卡洛模拟来稳定过程，对特征置换之前对观测值进行抽样，以及并行执行计算，这可以加快处理大数据集的运行时间，等等。
以下示例通过vip执行置换特征重要性计算。为了加快执行速度，我们对训练数据进行了30%的抽样，但重复模拟了3次以增加我们估计的稳定性

# 加载依赖包
library(dplyr)      #数据操纵
library(ggplot2)    # 可视化

# Modeling packages
library(h2o)       # H2O
library(recipes)   # 机器学习蓝图
library(rsample)   # 数据分割
library(xgboost)   # 拟合GBM模型

# 模型可解释性包
library(pdp)       # 偏依赖图及ICE曲线绘制
library(vip)       # 变量重要性VIP图
library(iml)       # 普遍IML相关函数
library(DALEX)     # 普遍IML相关函数
# devtools::install_github('thomasp85/lime')
library(lime)      # 局部可解释模型无关解释
load("inputdata.Rda")#加载示例数据
inputdata<-inputdata[,-1]
inputdata$Event<-factor(inputdata$Event,levels = c(0,1),labels = c("Alive","Death"))#结局变量因子化
set.seed(123)  # 设置随机种子保证可重复性
split <- initial_split(inputdata, strata = "Event")#数据分割
data_train <- training(split)
data_test <- testing(split)

# 分类变量ML蓝图
blueprint <- recipe(Event ~ ., data = data_train) %>%
  step_other(all_nominal(), threshold = 0.005)

#为h2o产生训练数据集和验证集
h2o.init()#启动h2o
train_h2o <- prep(blueprint, training = data_train, retain = TRUE) %>%
  juice() %>%
  as.h2o()
test_h2o <- prep(blueprint, training = data_train) %>%
  bake(new_data = data_test) %>%
  as.h2o()

#定义响应变量（目标变量）和特征变量
Y <- "Event"
X <- setdiff(names(data_train), Y)

#定义GBM网格搜索超参数
hyper_grid <- list(
  max_depth = c(1, 3, 5),
  min_rows = c(1, 5, 10),
  learn_rate = c(0.01, 0.05, 0.1),
  learn_rate_annealing = c(0.99, 1),
  sample_rate = c(0.5, 0.75, 1),
  col_sample_rate = c(0.8, 0.9, 1)
)

#定义随机网格搜索标准
search_criteria <- list(
  strategy = "RandomDiscrete",
  max_models = 25
)

#建立网格搜索
random_grid <- h2o.grid(
  algorithm = "gbm", grid_id = "gbm_grid", x = X, y = Y,
  training_frame = train_h2o, hyper_params = hyper_grid,
  search_criteria = search_criteria, ntrees = 5000, stopping_metric = "AUC",     
  stopping_rounds = 10, stopping_tolerance = 0, nfolds = 10, 
  fold_assignment = "Modulo", keep_cross_validation_predictions = TRUE,
  seed = 123
)

# 使用GBM算法拟合一个堆叠模型
ensemble <- h2o.stackedEnsemble(
  x = X, y = Y, training_frame = train_h2o, model_id = "ensemble_gbm_grid",
  base_models = random_grid@model_ids, metalearner_algorithm = "gbm"
)

#创建一个自定义函数，将预测值返回为一个向量
pred <- function(object, newdata)  {
  results <- as.vector(h2o.predict(object, as.h2o(newdata)))
  return(results)
}
#以下过程会花费较长时间
vip(
  ensemble,
  train = as.data.frame(train_h2o),
  method = "permute",
  target = "Event",
  metric = "AUC",
  reference_class="Alive",
  nsim =3,
  sample_frac = 0.3,#抽取30%的样本
  pred_wrapper = pred
)

通过VIP函数获得的堆叠模型的置换特征重要性解释

偏依赖图

在机器学习中，偏依赖图（Partial Dependence Plot，简称PDP）是一种可视化工具，用于分析模型中特征与预测结果之间的关系。它帮助我们了解单个特征对预测结果的影响，同时控制其他特征的平均效果。
偏依赖图的绘制过程如下：

首先，选择一个特征作为要分析的目标特征。
在这个特征的取值范围内，生成一组等间距的值，称为网格。
对于每个网格值，保持其他特征不变，使用训练好的模型进行预测，并记录预测结果。
绘制网格值和对应预测结果的散点图，并计算每个网格值的平均预测结果。
最后，绘制平均预测结果随目标特征取值的变化曲线，即偏依赖图。
通过偏依赖图，我们可以直观地观察到特征与预测结果之间的非线性关系，发现特征对预测结果的贡献程度，从而更好地理解模型的行为和预测结果的解释。它对于解释模型、调整特征、优化模型以及发现数据中的隐藏模式都具有重要的帮助作用。

偏依赖图的实现

pdq包（Brandon Greenwell 2018）是一个被广泛用于构建PDPs（Partial Dependence Plots）的R包。

pdp: An R Package for Constructing Partial Dependence Plotsby Brandon M. Greenwell

iml和DALEX包也提供PDP的功能。然而，pdq包内置了对许多模型的支持，但对于不支持的模型（例如h2o堆叠模型），我们需要封装一个自定义的预测函数，如下所示。首先，我们创建一个自定义预测函数，但这里我们返回预测值的平均值。

# 自定义函数，返回预测值的均值
pred.prob <- function(object, newdata) { 
  pred <- predict(object, as.h2o(newdata),type="prob")
  prob<-pred[, 2]
  mean(prob)
}

然后，我们使用pdp::partial()函数计算偏依赖值。

# 计算部分依赖值PDP
pd_values <- partial(
  ensemble,
  type="classification",
  train =as.data.frame(train_h2o), 
  pred.var ="ABHD2",
  pred.fun = pred.prob,
  grid.resolution = 20
)

我们可以查看下PDP值

head(pd_values) 
# ABHD2      yhat
# 1 1.249515 0.7794302
# 2 1.574986 0.7794302
# 3 1.900457 0.7794302
# 4 2.225927 0.7794302
# 5 2.551398 0.7794302
# 6 2.876869 0.7794302

我们可以使用autoplot()函数和ggplot2来查看PDPs。

# Partial dependence plot
autoplot(pd_values, rug = TRUE, train = as.data.frame(train_h2o))+
  theme_bw()+
  geom_line(linewidth=1,color="red")

ABHD2 部分依赖图

PDPs主要用于展示特征对预测响应值的边际效应。然而，Brandon M Greenwell、Boehmke和McCarthy（2018）提出了一种方法，使用偏依赖函数的相对“平坦程度”作为变量重要性的衡量指标。其思想是对预测响应值具有更大边际效应的特征更为重要。我们可以通过使用vip包，并设置method = "pdp"来实现基于PDP的特征重要性度量。此外，基于上述过程得到的变量重要性分数还保留了计算的偏依赖值。

个体条件期望（Individual Conditional Expectation，ICE）

个体条件期望（Individual Conditional Expectation，ICE）是一种用于解释机器学习模型预测结果的方法。它可以帮助理解某个样本在特征空间中的预测值是如何随着特征的变化而变化的。ICE图是通过在特征空间中对某个样本进行采样，然后固定其他特征的值，仅改变一个特征的值，来计算预测结果的期望。通过绘制每个特征的ICE图，我们可以直观地看到在不同特征取值下，模型的预测结果是如何变化的。这可以帮助我们发现模型中的非线性关系和交互作用，进一步理解模型的预测能力和局限性。
与PDPs类似，用于绘制ICE曲线的首选包是pdp包；但是，iml包也提供了ICE曲线绘制功能。使用pdp包创建ICE曲线的过程与PDPs基本相同，不同之处在于在自定义预测函数中不再对预测值进行平均化（不应用mean()）。默认情况下，autoplot()会绘制所有观测值，我们可以通过设置center = TRUE以使每个观测值所对应的曲线按照第一条曲线居中。

# Construct c-ICE curves
#先自定义函数用于返回预测值，与PDP不一样的是这里不取均值
pred.ice <- function(object, newdata)  {
  results <- as.vector(h2o.predict(object, as.h2o(newdata))[,2])
  return(results)
}
#绘制ICE图
partial(
  ensemble,
  train = as.data.frame(train_h2o), 
  pred.var = "ABHD2",
  pred.fun = pred,
  grid.resolution = 20,
  plot = TRUE,
  center = TRUE,
  plot.engine = "ggplot2"
)