评估提升模型

因果数据科学

如何比较和选择最佳的提升模型

因果推断在行业中最广泛应用之一是提升建模，即条件平均治疗效应的估计。

当估计一个治疗（药物、广告、产品等）对一个感兴趣的结果（疾病、公司收入、客户满意度等）的因果效应时，我们通常不仅仅关心治疗在平均情况下是否有效，而且希望知道对于哪些对象（患者、用户、客户等）它的效果更好或更差。

估计异质增量效应或提升是改进所关注政策的定位的重要中间步骤。例如，我们可能想要提醒某些人他们更有可能经历药物的副作用，或者仅向特定一组客户展示广告。

虽然存在许多用于建模提升的方法，但在特定应用中并不总是清楚应该使用哪种方法。关键是，由于因果推断的基本问题，我们永远无法观察到所关注的目标——提升，因此我们无法像使用机器学习预测算法那样验证我们的估计器。我们无法设置一个验证集并选择表现最佳的模型，因为我们没有基准真值，即使在验证集中也是如此，即使我们进行了随机实验也是如此。

那么我们能做什么呢？在本文中，我试图介绍用于评估提升模型的最流行方法。如果您对提升模型不熟悉，我建议先阅读我的入门文章。

理解元学习器

编辑描述

towardsdatascience.com

提升和促销邮件

想象一下，我们在一家对改进电子邮件营销活动感兴趣的产品公司的市场部门工作。从历史上看，我们主要向新客户发送电子邮件。然而，现在我们想采用数据驱动的方法，针对对收入影响最大的客户进行定位。这种影响也称为提升或增量。

让我们看看我们手头有哪些数据。我从src.dgp中导入数据生成过程dgp_promotional_email()。我还从src.utils导入一些绘图函数和库。

from src.utils import *from src.dgp import dgp_promotional_emaildgp = dgp_promotional_email(n=500)df = dgp.generate_data()df.head()

我们有关于500个客户的信息，我们观察到他们是否是新客户，他们的年龄，他们在邮件营销活动之前产生的销售额（sales_old），他们是否收到了邮件，以及邮件营销活动之后的销售额。

我们感兴趣的结果是销售额，我们用字母Y表示。我们希望改进的治疗或政策是邮件营销活动，我们用字母W表示。我们将所有其余变量称为混淆变量或控制变量，用X表示。

Y = 'sales'W = 'mail'X = ['age', 'sales_old', 'new']

表示变量之间因果关系的有向无环图（DAG）如下所示。感兴趣的因果关系用绿色表示。

从DAG中我们可以看到，new顾客指示器是一个混淆变量，需要进行控制以便识别mail对sales的影响。而age和sales_old则对估计来说不是必要的，但对于识别可能会有帮助。如果想了解更多关于DAG和控制变量的信息，可以查看我的入门文章。

DAG和控制变量

编辑描述

towardsdatascience.com

提升建模的目标是恢复个体处理效应（ITE） τᵢ，即发送促销mail对sales的增量效应。我们可以将ITE表示为两个假设数量之间的差异：如果顾客接收了邮件，他们的潜在结果Yᵢ⁽¹⁾，减去如果顾客没有接收邮件，他们的潜在结果Yᵢ⁽⁰⁾。

请注意，对于每个顾客，我们只观察到两个实现的结果中的一个，具体取决于他们是否实际收到了mail。因此，ITE本质上是不可观测的。相反，可以估计的是条件平均处理效应（CATE），即在协变量X条件下的预期个体处理效应τᵢ。例如，对于年长顾客（age > 50），mail对sales的平均效应。

为了能够恢复CATE，我们需要做三个假设。

1. 非混淆性：Y⁽⁰⁾，Y⁽¹⁾ ⊥ W | X
2. 重叠性：0 < e(X) < 1
3. 一致性：Y = W ⋅ Y⁽¹⁾ + (1−W) ⋅ Y⁽⁰⁾

其中e(X)是倾向得分，即在协变量X条件下的被处理的预期概率。

接下来，我们将使用机器学习方法来估计CATE τ(x)，倾向得分e(x)，以及结果的条件期望函数（CEF）μ(x)

我们使用随机森林回归算法来建模CATE和结果CEF，而使用逻辑回归来建模倾向得分。

from sklearn.ensemble import RandomForestRegressorfrom sklearn.linear_model import LogisticRegressionCVmodel_tau = RandomForestRegressor(max_depth=2)model_y = RandomForestRegressor(max_depth=2)model_e = LogisticRegressionCV()

在本文中，我们没有对底层的机器学习模型进行微调，但强烈建议进行微调以提高提升模型的准确性（例如，使用自动机器学习库，如FLAML）。

提升模型

有许多方法可以建模提升或者说估计条件平均处理效应（CATE）。由于本文的目标是比较评估提升模型的方法，我们不会详细解释这些方法。如果你想了解简介，可以查看我关于元学习器的入门文章。

我们将考虑以下学习器：

• S-learner或单学习器，由Kunzel、Sekhon、Bickel、Yu（2017）引入
• T-learner或双学习器，由Kunzel、Sekhon、Bickel、Yu（2017）引入
• X-learner或交叉学习器，由Kunzel、Sekhon、Bickel、Yu（2017）引入
• R-learner或Robinson学习器，由Nie、Wager（2017）引入
• DR-learner或双重稳健学习器，由Kennedy（2022）引入

我们从Microsoft的econml库中导入所有模型。

from src.learners_utils import *from econml.metalearners import SLearner, TLearner, XLearnerfrom econml.dml import NonParamDMLfrom econml.dr import DRLearnerS_learner = SLearner(overall_model=model_y)T_learner = TLearner(models=clone(model_y))X_learner = XLearner(models=model_y, propensity_model=model_e, cate_models=model_tau)R_learner = NonParamDML(model_y=model_y, model_t=model_e, model_final=model_tau, discrete_treatment=True)DR_learner = DRLearner(model_regression=model_y, model_propensity=model_e, model_final=model_tau)

我们在数据上使用fit()模型，指定结果变量Y，处理变量W和协变量X。

names = ['SL', 'TL', 'XL', 'RL', 'DRL']learners = [S_learner, T_learner, X_learner, R_learner, DR_learner]for learner in learners:    learner.fit(df[Y], df[W], X=df[X])

我们现在准备好评估模型了！我们应该选择哪个模型呢？

Oracle损失函数

评估提升模型的主要问题是，即使有验证集，即使有随机实验或AB测试，我们也无法观察到我们感兴趣的指标：个体处理效应。事实上，我们只能观察到未处理客户的实际结果Yᵢ⁽⁰⁾和处理客户的实际结果Yᵢ⁽¹⁾。因此，我们无法为任何客户在验证数据中计算个体处理效应τᵢ = Yᵢ⁽¹⁾ – Yᵢ⁽⁰⁾。

我们还能做些什么来评估我们的估计器吗？

答案是肯定的，但在提供更多细节之前，让我们首先了解一下如果我们能观察到个体处理效应τᵢ会做什么。

Oracle MSE损失

如果我们能观察到个体处理效应（但我们不能，因此称为“oracle”属性），我们可以尝试测量我们的估计值τ̂(Xᵢ)与真实值τᵢ之间的差距。这是我们在机器学习中通常做的评估预测方法的方式：我们留出一个验证数据集，在该数据集上比较预测值和真实值。有很多损失函数可以评估预测准确性，所以让我们集中讨论最流行的一个：均方误差（MSE）损失。

def loss_oracle_mse(data, learner):    tau = learner.effect(data[X])    return np.mean((tau - data['effect_on_sales'])**2)

函数compare_methods打印并绘制在单独的验证数据集上计算的评估指标。

def compare_methods(learners, names, loss, title=None, subtitle='较低为较好'):    data = dgp.generate_data(seed_data=1, seed_assignment=1, keep_po=True)    results = pd.DataFrame({        'learner': names,        'loss': [loss(data.copy(), learner) for learner in learners]    })    fig, ax = plt.subplots(1, 1, figsize=(6, 4))    sns.barplot(data=results, x="learner", y='loss').set(ylabel='')    plt.suptitle(title, y=1.02)    plt.title(subtitle, fontsize=12, fontweight=None, y=0.94)    return resultsresults = compare_methods(learners, names, loss_oracle_mse, title='Oracle MSE损失')

在此示例中，我们可以看到T-learner表现最差，S-learner紧随其后。另一方面，X-learner、R-learner和DR-learner表现明显更好，其中DR-learner获胜。

然而，这可能不是评估我们的提升模型的最佳损失函数。事实上，提升建模只是我们最终目标的一个中间步骤：提高收入。

Oracle策略增益

由于我们的最终目标是提高收入，我们可以通过它们在给定策略函数下增加的收入来评估估计器。例如，假设我们发送一封电子邮件的成本为0.01美元。那么，我们的策略将是对于每个预测的条件平均处理效应大于0.01美元的客户进行处理。

cost = 0.01

我们的收入实际上会增加多少？让我们定义d(τ̂)为我们的策略函数，使得如果τ ≥ 0.1，则 d=1 ，否则 d=0 。那么我们的增益（较高为较好）函数是：

同样，这是一个无法计算的“oracle”损失函数，因为我们无法观察到个体处理效应。

def gain_oracle_policy(data, learner):    tau_hat = learner.effect(data[X])    return np.sum((data['effect_on_sales'] - cost) * (tau_hat > cost))results = compare_methods(learners, names, gain_oracle_policy, title='Oracle策略增益', subtitle='较高为较好')

在这种情况下，S-learner明显是最差的表现者，没有对收入产生影响。T-learner导致了适度的收益，而X-learner、R-learner和DR-learner都带来了总体收益，其中X-learner略领先。

实用损失函数

在前一节中，我们已经看到了两个我们希望能够计算的损失函数的示例，如果我们能够观察到个体处理效应τᵢ。然而，在实践中，即使有随机实验和验证集，我们也无法观察到ITE，即我们感兴趣的对象。现在，我们将介绍一些基于这个实际约束条件试图评估提升模型的度量标准。

结果损失

第一种也是最简单的方法是切换到不同的损失变量。虽然我们无法观察到个体治疗效果τᵢ，但我们仍然可以观察到我们的结果Yᵢ。这虽然不是我们真正感兴趣的对象，但我们可能期望在预测y方面表现良好的提升模型也能产生良好的τ估计。

这样的损失函数可以是结果均方差损失，它是用于预测方法的常规均方差损失函数。

这里的问题是，并非所有模型都直接产生μ(x)的估计。因此，我们跳过此比较，转而使用可以评估任何提升模型的方法。

预测对预测损失

另一种非常简单的方法是将在训练集上训练的模型的预测与在验证集上训练的另一个模型的预测进行比较。虽然直观，但这种方法可能会非常误导人。

def loss_pred(data, learner):    tau = learner.effect(data[X])    learner2 = copy.deepcopy(learner).fit(data[Y], data[W], X=data[X])    tau2 = learner2.effect(data[X])    return np.mean((tau - tau2)**2)results = compare_methods(learners, names, loss_pred, '预测对预测损失')

毫不奇怪，这个度量表现非常糟糕，你应该永远不要使用它，因为它奖励一致的模型，而不管它们的质量如何。一个总是为每个观测预测随机常数CATE的模型将获得完美的分数。

分布损失

另一种方法是问：我们能有多好地匹配潜在结果的分布？我们可以针对已受治疗或未受治疗的潜在结果进行这个练习。让我们以最后一种情况为例。假设我们取未收到邮件的客户的观察到的销售额和收到邮件的客户的观察到的销售额减去估计的CATE τ̂(x)。根据无偏性假设，这两个未经处理的潜在结果的分布应该在给定协变量X的条件下相似。

因此，我们期望如果我们正确估计了治疗效应，这两个分布之间的距离会很近。

我们也可以对已接受治疗的潜在结果进行同样的练习。

我们使用能量距离作为距离度量。

from dcor import energy_distancedef loss_dist(data, learner):    tau = learner.effect(data[X])    data.loc[data.mail==1, 'sales'] -= tau[data.mail==1]    return energy_distance(data.loc[data.mail==0, [Y] + X], data.loc[data.mail==1, [Y] + X], exponent=2)results = compare_methods(learners, names, loss_dist, '分布损失')

这个度量非常嘈杂，并且奖励了S-learner和T-learner，实际上它们是表现最差的两个模型。

上下中位数差异

上下中位数损失试图回答这个问题：我们的提升模型是否检测到任何异质性？特别是，如果我们将验证集划分为预测提升τ̂(x)的上中位数和下中位数两个样本，使用均值差估计器估计平均效应的实际差异有多大？我们期望更好的估计器能够更好地将样本分成高效果和低效果。

from statsmodels.formula.api import ols def loss_ab(data, learner):    tau = learner.effect(data[X]) + np.random.normal(0, 1e-8, len(data))    data['above_median'] = tau >= np.median(tau)    param = ols('sales ~ mail * above_median', data=data).fit().params[-1]    return paramresults = compare_methods(learners, names, loss_ab, title='上下中位数差异', subtitle='越高越好')

不幸的是，上下中位数差异奖励了T-learner，它是最差的模型之一。

需要注意的是，两组（上中位数和下中位数τ̂(x)）的均值差估计器不能保证无偏，即使数据来自随机实验。事实上，我们将两组数据在高度内源性的变量τ̂(x)上进行了划分。因此，这个方法应该谨慎使用。

提升曲线

上下中位数测试的一个扩展是提升曲线。这个想法很简单：不再根据中位数（0.5分位数）将样本分成两组，而是将数据分成更多组（更多分位数）。

对于每一组，我们计算均值差估计，并将其累积总和绘制在相应的分位数上。结果称为提升曲线。解释很简单：曲线越高，我们越能将高效果和低效果的观测结果分离开来。然而，同样的免责声明也适用：均值差估计器不是无偏的。因此，它们应该谨慎使用。

def generate_uplift_curve(df):    Q = 20    df_q = pd.DataFrame()    data = dgp.generate_data(seed_data=1, seed_assignment=1, keep_po=True)    ate = np.mean(data[Y][data[W]==1]) - np.mean(data[Y][data[W]==0])    for learner, name in zip(learners, names):        data['tau_hat'] = learner.effect(data[X])        data['q'] = pd.qcut(-data.tau_hat + np.random.normal(0, 1e-8, len(data)), q=Q, labels=False)        for q in range(Q):            temp = data[data.q <= q]            uplift = (np.mean(temp[Y][temp[W]==1]) - np.mean(temp[Y][temp[W]==0])) * q / (Q-1)            df_q = pd.concat([df_q, pd.DataFrame({'q': [q], 'uplift': [uplift], 'learner': [name]})], ignore_index=True)        fig, ax = plt.subplots(1, 1, figsize=(8, 5))    sns.lineplot(x=range(Q), y=ate*range(Q)/(Q-1), color='k', ls='--', lw=3)    sns.lineplot(x='q', y='uplift', hue='learner', data=df_q);    plt.suptitle('提升曲线', y=1.02, fontsize=28, fontweight='bold')    plt.title('越高越好', fontsize=14, fontweight=None, y=0.96)generate_uplift_curve(df)

虽然可能不是评估提升模型的最佳方法，但提升曲线在理解和实施提升模型方面非常重要。事实上，对于每个模型，它告诉我们随着增加接受处理的人口比例（x轴），预期的平均治疗效应（y轴）是多少。

最近邻匹配

我们分析的最后几种方法使用聚合数据来了解这些方法是否适用于更大的群体。而最近邻匹配则试图了解提升模型预测个体治疗效应的能力。然而，由于个体治疗效应是不可观测的，它尝试通过匹配可观测特征X上的接受处理和对照观测来构建一个代理。

例如，如果我们取所有接受处理的观测（i: Wᵢ=1），并找到控制组中的最近邻（NN₀(Xᵢ)），相应的均方误差（MSE）损失函数为

from scipy.spatial import KDTreedef loss_nn(data, learner):    tau_hat = learner.effect(data[X])    nn0 = KDTree(data.loc[data[W]==0, X].values)    control_index = nn0.query(data.loc[data[W]==1, X], k=1)[-1]    tau_nn = data.loc[data[W]==1, Y].values - data.iloc[control_index, :][Y].values    return np.mean((tau_hat[data[W]==1] - tau_nn)**2)results = compare_methods(learners, names, loss_nn, title='最近邻损失')

在这种情况下，最近邻损失表现得相当好，识别出了两种效果最差的方法，即S-学习器和T-学习器。

IPW损失

逆概率加权（IPW）损失函数首次由Gutierrez、Gerardy（2017）提出，它是我们将要看到的三个使用伪结果Y*评估估计器的指标中的第一个。伪结果是其期望值是条件平均治疗效应的变量，但由于太不稳定而不能直接用作估计值。关于伪结果的更详细解释，请参阅我关于因果回归树的文章。与IPW损失相对应的伪结果为

因此，相应的损失函数为

def loss_ipw(data, learner):    tau_hat = learner.effect(data[X])    e_hat = clone(model_e).fit(data[X], data[W]).predict_proba(data[X])[:,1]    tau_gg = data[Y] * (data[W] - e_hat) / (e_hat * (1 - e_hat))    return np.mean((tau_hat - tau_gg)**2)results = compare_methods(learners, names, loss_ipw, title='IPW损失')

IPW损失函数非常嘈杂。解决方法是使用其更稳健的变体，即R损失或DR损失，我们将在下面介绍。

R损失

R损失是由Nie，Wager（2017）与R-learner一起引入的，它实质上是R-learner的目标函数。与IPW损失一样，其思想是尝试匹配一个伪结果，其期望值是条件平均处理效应。

相应的损失函数为

def loss_r(data, learner):    tau_hat = learner.effect(data[X])    y_hat = clone(model_y).fit(df[X + [W]], df[Y]).predict(data[X + [W]])    e_hat = clone(model_e).fit(df[X], df[W]).predict_proba(data[X])[:,1]    tau_nw = (data[Y] - y_hat) / (data[W] - e_hat)    return np.mean((tau_hat - tau_nw)**2)compare_methods(learners, names, loss_r, title='R损失')

R损失的噪声明显较小，而且明确地将S-learner与其他方法隔离开来。但是，它往往更偏向于其对应的学习器，即R-learner。

DR损失

DR损失是DR-learner的目标函数，由Saito，Yasui（2020）首次引入。与IPW损失和R损失一样，其思想是尝试匹配一个伪结果，其期望值是条件平均处理效应。DR伪结果与AIPW估计量密切相关，也被称为双重稳健估计量，因此得名DR。

相应的损失函数为

def loss_dr(data, learner):    tau_hat = learner.effect(data[X])    y_hat = clone(model_y).fit(df[X + [W]], df[Y]).predict(data[X + [W]])    mu1 = clone(model_y).fit(df[X + [W]], df[Y]).predict(data[X + [W]].assign(mail=1))    mu0 = clone(model_y).fit(df[X + [W]], df[Y]).predict(data[X + [W]].assign(mail=0))    e_hat = clone(model_e).fit(df[X], df[W]).predict_proba(data[X])[:,1]    tau_nw = mu1 - mu0 + (data[Y] - y_hat) * (data[W] - e_hat) / (e_hat * (1 - e_hat))    return np.mean((tau_hat - tau_nw)**2)results = compare_methods(learners, names, loss_dr, title='DR损失')

至于R-loss，DR-loss倾向于支持相应的学习器DR-learner。然而，从算法准确性的角度来看，它提供了更准确的排名。

经验策略增益

我们要分析的最后一个损失函数与迄今为止我们所见过的所有其他函数不同，因为它不关注我们能够如何准确地估计治疗效应，而是关注相应的最优治疗策略的表现。具体而言，Hitsch、Misra、Zhang（2023）提出了以下增益函数：

其中c是治疗成本，d是在估计的CATE τ̂(Xᵢ)下的最优治疗策略。在我们的案例中，我们假设个体治疗成本为c=0.01美元，因此最优策略是对每个估计的CATE大于0.01的客户进行治疗。

项Wᵢ⋅d(τ̂)和(1-Wᵢ)⋅(1-d(τ̂))意味着我们仅对实际治疗W与最优治疗策略d相对应的个体进行计算。

def gain_policy(data, learner):    tau_hat = learner.effect(data[X])    e_hat = clone(model_e).fit(data[X], data[W]).predict_proba(data[X])[:,1]    d = tau_hat > cost    return np.sum((d * data[W] * (data[Y] - cost)/ e_hat + (1-d) * (1-data[W]) * data[Y] / (1-e_hat)))results = compare_methods(learners, names, gain_policy, title='经验策略增益', subtitle='数值越大越好')

经验策略增益表现非常好，将S-learner和T-learner两种最差的方法排除在外。

元研究

在本文中，我们介绍了多种评估提升模型（即条件平均治疗效应估计器）的方法。我们还在我们的模拟数据集上进行了测试，这是一个非常特殊且有限的示例。这些指标在一般情况下的表现如何？

Schuler、Baiocchi、Tibshirani、Shah（2018）比较了S-loss、T-loss、R-loss在对应的估计器上的模拟数据。他们发现R-loss“是验证集度量值中，当优化时最一致地选择高性能模型的度量值”。作者还发现了所谓的亲和性偏差：R-loss或DR-loss等指标倾向于对应的学习器。

Curth、van der Schaar（2023）从理论角度研究了更广泛的学习器。他们发现“在我们考虑的所有实验条件下，没有现有的选择标准在全局范围内都是最好的”。

Mahajan、Mitliagkas、Neal、Syrgkanis（2023）在范围上进行了最全面的研究。作者比较了144个数据集和415个估计器上的许多指标。他们发现“没有任何指标明显优于其他指标”，但“使用DR元素的指标似乎总是候选赢家之一”。

结论

在本文中，我们探讨了多种评估提升模型的方法。最主要的挑战是无法观测到感兴趣的变量，即个体治疗效应。因此，不同的方法尝试使用其他变量、使用代理结果或近似隐含最优策略的效果来评估提升模型。

因为在理论和实证角度上都没有共识，所以很难推荐使用单一的方法来进行推荐。使用R-和DR-元素的损失函数往往表现得更好，但也对相应的学习器有偏向性。然而，了解这些度量如何工作可以帮助我们理解它们的偏见和局限性，从而根据具体情况做出最合适的决策。

参考文献

• Curth, van der Schaar (2023), “寻找洞见，而非魔法子弹：向异质处理效应估计的模型选择困境的揭秘”
• Gutierrez, Gerardy (2017), “因果推断和增益建模：文献综述”
• Hitsch, Misra, Zhang (2023), “异质处理效应和最优目标定向政策评估”
• Kennedy (2022), “朝着异质因果效应的最优双重鲁棒估计”
• Kunzel, Sekhon, Bickel, Yu (2017), “使用机器学习估计异质处理效应的元学习器”
• Mahajan, Mitliagkas, Neal, Syrgkanis (2023), “针对异质因果效应估计的模型选择的实证分析”
• Nie, Wager (2017), “伪正交估计异质处理效应”
• Saito, Yasui (2020), “反事实交叉验证：因果推断模型的稳定模型选择过程”
• Schuler, Baiocchi, Tibshirani, Shah (2018), “在估计个体处理效应时的方法比较”

相关文章

• 理解元学习器
• 理解AIPW，双重鲁棒估计器
• 理解因果树
• 从因果树到森林

代码

您可以在这里找到原始的Jupyter Notebook：

Blog-Posts/notebooks/evaluate_uplift.ipynb at main · matteocourthoud/Blog-Posts

用于VoAGI博客文章的代码和笔记本。通过创建一个来为matteocourthoud/Blog-Posts贡献代码和笔记本…

github.com

感谢您阅读！

我非常感谢！🤗 如果您喜欢这篇文章并希望看到更多，请关注我。我每周发布与因果推断和数据分析相关的主题。我尽量保持我的帖子简单而精确，始终提供代码，示例和模拟。

另外，一个小小的免责声明：我写作是为了学习，所以错误是常态，尽管我尽力而为。如果您发现错误，请告诉我。我也欢迎对新主题的建议！