Python-可解释的机器学习-四-Python 可解释的机器学习（四）第十章：可解释性特征选择与工程在前三章中，我

Python 可解释的机器学习（四）

原文：annas-archive.org/md5/9b99a159e8340372894ac9bde8bbd5d9

译者：飞龙

协议：CC BY-NC-SA 4.0

第十章：可解释性特征选择与工程

在前三章中，我们讨论了复杂性如何阻碍机器学习（ML）的可解释性。这里有一个权衡，因为你可能需要一些复杂性来最大化预测性能，但又不能达到无法依赖模型来满足可解释性原则：公平性、责任和透明度的程度。本章是四个专注于如何调整以实现可解释性的章节中的第一个。提高可解释性最简单的方法之一是通过特征选择。它有许多好处，例如加快训练速度并使模型更容易解释。但如果这两个原因不能说服你，也许另一个原因会。

一个常见的误解是，复杂的模型可以自行选择特征并仍然表现良好，那么为什么还要费心选择特征呢？是的，许多模型类别都有机制可以处理无用的特征，但它们并不完美。而且，随着每个剩余的机制的加入，过拟合的可能性也会增加。过拟合的模型是不可靠的，即使它们更准确。因此，虽然仍然强烈建议使用模型机制，如正则化，以避免过拟合，但特征选择仍然是有用的。

在本章中，我们将理解无关特征如何对模型的输出产生不利影响，从而了解特征选择对模型可解释性的重要性。然后，我们将回顾基于过滤器的特征选择方法，如斯皮尔曼相关系数，并了解嵌入式方法，如LASSO 和岭回归。然后，我们将发现包装方法，如顺序特征选择，以及混合方法，如递归特征消除（RFE）。最后，尽管特征工程通常在选择之前进行，但在特征选择完成后，探索特征工程仍有其价值。

这些是我们将在本章中讨论的主要主题：

理解无关特征的影响
回顾基于过滤器的特征选择方法
探索嵌入式特征选择方法
发现包装、混合和高级特征选择方法
考虑特征工程

让我们开始吧！

技术要求

本章的示例使用了mldatasets、pandas、numpy、scipy、mlxtend、sklearn-genetic-opt、xgboost、sklearn、matplotlib和seaborn库。有关如何安装所有这些库的说明见前言。

本章的 GitHub 代码位于此处：packt.link/1qP4P。

任务

据估计，全球有超过 1000 万个非营利组织，尽管其中很大一部分有公共资金，但大多数组织主要依赖私人捐赠者，包括企业和个人，以继续运营。因此，筹款是至关重要的任务，并且全年都在进行。

年复一年，捐款收入有所增长，但非营利组织面临一些问题：捐赠者的兴趣在变化，因此一年受欢迎的慈善机构可能在下一年被遗忘；非营利组织之间的竞争激烈，人口结构也在变化。在美国，平均捐赠者每年只捐赠两次慈善礼物，且年龄超过 64 岁。识别潜在捐赠者具有挑战性，而且吸引他们的活动可能成本高昂。

一个国家级退伍军人组织非营利分支拥有大约 190,000 名往届捐赠者的庞大邮件列表，并希望发送一份特别邮件请求捐款。然而，即使有特殊的批量折扣率，每地址的成本也高达 0.68 美元。这总计超过 130,000 美元。他们的市场预算只有 35,000 美元。鉴于他们已将此事列为高优先级，他们愿意扩展预算，但前提是投资回报率（ROI）足够高，以证明额外成本是合理的。

为了最大限度地减少使用他们有限的预算，他们希望尝试直接邮寄，目的是利用已知的信息来识别潜在捐赠者，例如过去的捐赠、地理位置和人口统计数据。他们将通过电子邮件联系其他捐赠者，这要便宜得多，整个列表的月成本不超过 1,000 美元。他们希望这种混合营销计划能产生更好的结果。他们还认识到，高价值捐赠者对个性化的纸质邮件响应更好，而较小的捐赠者无论如何对电子邮件的响应更好。

最多只有 6%的邮件列表捐赠者会对任何特定的活动进行捐赠。使用机器学习预测人类行为绝非易事，尤其是在数据类别不平衡的情况下。尽管如此，成功不是以最高的预测准确性来衡量的，而是以利润提升来衡量。换句话说，在测试数据集上评估的直接邮寄模型应该产生比如果他们向整个数据集进行群发邮件更多的利润。

他们寻求您的帮助，使用机器学习（ML）来生成一个模型，以识别最可能的捐赠者，但同时也保证一个高的 ROI。

您收到了非营利组织的数据集，该数据集大约平均分为训练数据和测试数据。如果您向测试数据集中的所有人发送邮件，您将获得 11,173 美元的利润，但如果您能够仅识别那些会捐赠的人，最大收益将达到 73,136 美元。您的目标是实现高利润提升和合理的 ROI。当活动进行时，它将识别整个邮件列表中最可能的捐赠者，非营利组织希望总支出不超过 35,000 美元。然而，数据集有 435 个列，一些简单的统计测试和建模练习表明，由于过度拟合，数据过于嘈杂，无法识别潜在捐赠者的可靠性。

方法

你决定首先使用所有特征拟合一个基础模型，并在不同的复杂度级别上评估它，以了解特征数量增加与预测模型过度拟合训练数据之间的关联。然后，你将采用一系列从简单的基于过滤的方法到最先进的方法的特征选择方法，以确定哪种方法实现了客户寻求的盈利性和可靠性目标。最后，一旦选定了最终特征列表，你就可以尝试特征工程。

由于问题的成本敏感性，阈值对于优化利润提升至关重要。我们将在稍后讨论阈值的作用，但一个显著的影响是，尽管这是一个分类问题，最好使用回归模型，然后使用预测来分类，这样只有一个阈值需要调整。也就是说，对于分类模型，你需要一个用于标签的阈值，比如那些捐赠超过 1 美元的，然后还需要另一个用于预测概率的阈值。另一方面，回归预测捐赠金额，阈值可以根据这个进行优化。

准备工作

此示例的代码可以在github.com/PacktPublishing/Interpretable-Machine-Learning-with-Python-2E/blob/main/10/Mailer.ipynb找到。

加载库

要运行此示例，我们需要安装以下库：

使用mldatasets加载数据集
使用pandas、numpy和scipy来操作它
使用mlxtend、sklearn-genetic-opt、xgboost和sklearn（scikit-learn）来拟合模型
使用matplotlib和seaborn创建和可视化解释

要加载库，请使用以下代码块：

import math
import os
import mldatasets
import pandas as pd
import numpy as np
import timeit
from tqdm.notebook import tqdm
from sklearn.feature_selection import VarianceThreshold,\
                                    mutual_info_classif, SelectKBest
from sklearn.feature_selection import SelectFromModel
from sklearn.linear_model import LogisticRegression,\
                                    LassoCV, LassoLarsCV, LassoLarsIC
from mlxtend.feature_selection import SequentialFeatureSelector
from sklearn.feature_selection import RFECV
from sklearn.decomposition import PCA import shap
from sklearn-genetic-opt import GAFeatureSelectionCV
from scipy.stats import rankdata
from sklearn.discriminant_analysis import
LinearDiscriminantAnalysis
from sklearn.ensemble import RandomForestRegressor
import xgboost as xgb
import matplotlib.pyplot as plt
import seaborn as sns

接下来，我们将加载并准备数据集。

理解和准备数据

我们将数据这样加载到两个 DataFrame（X_train和X_test）中，其中包含特征，以及两个相应的numpy数组标签（y_train和y_test）。请注意，这些 DataFrame 已经为我们预先准备，以删除稀疏或不必要的特征，处理缺失值，并对分类特征进行编码：

X_train, X_test, y_train, y_test = mldatasets.load(
    "nonprofit-mailer",
    prepare=True
)
y_train = y_train.squeeze()
y_test = y_test.squeeze()

所有特征都是数值型，没有缺失值，并且分类特征已经为我们进行了一热编码。在训练和测试邮件列表之间，应有超过 191,500 条记录和 435 个特征。你可以这样检查：

print(X_train.shape)
print(y_train.shape)
print(X_test.shape)
print(y_test.shape)

上一段代码应该输出以下内容：

(95485, 435)
(95485,)
(96017, 435)
(96017,)

接下来，我们可以使用变量成本 0.68（var_cost）验证测试标签是否有正确的捐赠者数量（test_donors）、捐赠金额（test_donations）和假设的利润范围（test_min_profit和test_max_profit）。我们可以打印这些信息，然后对训练数据集做同样的操作：

var_cost = 0.68
y_test_donors = y_test[y_test > 0]
test_donors = len(y_test_donors)
test_donations = sum(y_test_donors)
test_min_profit = test_donations - (len(y_test)*var_cost)
test_max_profit = test_donations - (test_donors*var_cost)
print(
    '%s test donors totaling $%.0f (min profit: $%.0f,\
    max profit: $%.0f)'
    %(test_donors, test_donations, test_min_profit,\
       test_max_profit))
y_train_donors = y_train[y_train > 0]
train_donors = len(y_train_donors)
train_donations = sum(y_train_donors)
train_min_profit = train_donations – (len(y_train)*var_cost)
train_max_profit = train_donations – (train_donors*var_cost)
print(
    '%s train donors totaling $%.0f (min profit: $%.0f,\
    max profit: $%.0f)'
    %(train_donors, train_donations, train_min_profit,\
    train_max_profit))

上一段代码应该输出以下内容：

4894 test donors totaling $76464 (min profit: $11173, max profit: $73136)
4812 train donors totaling $75113 (min profit: $10183, max profit: $71841)

事实上，如果非营利组织向测试邮件列表上的每个人大量邮寄，他们可能会获得大约 11,000 美元的利润，但为了实现这一目标，他们必须严重超支。非营利组织认识到，通过仅识别和针对捐赠者来获得最大利润几乎是一项不可能完成的任务。因此，他们宁愿生产一个能够可靠地产生超过最低利润但成本更低的模型，最好是低于预算。

理解无关特征的影响

特征选择也称为变量或属性选择。这是你可以自动或手动选择一组对构建机器学习模型有用的特定特征的方法。

并非更多的特征就一定能导致更好的模型。无关特征可能会影响学习过程，导致过拟合。因此，我们需要一些策略来移除可能对学习产生不利影响的任何特征。选择较小特征子集的一些优点包括以下内容：

理解简单的模型更容易：例如，对于使用 15 个变量的模型，其特征重要性比使用 150 个变量的模型更容易理解。
缩短训练时间：减少变量的数量可以降低计算成本，加快模型训练速度，而且最值得注意的是，简单的模型具有更快的推理时间。
通过减少过拟合来提高泛化能力：有时，预测价值很小，许多变量只是噪声。然而，机器学习模型却会从这些噪声中学习，并在最小化泛化的同时触发对训练数据的过拟合。通过移除这些无关或噪声特征，我们可以显著提高机器学习模型的泛化能力。
变量冗余：数据集中存在共线性特征是很常见的，这可能意味着某些特征是冗余的。在这些情况下，只要没有丢失显著信息，我们只需保留一个相关的特征，删除其他特征即可。

现在，我们将拟合一些模型来展示过多特征的影响。

创建基础模型

让我们为我们的邮件列表数据集创建一个基础模型，看看这将如何展开。但首先，让我们设置随机种子以确保可重复性：

rand = 9
os.environ['PYTHONHASHSEED']=str(rand)
np.random.seed(rand)

在本章中，我们将使用 XGBoost 的随机森林（RF）回归器（XGBRFRegressor）。它就像 scikit-learn 一样，但更快，因为它使用了目标函数的二阶近似。它还有更多选项，例如设置学习率和单调约束，这些在第十二章，单调约束和模型调优以提高可解释性中进行了考察。我们以保守的初始max_depth值4初始化XGBRFRegressor，并始终使用200估计量以确保一致性。然后，我们使用我们的训练数据对其进行拟合。我们将使用timeit来测量它需要多长时间，并将其保存在变量（baseline_time）中供以后参考：

stime = timeit.default_timer()
reg_mdl = xgb.XGBRFRegressor(max_depth=4, n_estimators=200, seed=rand)
fitted_mdl = reg_mdl.fit(X_train, y_train)
etime = timeit.default_timer()
baseline_time = etime-stime

现在我们已经有一个基础模型了，让我们来评估它。

评估模型

接下来，让我们创建一个字典(reg_mdls)来存放我们将在本章中拟合的所有模型，以测试哪些特征子集会产生最好的模型。在这里，我们可以使用evaluate_reg_mdl来评估具有所有特征和max_depth值为4的随机森林模型(rf_4_all)。它将生成一个总结和一个带有回归线的散点图：

reg_mdls = {}
reg_mdls['rf_4_all'] = mldatasets.evaluate_reg_mdl(
    fitted_mdl,
    X_train,
    X_test,
    y_train,
    y_test,
    plot_regplot=True,
    ret_eval_dict=True
)

之前的代码生成了图 10.1中显示的指标和图表：

图表描述自动生成

图 10.1：基础模型的预测性能

对于像图 10.1这样的图表，通常期望看到一条对角线，所以一眼看去就能判断出这个模型不具备预测性。此外，均方根误差（RMSEs）可能看起来并不糟糕，但在这种不平衡的问题背景下，它们却是令人沮丧的。考虑一下这个情况：只有 5%的人捐赠，而其中只有 20%的人捐赠额超过 20 美元，所以平均误差 4.3 美元至 4.6 美元是巨大的。

那么，这个模型有没有用呢？答案在于我们用它来分类所使用的阈值。让我们首先定义一个从 $0.40 到$ 25 的阈值数组(threshs)，我们首先以每 0.01 美元的间隔来设置这些阈值，直到达到 $1，然后以每 0.1 美元的间隔设置，直到达到$ 3，之后以每 1 美元的间隔设置：

threshs = np.hstack(
    [
      np.linspace(0.40,1,61),
      np.linspace(1.1,3,20),
      np.linspace(4,25,22)
    ]
)

mldatasets中有一个函数可以计算每个阈值下的利润(profits_by_thresh)。它只需要实际的(y_test)和预测标签，然后是阈值(threshs)、可变成本(var_costs)和所需的min_profit。只要利润高于min_profit，它就会生成一个包含每个阈值的收入、成本、利润和投资回报率的pandas DataFrame。记住，我们在本章开始时将这个最低值设置为$11,173，因为针对低于这个金额的捐赠者是没有意义的。在为测试和训练数据集生成这些利润 DataFrame 之后，我们可以将这些最大和最小金额放入模型的字典中，以供以后使用。然后，我们使用compare_df_plots来绘制每个阈值的测试和训练的成本、利润和投资回报率比率，只要它超过了利润最低值：

y_formatter = plt.FuncFormatter(
    lambda x, loc: "${:,}K".format(x/1000)
)
profits_test = mldatasets.profits_by_thresh(
    y_test,
    reg_mdls['rf_4_all']['preds_test'],
    threshs,
    var_costs=var_cost,
    min_profit=test_min_profit
)
profits_train = mldatasets.profits_by_thresh(
    y_train,
    reg_mdls['rf_4_all']['preds_train'],
    threshs,
    var_costs=var_cost,
    min_profit=train_min_profit
)
reg_mdls['rf_4_all']['max_profit_train'] =profits_train.profit.max()
reg_mdls['rf_4_all']['max_profit_test'] = profits_test.profit.max()
reg_mdls['rf_4_all']['max_roi'] = profits_test.roi.max()
reg_mdls['rf_4_all']['min_costs'] = profits_test.costs.min()
reg_mdls['rf_4_all']['profits_train'] = profits_train
reg_mdls['rf_4_all']['profits_test'] = profits_test
mldatasets.compare_df_plots(
    profits_test[['costs', 'profit', 'roi']],
    profits_train[['costs', 'profit', 'roi']],
    'Test',
    'Train',
    y_formatter=y_formatter,
    x_label='Threshold',\
    plot_args={'secondary_y':'roi'}
) 
Figure 10.2. You can tell that Test and Train are almost identical. Costs decrease steadily at a high rate and profit at a lower rate, while ROI increases steadily. However, some differences exist, such as ROI, which becomes a bit higher eventually, and although viable thresholds start at the same point, Train does end at a different threshold. It turns out the model can turn a profit, so despite the appearance of the plot in *Figure 10.1*, the model is far from useless:

图表，折线图描述自动生成

图 10.2：测试和训练数据集在阈值下基础模型的利润、成本和投资回报率比较

训练集和测试集的 RMSEs 差异是真实的。模型没有过拟合。主要原因是我们通过将max_depth值设置为4使用了相对较浅的树。我们可以通过计算有多少特征的feature_importances_值超过 0 来轻易地看到使用浅树的效果：

reg_mdls['rf_4_all']['total_feat'] =\
    reg_mdls['rf_4_all']['fitted'].feature_importances_.shape[0] reg_mdls['rf_4_all']['num_feat'] = sum(
    reg_mdls['rf_4_all']['fitted'].feature_importances_ > 0
)
print(reg_mdls['rf_4_all']['num_feat'])

之前的代码输出160。换句话说，只有 160 个在 435 个中使用了——在这样的浅树中只能容纳这么多的特征！自然地，这会导致降低过度拟合，但与此同时，在具有杂质度量的特征与随机选择特征之间的选择并不一定是最佳选择。

在不同的最大深度下训练基础模型

那么，如果我们使树更深会发生什么？让我们重复之前为浅层模型所做的所有步骤，但这次的最大深度在 5 到 12 之间：

for depth in tqdm(range(5, 13)):
mdlname = 'rf_'+str(depth)+'_all'
stime = timeit.default_timer()
reg_mdl = xgb.XGBRFRegressor(
    max_depth=depth,
    n_estimators=200,
    seed=rand
)
fitted_mdl = reg_mdl.fit(X_train, y_train)
etime = timeit.default_timer()
reg_mdls[mdlname] = mldatasets.evaluate_reg_mdl(
    fitted_mdl,
    X_train,
    X_test,
    y_train,
    y_test,
    plot_regplot=False,
    show_summary=False,
    ret_eval_dict=True
)
reg_mdls[mdlname]['speed'] = (etime - stime)/baseline_time
reg_mdls[mdlname]['depth'] = depth
reg_mdls[mdlname]['fs'] = 'all'
profits_test = mldatasets.profits_by_thresh(
    y_test,
    reg_mdls[mdlname]['preds_test'],
    threshs,
    var_costs=var_cost,
    min_profit=test_min_profit
)
profits_train = mldatasets.profits_by_thresh(
    y_train,
    reg_mdls[mdlname]['preds_train'],
    threshs,
    var_costs=var_cost,
    min_profit=train_min_profit
)
reg_mdls[mdlname]['max_profit_train'] = profits_train.profit.max()
reg_mdls[mdlname]['max_profit_test'] = profits_test.profit.max()
reg_mdls[mdlname]['max_roi'] = profits_test.roi.max()
reg_mdls[mdlname]['min_costs'] = profits_test.costs.min()
reg_mdls[mdlname]['profits_train'] = profits_train
reg_mdls[mdlname]['profits_test'] = profits_test
reg_mdls[mdlname]['total_feat'] =\
reg_mdls[mdlname]['fitted'].feature_importances_.shape[0]
reg_mdls[mdlname]['num_feat'] = sum(
    reg_mdls[mdlname]['fitted'].feature_importances_ > 0)

现在，让我们像之前使用 compare_df_plots 一样，绘制“最深”模型（最大深度为 12）的利润 DataFrame 的细节，生成图 10.3：

图表，折线图描述自动生成

图 10.3：比较测试和训练数据集对于“深”基础模型在阈值下的利润、成本和 ROI

看看这次图 10.3中不同的测试和训练。测试达到约 15,000 的最大值，而训练超过 20,000。训练的成本大幅下降，使得投资回报率比测试高几个数量级。此外，阈值范围也大不相同。你可能会问，这为什么会成为问题？如果我们必须猜测使用什么阈值来选择在下一封邮件中要针对的目标，训练的最佳阈值高于测试——这意味着使用过度拟合的模型，我们可能会错过目标，并在未见过的数据上表现不佳。

接下来，让我们将我们的模型字典（reg_mdls）转换为 DataFrame，并从中提取一些细节。然后，我们可以按深度排序它，格式化它，用颜色编码它，并输出它：

def display_mdl_metrics(reg_mdls, sort_by='depth', max_depth=None):
    reg_metrics_df = pd.DataFrame.from_dict( reg_mdls, 'index')\
                        [['depth', 'fs', 'rmse_train', 'rmse_test',\
                          'max_profit_train',\
                          'max_profit_test', 'max_roi',\
                          'min_costs', 'speed', 'num_feat']]
    pd.set_option('precision', 2) 
    html = reg_metrics_df.sort_values(
        by=sort_by, ascending=False).style.\
        format({'max_profit_train':'${0:,.0f}',\
        'max_profit_test':'${0:,.0f}', 'min_costs':'${0:,.0f}'}).\
        background_gradient(cmap='plasma', low=0.3, high=1,
                            subset=['rmse_train', 'rmse_test']).\
        background_gradient(cmap='viridis', low=1, high=0.3,
                            subset=[
                                'max_profit_train', 'max_profit_test'
                                ]
                            )
    return html
display_mdl_metrics(reg_mdls) 
display_mdl_metrics function to output the DataFrame shown in *Figure 10.4*. Something that should be immediately visible is how RMSE train and RMSE test are inverses. One decreases dramatically, and another increases slightly as the depth increases. The same can be said for profit. ROI tends to increase with depth and training speed and the number of features used as well:

表格描述自动生成

图 10.4：比较所有基础 RF 模型在不同深度下的指标

我们可能会倾向于使用具有最高盈利能力的 rf_11_all，但使用它是有风险的！一个常见的误解是，黑盒模型可以有效地消除任何数量的无关特征。虽然它们通常能够找到有价值的东西并充分利用它，但过多的特征可能会因为过度拟合训练数据集中的噪声而降低它们的可靠性。幸运的是，存在一个甜蜜点，你可以以最小的过度拟合达到高盈利能力，但为了达到这一点，我们首先必须减少特征的数量！

检查基于过滤的特征选择方法

基于过滤的方法独立地从数据集中选择特征，而不使用任何机器学习。这些方法仅依赖于变量的特征，并且相对有效、计算成本低、执行速度快。因此，作为特征选择方法的低垂之果，它们通常是任何特征选择流程的第一步。

基于过滤的方法可以分为：

单变量：它们独立于特征空间，一次评估和评级一个特征。单变量方法可能存在的问题是，由于它们没有考虑特征之间的关系，可能会过滤掉太多信息。
多元性：这些方法考虑整个特征空间以及特征之间的相互作用。

总体而言，对于移除过时、冗余、常数、重复和不相关的特征，过滤方法非常有效。然而，由于它们没有考虑到只有机器学习模型才能发现的复杂、非线性、非单调的相关性和相互作用，当这些关系在数据中突出时，它们并不有效。

我们将回顾三种基于过滤的方法：

基础
相关性
排序

我们将在各自的章节中进一步解释它们。

基础过滤方法

在数据准备阶段，我们采用基本过滤方法，特别是在任何建模之前的数据清洗阶段。这样做的原因是，做出可能对模型产生不利影响的特征选择决策的风险很低。这些方法涉及常识性操作，例如移除不携带信息或重复信息的特征。

带有方差阈值的常数特征

常数特征在训练数据集中不发生变化，因此不携带任何信息，模型无法从中学习。我们可以使用一个名为VarianceThreshold的单变量方法，它移除低方差特征。我们将使用零作为阈值，因为我们只想过滤掉具有零方差的特征——换句话说，就是常数特征。它仅适用于数值特征，因此我们必须首先确定哪些特征是数值的，哪些是分类的。一旦我们将方法拟合到数值列上，get_support()返回的不是常数特征的列表，我们可以使用集合代数来返回仅包含常数特征的集合（num_const_cols）：

num_cols_l = X_train.select_dtypes([np.number]).columns
cat_cols_l = X_train.select_dtypes([np.bool, np.object]).columns
num_const = VarianceThreshold(threshold=0)
num_const.fit(X_train[num_cols_l])
num_const_cols = list(
    set(X_train[num_cols_l].columns) -
    set(num_cols_l[num_const.get_support()])
) 
nunique() function on categorical features. It will return a pandas series, and then a lambda function can filter out only those with one unique value. Then, .index.tolist() returns the name of the features as a list. Now, we just join both lists of constant features, and voilà! We have all constants (all_const_cols). We can print them; there should be three:

cat_const_cols = X_train[cat_cols_l].nunique()[lambda x:\
                                               x<2].index.tolist()
all_const_cols = num_const_cols + cat_const_cols
print(all_const_cols)

在大多数情况下，仅移除常数特征是不够的。一个冗余特征可能几乎是常数或准常数。

带有 value_counts 的准常数特征

准常数特征几乎都是相同的值。与常数过滤不同，使用方差阈值不会起作用，因为高方差和准常数性不是互斥的。相反，我们将迭代所有特征并获取 value_counts()，它返回每个值的行数。然后，将这些计数除以总行数以获得百分比，并按最高百分比排序。如果最高值高于预先设定的阈值（thresh），则将其追加到准常数列列表（quasi_const_cols）中。请注意，选择此阈值必须非常谨慎，并且需要对问题有深入的理解。例如，在这种情况下，我们知道这是不平衡的，因为只有 5% 的人捐赠，其中大多数人捐赠的金额很低，所以即使是特征的一小部分也可能产生影响，这就是为什么我们的阈值如此之高，达到 99.9%：

thresh = 0.999
quasi_const_cols = []
num_rows = X_train.shape[0]
for col in tqdm(X_train.columns):
    top_val = (
        X_train[col].value_counts() / num_rows
        ).sort_values(ascending=False).values[0]
    if top_val >= thresh:
        quasi_const_cols.append(col)
print(quasi_const_cols)

前面的代码应该已经打印出五个特征，其中包括之前获得的三个。接下来，我们将处理另一种形式的不相关特征：重复项！

重复特征

通常，当我们讨论数据中的重复项时，我们首先想到的是重复的行，但重复的列也是问题所在。我们可以像查找重复行一样找到它们，使用 pandas duplicated() 函数，但首先需要将 DataFrame 转置，反转列和行：

X_train_transposed = X_train.T
dup_cols = X_train_transposed[
    X_train_transposed.duplicated()].index.tolist()
print(dup_cols)

前面的代码片段输出了一个包含两个重复列的列表。

移除不必要的特征

与其他特征选择方法不同，您应该用模型测试这些方法，您可以直接通过移除您认为无用的特征来应用基于基本过滤的特征选择方法。但以防万一，制作原始数据的副本是一个好习惯。请注意，我们不将常数列（all_constant_cols）包括在我们打算删除的列（drop_cols）中，因为准常数列已经包含它们：

X_train_orig = X_train.copy()
X_test_orig = X_test.copy()
drop_cols = quasi_const_cols + dup_cols
X_train.drop(labels=drop_cols, axis=1, inplace=True)
X_test.drop(labels=drop_cols, axis=1, inplace=True)

接下来，我们将探索剩余特征上的多变量过滤方法。

基于相关性的过滤方法

基于相关性的过滤方法量化两个特征之间关系的强度。这对于特征选择很有用，因为我们可能想要过滤掉高度相关的特征或那些与其他特征完全不相关的特征。无论如何，它是一种多变量特征选择方法——更确切地说，是双变量特征选择方法。

但首先，我们应该选择一个相关性方法：

皮尔逊相关系数：衡量两个特征之间的线性相关性。它输出一个介于 -1（负相关）和 1（正相关）之间的系数，0 表示没有线性相关性。与线性回归类似，它假设线性、正态性和同方差性——也就是说，线性回归线周围的误差项在所有值中大小相似。
斯皮尔曼秩相关系数：衡量两个特征单调性的强度，无论它们是否线性相关。单调性是指一个特征增加时，另一个特征持续增加或减少的程度。它在-1 和 1 之间衡量，0 表示没有单调相关性。它不做分布假设，可以与连续和离散特征一起使用。然而，它的弱点在于非单调关系。
肯德尔 tau 相关系数：衡量特征之间的序数关联——也就是说，它计算有序数字列表之间的相似性。它也介于-1 和 1 之间，但分别代表低和高。对于离散特征来说，它很有用。

数据集是连续和离散的混合，我们不能对其做出任何线性假设，因此spearman是正确的选择。尽管如此，所有三个都可以与pandas的corr函数一起使用：

corrs = X_train.corr(method='spearman')
print(corrs.shape)

前面的代码应该输出相关矩阵的形状，即(428, 428)。这个维度是有意义的，因为还剩下 428 个特征，每个特征都与 428 个特征有关，包括它自己。

我们现在可以在相关矩阵（corrs）中寻找要删除的特征。请注意，为了做到这一点，我们必须建立阈值。例如，我们可以说一个高度相关的特征具有超过 0.99 的绝对值系数，而对于一个不相关的特征则小于 0.15。有了这些阈值，我们可以找到只与一个特征相关并且与多个特征高度相关的特征。为什么是一个特征？因为在相关矩阵的对角线总是 1，因为一个特征总是与自己完美相关。以下代码中的lambda函数确保我们考虑到这一点：

extcorr_cols = (abs(corrs) > 0.99).sum(axis=1)[lambda x: x>1]\
.index.tolist()
print(extcorr_cols)
uncorr_cols = (abs(corrs) > 0.15).sum(axis=1)[lambda x: x==1]\
.index.tolist()
print(uncorr_cols)

前面的代码以如下方式输出两个列表：

['MAJOR', 'HHAGE1', 'HHAGE3', 'HHN3', 'HHP1', 'HV1', 'HV2', 'MDMAUD_R', 'MDMAUD_F', 'MDMAUD_A']
['TCODE', 'MAILCODE', 'NOEXCH', 'CHILD03', 'CHILD07', 'CHILD12', 'CHILD18', 'HC15', 'MAXADATE']

第一个列表包含与除自身以外的其他特征高度相关的特征。虽然了解这一点很有用，但你不应在没有理解它们与哪些特征以及如何相关，以及与目标相关的情况下从该列表中删除特征。然后，只有在发现冗余的情况下，确保只删除其中一个。第二个列表包含与除自身以外的任何其他特征都不相关的特征，鉴于特征的数量众多，这在当前情况下是可疑的。话虽如此，我们也应该逐个检查它们，特别是要衡量它们与目标的相关性，看看它们是否冗余。然而，我们将冒险排除不相关的特征，创建一个特征子集（corr_cols）：

corr_cols = X_train.columns[
    ~X_train.columns.isin(uncorr_cols)
].tolist()
print(len(corr_cols))

前面的代码应该输出419。现在让我们只使用这些特征来拟合 RF 模型。鉴于仍有超过 400 个特征，我们将使用max_depth值为11。除了这一点和一个不同的模型名称（mdlname）之外，代码与之前相同：

mdlname = 'rf_11_f-corr'
stime = timeit.default_timer()
reg_mdl = xgb.XGBRFRegressor(
    max_depth=11,
    n_estimators=200,
    seed=rand
)
fitted_mdl = reg_mdl.fit(X_train[corr_cols], y_train)
reg_mdls[mdlname]['num_feat'] = sum(
    reg_mdls[mdlname]['fitted'].feature_importances_ > 0
)

在比较前面模型的输出结果之前，让我们了解一下排名滤波方法。

基于排序过滤的方法

基于排序过滤的方法基于统计单变量排序测试，这些测试评估特征与目标之间的依赖强度。这些是一些最受欢迎的方法：

方差分析 F 检验（ANOVA）F 检验衡量特征与目标之间的线性依赖性。正如其名所示，它是通过分解方差来做到这一点的。它做出了与线性回归类似的假设，例如正态性、独立性和同方差性。在 scikit-learn 中，您可以使用 f_regression 和 f_classification 分别对回归和分类进行排序，以 F 检验产生的 F 分数来排序特征。
卡方检验独立性：这个测试衡量非负分类变量与二元目标之间的关联性，因此它只适用于分类问题。在 scikit-learn 中，您可以使用 chi2。
互信息（MI）：与前面两种方法不同，这种方法是从信息理论而不是经典统计假设检验中推导出来的。虽然名称不同，但这个概念我们在本书中已经讨论过，称为库尔巴克-莱布勒（KL）散度，因为它是对特征 X 和目标 Y 的 KL。scikit-learn 中的 Python 实现使用了一个数值稳定的对称 KL 衍生品，称为Jensen-Shannon（JS）散度，并利用 k-最近邻来计算距离。可以使用 mutual_info_regression 和 mutual_info_classif 分别对回归和分类进行特征排序。

在提到的三种选项中，最适合这个数据集的是 MI，因为我们不能假设特征之间存在线性关系，而且其中大部分也不是分类数据。我们可以尝试使用阈值为 $0.68 的分类，这至少可以覆盖发送邮件的成本。为此，我们必须首先使用该阈值创建一个二元分类目标（y_train_class）：

y_train_class = np.where(y_train > 0.68, 1, 0)

接下来，我们可以使用 SelectKBest 根据互信息分类（MIC）获取前 160 个特征。然后我们使用 get_support() 获取一个布尔向量（或掩码），它告诉我们哪些特征在前 160 个中，并使用这个掩码对特征列表进行子集化：

mic_selection = SelectKBest(
    mutual_info_classif, k=160).fit(X_train, y_train_class)
mic_cols = X_train.columns[mic_selection.get_support()].tolist()
print(len(mic_cols))

前面的代码应该确认 mic_cols 列表中确实有 160 个特征。顺便说一下，这是一个任意数字。理想情况下，我们可以测试分类目标的不同阈值和 MI 的 k 值，寻找在最小过拟合的同时实现最高利润提升的模型。接下来，我们将使用与之前相同的 MIC 特征拟合 RF 模型。这次，我们将使用最大深度为 5，因为特征数量显著减少：

mdlname = 'rf_5_f-mic'
stime = timeit.default_timer()
reg_mdl = xgb.XGBRFRegressor(max_depth=5, n_estimators=200, seed=rand)
fitted_mdl = reg_mdl.fit(X_train[mic_cols], y_train)
reg_mdls[mdlname]['num_feat'] = sum(
    reg_mdls[mdlname]['fitted'].feature_importances_ > 0
)

现在，让我们像在 图 10.3 中所做的那样绘制测试和训练的利润，但这次是针对 MIC 模型。它将产生 图 10.5 中所示的内容：

图表，折线图自动生成的描述

图 10.5：具有 MIC 特征的模型在阈值之间的利润、成本和 ROI 测试和训练数据集的比较

在图 10.5中，你可以看出测试和训练之间存在相当大的差异，但相似之处表明过拟合最小。例如，训练的最高盈利性可以在 0.66 和 0.75 之间找到，而测试主要在 0.66 和 0.7 之间，之后逐渐下降。

尽管我们已经视觉检查了 MIC 模型，但查看原始指标也是一种令人放心的方式。接下来，我们将使用一致的指标比较我们迄今为止训练的所有模型。

比较基于滤波的方法

我们已经将指标保存到一个字典（reg_mdls）中，我们很容易将其转换为 DataFrame，并像之前那样输出，但这次我们按max_profit_test排序：

display_mdl_metrics(reg_mdls, 'max_profit_test') 
Figure 10.6. It is evident that the filter MIC model is the least overfitted of all. It ranked higher than more complex models with more features and took less time to train than any model. Its speed is an advantage for hyperparameter tuning. What if we wanted to find the best classification target thresholds or MIC *k*? We won’t do this now, but we would likely get a better model if we ran every combination, but it would take time to do and even more with more features:

应用，表格，自动生成中等置信度的描述

图 10.6：比较所有基础模型和基于滤波的特征选择模型的指标

在图 10.6中，我们可以看出，与具有更多特征和相同max_depth量的模型（rf_11_all）相比，相关滤波模型（rf_11_f-corr）的表现更差，这表明我们可能移除了一个重要的特征。正如该部分所警告的，盲目设置阈值并移除其上所有内容的问题在于你可能会无意中移除有用的东西。并非所有高度相关和无关的特征都是无用的，因此需要进一步检查。接下来，我们将探索一些嵌入方法，当与交叉验证结合使用时，需要更少的监督。

探索嵌入特征选择方法

嵌入方法存在于模型本身中，通过训练过程中自然选择特征。你可以利用具有这些特性的任何模型的内在属性来捕获所选特征：

基于树的模型：例如，我们多次使用以下代码来计算 RF 模型使用的特征数量，这是学习过程中自然发生特征选择的证据：
```
sum(reg_mdls[mdlname]['fitted'].feature_importances_ > 0) 
```
XGBoost 的 RF 默认使用增益，这是在所有使用该特征进行特征重要性计算的分割中平均错误减少。我们可以将阈值提高到 0 以上，根据它们的相对贡献选择更少的特征。然而，通过限制树的深度，我们迫使模型选择更少的特征。
具有系数的正则化模型：我们将在第十二章，单调约束和模型调优以提高可解释性中进一步研究这个问题，但许多模型类可以采用基于惩罚的正则化，如 L1、L2 和弹性网络。然而，并非所有这些模型都具有可以提取以确定哪些特征被惩罚的内在参数，如系数。

本节将仅涵盖正则化模型，因为我们已经使用了一个基于树的模型。最好利用不同的模型类别来获得对哪些特征最重要的不同视角。

我们在第三章，解释挑战中介绍了一些这些模型，但这些都是一些结合基于惩罚的正则化和输出特征特定系数的模型类别：

最小绝对收缩和选择算子（LASSO）：因为它在损失函数中使用 L1 惩罚，所以 LASSO 可以将系数设置为 0。
最小角度回归（LARS）：类似于 LASSO，但基于向量，更适合高维数据。它也对等相关的特征更加公平。
岭回归：在损失函数中使用 L2 惩罚，因此只能将不相关的系数缩小到接近 0，但不能缩小到 0。
弹性网络回归：使用 L1 和 L2 范数的混合作为惩罚。
逻辑回归：根据求解器，它可以处理 L1、L2 或弹性网络惩罚。

之前提到的模型也有一些变体，例如LASSO LARS，它使用 LARS 算法进行 LASSO 拟合，或者甚至是LASSO LARS IC，它与前者相同，但在模型部分使用 AIC 或 BIC 准则：

赤池信息准则（AIC）：基于信息理论的一种相对拟合优度度量
贝叶斯信息准则（BIC）：与 AIC 具有相似的公式，但具有不同的惩罚项

好的，现在让我们使用SelectFromModel从 LASSO 模型中提取顶级特征。我们将使用LassoCV，因为它可以自动进行交叉验证以找到最优的惩罚强度。一旦拟合，我们就可以使用get_support()获取特征掩码。然后我们可以打印特征数量和特征列表：

lasso_selection = SelectFromModel(
    LassoCV(n_jobs=-1, random_state=rand)
)
lasso_selection.fit(X_train, y_train)
lasso_cols = X_train.columns[lasso_selection.get_support()].tolist()
print(len(lasso_cols))
print(lasso_cols)

上一段代码输出以下内容：

7
['ODATEDW', 'TCODE', 'POP901', 'POP902', 'HV2', 'RAMNTALL', 'MAXRDATE']

现在，让我们尝试使用LassoLarsCV进行相同的操作：

llars_selection = SelectFromModel(LassoLarsCV(n_jobs=-1))
llars_selection.fit(X_train, y_train)
llars_cols = X_train.columns[llars_selection.get_support()].tolist()
print(len(llars_cols))
print(llars_cols)

上一段代码生成以下输出：

8
['RECPGVG', 'MDMAUD', 'HVP3', 'RAMNTALL', 'LASTGIFT', 'AVGGIFT', 'MDMAUD_A', 'DOMAIN_SOCIALCLS']

LASSO 将除七个特征外的所有系数缩小到 0，而 LASSO LARS 也将八个系数缩小到 0。然而，请注意这两个列表之间没有重叠！好的，那么让我们尝试将 AIC 模型选择与 LASSO LARS 结合使用LassoLarsIC：

llarsic_selection = SelectFromModel(LassoLarsIC(criterion='aic'))
llarsic_selection.fit(X_train, y_train)
llarsic_cols = X_train.columns[
    llarsic_selection.get_support()
].tolist()
print(len(llarsic_cols))
print(llarsic_cols)

上一段代码生成以下输出：

111
['TCODE', 'STATE', 'MAILCODE', 'RECINHSE', 'RECP3', 'RECPGVG', 'RECSWEEP',..., 'DOMAIN_URBANICITY', 'DOMAIN_SOCIALCLS', 'ZIP_LON']

这是一种相同的算法，但采用了不同的方法来选择正则化参数的值。注意这种不那么保守的方法将特征数量扩展到 111 个。到目前为止，我们使用的方法都具有 L1 范数。让我们尝试一个使用 L2 的——更具体地说，是 L2 惩罚逻辑回归。我们做的是之前所做的，但这次，我们使用二元分类目标（y_train_class）进行拟合：

log_selection = SelectFromModel(
    LogisticRegression(
        C=0.0001,
        solver='sag',
        penalty='l2',
        n_jobs=-1,
        random_state=rand
    )
)
log_selection.fit(X_train, y_train_class)
log_cols = X_train.columns[log_selection.get_support()].tolist()
print(len(log_cols))
print(log_cols)

上一段代码生成以下输出：

87
['ODATEDW', 'TCODE', 'STATE', 'POP901', 'POP902', 'POP903', 'ETH1', 'ETH2', 'ETH5', 'CHIL1', 'HHN2',..., 'AMT_7', 'ZIP_LON']

现在我们有几个特征子集要测试，我们可以将它们的名称放入一个列表（fsnames）中，将特征子集列表放入另一个列表（fscols）中：

fsnames = ['e-lasso', 'e-llars', 'e-llarsic', 'e-logl2']
fscols = [lasso_cols, llars_cols, llarsic_cols, log_cols]

然后，我们可以遍历所有列表名称，并在每次迭代中增加max_depth，就像我们之前做的那样来拟合和评估我们的XGBRFRegressor模型：

def train_mdls_with_fs(reg_mdls, fsnames, fscols, depths):
    for i, fsname in tqdm(enumerate(fsnames), total=len(fsnames)):
       depth = depths[i]
       cols = fscols[i]
       mdlname = 'rf_'+str(depth)+'_'+fsname
       stime = timeit.default_timer()
       reg_mdl = xgb.XGBRFRegressor(
           max_depth=depth, n_estimators=200, seed=rand
       )
       fitted_mdl = reg_mdl.fit(X_train[cols], y_train)
       reg_mdls[mdlname]['num_feat'] = sum(
           reg_mdls[mdlname]['fitted'].feature_importances_ > 0
       )
train_mdls_with_fs(reg_mdls, fsnames, fscols, [3, 4, 5, 6])

现在，让我们看看我们的嵌入式特征选择模型与过滤模型相比的表现如何。我们将重新运行之前运行的代码，输出图 10.6中显示的内容。这次，我们将得到图 10.7中显示的内容：

表描述自动生成，置信度中等

图 10.7：比较所有基础模型和基于过滤和嵌入式特征选择模型的指标

根据图 10.7，我们尝试的四种嵌入式方法中有三种产生了具有最低测试 RMSE（rf_5_e-llarsic、rf_e-lasso和rf_4_e-llars）的模型。它们也都比其他模型训练得快得多，并且比任何同等复杂性的模型都更有利可图。其中之一（rf_5_e-llarsic）甚至非常有利可图。与具有相似测试盈利能力的rf_9_all进行比较，看看性能如何从训练数据中偏离。

发现包装、混合和高级特征选择方法

到目前为止研究的特征选择方法在计算上成本较低，因为它们不需要模型拟合或拟合更简单的白盒模型。在本节中，我们将了解其他更全面的方法，这些方法具有许多可能的调整选项。这里包括的方法类别如下：

包装：通过使用测量指标改进的搜索策略来拟合机器学习模型，彻底搜索最佳特征子集。
混合：一种结合嵌入式和过滤方法以及包装方法的方法。
高级：一种不属于之前讨论的任何类别的的方法。例如包括降维、模型无关特征重要性和遗传算法（GAs）。

现在，让我们开始包装方法吧！

包装方法

包装方法背后的概念相当简单：评估特征的不同子集在机器学习模型上的表现，并选择在预定的目标函数上实现最佳得分的那个。这里变化的是搜索策略：

顺序正向选择（SFS）：这种方法开始时没有特征，然后每次添加一个。
顺序正向浮点选择（SFFS）：与之前相同，除了每次添加一个特征时，它可以移除一个特征，只要目标函数增加。
顺序向后选择（SBS）：这个过程从所有特征都存在开始，每次消除一个特征。
顺序浮点向后选择（SFBS）：与之前相同，除了每次移除一个特征时，它还可以添加一个特征，只要目标函数增加。
穷举特征选择（EFS）：这种方法寻求所有可能的特征组合。
双向搜索（BDS）：这个方法同时允许向前和向后进行函数选择，以获得一个独特的解决方案。

这些方法是贪婪算法，因为它们逐个解决问题，根据它们的即时利益选择部分。尽管它们可能达到全局最大值，但它们采取的方法更适合寻找局部最大值。根据特征的数量，它们可能过于计算密集，以至于不实用，特别是 EFS，它呈指数增长。另一个重要的区别是，向前方法随着特征的添加而提高准确性，而向后方法则随着特征的移除而监控准确性下降。为了缩短搜索时间，我们将做两件事：

我们从其他方法共同选出的特征开始搜索，以拥有更小的特征空间进行选择。为此，我们将来自几种方法的特征列表合并成一个单一的top_cols列表：
```
top_cols = list(set(mic_cols).union(set(llarsic_cols)\
).union(set(log_cols)))
len(top_cols) 
```

样本我们的数据集，以便机器学习模型加速。我们可以使用np.random.choice进行随机选择行索引，而不进行替换：

sample_size = 0.1
sample_train_idx = np.random.choice(
    X_train.shape[0],
    math.ceil(X_train.shape[0]*sample_size),
    replace=False
)
sample_test_idx = np.random.choice(
    X_test.shape[0],
    math.ceil(X_test.shape[0]*sample_size),
    replace=False
)

在所提出的包装方法中，我们只执行 SFS，因为它们非常耗时。然而，对于更小的数据集，你可以尝试其他选项，这些选项mlextend库也支持。

顺序前向选择（SFS）

包装方法的第一参数是一个未拟合的估计器（一个模型）。在SequentialFeatureSelector中，我们放置了一个LinearDiscriminantAnalysis模型。其他参数包括方向（forward=true），是否浮动（floating=False），这意味着它可能会撤销之前对特征的排除或包含，我们希望选择的特征数量（k_features=27），交叉验证的数量（cv=3），以及要使用的损失函数（scoring=f1）。一些推荐的可选参数包括详细程度（verbose=2）和并行运行的工作数量（n_jobs=-1）。由于它可能需要一段时间，我们肯定希望它输出一些内容，并尽可能多地使用处理器：

sfs_lda = SequentialFeatureSelector(
    LinearDiscriminantAnalysis(n_components=1),
    forward=True,
    floating=False,
    k_features=100,
    cv=3,
    scoring='f1',
    verbose=2,
    n_jobs=-1
)
sfs_lda = sfs_lda.fit(X_train.iloc[sample_train_idx][top_cols],\
                      y_train_class[sample_train_idx])
sfs_lda_cols = X_train.columns[list(sfs_lda.k_feature_idx_)].tolist()

一旦我们拟合了 SFS，它将返回使用k_feature_idx_选定的特征的索引，我们可以使用这些索引来子集列并获取特征名称列表。

混合方法

从 435 个特征开始，仅 27 个特征子集的组合就有超过 10⁴²种！所以，你可以看到在如此大的特征空间中 EFS 是如何不切实际的。因此，除了在整个数据集上使用 EFS 之外，包装方法不可避免地会采取一些捷径来选择特征。无论你是向前、向后还是两者都进行，只要你不评估每个特征的组合，你就很容易错过最佳选择。

然而，我们可以利用包装方法的更严格、更全面的搜索方法，同时结合筛选和嵌入方法的效率。这种方法的结果是混合方法。例如，你可以使用筛选或嵌入方法仅提取前 10 个特征，并在这些特征上仅执行 EFS 或 SBS。

递归特征消除（RFE）

另一种更常见的方法是 SBS，但它不是仅基于改进一个指标来删除特征，而是使用模型的内在参数来对特征进行排序，并仅删除排名最低的特征。这种方法被称为递归特征消除（RFE），它是嵌入和包装方法之间的混合。我们只能使用具有feature_importances_或系数（coef_）的模型，因为这是该方法知道要删除哪些特征的方式。具有这些属性的 scikit-learn 模型类别被归类为linear_model、tree和ensemble。此外，XGBoost、LightGBM 和 CatBoost 的 scikit-learn 兼容版本也具有feature_importances_。

我们将使用交叉验证版本的递归特征消除（RFE），因为它更可靠。RFECV首先采用估计器（LinearDiscriminantAnalysis）。然后我们可以定义step，它设置每次迭代应删除多少特征，交叉验证的次数（cv），以及用于评估的指标（scoring）。最后，建议设置详细程度（verbose=2）并尽可能利用更多处理器（n_jobs=-1）。为了加快速度，我们将再次使用样本进行训练，并从top_cols的 267 开始：

rfe_lda = RFECV(
    LinearDiscriminantAnalysis(n_components=1),
    step=2, cv=3, scoring='f1', verbose=2, n_jobs=-1
)
rfe_lda.fit(
    X_train.iloc[sample_train_idx][top_cols],
    y_train_class[sample_train_idx]
)
rfe_lda_cols = np.array(top_cols)[rfe_lda.support_].tolist()

接下来，我们将尝试与主要三个特征选择类别（筛选、嵌入和包装）无关的不同方法。

高级方法

许多方法可以归类为高级特征选择方法，包括以下子类别：

模型无关特征重要性：任何在第四章、全局模型无关解释方法中提到的特征重要性方法都可以用于获取模型的特征选择中的顶级特征。
遗传算法：这是一种包装方法，因为它“包装”了一个评估多个特征子集预测性能的模型。然而，与我们所检查的包装方法不同，它并不总是做出最局部最优的选择。它更适合与大型特征空间一起工作。它被称为遗传算法，因为它受到了生物学的启发——自然选择，特别是。
降维：一些降维方法，如主成分分析（PCA），可以在特征基础上返回解释方差。对于其他方法，如因子分析，它可以从其他输出中推导出来。解释方差可以用于对特征进行排序。
自动编码器：我们不会深入探讨这一点，但深度学习可以利用自动编码器进行特征选择。这种方法在 Google Scholar 上有许多变体，但在工业界并不广泛采用。

在本节中，我们将简要介绍前两种方法，以便您了解它们如何实现。让我们直接进入正题！

模型无关特征重要性

在这本书的整个过程中，我们使用的一个流行的模型无关特征重要性方法是 SHAP，它有许多属性使其比其他方法更可靠。在下面的代码中，我们可以使用TreeExplainer提取我们最佳模型的shap_values：

fitted_rf_mdl = reg_mdls['rf_11_all']['fitted']
shap_rf_explainer = shap.TreeExplainer(fitted_rf_mdl)
shap_rf_values = shap_rf_explainer.shap_values(
    X_test_orig.iloc[sample_test_idx]
)
shap_imps = pd.DataFrame(
    {'col':X_train_orig.columns, 'imp':np.abs(shap_rf_values).mean(0)}
).sort_values(by='imp',ascending=False)
shap_cols = shap_imps.head(120).col.tolist()

然后，SHAP 值绝对值的平均值在第一维上为我们提供了每个特征的排名。我们将这个值放入一个 DataFrame 中，并按我们为 PCA 所做的方式对其进行排序。最后，也将前 120 个放入一个列表（shap_cols）。

遗传算法

算法遗传学（GAs）是一种受自然选择启发的随机全局优化技术，它像包装方法一样包装一个模型。然而，它们不是基于一步一步的序列。GAs 没有迭代，但有代，包括染色体的种群。每个染色体是特征空间的二进制表示，其中 1 表示选择一个特征，0 表示不选择。每一代都是通过以下操作产生的：

选择：就像自然选择一样，这部分是随机的（探索）和部分是基于已经有效的东西（利用）。有效的是其适应性。适应性是通过一个“scorer”来评估的，就像包装方法一样。适应性差的染色体被移除，而好的染色体则通过“交叉”繁殖。
交叉：随机地，一些好的位（或特征）从每个父代传递给子代。
变异：即使染色体已经证明有效，给定一个低的突变率，它偶尔也会突变或翻转其位之一，换句话说，特征。

我们将要使用的 Python 实现有许多选项。在这里我们不会解释所有这些选项，但如果您感兴趣，它们在代码中都有很好的文档说明。第一个属性是估计器。我们还可以定义交叉验证迭代次数（cv=3）和scoring来决定染色体是否适合。有一些重要的概率属性，例如突变位（mutation_probability）的概率和位交换（crossover_probability）的概率。在每一代中，n_gen_no_change提供了一种在代数没有改进时提前停止的手段，默认的generations是 40，但我们将使用 5。我们可以像任何模型一样拟合GeneticSelectionCV。这可能需要一些时间，因此最好定义详细程度并允许它使用所有处理能力。一旦完成，我们可以使用布尔掩码（support_）来子集特征：

ga_rf = GAFeatureSelectionCV(
    RandomForestRegressor(random_state=rand, max_depth=3),
    cv=3,
    scoring='neg_root_mean_squared_error',
    crossover_probability=0.8,
    mutation_probability=0.1,
    generations=5, n_jobs=-1
)
ga_rf = ga_rf.fit(
    X_train.iloc[sample_train_idx][top_cols].values,
    y_train[sample_train_idx]
)
ga_rf_cols = np.array(top_cols)[ga_rf.best_features_].tolist()

好的，现在我们已经在本节中介绍了各种包装、混合和高级特征选择方法，让我们一次性评估它们并比较结果。

评估所有特征选择模型

正如我们对待嵌入方法一样，我们可以将特征子集名称 (fsnames)、列表 (fscols) 和相应的 depths 放入列表中：

fsnames = ['w-sfs-lda', 'h-rfe-lda', 'a-shap', 'a-ga-rf']
fscols = [sfs_lda_cols, rfe_lda_cols, shap_cols, ga_rf_cols]
depths = [5, 6, 5, 6]

然后，我们可以使用我们创建的两个函数，首先遍历所有特征子集，用它们训练和评估一个模型。然后第二个函数输出评估结果，以 DataFrame 的形式包含先前训练的模型：

train_mdls_with_fs(reg_mdls, fsnames, fscols, depths) 
display_mdl_metrics(reg_mdls, 'max_profit_test', max_depth=7) 
Figure 10.8:

图形用户界面，应用程序描述自动生成

图 10.8：比较所有特征选择模型的指标

图 10.8 展示了与包含所有特征相比，特征选择模型在相同深度下的盈利能力更强。此外，嵌入的 LASSO LARS 与 AIC (e-llarsic) 方法和 MIC (f-mic) 过滤方法在相同深度下优于所有包装、混合和高级方法。尽管如此，我们还是通过使用训练数据集的一个样本来阻碍了这些方法，这是加快过程所必需的。也许在其他情况下，它们会优于最顶尖的模型。然而，接下来的三种特征选择方法竞争力相当强：

基于 LDA 的 RFE：混合方法 (h-rfe-lda)
带有 L2 正则化的逻辑回归：嵌入方法 (e-logl2)
基于 RF 的 GAs：高级方法 (a-ga-rf)

在这本书中回顾的方法有很多变体，花费很多天去运行这些变体是有意义的。例如，也许 RFE 与 L1 正则化的逻辑回归或 GA 与支持向量机以及额外的突变会产生最佳模型。有如此多的不同可能性！然而，如果你被迫仅基于 图 10.8 中的利润来做出推荐，那么 111 特征的 e-llarsic 是最佳选择，但它也有比任何顶级模型更高的最低成本和更低的最高回报率。这是一个权衡。尽管它的测试 RMSE 值最高，但 160 特征的模型 (f-mic) 在最大利润训练和测试之间的差异相似，并且在最大回报率和最低成本方面超过了它。因此，这两个选项是合理的。但在做出最终决定之前，必须将不同阈值下的盈利能力进行比较，以评估每个模型在什么成本和回报率下可以做出最可靠的预测。

考虑特征工程

假设非营利组织选择了使用具有 LASSO LARS 与 AIC (e-llarsic) 选择特征的模型，但想评估你是否可以进一步改进它。现在你已经移除了可能只略微提高预测性能但主要增加噪声的 300 多个特征，你剩下的是更相关的特征。然而，你也知道，e-llars 选出的 8 个特征产生了与 111 个特征相同的 RMSE。这意味着虽然那些额外特征中有些东西可以提高盈利能力，但它并没有提高 RMSE。

从特征选择的角度来看，可以采取许多方法来解决这个问题。例如，检查e-llarsic和e-llars之间特征的交集和差异，并在这些特征上严格进行特征选择，以查看 RMSE 是否在任何组合中下降，同时保持或提高当前的盈利能力。然而，还有一种可能性，那就是特征工程。在这个阶段进行特征工程有几个重要的原因：

使模型解释更容易理解：例如，有时特征有一个不直观的尺度，或者尺度是直观的，但分布使得理解变得困难。只要对这些特征的转换不会降低模型性能，转换特征以更好地理解解释方法的输出是有价值的。随着你在更多工程化特征上训练模型，你会意识到什么有效以及为什么有效。这将帮助你理解模型，更重要的是，理解数据。
对单个特征设置护栏：有时，特征分布不均匀，模型倾向于在特征直方图的稀疏区域或存在重要异常值的地方过拟合。
清理反直觉的交互：一些模型发现的不合逻辑的交互，仅因为特征相关，但并非出于正确的原因而存在。它们可能是混淆变量，甚至可能是冗余的（例如我们在第四章，全局模型无关解释方法中找到的）。你可以决定设计一个交互特征或删除一个冗余的特征。

关于最后两个原因，我们将在第十二章，单调约束和模型调优以实现可解释性中更详细地研究特征工程策略。本节将重点介绍第一个原因，尤其是因为它是一个很好的起点，因为它将允许你更好地理解数据，直到你足够了解它，可以做出更转型的改变。

因此，我们剩下 111 个特征，但不知道它们如何与目标或彼此相关。我们首先应该做的是运行一个特征重要性方法。我们可以在e-llarsic模型上使用 SHAP 的TreeExplainer。TreeExplainer的一个优点是它可以计算 SHAP 交互值，shap_interaction_values。与shap_values输出一个(N, 111)维度的数组不同，其中N是观察数量，它将输出(N, 111, 111)。我们可以用它生成一个summary_plot图，该图对单个特征和交互进行排名。交互值唯一的区别是您使用plot_type="compact_dot"：

winning_mdl = 'rf_5_e-llarsic'
fitted_rf_mdl = reg_mdls[winning_mdl]['fitted']
shap_rf_explainer = shap.TreeExplainer(fitted_rf_mdl)
shap_rf_interact_values = \
    shap_rf_explainer.shap_interaction_values(
        X_test.iloc[sample_test_idx][llarsic_cols]
    )
shap.summary_plot(
    shap_rf_interact_values,
    X_test.iloc[sample_test_idx][llarsic_cols],
    plot_type="compact_dot",
    sort=True
) 
Figure 10.9:

图形用户界面，应用程序，表格自动生成的描述

图 10.9：SHAP 交互总结图

我们可以像阅读任何总结图一样阅读图 10.9，除了它包含了两次双变量交互——首先是一个特征，然后是另一个。例如，MDMAUD_A* - CLUSTER是从MDMAUD_A的角度来看该交互的交互 SHAP 值，因此特征值对应于该特征本身，但 SHAP 值是针对交互的。我们在这里可以达成一致的是，考虑到重要性值的规模和比较无序的双变量交互的复杂性，这个图很难阅读。我们将在稍后解决这个问题。

在这本书中，带有表格数据的章节通常以数据字典开始。这个例外是因为一开始有 435 个特征。现在，至少了解哪些是顶级特征是有意义的。完整的数据字典可以在kdd.ics.uci.edu/databases/kddcup98/epsilon_mirror/cup98dic.txt找到，但由于分类编码，一些特征已经发生了变化，因此我们将在这里更详细地解释它们：

MAXRAMNT: 连续型，迄今为止最大赠礼的美元金额
HVP2: 离散型，捐赠者社区中价值>= $150,000 的房屋比例（值在 0 到 100 之间）
LASTGIFT: 连续型，最近一次赠礼的美元金额
RAMNTALL: 连续型，迄今为止终身赠礼的美元金额
AVGGIFT: 连续型，迄今为止赠礼的平均美元金额
MDMAUD_A: 序数型，对于在其捐赠历史中任何时间点都捐赠了 $100+赠礼的捐赠者的捐赠金额代码（值在 0 到 3 之间，对于从未超过$ 100 的捐赠者为-1）。金额代码是RFA（最近/频率/金额）主要客户矩阵代码的第三个字节，即捐赠的金额。类别如下：

0：少于$100（低金额）

1：$100 – 499（核心）

2: $500 – 999 (major)

3: $1,000 + (top)

NGIFTALL: 离散型，迄今为止终身赠礼的数量
AMT_14: 序数型，14 次之前推广的 RFA 捐赠金额代码，这对应于当时最后一次捐赠的金额：

0: $0.01 – 1.99

1: $2.00 – 2.99

2: $3.00 – 4.99

3: $5.00 – 9.99

4: $10.00 – 14.99

5: $15.00 – 24.99

6: $25.00 及以上

DOMAIN_SOCIALCLS: 名义型，社会经济地位（SES）的社区，它与DOMAIN_URBANICITY（0：城市，1：城市，2：郊区，3：镇，4：农村）结合，意味着以下：

1: 最高社会经济地位

2: 平均社会经济地位，但城市社区的平均水平以上

3: 最低社会经济地位，但城市社区的平均水平以下

4: 仅城市社区最低社会经济地位

CLUSTER: 名义型，表示捐赠者所属的集群组的代码
MINRAMNT: 连续型，迄今为止最小赠礼的美元金额
LSC2: 离散型，捐赠者社区中西班牙语家庭的比例（值在 0 到 100 之间）
IC15：离散值，捐赠者所在地区家庭收入低于$15,000 的家庭百分比（值在 0 到 100 之间）

可以从前面的字典和图 10.9中提炼出以下见解：

赠款金额优先：其中七个顶级功能与赠款金额相关，无论是总额、最小值、最大值、平均值还是最后值。如果你包括赠款总数（NGIFTALL），则有八个特征涉及捐赠历史，这完全合理。那么，这有什么相关性呢？因为这些特征很可能高度相关，理解它们可能是提高模型的关键。也许可以创建其他特征，更好地提炼这些关系。
连续赠款金额特征的值较高具有高 SHAP 值：像这样绘制任何这些特征的箱线图，例如plt.boxplot(X_test.MAXRAMNT)，你会看到这些特征是如何右偏斜的。也许通过将它们分成区间——称为“离散化”——或使用不同的尺度，如对数尺度（尝试plt.boxplot(np.log(X_test.MAXRAMNT))），可以帮助解释这些特征，同时也有助于找到捐赠可能性显著增加的区域。
与第十四次促销的关系：他们在两年前进行的促销与数据集标签中标记的促销之间发生了什么？促销材料是否相似？是否每两年发生一次季节性因素？也许你可以设计一个特征来更好地识别这种现象。
分类不一致：DOMAIN_SOCIALCLS根据DOMAIN_URBANITY值的不同而具有不同的类别。我们可以通过使用量表中的所有五个类别（最高、高于平均水平、平均水平、低于平均水平、最低）来使这一分类一致，即使这意味着非城市捐赠者只会使用三个类别。这样做的好处是更容易解释，而且不太可能对模型的性能产生不利影响。

SHAP 交互摘要图可以用来识别特征和交互排名以及它们之间的某些共同点，但在这种情况下（见图 10.9），阅读起来很困难。但要深入挖掘交互，你首先需要量化它们的影响。为此，让我们创建一个热图，只包含按其平均绝对 SHAP 值（shap_rf_interact_avgs）测量的顶级交互。然后，我们应该将所有对角线值设置为 0（shap_rf_interact_avgs_nodiag），因为这些不是交互，而是特征 SHAP 值，没有它们更容易观察交互。我们可以将这个矩阵放入 DataFrame 中，但它是一个有 111 列和 111 行的 DataFrame，所以为了过滤出具有最多交互的特征，我们使用scipy的rankdata对它们求和并排名。然后，我们使用排名来识别 12 个最具交互性的特征（most_interact_cols），并按这些特征子集 DataFrame。最后，我们将 DataFrame 绘制成热图：

shap_rf_interact_avgs = np.abs(shap_rf_interact_values).mean(0)
shap_rf_interact_avgs_nodiag = shap_rf_interact_avgs.copy()
np.fill_diagonal(shap_rf_interact_avgs_nodiag, 0)
shap_rf_interact_df = pd.DataFrame(shap_rf_interact_avgs_nodiag)
shap_rf_interact_df.columns = X_test[llarsic_cols].columns
shap_rf_interact_df.index = X_test[llarsic_cols].columns
shap_rf_interact_ranks = 112 -rankdata(np.sum(
     shap_rf_interact_avgs_nodiag, axis=0)
)
most_interact_cols = shap_rf_interact_df.columns[
    shap_rf_interact_ranks < 13
]
shap_rf_interact_df = shap_rf_interact_df.loc[
most_interact_cols,most_interact_cols
]
sns.heatmap(
    shap_rf_interact_df,
    cmap='Blues',
    annot=True,
    annot_kws={'size':10},
    fmt='.3f',
    linewidths=.5
) 
Figure 10.10. It depicts the most salient feature interactions according to SHAP interaction absolute mean values. Note that these are averages, so given how right-skewed most of these features are, it is likely much higher for many observations. However, it’s still a good indication of relative impact:

图表描述自动生成

图 10.10：SHAP 交互热图

我们可以通过 SHAP 的dependence_plot逐个理解特征交互。例如，我们可以选择我们的顶级特征MAXRAMNT，并将其与RAMNTALL、LSC4、HVP2和AVGGIFT等特征进行颜色编码的交互绘图。但首先，我们需要计算shap_values。然而，还有一些问题需要解决，我们之前已经提到了。这些问题与以下内容有关：

异常值的普遍性：我们可以通过使用特征和 SHAP 值的百分位数来限制x轴和y轴，分别用plt.xlim和plt.ylim来将这些异常值从图中剔除。这本质上是在 1st 和 99th 百分位数之间的案例上进行放大。
金额特征的偏斜分布：在涉及金钱的任何特征中，它通常是右偏斜的。有许多方法可以简化它，例如使用百分位数对特征进行分箱，但一个快速的方法是使用对数刻度。在matplotlib中，您可以通过plt.xscale('log')来实现这一点，而无需转换特征。

以下代码考虑了两个问题。您可以尝试取消注释xlim、ylim或xscale，以查看它们各自在理解dependence_plot时产生的巨大差异：

shap_rf_values = shap_rf_explainer.shap_values(
    X_test.iloc[sample_test_idx] [llarsic_cols]
)
maxramt_shap = shap_rf_values[:,llarsic_cols.index("MAXRAMNT")]
shap.dependence_plot(
    "MAXRAMNT",
    shap_rf_values,
    X_test.iloc[sample_test_idx][llarsic_cols],
    interaction_index="AVGGIFT",
    show=False, alpha=0.1
)
plt.xlim(xmin=np.percentile(X_test.MAXRAMNT, 1),\
         xmax=np.percentile(X_test.MAXRAMNT, 99))
plt.ylim(ymin=np.percentile(maxramt_shap, 1),\
         ymax=np.percentile(maxramt_shap, 99))
plt.xscale('log')

上一段代码生成了图 10.11中所示的内容。它显示了MAXRAMNT在 10 到 100 之间有一个转折点，模型输出的平均影响开始逐渐增加，这些与更高的AVGGIFT值相关：

图表描述自动生成，中等置信度

图 10.11：MAXRAMNT 和 AVGGIFT 之间的 SHAP 交互图

从图 10.11中可以得到的教训是，这些特征的一定值以及可能的一些其他值可以增加捐赠的可能性，从而形成一个簇。从特征工程的角度来看，我们可以采用无监督方法，仅基于您已识别为相关的少数特征来创建特殊的簇特征。或者，我们可以采取更手动的方法，通过比较不同的图表来了解如何最好地识别簇。我们可以从这个过程中推导出二元特征，甚至可以推导出特征之间的比率，这些比率可以更清楚地描述交互或簇归属。

这里的想法不是试图重新发明轮子，去做模型已经做得很好的事情，而是首先追求一个更直观的模型解释。希望这甚至可以通过整理特征对预测性能产生积极影响，因为如果您更好地理解它们，也许模型也会！这就像平滑一个颗粒感强的图像；它可能会让您和模型都少一些困惑（有关更多信息，请参阅第十三章，对抗鲁棒性）！但通过模型更好地理解数据还有其他积极的影响。

事实上，课程不仅仅是关于特征工程或建模，还可以直接应用于促销活动。如果能够识别出转折点，能否用来鼓励捐款呢？或许如果你捐款超过X美元，就可以获得一个免费的杯子？或者设置一个每月捐款X美元的定期捐款，并成为“银牌”赞助者的专属名单之一？

我们将以这个好奇的笔记结束这个话题，但希望这能激发你去欣赏我们如何将模型解释的教训应用到特征选择、工程以及更多方面。

任务完成

为了完成这个任务，你主要使用特征选择工具集来减少过拟合。非营利组织对大约 30%的利润提升感到满意，总成本为 35,601 美元，比向测试数据集中的每个人发送邮件的成本低 30,000 美元。然而，他们仍然希望确保他们可以安全地使用这个模型，而不用担心会亏损。

在本章中，我们探讨了过拟合如何导致盈利曲线不一致。不一致性是关键的，因为它可能意味着基于训练数据选择的阈值在样本外数据上不可靠。因此，你使用compare_df_plots来比较测试集和训练集之间的盈利，就像你之前做的那样，但这次是为了选定的模型（rf_5_e-llarsic）：

profits_test = reg_mdls['rf_5_e-llarsic']['profits_test']
profits_train = reg_mdls['rf_5_e-llarsic']['profits_train']
mldatasets.compare_df_plots(
    profits_test[['costs', 'profit', 'roi']],
    profits_train[['costs', 'profit', 'roi']],
    'Test',
    'Train',
    x_label='Threshold',
    y_formatter=y_formatter,
    plot_args={'secondary_y':'roi'}
)

上述代码生成了图 10.12中所示的内容。你可以向非营利组织展示，以证明在测试中，0.68 美元是一个甜点，是可获得的第二高利润。它也在他们的预算范围内，实现了 41%的投资回报率。更重要的是，这些数字与训练数据非常接近。另一个令人高兴的是，训练和测试的利润曲线缓慢下降，而不是突然跌落悬崖。非营利组织可以确信，如果他们选择提高阈值，运营仍然会盈利。毕竟，他们希望针对整个邮件列表的捐赠者，为了使这从财务上可行，他们必须更加专属。比如说，他们在整个邮件列表上使用 0.77 美元的阈值，活动成本约为 46,000 美元，但利润超过 24,000 美元：

图表，折线图自动生成的描述

图 10.12：通过 AIC 特征在不同阈值下，模型使用 LASSO LARS 的测试集和训练集的盈利、成本和投资回报率比较

恭喜！你已经完成了这个任务！

但有一个关键细节，如果我们不提出来，我们可能会疏忽。

尽管我们考虑到下一场活动来训练这个模型，但这个模型很可能会在未来直接营销活动中使用，而无需重新训练。这种模型的重用带来一个问题。有一个概念叫做数据漂移，也称为特征漂移，即随着时间的推移，模型关于目标变量特征的所学内容不再成立。另一个概念，概念漂移，是关于目标特征定义随时间变化的情况。例如，构成有利捐赠者的条件可能会改变。这两种漂移可能同时发生，并且涉及人类行为的问题，这是可以预料的。行为受到文化、习惯、态度、技术和时尚的影响，这些总是在不断发展。您可以警告非营利组织，您只能保证模型在下一场活动中是可靠的，但他们无法承担每次都雇佣您进行模型重新训练的费用！

您可以向客户提议创建一个脚本，直接监控他们的邮件列表数据库中的漂移情况。如果它检测到模型使用的特征有显著变化，它将向他们和您发出警报。在这种情况下，您可以触发模型的自动重新训练。然而，如果漂移是由于数据损坏造成的，您将没有机会解决这个问题。即使进行了自动重新训练，如果性能指标没有达到预定的标准，也无法部署。无论如何，您都应该密切关注预测性能，以确保可靠性。可靠性是模型可解释性的一个基本主题，因为它与问责制密切相关。本书不会涵盖漂移检测，但未来的章节将讨论数据增强（第十一章，偏差缓解和因果推断方法）和对抗鲁棒性（第十三章，对抗鲁棒性），这些都关乎可靠性。

摘要

在本章中，我们学习了无关特征如何影响模型结果，以及特征选择如何提供一套工具来解决此问题。然后，我们探讨了这套工具中的许多不同方法，从最基本过滤器方法到最先进的方法。最后，我们讨论了特征工程的可解释性问题。特征工程可以使模型更具可解释性，从而表现更好。我们将在第十二章，单调约束和模型调优以实现可解释性中更详细地介绍这个主题。

在下一章中，我们将讨论偏差缓解和因果推断的方法。

数据集来源

Ling, C. 和 Li, C.，1998 年，《直接营销的数据挖掘：问题和解决方案》。在第四届国际知识发现和数据挖掘会议（KDD’98）论文集中。AAAI 出版社，第 73-79 页：dl.acm.org/doi/10.5555/3000292.3000304
UCI 机器学习仓库，1998，KDD Cup 1998 数据集：archive.ics.uci.edu/ml/datasets/KDD+Cup+1998+Data

进一步阅读

Ross, B.C.，2014，离散和连续数据集之间的互信息。PLoS ONE，9：journals.plos.org/plosone/article?id=10.1371/journal.pone.0087357
Geurts, P., Ernst, D., 和 Wehenkel, L.，2006，极端随机树。Machine Learning，63(1)，3-42：link.springer.com/article/10.1007/s10994-006-6226-1
Abid, A.，Balin, M.F.，和 Zou, J.，2019，用于可微分特征选择和重建的混凝土自编码器。ICML：arxiv.org/abs/1901.09346
Tan, F., Fu, X., Zhang, Y., 和 Bourgeois, A.G.，2008，基于遗传算法的特征子集选择方法。Soft Computing，12，111-120：link.springer.com/article/10.1007/s00500-007-0193-8
Calzolari, M.，2020，10 月 12 日，manuel-calzolari/sklearn-genetic：sklearn-genetic 0.3.0（版本 0.3.0）。Zenodo：doi.org/10.5281/zenodo.4081754

在 Discord 上了解更多

要加入这本书的 Discord 社区——在那里您可以分享反馈、向作者提问，并了解新书发布——请扫描下面的二维码：

packt.link/inml

二维码

第十一章：偏差缓解和因果推断方法

在第六章“锚点和反事实解释”中，我们探讨了公平及其与决策的关系，但仅限于事后模型解释方法。在第十章“用于可解释性的特征选择和工程”中，我们提出了成本敏感性的问题，这通常与平衡或公平相关。在本章中，我们将探讨平衡数据和调整模型以实现公平的方法。

使用信用卡违约数据集，我们将学习如何利用目标可视化工具，如类平衡，来检测不希望的偏差，然后通过预处理方法，如重新加权和不平等影响移除（用于处理过程中）和等概率（用于后处理）来减少它。从第六章的“锚点和反事实解释”和第十章的“用于可解释性的特征选择和工程”等主题扩展，我们还将研究政策决策可能产生意外、反直觉或有害的影响。在假设检验的背景下，一个决策被称为治疗。对于许多决策场景，估计其效果并确保这个估计是可靠的至关重要。

因此，我们将假设针对最脆弱的信用卡违约人群的治疗方法，并利用因果模型来确定其平均治疗效果（ATE）和条件平均治疗效果（CATE）。最后，我们将使用各种方法测试因果假设和估计的稳健性。

这些是我们将要涵盖的主要主题：

检测偏差
缓解偏差
创建因果模型
理解异质治疗效果
测试估计的稳健性

技术要求

本章的示例使用了mldatasets、pandas、numpy、sklearn、lightgbm、xgboost、matplotlib、seaborn、xai、aif360、econml和dowhy库。有关如何安装所有这些库的说明见前言。

本章的代码位于此处：

packt.link/xe6ie

任务

全球流通的信用卡超过 28 亿张，我们每年在它们上的总消费超过 25 万亿美元（美元）（www.ft.com/content/ad826e32-2ee8-11e9-ba00-0251022932c8）。这无疑是一个天文数字，但衡量信用卡行业的规模，不应仅看消费额，而应看债务总额。银行等发卡机构主要通过收取利息来赚取大部分收入。因此，消费者（2022 年）欠下的超过 60 万亿美元的债务，其中信用卡债务占很大一部分，为贷方提供了稳定的利息收入。这可能有利于商业，但也带来了充足的风险，因为如果借款人在本金加运营成本偿还之前违约，贷方可能会亏损，尤其是当他们已经用尽法律途径追讨债务时。

当出现信用泡沫时，这个问题会加剧，因为不健康的债务水平可能会损害贷方的财务状况，并在泡沫破裂时将他们的股东拖下水。2008 年的住房泡沫，也称为次贷危机，就是这种情况。这些泡沫通常始于对增长的投机和对无资格需求的寻求，以推动这种增长。在次贷危机的情况下，银行向那些没有证明还款能力的个人提供了抵押贷款。遗憾的是，他们也针对了少数族裔，一旦泡沫破裂，他们的全部净资产就会被清零。金融危机、萧条以及介于两者之间的每一次灾难，往往以更高的比率影响最脆弱的人群。

信用卡也涉及到了灾难性的泡沫，特别是在 2003 年的韩国(www.bis.org/repofficepubl/arpresearch_fs_200806.10.pdf)和 2006 年的台湾。本章将考察 2005 年的数据，导致台湾信用卡危机。到 2006 年，逾期信用卡债务达到 2680 亿美元，由超过 70 万人欠下。超过 3%的台湾人口甚至无法支付信用卡的最低还款额，俗称为信用卡奴隶。随之而来的是重大的社会影响，如无家可归者数量的急剧增加、毒品走私/滥用，甚至自杀。在 1997 年亚洲金融危机之后，该地区的自杀率稳步上升。2005 年至 2006 年间的 23%的增幅将台湾的自杀率推到了世界第二高。

如果我们将危机追溯到其根本原因，那是因为新发卡银行已经耗尽了一个饱和的房地产市场，削减了获取信用卡的要求，而当时这些信用卡的监管由当局执行得并不好。

这对年轻人影响最大，因为他们通常收入较低，管理资金的经验也较少。2005 年，台湾金融监督管理委员会发布了新规定，提高了信用卡申请人的要求，防止出现新的信用卡奴隶。然而，还需要更多的政策来处理系统中已经存在的债务和债务人。

当局开始讨论创建资产管理公司（AMCs）以从银行的资产负债表中剥离不良债务。他们还希望通过一项债务人还款规定，为谈判合理的还款计划提供一个框架。这两项政策直到 2006 年才被纳入法律。

假设一下，现在是 2005 年 8 月，你从未来带着新颖的机器学习和因果推理方法来到这里！一家台湾银行希望创建一个分类模型来预测将违约的客户。他们为你提供了一份包含 30,000 名信用卡客户的数据库。监管机构仍在起草法律，因此有机会提出既有利于银行又有利于债务人的政策。当法律通过后，他们可以使用分类模型预测哪些债务应该卖给资产管理公司（AMCs），并通过因果模型估计哪些政策将有利于其他客户和银行，但他们希望公平且稳健地完成这项任务——这是你的使命！

方法

银行已经强调，将公平性嵌入到你的方法中是多么重要，因为监管机构和公众普遍希望确保银行不会造成更多的伤害。他们的声誉也依赖于这一点，因为在最近几个月里，媒体一直在无情地指责他们进行不诚实和掠夺性贷款行为，导致消费者失去信任。因此，他们希望使用最先进的稳健性测试来证明规定的政策将减轻问题。你提出的方法包括以下要点：

据报道，年轻放贷人更容易违约还款，因此你预计会发现年龄偏差，但你也会寻找其他受保护群体，如性别，的偏差。
一旦检测到偏差，你可以使用AI Fairness 360（AIF360）库中的预处理、处理和后处理算法来减轻偏差。在这个过程中，你将使用每个算法训练不同的模型，评估它们的公平性，并选择最公平的模型。
为了能够理解政策的影响，银行对一小部分客户进行了一项实验。通过实验结果，你可以通过dowhy库拟合一个因果模型，这将识别出因果效应。这些效应进一步由因果模型分解，以揭示异质的治疗效应。
然后，你可以评估异质处理效果来理解它们并决定哪种处理最有效。
最后，为了确保你的结论是稳健的，你需要用几种方法来反驳结果，看看效果是否仍然存在。

让我们深入探讨！

准备工作

您可以在以下链接找到本例的代码：github.com/PacktPublishing/Interpretable-Machine-Learning-with-Python/blob/master/Chapter11/CreditCardDefaults.ipynb.

加载库

要运行此示例，您需要安装以下库：

mldatasets用于加载数据集
pandas和numpy用于操作数据
sklearn（scikit-learn）、xgboost、aif360和lightgbm用于分割数据和拟合模型
matplotlib、seaborn和xai用于可视化解释
econml和dowhy用于因果推断

您应该首先加载所有这些库，如下所示：

import math
import os
import mldatasets
import pandas as pd
import numpy as np
from tqdm.notebook import tqdm
from sklearn import model_selection, tree
import lightgbm as lgb
import xgboost as xgb
from aif360.datasets import BinaryLabelDataset
from aif360.metrics import BinaryLabelDatasetMetric,\
                           ClassificationMetric
from aif360.algorithms.preprocessing import Reweighing,\
                                           DisparateImpactRemover
from aif360.algorithms.inprocessing import ExponentiatedGradientReduction, GerryFairClassifier
from aif360.algorithms.postprocessing.\
                      calibrated_eq_odds_postprocessing \
                            import CalibratedEqOddsPostprocessing
from aif360.algorithms.postprocessing.eq_odds_postprocessing\
                            import EqOddsPostprocessing
from econml.dr import LinearDRLearner
import dowhy
from dowhy import CausalModel
import xai
from networkx.drawing.nx_pydot import to_pydot
from IPython.display import Image, display
import matplotlib.pyplot as plt
import seaborn as sns

理解和准备数据

我们将数据如下加载到名为ccdefault_all_df的 DataFrame 中：

ccdefault_all_df = mldatasets.load("cc-default", prepare=True)

应该有 30,000 条记录和 31 列。我们可以使用info()来验证这一点，如下所示：

ccdefault_all_df.info()

上述代码输出以下内容：

Int64Index: 30000 entries, 1 to 30000
Data columns (total 31 columns):
#   Column            Non-Null Count  Dtype  
---  ------            --------------  -----  
0   CC_LIMIT_CAT      30000 non-null  int8   
1   EDUCATION         30000 non-null  int8   
2   MARITAL_STATUS    30000 non-null  int8   
3   GENDER            30000 non-null  int8   
4   AGE_GROUP         30000 non-null  int8   
5   pay_status_1      30000 non-null  int8   
6   pay_status_2      30000 non-null  int8   
7   pay_status_3      30000 non-null  int8   
8   pay_status_4      30000 non-null  int8   
9   pay_status_5      30000 non-null  int8   
10  pay_status_6      30000 non-null  int8   
11  paid_pct_1        30000 non-null  float64
12  paid_pct_2        30000 non-null  float64
13  paid_pct_3        30000 non-null  float64
14  paid_pct_4        30000 non-null  float64
15  paid_pct_5        30000 non-null  float64
16  paid_pct_6        30000 non-null  float64
17  bill1_over_limit  30000 non-null  float64
18  IS_DEFAULT        30000 non-null  int8   
19  _AGE              30000 non-null  int16  
20  _spend            30000 non-null  int32  
21  _tpm              30000 non-null  int16  
22  _ppm              30000 non-null  int16  
23  _RETAIL           30000 non-null  int8   
24  _URBAN            30000 non-null  int8   
25  _RURAL            30000 non-null  int8   
26  _PREMIUM          30000 non-null  int8   
27  _TREATMENT        30000 non-null  int8   
28  _LTV              30000 non-null  float64
29  _CC_LIMIT         30000 non-null  int32  
30  _risk_score       30000 non-null  float64

输出检查无误。所有特征都是数值型，没有缺失值，因为我们使用了prepare=True，这确保了所有空值都被填补。分类特征都是int8，因为它们已经被编码。

数据字典

有 30 个特征，但我们不会一起使用它们，因为其中 18 个用于偏差缓解练习，剩下的 12 个以下划线（_）开头的用于因果推断练习。很快，我们将把数据分割成每个练习对应的相应数据集。重要的是要注意，小写特征与每个客户的交易历史有关，而客户账户或目标特征是大写的。

我们将在以下偏差缓解练习中使用以下特征：

CC_LIMIT_CAT: 序数；信用卡额度（_CC_LIMIT）分为八个大致均匀分布的四分位数
EDUCATION: 序数；客户的受教育程度（0: 其他，1: 高中，2: 本科，3: 研究生）
MARITAL_STATUS: 名义变量；客户的婚姻状况（0: 其他，1: 单身，2: 已婚）
GENDER: 名义变量；客户的性别（1: 男，2: 女）
AGE GROUP: 二元变量；表示客户是否属于特权年龄组（1: 特权组（26-47 岁），0: 非特权组（其他所有年龄））
pay_status_1 pay_status_6: 序数；从 2005 年 4 月的pay_status_6到 8 月的还款状态（-1: 按时还款，1: 延迟 1 个月还款，2: 延迟 2 个月还款，8: 延迟 8 个月，9: 延迟 9 个月及以上）
paid_pct_1 paid_pct_6: 连续型；从 4 月到 2005 年 8 月，每月应付款项的百分比，paid_pct_6为 8 月，paid_pct_1为 4 月。
bill1_over_limit: 连续型；2005 年 8 月最后一张账单与相应信用额度的比率
IS_DEFAULT: 二元型；目标变量；客户是否违约

这些是我们将在因果推断练习中使用的特征：

_AGE: 连续型；客户的年龄，单位为年。
_spend: 连续型；每位客户在新台币（NT$）中的消费金额。
_tpm: 连续型；客户在前 6 个月内使用信用卡的月均交易量。
_ppm: 连续型；客户在前 6 个月内使用信用卡的月均购买量。
_RETAIL: 二元型；如果客户是零售客户，而不是通过雇主获得的客户。
_URBAN: 二元型；如果客户是城市客户。
_RURAL: 二元型；如果客户是农村客户。
_PREMIUM: 二元型；如果客户是“高级”客户。高级客户会获得现金返还和其他消费激励。
_TREATMENT: 名义型；针对每位客户指定的干预或政策（-1：非实验部分，0：对照组，1：降低信用额度，2：支付计划，3：支付计划和信用额度）。
_LTV: 连续型；干预的结果，即在过去的 6 个月信用支付行为的基础上，估计的新台币（NT$）终身价值。
_CC_LIMIT: 连续型；客户在治疗前的原始信用卡额度，单位为新台币（NT$）。银行家们预期治疗结果将受到这一特征的极大影响。
_risk_score: 连续型；银行在 6 个月前根据信用卡账单与信用额度比率计算出的每位客户的风险评分。它与bill1_over_limit类似，但它是对 6 个月支付历史记录的加权平均值，且在选择治疗措施前 5 个月产生。

我们将在接下来的几节中更详细地解释因果推断特征及其目的。同时，让我们通过value_counts()按其值分解_TREATMENT特征，以了解我们将如何分割这个数据集，如下所示：

ccdefault_all_df._TREATMENT.value_counts()

上述代码输出了以下内容：

-1    28904
3      274
2      274
1      274
0      274

大多数观测值是治疗-1，因此它们不是因果推断的一部分。其余部分在三种治疗（1-3）和对照组（0）之间平均分配。自然地，我们将使用这四个组进行因果推断练习。然而，由于对照组没有指定治疗措施，我们可以将其与-1治疗措施一起用于我们的偏差缓解练习。我们必须小心排除在偏差缓解练习中行为被操纵的客户。整个目的是在尝试减少偏差的同时，预测在“照常营业”的情况下，哪些客户最有可能违约。

数据准备

目前，我们的单一数据准备步骤是将数据集拆分，这可以通过使用_TREATMENT列对pandas DataFrame 进行子集化轻松完成。我们将为每个练习创建一个 DataFrame，使用这种子集化：偏差缓解（ccdefault_bias_df）和因果推断（ccdefault_causal_df）。这些可以在以下代码片段中看到：

ccdefault_bias_df = ccdefault_all_df[
    ccdefault_all_df._TREATMENT < 1
]
ccdefault_causal_df =ccdefault_all_df[
    ccdefault_all_df._TREATMENT >= 0
]

我们将在深入部分进行一些其他数据准备步骤，但现在我们可以开始着手了！

检测偏差

机器学习中存在许多偏差来源。如第一章中所述，“解释、可解释性和可解释性；以及这一切为什么都重要？”，存在大量的偏差来源。那些根植于数据所代表真相的，如系统性和结构性偏差，导致数据中的偏见偏差。还有根植于数据的偏差，如样本、排除、关联和测量偏差。最后，还有我们从数据或模型中得出的见解中的偏差，我们必须小心处理，如保守主义偏差、显著性偏差和基本归因错误。

对于这个例子，为了正确地解开这么多偏差水平，我们应该将我们的数据与 2005 年台湾的普查数据和按人口统计划分的历史贷款数据联系起来。然后，使用这些外部数据集，控制信用卡合同条件，以及性别、收入和其他人口统计数据，以确定年轻人是否特别被针对，获得他们不应有资格的高利率信用卡。我们还需要追踪数据集到作者那里，并与他们以及领域专家协商，检查数据集中与偏差相关的数据质量问题。理想情况下，这些步骤对于验证假设是必要的，但这将是一个需要几章解释的巨大任务。

因此，本着简便的原则，我们直接接受本章的前提。也就是说，由于掠夺性贷款行为，某些年龄群体更容易受到信用卡违约的影响，这不是他们自己的任何过错。我们还将直接接受数据集的质量。有了这些保留，这意味着如果我们发现数据或由此数据派生的任何模型中年龄群体之间存在差异，这可以归因于掠夺性贷款行为。

这里还概述了两种公平性类型：

程序公平性：这是关于公平或平等对待。在法律上很难定义这个术语，因为它在很大程度上取决于上下文。
结果公平性：这完全是关于衡量公平的结果。

这两个概念并不是相互排斥的，因为程序可能是公平的，但结果可能不公平，反之亦然。在这个例子中，不公平的程序是向不合格的客户提供高利率信用卡。尽管如此，我们将在本章中关注结果公平性。

当我们讨论机器学习中的偏差时，它将影响受保护的群体，并且在这些群体中，将存在特权和弱势群体。后者是一个受到偏差负面影响的群体。偏差的表现方式也有很多，以下是如何应对偏差的：

代表性：可能存在代表性不足或弱势群体过度代表的情况。与其它群体相比，模型将学习关于这个群体的信息要么太少要么太多。
分布：特征在群体之间的分布差异可能导致模型做出有偏的关联，这些关联可能直接或间接地影响模型的结果。
概率：对于分类问题，如第六章中讨论的，群体之间的类别平衡差异可能导致模型学习到某个群体有更高的概率属于某一类或另一类。这些可以通过混淆矩阵或比较它们的分类指标（如假阳性或假阴性率）来轻松观察到。
混合：上述任何表现形式的组合。

在缓解偏差部分讨论了任何偏差表现策略，但我们在本章中讨论的偏差类型与我们的主要受保护属性（_AGE）的概率差异有关。我们将通过以下方式观察这一点：

可视化数据集偏差：通过可视化观察受保护特征的差异。
量化数据集偏差：使用公平性指标来衡量偏差。
量化模型偏差：我们将训练一个分类模型并使用为模型设计的其他公平性指标。

模型偏差可以像我们在第六章，锚点和反事实解释中已经做的那样可视化，或者像我们在第十二章，单调约束和模型调优以提高可解释性中将要做的那样可视化。我们将在本章稍后的一个子节将所有内容结合起来中快速探索一些其他可视化。现在，让我们不耽搁，继续本节的实际部分。

可视化数据集偏差

数据本身讲述了某个群体属于正类与另一类相比的可能性有多大。如果是一个分类特征，可以通过将正类的value_counts()函数除以所有类别的总和来获得这些概率。例如，对于性别，我们可以这样做：

ccdefault_bias_df[
    ccdefault_bias_df.**IS_DEFAULT==1**
].**GENDER.value_counts()**/ccdefault_bias_df.**GENDER.value_counts()**

前面的代码片段产生了以下输出，显示了男性平均有更高的概率违约他们的信用卡：

2    0.206529
1    0.241633

对于连续特征的这种操作代码要复杂一些。建议您首先使用pandas的qcut将特征划分为四分位数，然后使用与分类特征相同的方法。幸运的是，plot_prob_progression函数为您完成了这项工作，并绘制了每个四分位数的概率进展。第一个属性是pandas系列，一个包含受保护特征（_AGE）的数组或列表，第二个是相同的，但用于目标特征（IS_DEFAULT）。然后我们选择要设置为四分位数的区间数（x_intervals）（use_quantiles=True）。

其余的属性是美学上的，例如标签、标题和添加mean_line。代码可以在下面的片段中看到：

mldatasets.**plot_prob_progression(**
    ccdefault_bias_df.**_AGE**,
    ccdefault_bias_df.**IS_DEFAULT**,
    **x_intervals**=8,
    use_quantiles=True,
    xlabel='Age',
    **mean_line**=True,
    title='Probability of Default by Age'
)

上述代码生成了以下输出，展示了最年轻（21-25）和最年长（47-79）的人群最有可能违约。其他所有群体仅代表超过一个标准差：

图表，折线图自动生成的描述

图 11.1：按年龄划分的 CC 违约概率

我们可以将最年轻和最年长的四分位数称为弱势群体，而其他所有人称为特权群体。为了检测和减轻不公平性，最好将它们编码为二元特征——我们正是这样使用AGE_GROUP做到的。我们可以再次利用plot_prob_progression，但这次用AGE_GROUP代替AGE，并将数字替换为我们更容易理解的标签。代码可以在下面的片段中看到：

mldatasets.**plot_prob_progression(**
    ccdefault_bias_df.**AGE_GROUP**.**replace**({0:'21-25,48+',1:'26-47'}),
    ccdefault_bias_df.**IS_DEFAULT**,
    xlabel='Age Group',
    title='Probability of Default by Age Group',
    **mean_line**=True
)

上述片段生成了以下输出，其中两组之间的差异相当明显：

图表，折线图自动生成的描述

图 11.2：按年龄组划分的 CC 违约概率

接下来，让我们将GENDER重新引入画面。我们可以使用plot_prob_contour_map，它类似于plot_prob_progression，但在二维空间中，用颜色编码概率而不是绘制线条。因此，前两个属性是我们希望在x轴（GENDER）和y轴（AGE_GROUP）上的特征，第三个是目标（IS_DEFAULT）。由于我们的特征都是二元的，最好使用plot_type='grid'而不是contour。代码可以在下面的片段中看到：

mldatasets.plot_prob_contour_map(
    ccdefault_bias_df.**GENDER**.replace({1:'Male',2:'Female'}),
    ccdefault_bias_df.**AGE_GROUP**.replace({0:'21-25,48+',1:'26-47'}),
    ccdefault_bias_df.**IS_DEFAULT**,
    xlabel='Gender',
    ylabel='Age Group',
    annotate=True,
    plot_type='grid',
    title='Probability of Default by Gender/Age Group'
) 
group is 26- 47-year-old females, followed by their male counterparts at about 3-4% apart. The same happens with the underprivileged age group:

图表，树状图自动生成的描述

图 11.3：按性别和年龄组划分的 CC 默认概率网格

性别差异是一个有趣的观察结果，我们可以提出许多假设来解释为什么女性违约较少。她们是否只是简单地更擅长管理债务？这是否与她们的婚姻状况或教育有关？我们不会深入探讨这些问题。鉴于我们只知道基于年龄的歧视，我们将在特权组中仅使用AGE_GROUP，但将GENDER保持为受保护属性，这将在我们监控的一些公平性指标中考虑。说到指标，我们将量化数据集偏差。

量化数据集偏差

公平性指标分为三类，如下概述：

个体公平性：个体观察值在数据中与同龄人的接近程度。例如，欧几里得距离和曼哈顿距离等距离度量可以用于此目的。
组公平性：组与组之间标签或结果的平均距离。这可以在数据或模型中进行衡量。
两者都：一些指标通过在组内和组间同时考虑不平等来衡量熵或方差，例如Theil 指数和变异系数。

在本章中，我们将专注于组公平性指标。

在我们计算公平性指标之前，有一些待处理的数据准备步骤。让我们确保我们将用于偏差缓解练习的数据集（ccdefault_bias_df）只包含相关的列，这些列不以下划线（"_"）开头。另一方面，因果推断练习将包括以下划线开头的列以及AGE_GROUP和IS_DEFAULT。代码可以在下面的代码片段中看到：

cols_bias_l = ccdefault_all_df.columns[
    **~ccdefault_all_df.columns.str.startswith('_')**
].tolist()
cols_causal_l = [**'AGE_GROUP'**,**'IS_DEFAULT'**] +\
    ccdefault_all_df.columns[
        **ccdefault_all_df.columns.str.startswith('_')**
    ].tolist()
ccdefault_bias_df = ccdefault_bias_df[**cols_bias_l**]
ccdefault_causal_df = ccdefault_causal_df[**cols_causal_l**]

此外，量化训练数据集中的数据集偏差更为重要，因为这是模型将从中学习的数据，所以让我们继续将数据分成训练集和测试集X和y对。我们在初始化随机种子以实现某些可重复性之后进行此操作。代码可以在下面的代码片段中看到：

rand = 9
os.environ['PYTHONHASHSEED']=str(rand)
np.random.seed(rand)
y = ccdefault_bias_df[**'IS_DEFAULT'**]
X = ccdefault_bias_df.drop([**'IS_DEFAULT'**], axis=1).copy()
X_train, X_test, y_train, y_test = model_selection.**train_test_split**(
    X, y, test_size=0.25, random_state=rand
)

尽管我们将使用我们刚刚分割的pandas数据用于训练和性能评估，但我们将在这次练习中使用的库，称为 AIF360，将数据集抽象为基类。这些类包括转换为numpy数组的数据，并存储与公平性相关的属性。

对于回归，AIF360 有RegressionDataset，但对于这个二元分类示例，我们将使用BinaryLabelDataset。你可以使用包含特征和标签的pandas DataFrame 来初始化它（X_train.join(y_train)）。然后，你指定标签的名称（label_names）和受保护的属性（protected_attribute_names），并且建议你为favorable_label和unfavorable_label输入一个值，这样 AIF360 就可以将其纳入评估公平性的考量。尽管这可能听起来很复杂，但在二元分类中，正数和负数仅与我们要预测的内容相关——正类——而不是它是否是一个有利的结果。代码可以在下面的片段中看到：

train_ds = **BinaryLabelDataset**(
    df=**X_train.join(y_train)**,
    label_names=['IS_DEFAULT'],
    protected_attribute_names=['AGE_GROUP', 'GENDER'],
    favorable_label=0,
    unfavorable_label=1
)
test_ds = **BinaryLabelDataset**(
    **df=X_test.join(y_test)**,
    label_names=['IS_DEFAULT'],
    protected_attribute_names=['AGE_GROUP', 'GENDER'],
    favorable_label=0, unfavorable_label=1
)

接下来，我们为underprivileged groups和privileged_groups创建数组。在AGE_GROUP=1中的成员有较低的违约概率，因此他们是特权组，反之亦然。然后，使用这些数组以及用于训练的抽象数据集(train_ds)，我们可以通过BinaryLabelDatasetMetric初始化一个度量类。这个类有计算几个群体公平度度量的函数，仅凭数据本身进行判断。我们将输出其中的三个，并解释它们的含义。代码可以在下面的片段中看到：

underprivileged_groups=[**{'AGE_GROUP': 0}**]
privileged_groups=[**{'AGE_GROUP': 1}**]
metrics_train_ds = **BinaryLabelDatasetMetric**(
    train_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
print('Statistical Parity Difference (SPD): %.4f' %
      metrics_train_ds.**statistical_parity_difference**())
print('Disparate Impact (DI): %.4f' % 
      metrics_train_ds.**disparate_impact**())
print('Smoothed Empirical Differential Fairness (SEDF): %.4f' %\
      metrics_train_ds.**smoothed_empirical_differential_fairness**())

前面的代码片段生成了以下输出：

Statistical Parity Difference (SPD):               -0.0437
Disparate Impact (DI):                              0.9447
Smoothed Empirical Differential Fairness (SEDF):    0.3514

现在，让我们解释每个度量分别代表什么，如下：

统计差异差异（SPD）：也称为平均差异，这是弱势群体和特权群体之间有利结果平均概率的差异。负数表示对弱势群体的不公平，正数表示更好，但接近零的数字表示公平的结果，特权群体和弱势群体之间没有显著差异。它使用以下公式计算，其中f是有利类的值，D是客户组，Y是客户是否会违约：

差异影响（DI）：DI 与 SPD 完全相同，只是它是比率而不是差异。在比率方面，越接近一，对弱势群体来说越好。换句话说，一代表群体之间公平的结果，没有差异，低于一表示与特权群体相比，弱势群体不利的结果，而高于一表示与特权群体相比，弱势群体有利的结果。公式如下：

平滑经验差异公平性（SEDF）：这个公平性指标是从一篇名为*“交叉性公平性的定义。”的论文中提出的许多新指标之一。与前面两个指标不同，它不仅限于预定的特权群体和弱势群体，而是扩展到包括受保护属性中的所有类别——在本例中是图 11.3*中的四个。论文的作者们认为，当有受保护属性的交叉表时，公平性尤其棘手。这是因为辛普森悖论，即一个群体在总体上可能是有利的或是不利的，但在细分到交叉表时则不是。我们不会深入数学，但他们的方法在测量交叉性场景中的合理公平性水平时考虑到这种可能性。为了解释它，零代表绝对公平，越远离零，公平性越低。

接下来，我们将量化模型的群体公平性指标。

量化模型偏差

在我们计算指标之前，我们需要训练一个模型。为此，我们将使用最佳超参数(lgb_params)初始化一个 LightGBM 分类器(LGBMClassifier)。这些参数已经为我们进行了超参数调整（更多关于如何做这件事的细节在第十二章，单调约束和模型调优以提高可解释性）。

请注意，这些参数包括scale_pos_weight，这是用于类别加权的。由于这是一个不平衡的分类任务，这是一个重要的参数，以便使分类器进行成本敏感训练，对一种误分类形式进行惩罚，而不是另一种。一旦分类器初始化，它将通过evaluate_class_mdl进行fit和评估，该函数返回一个包含预测性能指标的字典，我们可以将其存储在模型字典(cls_mdls)中。代码可以在下面的代码片段中看到：

cls_mdls = {}
lgb_params = {
    'learning_rate': 0.4,
    'reg_alpha': 21,
    'reg_lambda': 1,
    **'scale_pos_weight'**: 1.8
}
lgb_base_mdl = lgb.LGBMClassifier(
    random_seed=rand,
    max_depth=6,
    num_leaves=33,
    **lgb_params
)
lgb_base_mdl.fit(X_train, y_train)
**cls_mdls**['lgb_0_base'] = mldatasets.**evaluate_class_mdl**(
    lgb_base_mdl,
    X_train,
    X_test,
    y_train,
    y_test,
    plot_roc=False,
    plot_conf_matrix=True,
    show_summary=True,
    ret_eval_dict=True
) 
Figure 11.4. The scale_pos_weight parameter ensures a healthier balance between false positives in the top-right corner and false negatives at the bottom left. As a result, precision and recall aren’t too far off from each other. We favor high precision for a problem such as this one because we want to maximize true positives, but, not at the great expense of recall, so a balance between both is critical. While hyperparameter tuning, the F1 score, and the Matthews correlation coefficient (MCC) are useful metrics to use to this end. The evaluation of the LightGBM base model is shown here:

图 11.4：LightGBM 基础模型的评估

接下来，让我们计算模型的公平性指标。为此，我们需要对 AIF360 数据集进行“深度”复制（deepcopy=True），但我们将labels和scores更改为我们的模型预测的值。compute_aif_metrics函数使用 AIF360 的ClassificationMetric类为模型执行BinaryLabelDatasetMetric为数据集所执行的操作。然而，它并不直接与模型交互。它使用原始数据集（test_ds）和修改后的数据集（test_pred_ds）包含模型的预测来计算公平性。compute_aif_metrics函数创建一个包含几个预先计算的指标（metrics_test_dict）和指标类（metrics_test_cls）的字典，可以用来逐个获取指标。代码可以在下面的代码片段中看到：

test_pred_ds = test_ds.copy(deepcopy=True)
test_pred_ds.labels =\
    cls_mdls['lgb_0_base']['preds_test'].reshape(-1,1)
test_pred_ds.scores = \
    cls_mdls['lgb_0_base']['probs_test'].reshape(-1,1)
metrics_test_dict, metrics_test_cls = \
    mldatasets.**compute_aif_metrics**(
        test_ds,
        test_pred_ds,
        unprivileged_groups=underprivileged_groups,
        privileged_groups=privileged_groups
    )
cls_mdls['lgb_0_base'].**update**(metrics_test_dict)
print('Statistical Parity Difference (SPD): %.4f' %
      metrics_test_cls.**statistical_parity_difference**())
print('Disparate Impact (DI): %.4f' %
      metrics_test_cls.**disparate_impact**())
print('Average Odds Difference (AOD): %.4f' %
      metrics_test_cls.**average_odds_difference**())
print('Equal Opportunity Difference (EOD): %.4f' %
      metrics_test_cls.**equal_opportunity_difference**())
print('Differential Fairness Bias Amplification(DFBA): %.4f' % \
    metrics_test_cls.**differential_fairness_bias_amplification**())

前面的代码片段生成了以下输出：

Statistical Parity Difference (SPD):               -0.0679
Disparate Impact (DI):                                      0.9193
Average Odds Difference (AOD):                        -0.0550
Equal Opportunity Difference (EOD):                 -0.0265
Differential Fairness Bias Amplification (DFBA):    0.2328

现在，把我们已经解释过的指标放在一边，让我们解释其他指标的含义，如下所示：

平均机会差异（AOD）：这是对特权组和弱势群体假阳性率（FPR）的平均值与假阴性率（FNR）之间的差异。负值表示弱势群体存在不利，越接近零越好。公式如下：

平等机会差异（EOD）：它只是 AOD 的真正阳性率（TPR）差异，因此它只用于测量 TPR 的机会。与 AOD 一样，负值确认了弱势群体存在不利，值越接近零意味着组间没有显著差异。公式如下：

差异公平性偏差放大（DFBA）：这个指标与 SEDF 在同一篇论文中定义，同样以零作为公平性的基准，并且也是交叉的。然而，它只测量了在称为偏差放大的现象中，模型和数据在不公平性比例上的差异。换句话说，这个值表示模型相对于原始数据增加了多少不公平性。

如果你将模型的SPD和DI指标与数据相比，它们确实更差。这并不奇怪，因为这是预期的，因为模型学习到的表示往往会放大偏差。你可以用DFBA指标来证实这一点。至于AOD和EOD，它们通常与SPD指标在同一区域，但理想情况下，EOD指标应该比AOD指标更接近零，因为我们在这个例子中更关心 TPR。

接下来，我们将介绍减轻模型偏差的方法。

减轻偏差

我们可以通过在以下三个不同层面上操作的方法来在三个不同层面上减轻偏差：

预处理：这些是在训练模型之前检测和消除训练数据中偏差的干预措施。利用预处理的方法的优点是它们在源头解决偏差。另一方面，任何未检测到的偏差仍可能被模型放大。
内处理：这些方法在模型训练期间减轻偏差，因此高度依赖于模型，并且通常不像预处理和后处理方法那样不依赖于模型。它们还需要调整超参数来校准公平性指标。
后处理：这些方法在模型推理期间缓解偏差。在第六章，锚点和反事实解释中，我们提到了使用 What-If 工具来选择正确的阈值（参见该章节中的图 6.13），并且我们手动调整它们以达到与假阳性相同的效果。就像那时一样，后处理方法旨在直接在结果中检测和纠正公平性，但需要进行的调整将取决于哪些指标对你的问题最重要。后处理方法的优势在于它们可以解决结果不公平性，这在可以产生最大影响的地方，但由于它与模型开发的其余部分脱节，可能会扭曲事物。

请注意，偏差缓解方法可能会损害预测性能，因此通常存在权衡。可能会有相反的目标，尤其是在数据反映了有偏见的真相的情况下。我们可以选择追求更好的真相：一个正义的真相——我们想要的，而不是我们拥有的那个。

本节将解释每个级别的几种方法，但只为每种方法实现和评估两个。此外，我们不会在本章中这样做，但你可以将不同类型的方法结合起来以最大化缓解——例如，你可以使用预处理方法来去偏数据，然后用它来训练模型，最后使用后处理方法来移除模型添加的偏差。

预处理偏差缓解方法

这些是一些最重要的预处理或数据特定偏差缓解方法：

无意识：也称为压制。移除偏差最直接的方法是排除数据集中的有偏特征，但这是一种天真方法，因为你假设偏差严格包含在这些特征中。
特征工程：有时，连续特征会捕捉到偏差，因为存在许多稀疏区域，模型可以用假设来填补空白或从异常值中学习。它也可以与交互做同样的事情。特征工程可以设置护栏。我们将在第十二章，单调约束和模型调优以实现可解释性中讨论这个话题。
平衡：也称为重采样。单独来看，通过平衡数据集可以相对容易地修复表示问题。XAI 库（github.com/EthicalML/xai）有一个balance函数，通过随机下采样和上采样组表示来实现这一点。下采样，或欠采样，就是我们通常所说的采样，即只取一定比例的观察结果，而上采样，或过采样，则是创建一定比例的随机重复。一些策略会合成上采样而不是重复，例如合成少数过采样技术（SMOTE）。然而，我们必须警告，如果你有足够的数据，总是优先下采样而不是上采样。如果有其他可能的偏见问题，最好不要只使用平衡策略。
重新标记：也称为调整，这是一种算法改变最可能存在偏见的观察结果的标签，通过排名来产生调整后的数据。通常，这使用朴素贝叶斯分类器执行，为了保持类别分布，它不仅提升了一些观察结果，还降低了一样多的数量。
重新加权：这种方法与重新标记类似，但不是翻转它们的标签，而是为每个观察结果推导出一个权重，我们可以在学习过程中实现它。就像类别权重应用于每个类别一样，样本权重应用于每个观察结果或样本。许多回归器和分类器，包括LGBMClassifier，都支持样本权重。尽管技术上重新加权不接触数据和模型应用到的解决方案，但它是一种预处理方法，因为我们检测到数据中的偏见。
差异影响消除器：这种方法的设计者非常小心，遵守法律对偏见的定义，并在不改变标签或受保护属性的情况下保持数据的完整性。它实施了一个修复过程，试图从剩余的特征中消除偏见。当我们怀疑偏见主要存在于那里时，这是一个非常好的过程——也就是说，特征与受保护属性高度相关，但它不解决其他地方的偏见。在任何情况下，它都是一个很好的基线，用于了解有多少偏见是非受保护特征。
学习公平表示：这种方法利用了对抗性学习框架。有一个生成器（自动编码器）创建排除受保护属性的数据表示，还有一个评论家，其目标是使特权组和弱势群体中学习的表示尽可能接近。
用于歧视预防的优化预处理：这种方法通过数学优化数据，以保持整体概率分布。同时，保护属性与目标之间的相关性被消除。这个过程的结果是数据略微扭曲，以消除偏差。

由于存在许多预处理方法，我们将在本章中仅使用其中两种。尽管如此，如果你对使用我们未涉及的方法感兴趣，它们在 AIF360 库中可用，你可以在其文档中了解它们（aif360.res.ibm.com/)。

重新加权方法

重新加权方法相对简单易行。你通过指定组来初始化它，然后像使用任何 scikit-learn 编码器或缩放器一样fit和transform数据。对于那些不熟悉fit的人来说，算法学习如何转换提供的数据，而transform使用学到的知识来转换它。以下代码片段中可以看到代码：

reweighter= **Reweighing**(
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
reweighter.**fit**(train_ds)
train_rw_ds = reweighter.**transform**(train_ds)

从这个过程得到的转换不会改变数据，但为每个观测值创建权重。AIF360 库能够将这些权重因素纳入公平性的计算中，因此我们可以使用之前使用的BinaryLabelDatasetMetric来计算不同的指标。以下代码片段中可以看到代码：

metrics_train_rw_ds = **BinaryLabelDatasetMetric**(
    train_rw_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
print('Statistical Parity Difference (SPD): %.4f' %
      metrics_train_rw_ds.**statistical_parity_difference**())
print('Disparate Impact (DI): %.4f' %
       metrics_train_rw_ds.**disparate_impact**())
print('Smoothed Empirical Differential Fairness(SEDF): %.4f'%
metrics_train_rw_ds.**smoothed_empirical_differential_fairness**())

上述代码输出以下内容：

Statistical Parity Difference (SPD):                    -0.0000
Disparate Impact (DI):                                   1.0000
Smoothed Empirical Differential Fairness (SEDF):    0.1942

权重对 SPD 和 DI 有完美的影响，使它们在这些指标方面绝对公平。然而，请注意，SEDF 比以前更好，但不是零。这是因为特权群体和弱势群体仅与AGE_GROUP保护属性相关，但不与GENDER相关。SEDF 是交叉公平性的度量，重新加权没有涉及。

你可能会认为给观测值添加权重会对预测性能产生不利影响。然而，这种方法被设计用来保持平衡。在未加权的数据集中，所有观测值都有一个权重为 1，因此所有权重的平均值是 1。在重新加权时，会改变观测值的权重，但平均值仍然大约是 1。你可以通过比较原始数据集和重新加权的数据集中instance_weights的平均值的绝对差异来检查这一点。它应该是微不足道的。以下代码片段中可以看到代码：

np.**abs**(train_ds.**instance_weights.mean**() -\
       train_rw_ds.**instance_weights.mean**()) < 1e-6

那么，你可能会问，如何应用instance_weights？许多模型类在fit方法中有一个不太为人所知的属性，称为sample_weight。你只需将其插入其中，在训练过程中，它将根据相应的权重从观测值中学习。以下代码片段展示了这种方法：

lgb_rw_mdl = lgb.LGBMClassifier(
    random_seed=rand,
    max_depth=6,
    num_leaves=33,
    **lgb_params
)
lgb_rw_mdl.fit(
    X_train,
    y_train,
    **sample_weight**=train_rw_ds.instance_weights
)

我们可以使用与基础模型相同的方法评估此模型，使用evaluate_class_mdl。然而，当我们使用compute_aif_metrics计算公平性指标时，我们将它们保存在模型字典中。我们不会逐个查看每种方法的输出，而是在本节结束时进行比较。以下是一个代码片段：

cls_mdls['lgb_1_rw'] = mldatasets.**evaluate_class_mdl**(
    lgb_rw_mdl,
    train_rw_ds.features,
    X_test,
    train_rw_ds.labels,
    y_test,
    plot_roc=False,
    plot_conf_matrix=True,
    show_summary=True,
    ret_eval_dict=True
)
test_pred_rw_ds = test_ds.copy(deepcopy=True)
test_pred_rw_ds.labels = cls_mdls['lgb_1_rw']['preds_test'
    ].reshape(-1,1)
test_pred_rw_ds.scores = cls_mdls['lgb_1_rw']['probs_test'
    ].reshape(-1,1)
metrics_test_rw_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_rw_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['lgb_1_rw'].update(metrics_test_rw_dict) 
Figure 11.5:

图表，瀑布图描述自动生成

图 11.5：LightGBM 重新加权模型的评估

如果你将图 11.5与图 11.4进行比较，你可以得出结论，重新加权模型和基础模型之间的预测性能没有太大差异。这个结果是可以预料的，但仍然值得验证。一些偏差缓解方法可能会对预测性能产生不利影响，但重新加权并没有。同样，DI 消除器（差异影响，DI）也不应该如此，我们将在下一节中讨论！

差异影响消除方法

此方法专注于不在受保护属性（AGE_GROUP）中的偏差，因此我们将在过程中删除此特征。为此，我们需要它的索引——换句话说，它在列列表中的位置。我们可以将此位置（protected_index）保存为变量，如下所示：

protected_index = train_ds.feature_names.index('AGE_GROUP')

DI 消除器是参数化的。它需要一个介于零和一之间的修复水平，因此我们需要找到最佳值。为此，我们可以遍历一个具有不同修复水平值的数组（levels），使用每个level初始化DisparateImpactRemover，并对数据进行fit_transform，这将消除数据中的偏差。然而，我们随后在不包含受保护属性的情况下训练模型，并使用BinaryLabelDatasetMetric评估disparate_impact。记住，DI 是一个比率，因此它是一个可以在超过和低于一之间的指标，最佳 DI 是最接近一的。因此，当我们遍历不同的修复水平时，我们将持续保存 DI 最接近一的模型。我们还将 DI 追加到数组中，以供以后使用。以下是一个代码片段：

di = np.array([])
train_dir_ds = None
test_dir_ds = None
lgb_dir_mdl = None
X_train_dir = None
X_test_dir = None
levels = np.hstack(
    [np.linspace(0., 0.1, 41), np.linspace(0.2, 1, 9)]
)
for level in tqdm(levels):
    di_remover = **DisparateImpactRemover**(repair_level=level)
    train_dir_ds_i = di_remover.**fit_transform**(train_ds)
    test_dir_ds_i = di_remover.**fit_transform**(test_ds)
    X_train_dir_i = np.**delete**(
        train_dir_ds_i.features,
        protected_index,
        axis=1
    )
    X_test_dir_i = np.**delete**(
        test_dir_ds_i.features,
        protected_index,
        axis=1
    )
    lgb_dir_mdl_i = lgb.**LGBMClassifier**(
        random_seed=rand,
        max_depth=5,
        num_leaves=33,
        **lgb_params
    )
    lgb_dir_mdl_i.**fit**(X_train_dir_i, train_dir_ds_i.labels)
    test_dir_ds_pred_i = test_dir_ds_i.copy()
    test_dir_ds_pred_i.labels = lgb_dir_mdl_i.predict(
        X_test_dir_i
    )
    metrics_test_dir_ds = **BinaryLabelDatasetMetric**(
        test_dir_ds_pred_i,
        unprivileged_groups=underprivileged_groups,
        privileged_groups=privileged_groups
    )
    di_i = metrics_test_dir_ds.disparate_impact()
    if (di.shape[0]==0) or (np.min(np.abs(di-1)) >= abs(di_i-1)):
        print(abs(di_i-1))
        train_dir_ds = train_dir_ds_i
        test_dir_ds = test_dir_ds_i
        X_train_dir = X_train_dir_i
        X_test_dir = X_test_dir_i
        lgb_dir_mdl = lgb_dir_mdl_i
    di = np.append(np.array(di), di_i)

为了观察不同修复水平下的 DI，我们可以使用以下代码，如果你想要放大最佳 DI 所在区域，只需取消注释xlim行：

plt.plot(**levels**, **di**, marker='o')

上述代码生成以下输出。正如你所看到的，最佳修复水平位于 0 和 0.1 之间，因为那里的值最接近 1：

图表，折线图描述自动生成

图 11.6：不同 DI 消除修复水平下的 DI

现在，让我们使用evaluate_class_mdl评估最佳的 DI 修复模型，并计算公平性指标（compute_aif_metrics）。这次我们甚至不会绘制混淆矩阵，但我们会将所有结果保存到cls_mdls字典中，以供后续检查。代码如下所示：

cls_mdls['lgb_1_dir'] = mldatasets.**evaluate_class_mdl**(
    lgb_dir_mdl,
    X_train_dir,
    X_test_dir,
    train_dir_ds.labels,
    test_dir_ds.labels,
    plot_roc=False,
    plot_conf_matrix=False,
    show_summary=False,
    ret_eval_dict=True
)
test_pred_dir_ds = test_ds.copy(deepcopy=True)
test_pred_dir_ds.labels = cls_mdls['lgb_1_dir']['preds_test'
].reshape(-1,1)
metrics_test_dir_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_dir_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['lgb_1_dir'].**update**(metrics_test_dir_dict)

数据链中的下一个链接是模型，因此即使我们去除了数据中的偏差，模型本身也会引入偏差，因此训练能够处理这种偏差的模型是有意义的，这正是我们接下来将要学习如何做的！

处理中偏差减轻方法

这些是一些最重要的处理中或模型特定的偏差减轻方法：

成本敏感训练：我们已经在本章训练的每个 LightGBM 模型中通过scale_pos_weight参数整合了这种方法。它通常用于不平衡分类问题，并被简单地视为提高少数类准确率的一种手段。然而，鉴于类别不平衡往往倾向于使某些群体优于其他群体，这种方法也可以用来减轻偏差，但并不能保证一定会这样做。它可以作为类权重或通过创建自定义损失函数来整合。实现方式将根据模型类别和与偏差相关的成本而有所不同。如果它们与误分类线性增长，则类权重就足够了，否则建议使用自定义损失函数。
约束：许多模型类别支持单调性和交互约束，TensorFlow Lattice（TFL）提供了更高级的自定义形状约束。这些确保了特征和目标之间的关系被限制在某种模式中，在模型级别上设置了护栏。你会有很多理由想要使用它们，但其中最重要的是减轻偏差。我们将在第十二章单调约束和模型调优以实现可解释性中讨论这个话题。
偏见消除正则化器：这种方法将偏见定义为敏感变量和目标变量之间的统计依赖性。然而，这种方法的目标是最大限度地减少间接偏见，排除可以通过简单地删除敏感变量来避免的偏见。因此，该方法首先通过偏见指数（PI）对其进行量化，这是目标和敏感变量之间的互信息。顺便提一下，我们在第十章可解释性特征选择和工程中介绍了互信息。然后，与 L2 一起，PI 被整合到一个自定义正则化项中。从理论上讲，任何模型分类器都可以使用基于 PI 的正则化器进行正则化，但到目前为止，唯一实现的例子是逻辑回归。
Gerry fair 分类器：这是受公平性划分概念的启发，它在某一群体中看似公平，但在细分到子群体时却缺乏公平性。该算法利用一种基于博弈论的虚构博弈方法，其中你有一个学习者和审计员之间的零和游戏。学习者最小化预测误差和基于公平性的总惩罚项。审计员通过基于在最不公平对待的子群体中观察到的最坏结果来进一步惩罚学习者。

游戏的目标是达到纳什均衡，这是在两个可能具有矛盾目标的非合作玩家达成部分满足双方的解决方案时实现的。在这种情况下，学习者获得最小的预测误差和总体不公平性，审计员获得最小的子群体不公平性。该方法的实现是模型无关的。

对抗性去偏：与 gerry fair 分类器类似，对抗性去偏利用两个对立的演员，但这次是两个神经网络：预测器和对手。我们最大化预测器预测目标的能力，同时最小化对手预测受保护特征的能力，从而增加特权群体和弱势群体之间机会的平等性。
指数梯度下降法：这种方法通过将其简化为一系列此类问题来自动化成本敏感优化，并使用关于受保护属性（如人口统计学平等或均衡机会）的公平性约束。它是模型无关的，但仅限于与 scikit-learn 兼容的二分类器。

由于存在如此多的预处理方法，我们将在本章中仅使用其中两种。尽管如此，如果你对我们将不涉及的方法感兴趣，它们可以在 AIF360 库和文档中找到。

指数梯度下降法

ExponentiatedGradientReduction方法是对具有约束的成本敏感训练的实现。我们用基估计器初始化它，指定要执行的迭代次数的最大值（max_iter），并指定要使用的差异constraints。然后，我们fit它。这种方法可以在下面的代码片段中看到：

lgb_egr_mdl = ExponentiatedGradientReduction(
    estimator=lgb_base_mdl,
    max_iter=50,
    constraints='DemographicParity'
)
lgb_egr_mdl.fit(train_ds)

我们可以使用predict函数来获取训练和测试预测，然后使用evaluate_class_metrics_mdl和compute_aif_metrics分别获取预测性能和公平性指标。我们将它们都放入cls_mdls字典中，如下面的代码片段所示：

train_pred_egr_ds = lgb_egr_mdl.**predict**(train_ds)
test_pred_egr_ds = lgb_egr_mdl.**predict**(test_ds)
cls_mdls['lgb_2_egr'] = mldatasets.**evaluate_class_metrics_mdl**(
    lgb_egr_mdl,
    train_pred_egr_ds.labels,
    test_pred_egr_ds.scores,
    test_pred_egr_ds.labels,
    y_train,
    y_test
)
metrics_test_egr_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_egr_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['lgb_2_egr'].**update**(metrics_test_egr_dict)

接下来，我们将了解一种部分模型无关的预处理方法，它考虑了交叉性。

gerry fair 分类器方法

Gerry fair 分类器部分是模型无关的。它只支持线性模型、支持向量机（SVMs）、核回归和决策树。我们通过定义正则化强度（C）、用于早期停止的公平近似（gamma）、是否详细输出（printflag）、最大迭代次数（max_iters）、模型（predictor）以及要采用的公平概念（fairness_def）来初始化GerryFairClassifier。我们将使用错误的负例（"FN"）的公平概念来计算公平违规的加权差异。一旦初始化完成，我们只需要调用fit方法并启用early_termination，如果它在五次迭代中没有改进，则停止。以下代码片段展示了代码：

dt_gf_mdl = **GerryFairClassifier**(
    C=100,
    gamma=.005,
    max_iters=50,
    **fairness_def**='FN',
    printflag=True,
    **predictor**=tree.DecisionTreeRegressor(max_depth=3)
)
dt_gf_mdl.**fit**(train_ds, early_termination=True)

我们可以使用predict函数来获取训练和测试预测，然后使用evaluate_class_metrics_mdl和compute_aif_metrics来分别获得预测性能和公平性指标。我们将它们都放入cl_smdls字典中，如下面的代码片段所示：

train_pred_gf_ds = dt_gf_mdl.**predict**(train_ds, threshold=False)
test_pred_gf_ds = dt_gf_mdl.**predict**(test_ds, threshold=False)
cls_mdls['dt_2_gf'] = mldatasets.evaluate_class_metrics_mdl(
    dt_gf_mdl,
    train_pred_gf_ds.labels,
    None,
    test_pred_gf_ds.labels,
    y_train,
    y_test
)
metrics_test_gf_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_gf_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['dt_2_gf'].**update**(metrics_test_gf_dict)

在模型推理之后的链中的下一个和最后一个链接，因此即使你去除了数据和模型的偏差，也可能还剩下一些偏差，因此在这个阶段处理它也是有意义的，这正是我们接下来将要学习如何做的！

后处理偏差缓解方法

这些是一些最重要的后处理或推理特定偏差缓解方法：

预测弃权：这有许多潜在的好处，如公平性、安全性和控制成本，但具体应用取决于你的问题。通常，模型会返回所有预测，即使是低置信度的预测——也就是说，接近分类阈值的预测，或者当模型返回的置信区间超出预定阈值时。当涉及公平性时，如果我们将在低置信度区域将预测改为我不知道（IDK），那么在评估公平性指标时，仅针对所做的预测，模型可能会因为副作用而变得更加公平。还可能将预测弃权作为一个内部处理方法。一篇名为《负责任地预测：通过学习推迟来提高公平性》的论文讨论了两种方法，通过训练模型来回避（学习预测 IDK）或推迟（当正确率低于专家意见时预测 IDK）。另一篇名为《在二元分类中弃权的效用》的论文采用了一个名为Knows What It Knows（KWIK）的强化学习框架，它对自己的错误有自我意识，但允许弃权。
均衡机会后处理：也称为不同对待，这确保了特权群体和弱势群体在错误分类方面得到平等对待，无论是假阳性还是假阴性。它找到最佳概率阈值，通过改变标签来平衡组之间的机会。
校准的均等机会后处理：这种方法不是改变标签，而是修改概率估计，使它们平均相等。它称之为校准。然而，这个约束不能同时满足假阳性和假阴性，因此你被迫在两者之间做出选择。因此，在召回率远比精确度更重要或反之亦然的情况下，校准估计的概率是有利的。
拒绝选项分类法：这种方法利用了直觉，即决策边界周围的预测往往是最不公平的。然后，它找到决策边界周围的一个最优带，在这个带中，翻转弱势和优势群体的标签可以产生最公平的结果。

在本章中，我们只会使用这两种后处理方法。拒绝选项分类法在 AIF360 库和文档中可用。

均等机会后处理方法

均等机会后处理方法（EqOddsPostprocessing）初始化时，需要指定我们想要均等机会的群体和随机种子。然后，我们fit它。请注意，拟合需要两个数据集：原始数据集（test_ds）以及为我们基础模型提供预测的数据集（test_pred_ds）。fit所做的就是计算最优概率阈值。然后，predict创建一个新的数据集，其中这些阈值已经改变了labels。代码可以在下面的片段中看到：

epp = **EqOddsPostprocessing**(
    privileged_groups=privileged_groups,
    unprivileged_groups=underprivileged_groups,
    seed=rand
)
epp = epp.**fit**(test_ds, test_pred_ds)
test_pred_epp_ds = epp.**predict**(test_pred_ds)

我们可以使用evaluate_class_metrics_mdl和compute_aif_metrics来分别获得等比例概率（EPP）的预测性能和公平性指标。我们将它们都放入cls_mdls字典中。代码可以在下面的片段中看到：

cls_mdls['lgb_3_epp'] = mldatasets.**evaluate_class_metrics_mdl**(
    lgb_base_mdl,
    cls_mdls['lgb_0_base']['preds_train'],
    test_pred_epp_ds.scores,
    test_pred_epp_ds.labels,
    y_train,
    y_test
)
metrics_test_epp_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_epp_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['lgb_3_epp'].**update**(metrics_test_epp_dict)

接下来，我们将了解另一种后处理方法。主要区别在于它校准概率分数，而不仅仅是改变预测标签。

校准的均等机会后处理方法

校准的均等机会（CalibratedEqOddsPostprocessing）的实现方式与均等机会完全相同，但它有一个更关键的属性（cost_constraint）。这个属性定义了要满足哪个约束，因为它不能同时使分数对 FPRs 和 FNRs 都是公平的。我们选择 FPR，然后fit、predict和evaluate，就像我们对均等机会所做的那样。代码可以在下面的片段中看到：

cpp = **CalibratedEqOddsPostprocessing**(
    privileged_groups=privileged_groups,
    unprivileged_groups=underprivileged_groups,
    **cost_constraint**="fpr",
    seed=rand
)
cpp = cpp.**fit**(test_ds, test_pred_ds)
test_pred_cpp_ds = cpp.**predict**(test_pred_ds)
cls_mdls['lgb_3_cpp'] = mldatasets.**evaluate_class_metrics_mdl**(
    lgb_base_mdl,
    cls_mdls['lgb_0_base']['preds_train'],
    test_pred_cpp_ds.scores,
    test_pred_cpp_ds.labels,
    y_train,
    y_test
)
metrics_test_cpp_dict, _ = mldatasets.**compute_aif_metrics**(
    test_ds,
    test_pred_cpp_ds,
    unprivileged_groups=underprivileged_groups,
    privileged_groups=privileged_groups
)
cls_mdls['lgb_3_cpp'].**update**(metrics_test_cpp_dict)

现在我们已经尝试了六种偏差缓解方法，每个级别两种，我们可以将它们相互比较，并与基础模型进行比较！

将所有这些结合起来！

为了比较所有方法的指标，我们可以将字典（cls_mdls）放入 DataFrame（cls_metrics_df）中。我们只对一些性能指标和记录的大多数公平性指标感兴趣。然后，我们输出按测试准确率排序的 DataFrame，并使用所有公平性指标进行颜色编码。代码可以在下面的片段中看到：

cls_metrics_df = pd.DataFrame.from_dict(cls_mdls, 'index')[
    [
        'accuracy_train',
        'accuracy_test',
        'f1_test',
        'mcc_test',
        'SPD',
        'DI',
        'AOD',
        'EOD',
        'DFBA'
    ]
]
metrics_fmt_dict = dict(
    zip(cls_metrics_df.columns,['{:.1%}']*3+ ['{:.3f}']*6)
)
cls_metrics_df.sort_values(
    by='accuracy_test',
    ascending=False
).style.format(metrics_fmt_dict)

前面的代码片段输出了以下 DataFrame：

图 11.7：所有偏差缓解方法与不同公平性指标的对比

图 11.7 显示，大多数方法在 SPD、DI、AOD 和 EOD 方面产生的模型比基础模型更公平。校准等概率后处理 (lgb_3_cpp) 是一个例外，但它具有最佳的 DFBAs 之一，但由于校准的不平衡性质，它产生了次优的 DI。请注意，这种方法在校准分数时特别擅长实现 FPR 或 FNR 的平衡，但所有这些公平性指标对于捕捉这一点都没有用。相反，你可以创建一个指标，它是 FPRs 的比率，就像我们在 第六章，锚点和反事实解释 中所做的那样。偶然的是，这将是一个完美的 校准等概率 （CPP）的用例。

获得最佳 SPD、DI、AOD 和 DFBA 以及次优 EOD 的方法是等概率后处理 (lgb_3_epp)，因此让我们使用 XAI 的图表来可视化其公平性。为此，我们首先创建一个包含测试示例的 DataFrame (test_df)，然后使用 replace 将 AGE_GROUP 转换为分类变量，并获取分类列的列表 (cat_cols_l)。然后，我们可以使用真实标签 (y_test)、EPP 模型的预测概率分数、DataFrame (test_df)、受保护属性 (cross_cols) 和分类列来比较不同的指标 (metrics_plot)。我们也可以为 受试者工作特征 （ROC）图表 (roc_plot) 和 精确率-召回率 （PR）图表 (pr_plot) 做同样的事情。代码可以在下面的代码片段中看到：

test_df = ccdefault_bias_df.loc[**X_test.index**]
test_df['AGE_GROUP'] = test_df.AGE_GROUP.**replace**(
    {0:'underprivileged', 1:'privileged'}
)
cat_cols_l = ccdefault_bias_df.dtypes[**lambda x: x==np.int8**
                                      ].index.tolist()
_ = xai.**metrics_plot**(
    y_test,cls_mdls['lgb_3_epp']['probs_test'],
    df=test_df, cross_cols=['AGE_GROUP'],
    categorical_cols=cat_cols_l
)
_ = xai.**roc_plot**(
    y_test, cls_mdls['lgb_3_epp']['probs_test'],
    df=test_df, cross_cols=['AGE_GROUP'],
    categorical_cols=cat_cols_l
)
_ = xai.**pr_plot**(
    y_test,
    cls_mdls['lgb_3_epp']['probs_test'],
    df=test_df, cross_cols=['AGE_GROUP'],
    categorical_cols=cat_cols_l
) 
Figure 11.8. The first one shows that even the fairest model still has some disparities between both groups, especially between precision and recall and, by extension, F1 score, which is their average. However, the ROC curve shows how close both groups are from an FPR versus a TPR standpoint. The third plot is where the disparities in precision and recall become even more evident. This all demonstrates how hard it is to keep a fair balance on all fronts! Some methods are best for making one aspect perfect but nothing else, while others are pretty good on a handful of aspects but nothing else. Despite the shortcomings of the methods, most of them achieved a sizable improvement. Ultimately, choosing methods will depend on what you most care about, and combining them is also recommended for maximum effect! The output is shown here:

图 11.8：展示最公平模型的公平性图表

我们已经完成了偏差缓解练习，并将继续进行因果推断练习，我们将讨论如何确保公平和稳健的政策。

创建因果模型

决策通常需要理解因果关系。如果效果是可取的，你可以决定复制其原因，或者避免它。你可以故意改变某些东西来观察它如何改变结果，或者将意外效应追溯到其原因，或者模拟哪种改变会产生最大的积极影响。因果推断可以通过创建因果图和模型来帮助我们完成所有这些，这些图将所有变量联系起来并估计效应，以便做出更原则性的决策。然而，为了正确评估原因的影响，无论是设计还是意外，你需要将其效应与混杂变量分开。

因果推断与本章相关的原因是银行的决策具有显著影响持卡人生计的力量，鉴于自杀率的上升，甚至关系到生死。因此，有必要极其谨慎地评估政策决策。

台湾银行进行了一项为期 6 个月的贷款政策实验。银行看到了形势的严峻，知道那些最高风险的违约客户将 somehow 从资产负债表中注销，从而减轻了这些客户的财务义务。因此，实验的重点仅涉及银行认为可以挽救的部分，即低至中等风险的违约客户，现在实验已经结束，他们想了解以下政策如何影响了客户行为：

降低信用额度：一些客户的信用额度降低了 25%。
付款计划：他们被给予 6 个月的时间来偿还当前的信用卡债务。换句话说，债务被分成六部分，每个月他们必须偿还一部分。
两项措施：降低信用额度和付款计划。

此外，2005 年台湾普遍的信用卡利率约为 16-20%，但银行得知这些利率将被台湾金融监督管理委员会限制在 4%。因此，他们确保所有参与实验的客户都能自动获得该水平的利率。一些银行高管认为这只会加剧债务负担，并在这个过程中创造更多的“信用卡奴隶”。这些担忧促使提出以较低的信用卡额度作为对策进行实验的建议。另一方面，制定付款计划是为了了解债务减免是否给了客户使用信用卡而不必担心的空间。

在业务方面，理由是必须鼓励健康水平的消费，因为随着利率的降低，大部分利润将来自支付处理、现金返还合作伙伴关系和其他与消费相关的来源，从而增加客户的使用寿命。这对客户也有好处，因为如果他们在作为消费者比作为债务人更有利可图，这意味着激励措施已经到位，以防止他们成为后者。所有这些都证明了使用估计的终身价值（_LTV）作为代理指标来衡量实验结果如何使银行和客户受益的合理性。多年来，银行一直在使用一种相当准确的计算方法来估计信用卡持卡人根据他们的消费和支付历史以及诸如额度、利率等参数将为银行提供多少价值。

在实验设计的术语中，选择的政策被称为治疗，除了三个接受治疗的组别外，还有一个未接受治疗的对照组——即政策没有任何变化，甚至没有降低利率。在我们继续前进之前，让我们首先初始化一个包含治疗名称的列表（treatment_names）和一个包含甚至对照组的列表（all_treatment_names），如下所示：

treatment_names = [
    'Lower Credit Limit',
    'Payment Plan',
    'Payment Plan &Credit Limit'
]
all_treatment_names = np.array(["None"] + treatment_names)

现在，让我们检查实验的结果，以帮助我们设计一个最优的因果模型。

理解实验结果

评估治疗有效性的一个相当直观的方法是通过比较它们的成果。我们想知道以下两个简单问题的答案：

相比对照组，治疗是否降低了违约率？
支出行为是否有利于提高终身价值估计？

我们可以在一个图表中可视化这两个因素。为此，我们获得一个包含每个组违约百分比的pandas系列（pct_s），然后另一个包含每个组终身价值总和的系列（ltv_s），单位为千新台币（NTD）（K$）。我们将这两个系列放入pandas DataFrame 中，并绘制它，如下面的代码片段所示：

pct_s = ccdefault_causal_df[
    ccdefault_causal_df.IS_DEFAULT==1]
    .groupby(['_TREATMENT'])
    .size()
    /ccdefault_causal_df.groupby(['_TREATMENT']).size()
ltv_s = ccdefault_causal_df.groupby(
    ['_TREATMENT'])['_LTV'].sum()/1000
plot_df = pd.DataFrame(
    {'% Defaulted':pct_s,
     'Total LTV, K$':ltv_s}
)
plot_df.index = all_treatment_names
ax = plot_df.plot(secondary_y=['Total LTV, K$'], figsize=(8,5))
ax.get_legend().set_bbox_to_anchor((0.7, 0.99))
plt.grid(False) 
Figure 11.9. It can be inferred that all treatments fare better than the control group. The lowering of the credit limit on its own decreases the default rate by over 12% and more than doubles the estimated LTV, while the payment plan only decreases the defaults by 3% and increases the LTV by about 85%. However, both policies combined quadrupled the control group’s LTV and reduced the default rate by nearly 15%! The output can be seen here:

图表，折线图描述自动生成

图 11.9：不同信用政策的治疗实验结果

在银行高管们为找到了获胜政策而欢欣鼓舞之前，我们必须检查他们是如何在实验中的信用卡持卡人之间分配它的。我们了解到，他们根据风险因素（由_risk_score变量衡量）选择治疗。然而，终身价值在很大程度上受到可用信用额（_CC_LIMIT）的影响，因此我们必须考虑这一点。理解分布的一种方法是通过将两个变量以散点图的形式相互对比，并按_TREATMENT进行颜色编码。以下代码片段展示了如何实现这一点：

sns.**scatterplot**(
    x=ccdefault_causal_df['_CC_LIMIT'].values,
    y=ccdefault_causal_df['_risk_score'].values,
    hue=all_treatment_names[ccdefault_causal_df['_TREATMENT'].values],
    hue_order = all_treatment_names
)

上述代码生成了图 11.10中的图表。它显示三种治疗对应不同的风险水平，而对照组（None）在垂直方向上分布得更广。基于风险水平分配治疗的选择也意味着他们基于_CC_LIMIT不均匀地分配了治疗。我们应该问自己，这个实验的偏见条件是否使得结果解释甚至变得可行。请看以下输出：

图表，散点图描述自动生成

图 11.10：风险因素与原始信用额的比较

图 11.10中的散点图展示了治疗在风险因素上的分层。然而，散点图在理解分布时可能具有挑战性。为此，最好使用核密度估计（KDE）图。因此，让我们看看_CC_LIMIT和终身价值（_LTV）在所有治疗中的分布情况，使用 Seaborn 的displot。请看以下代码片段：

sns.**displot**(
    ccdefault_causal_df,
    x="_CC_LIMIT",
    hue="_TREATMENT",
    kind="kde",
    fill=True
)
sns.**displot**(
    ccdefault_causal_df,
    x="_LTV",
    hue="_TREATMENT",
    kind="kde", fill=True
) 
Figure 11.11. We can easily tell how far apart all four distributions are for both plots, mostly regarding treatment #3 (Payment Plan & Lower Credit Limit), which tends to be centered significantly more to the right and has a longer and fatter right tail. You can view the output here:

图表描述自动生成

图 11.11：根据 _TREATMENT 的 _CC_LIMIT 和 _LTV 的 KDE 分布

理想情况下，当你设计此类实验时，你应该根据可能改变结果的相关因素，在所有组别中追求平等分布。然而，这并不总是可行的，可能因为物流或战略限制。在这种情况下，结果（_LTV）根据客户信用卡额度（_CC_LIMIT）、异质性特征——换句话说，直接影响处理效果的变量，也称为异质性处理效应调节因子而变化。我们可以创建一个包含_TREATMENT特征和效应调节因子（_CC_LIMIT）的因果模型。

理解因果模型

我们将要构建的因果模型可以分为以下四个部分：

结果 (Y): 因果模型的结果变量。
处理 (T): 影响结果的处理变量。
效应调节因子 (X): 影响效应异质性的变量，它位于处理和结果之间。
控制变量 (W): 也称为共同原因或混杂因素。它们是影响结果和处理的特征。

我们将首先将这些组件在数据中识别为单独的pandas数据框，如下所示：

**W** = ccdefault_causal_df[
    [
      '_spend','_tpm', '_ppm', '_RETAIL','_URBAN', '_RURAL',
      '_PREMIUM'
    ]
]
**X** = ccdefault_causal_df[['_CC_LIMIT']]
**T** = ccdefault_causal_df[['_TREATMENT']]
**Y** = ccdefault_causal_df[['_LTV']]

我们将使用双重稳健学习（DRL）方法来估计处理效应。它被称为“双重”，因为它利用了两个模型，如下所示：

它使用回归模型预测结果，如图所示：

它使用倾向模型预测处理，如图所示：

由于最终阶段结合了两种模型，同时保持了多个理想的统计特性，如置信区间和渐近正态性，因此它也是稳健的。更正式地说，估计利用了条件在处理t上的回归模型g和倾向模型p，如下所示：

它还做了以下操作：

目标是推导出与每个处理t相关的异质效应X的条件平均处理效应（CATE），表示为。首先，DRL 方法通过应用逆倾向来去偏回归模型，如下所示：

如何精确地估计模型中的系数将取决于所采用的 DRL 变体。我们将使用线性变体（LinearDRLearner），以便它返回系数和截距，这些可以很容易地解释。它通过在处理组t和控制组在x[t]上的结果差异中运行普通线性回归（OLS）来推导。这种直观的做法是有意义的，因为处理的估计效应减去没有处理的估计效应（t = 0）是这种处理的净效应。

现在，所有理论都已经讲完，让我们深入挖掘吧！

初始化线性双重稳健学习器

我们可以通过指定任何与 scikit-learn 兼容的回归器（model_regression）和分类器（model_propensity）来从econml库初始化一个LinearDRLearner，我们称之为drlearner。我们将使用 XGBoost 来处理这两个，但请注意，分类器有一个objective=multi:softmax属性。记住，我们有多个处理，所以这是一个多类分类问题。代码可以在下面的片段中看到：

drlearner = **LinearDRLearner**(
    model_regression=xgb.XGBRegressor(learning_rate=0.1),
    model_propensity=xgb.XGBClassifier(learning_rate=0.1,
    max_depth=2,
    objective="multi:softmax"),
    random_state=rand
)

如果你想了解回归模型和倾向性模型都在做什么，你可以轻松地拟合xgb.XGBRegressor().fit(W.join(X),Y)和xgb.XGBClassifier(objective="multi:softmax").fit(W.join(X), T)模型。我们现在不会这样做，但如果你好奇，你可以评估它们的性能，甚至运行特征重要性方法来了解它们各自预测的影响。因果模型将它们与 DRL 框架结合在一起，导致不同的结论。

拟合因果模型

我们可以使用drlearner中的fit来拟合因果模型，利用econml的dowhy包装器。首先的属性是Y、T、X和Y组件：pandas数据框。可选地，你可以为这些组件提供变量名称：每个pandas数据框的列名。最后，我们希望估计处理效应。可选地，我们可以提供用于此的效果修饰符（X），我们将使用其中的一半数据来这样做，如下面的代码片段所示：

causal_mdl = drlearner.dowhy.fit(
    Y,
    T,
    X=X,
    W=W,
    outcome_names=Y.columns.to_list(),
    treatment_names=T.columns.to_list(),
    feature_names=X.columns.to_list(),
    confounder_names=W.columns.to_list(),
    target_units=X.iloc[:550].values
)

在因果模型初始化后，我们可以可视化它。pydot库与pygraphviz可以为我们完成这项工作。请注意，这个库在某些环境中配置困难，所以它可能无法加载并显示view_model的默认图形。如果发生这种情况，请不要担心。看看下面的代码片段：

try:
    display(**Image**(to_pydot(causal_mdl._graph._graph).create_png()))
except:
    causal_mdl.**view_model**()

前一个代码片段中的代码输出了此处显示的模型图。有了它，你可以欣赏到所有变量是如何相互连接的：

图描述自动生成

图 11.12：因果模型图

因果模型已经拟合，那么让我们检查和解释结果，好吗？

理解异质处理效应

首先，需要注意的是，econml 的 dowhy 包装器通过 dowhy.fit 方法简化了一些步骤。通常，当你直接使用 dowhy 构建 CausalModel（如本例所示）时，它有一个名为 identify_effect 的方法，该方法推导出要估计的效果的概率表达式（即 识别估计量）。在这种情况下，这被称为 平均处理效应（ATE）。然后，另一个名为 estimate_effect 的方法接受这个表达式以及它应该与之关联的模型（回归和倾向）。有了它们，它为每个结果 i 和处理 t 计算 ATE，，和 CATE，。然而，由于我们使用了包装器来 fit 因果模型，它自动处理了识别和估计步骤。

您可以通过因果模型的 identified_estimand_ 属性访问识别的 ATE，并通过 estimate_ 属性访问估计结果。以下代码片段显示了代码：

identified_ate = causal_mdl.identified_estimand_
print(identified_ate)
drlearner_estimate = causal_mdl.estimate_
print(drlearner_estimate) 
 the estimand expression for identified_estimand_, which is a derivation of the expected value for , with some assumptions. Then, the causal-realized estimate_ returns the ATE for treatment #1, as illustrated in the following code snippet:

Estimand type: nonparametric-ate
### Estimand : 1
Estimand name: backdoor1 (Default)
Estimand expression:
      d 
─────────────(E[_LTV|_RETAIL,_URBAN,_PREMIUM,_RURAL,_CC_LIMIT,
d[_TREATMENT]
])
Estimand assumption 1, Unconfoundedness: If U→{_TREATMENT} and U→_LTV then \ P(_LTV|_TREATMENT,_RETAIL,_URBAN,_PREMIUM,_RURAL,_CC_LIMIT,_spend,_ppm,_tpm,U) = \ P(_LTV|_TREATMENT,_RETAIL,_URBAN,_PREMIUM,_RURAL,_CC_LIMIT,_spend,_ppm,_tpm)
*** Causal Estimate ***
## Identified estimand
Estimand type: nonparametric-ate
## Realized estimand
b:_LTV ~ _TREATMENT + _RETAIL + _URBAN + _PREMIUM + _RURAL + \ _CC_LIMIT + _spend + _ppm + _tpm | _CC_LIMIT
Target units:
## Estimate
Mean value: 7227.904763676559
Effect estimates: [6766.07978487 7337.39526574 7363.36013004
                   7224.20893104 7500.84310705 7221.40328496]

接下来，我们可以遍历因果模型中的所有处理，并为每个处理返回一个总结，如下所示：

for i in range(causal_mdl._d_t[0]):
    print("Treatment: %s" % treatment_names[i])
     display(econml_mdl.**summary**(T=i+1))

前面的代码输出了三个线性回归总结。第一个看起来像这样：

图形用户界面，表格自动生成的描述

图 11.13：某处理总结

为了更好地理解系数和截距，我们可以用它们各自的置信区间来绘制它们。为此，我们首先创建一个处理索引（idxs）。有三个处理，所以这是一个介于 0 和 2 之间的数字数组。然后，使用列表推导将所有系数（coef_）和截距（intercept_）放入一个数组中。然而，对于系数和截距的 90%置信区间来说，这要复杂一些，因为 coef__interval 和 intercept__interval 返回这些区间的下限和上限。我们需要两个方向的误差范围的长度，而不是界限。我们从这些界限中减去系数和截距，以获得它们各自的误差范围，如下面的代码片段所示：

idxs = np.arange(0, causal_mdl._d_t[0])
coefs = np.hstack([causal_mdl.**coef_**(T=i+1) for i in idxs])
intercepts = np.hstack(
    [causal_mdl.**intercept_**(T=i+1)for i in idxs]
)
coefs_err = np.hstack(
    [causal_mdl.**coef__interval**(T=i+1) for i in idxs]
)
coefs_err[0, :] = coefs - coefs_err[0, :]
coefs_err[1, :] = coefs_err[1, :] - coefs
intercepts_err = np.vstack(
    [causal_mdl.**intercept__interval**(T=i+1) for i in idxs]
).Tintercepts_err[0, :] = intercepts - intercepts_err[0, :]
intercepts_err[1, :] = intercepts_err[1, :] - intercepts

接下来，我们使用 errorbar 函数绘制每个处理及其相应误差的系数。我们还可以将截距作为另一个子图进行相同的操作，如下面的代码片段所示：

ax1 = plt.subplot(2, 1, 1)
plt.errorbar(**idxs**, **coefs**, **coefs_err**, fmt="o")
plt.xticks(idxs, treatment_names)
plt.setp(ax1.get_xticklabels(), visible=False)
plt.title("Coefficients")
plt.subplot(2, 1, 2)
plt.errorbar(**idxs**, **intercepts**, **intercepts_err**, fmt="o")
plt.xticks(idxs, treatment_names)
plt.title("Intercepts")

前面的代码片段输出以下内容：

包含图表的图片自动生成的描述

图 11.14：所有处理的系数和截距

通过图 11.14，你可以欣赏到所有截距和系数的相对误差范围有多大。尽管如此，很明显，仅从系数来看，从左到右读取时，治疗的效果会逐渐变好。但在我们得出支付计划 & 降低信用额度是最佳政策的结论之前，我们必须考虑截距，这个截距对于这种治疗比第一个要低。本质上，这意味着具有最低信用卡额度的客户更有可能通过第一种政策提高终身价值，因为系数是乘以限制的，而截距是起点。鉴于没有一种最佳政策适用于所有客户，让我们来看看如何使用因果模型为每个客户选择政策。

选择政策

我们可以使用const_marginal_effect方法根据客户基础制定信用政策，该方法考虑了X效果修正器（_CC_LIMIT）并计算反事实 CATE，。换句话说，它返回了所有观察到的X中所有治疗的估计_LTV。

然而，它们并不都花费相同。制定支付计划需要每份合同约NT $*1,000 的行政和法律费用，根据银行的精算部门，降低信用额度 25 的机遇成本估计为每月平均支付*NT$ 72（_ppm），在整个客户生命周期内。为了考虑这些成本，我们可以设置一个简单的lambda函数，该函数接受所有治疗的支付计划成本并将它们添加到变量信用额度成本中，这自然地乘以_ppm。给定一个长度为n的信用卡额度数组，成本函数返回一个(n, 3)维度的数组，其中包含每个治疗的成本。然后，我们获得反事实 CATE 并扣除成本（treatment_effect_minus_costs）。然后，我们将数组扩展以包括一列表示无治疗的零，并使用argmax返回每个客户的推荐治疗索引（recommended_T），如下面的代码片段所示：

**cost_fn** = lambda X: np.repeat(
    np.array([[0, 1000, 1000]]),
    X.shape[0], axis=0) + (np.repeat(np.array([[72, 0, 72]]),
    X.shape[0], axis=0)
    *X._ppm.values.reshape(-1,1)
)
**treatment_effect_minus_costs** = causal_mdl.const_marginal_effect(
    X=X.values) - **cost_fn**(ccdefault_causal_df)
treatment_effect_minus_costs = np.hstack(
    [
        np.zeros(X.shape),
        **treatment_effect_minus_costs**
    ]
)
recommended_T = np.**argmax**(treatment_effect_minus_costs, axis=1)

我们可以使用scatterplot _CC_LIMIT和_ppm，按推荐治疗进行颜色编码，以观察客户的最佳信用政策，如下所示：

sns.scatterplot(
    x=ccdefault_causal_df['_CC_LIMIT'].values,
    y=ccdefault_causal_df['_ppm'].values,
    hue=all_treatment_names[recommended_T],
    hue_order=all_treatment_names
)
plt.title("Optimal Credit Policy by Customer")
plt.xlabel("Original Credit Limit")
plt.ylabel("Payments/month")

前面的代码片段输出以下散点图：

图表，散点图自动生成的描述图 11.15：根据原始信用额度和卡片使用情况，客户最优信用政策

在 图 11.15 中很明显，“None”（无治疗）永远不会被推荐给任何客户。即使不扣除成本，这一事实也成立——你可以从 treatment_effect_minus_costs 中移除 cost_fn 并重新运行输出图表的代码来验证，无论成本如何，治疗总是被推荐的。你可以推断出所有治疗对客户都有益，其中一些比其他更多。当然，根据客户的不同，一些治疗比其他治疗对银行更有利。这里有一条很细的界限。

最大的担忧之一是客户的公平性，特别是那些银行伤害最严重的客户：弱势年龄群体。仅仅因为一项政策对银行的成本比另一项更高，并不意味着应该排除访问其他政策的机会。评估这一点的一种方法可以使用所有推荐政策的百分比堆叠条形图。这样，我们可以观察推荐政策在特权群体和弱势群体之间的分配情况。看看下面的代码片段：

ccdefault_causal_df['recommended_T'] = recommended_T
plot_df = ccdefault_causal_df.groupby(
    ['recommended_T','AGE_GROUP']).size().reset_index()
plot_df['AGE_GROUP'] = plot_df.AGE_GROUP.**replace**(
    {0:'underprivileged', 1:'privileged'}
)
plot_df = plot_df.pivot(
    columns='AGE_GROUP',
    index='recommended_T',
    values=0
)
plot_df.index = treatment_names
plot_df = plot_df.apply(lambda r: **r/r.sum()*100**, axis=1)
plot_df.plot.bar(stacked=True, rot=0)
plt.xlabel('Optimal Policy')
plt.ylabel('%')

前一个代码片段中的代码输出如下：

图 11.16：最优策略分布的公平性

图 11.16 展示了特权群体被分配到具有支付计划的政策的比例更高。这种差异主要是由于银行的成本是一个因素，所以如果银行能够承担一些这些成本，那么它可能会更加公平。但什么是公平的解决方案呢？选择信贷政策是程序公平性的一个例子，并且有许多可能的定义。平等对待是否字面意义上的平等对待或比例对待？它是否包括选择自由的概念？如果客户更喜欢一项政策而不是另一项，他们应该被允许切换吗？无论定义如何，都可以通过因果模型的帮助来解决。我们可以将相同的政策分配给所有客户，或者调整推荐政策的分布，使得比例相等，或者每个客户都可以在第一和第二最优政策之间进行选择。有如此多的方法可以这样做！

测试估计的鲁棒性

dowhy 库提供了四种方法来测试估计因果效应的鲁棒性，具体如下：

随机共同原因：添加一个随机生成的混杂因素。如果估计是鲁棒的，ATE（平均处理效应）不应该变化太多。
安慰剂治疗反驳者：用随机变量（安慰剂）替换治疗。如果估计是鲁棒的，ATE 应该接近零。
数据子集反驳者：移除数据的一个随机子集。如果估计器泛化良好，ATE 不应该变化太多。
添加未观察到的共同原因：添加一个与处理和结果都相关的未观察到的混杂因素。估计量假设存在一定程度的未混杂性，但添加更多应该会偏误估计。根据混杂因素效应的强度，它应该对 ATE 有相同的影响。

我们将用前两个来测试稳健性。

添加随机共同原因

此方法通过调用refute_estimate并指定method_name="random_common_cause"来实现，这是最简单的实现方式。这将返回一个可以打印的摘要。请看以下代码片段：

ref_random = causal_mdl.refute_estimate(
    method_name="random_common_cause"
)
print(ref_random)

前述代码片段输出如下：

Refute: Add a Random Common Cause
Estimated effect:7227.904763676559
New effect:7241.433599647397

前面的输出告诉我们，一个新的共同原因，或称 W 变量，对平均处理效应（ATE）没有显著影响。

用随机变量替换处理变量

使用此方法，我们将用噪声替换处理变量。如果处理与结果有稳健的相关性，这应该将平均效应降至零。为了实现它，我们同样调用refute_estimate函数，但使用placebo_treatment_refuter作为方法。我们还必须指定placebo_type和模拟次数（num_simulations）。我们将使用的安慰剂类型是permute，模拟次数越多越好，但这也会花费更长的时间。代码可以在以下片段中看到：

ref_placebo = causal_mdl.refute_estimate(
    method_name="placebo_treatment_refuter",
    placebo_type="permute", num_simulations=20
)
print(ref_placebo)

前面的代码输出如下：

Refute: Use a Placebo Treatment
Estimated effect:7227.904763676559
New effect:381.05420029741083
p value:0.32491556283289624

如前述输出所示，新的效应接近于零。然而，鉴于 p 值高于 0.05，我们不能拒绝 ATE 大于零的零假设。这告诉我们，估计的因果效应并不非常稳健。我们可能通过添加相关的混杂因素或使用不同的因果模型来改进它，但同样，实验设计存在我们无法修复的缺陷，例如银行根据风险因素偏袒地指定治疗方式。

任务完成

本章的任务有两个，如下所述：

创建一个公平的预测模型来预测哪些客户最有可能违约。
创建一个稳健的因果模型来估计哪些政策对客户和银行最有益。

关于第一个目标，我们已经根据四个公平性指标（SPD、DI、AOD、EOD）——在比较特权群体和弱势群体年龄组时——产生了四个具有偏差缓解方法的模型，这些模型在客观上比基础模型更公平。然而，根据 DFBA（参见图 11.7），只有其中两个模型在同时使用年龄组和性别时具有交叉公平性。通过结合方法，我们仍然可以显著提高公平性，但任何一种模型都能改进基础模型。

对于第二个目标，因果推断框架确定，所测试的任何政策对于双方来说都比没有政策要好。太好了！然而，它得出的估计并没有确立一个单一的获胜者。尽管如此，正如预期的那样，推荐的政策会根据客户的信用额度而变化——另一方面，如果我们旨在最大化银行利润，我们必须考虑信用卡的平均使用情况。盈利性的问题提出了我们必须协调的两个目标：制定对客户或银行最有利的推荐政策。

因此，如何程序上公平是一个复杂的问题，有许多可能的答案，任何解决方案都可能导致银行吸收与实施政策相关的部分成本。至于鲁棒性，尽管实验存在缺陷，但我们可以说我们的估计具有中等水平的鲁棒性，通过了一个鲁棒性测试但没有通过另一个。话虽如此，这完全取决于我们认为足够鲁棒以验证我们的发现。理想情况下，我们会要求银行开始一个新的无偏实验，但等待另外 6 个月可能不可行。

在数据科学中，我们经常发现自己在与有缺陷的实验和有偏差的数据打交道，并必须充分利用它们。因果推断通过分离原因和效果，包括估计及其相应的置信区间，提供了一种这样做的方法。然后我们可以提供带有所有免责声明的发现，以便决策者可以做出明智的决策。有偏差的决策会导致有偏差的结果，因此解决偏差的道德必要性可以从塑造决策开始。

摘要

在阅读本章之后，你应该了解如何通过视觉和指标在数据和模型中检测偏差，然后通过预处理、处理和后处理方法来减轻偏差。我们还通过估计异质处理效应、用它们做出公平的政策决策以及测试它们的鲁棒性来了解因果推断。在下一章中，我们也将讨论偏差，但学习如何调整模型以满足多个目标，包括公平性。

数据集来源

Yeh, I. C., & Lien, C. H. (2009). 比较数据挖掘技术在预测信用卡客户违约概率方面的准确性. 《专家系统与应用》，36(2)，2473-2480: dl.acm.org/doi/abs/10.1016/j.eswa.2007.12.020

进一步阅读

Chang, C., Chang, H.H., and Tien, J., 2017, 关于信息不对称下金融监管机构应对策略的研究：台湾信用卡市场案例研究. 《通用管理杂志》，5，429-436: doi.org/10.13189/ujm.2017.050903
Foulds, J., and Pan, S., 2020, An Intersectional Definition of Fairness. 2020 IEEE 36th International Conference on Data Engineering (ICDE), 1918-1921: arxiv.org/abs/1807.08362
Kamiran, F., and Calders, T., 2011, Data preprocessing techniques for classification without discrimination. Knowledge and Information Systems, 33, 1-33: link.springer.com/article/10.1007/s10115-011-0463-8
Feldman, M., Friedler, S., Moeller, J., Scheidegger, C., and Venkatasubramanian, S., 2015, Certifying and Removing DI. Proceedings of the 21st ACM SIGKDD International Conference on Knowledge Discovery and Data Mining: arxiv.org/abs/1412.3756
Kamishima, T., Akaho, S., Asoh, H., and Sakuma, J., 2012, Fairness-Aware Classifier with Prejudice Remover Regularizer. ECML/PKDD: dl.acm.org/doi/10.5555/3120007.3120011
A. Agarwal, A. Beygelzimer, M. Dudik, J. Langford, and H. Wallach, A Reductions Approach to Fair Classification, International Conference on Machine Learning, 2018. arxiv.org/pdf/1803.02453.pdf
Kearns, M., Neel, S., Roth, A., and Wu, Z., 2018, Preventing Fairness Gerrymandering: Auditing and Learning for Subgroup Fairness. ICML: arxiv.org/pdf/1711.05144.pdf
Pleiss, G., Raghavan, M., Wu, F., Kleinberg, J., and Weinberger, K.Q., 2017, On Fairness and Calibration. NIPS: arxiv.org/abs/1709.02012
Foster, D. and Syrgkanis, V., 2019, Orthogonal Statistical Learning. ICML: arxiv.org/abs/1901.09036

在 Discord 上了解更多

要加入这本书的 Discord 社区——在那里您可以分享反馈、向作者提问，并了解新版本——请扫描下面的二维码：

packt.link/inml