本章介绍了因果推断应用于行业时可能最有趣的发展:效应异质性。在此之前,你了解的是一种处理方式的总体影响。现在,你将聚焦于探究它如何对不同的人产生不同影响。“处理效应并非恒定不变”这个理念很简单,但却极具力量。知晓哪些个体对特定处理反应更好,是决定谁能接受该处理的关键。效应异质性为备受推崇的个性化理念提供了一种因果推断方法。

你将首先从理论层面理解效应异质性,了解估计效应异质性时存在的挑战,以及如何拓展已学知识来应对这些挑战。接下来,你会发现估计异质效应与预测问题密切相关,而预测问题对数据科学家来说已经非常熟悉。因此,你会看到交叉验证和模型选择的理念在处理异质性模型中仍然适用。不过,验证你的效应估计要比评估简单的预测模型困难得多,这就是为什么你会看到一些关于如何进行验证的新颖思路。

本章以一些关于如何利用效应异质性来指导决策的准则和例子收尾。尽管这些例子并非详尽无遗,但我希望它们能让你了解如何将这些理念应用到你自己的业务问题中。

从平均处理效应到条件平均处理效应

到目前为止,每次你估计一种处理的因果影响时,大多是平均处理效应:

\(\tau = E[Y_1 - Y_0]\)

或者是连续处理的等价形式:

\(\tau = E[y^\prime(t)]\)

其中 \(y^\prime(t)\) 是处理响应函数的导数。也就是说,你已经学习了用于揭示一种处理总体有效性的技术。平均处理效应估计是因果推断的基石。对于被称为项目评估的决策问题来说,它是一个极其有用的工具:当你想知道是否应该向全体人群推行一种处理时,就会用到它。

现在,是时候学习如何为另一种类型的决策提供依据了:你应该对哪个单元进行处理?为此,你需要允许决策因单元不同而变化。例如,给一位顾客发放折扣券可能是有益的,但对另一位顾客可能并非如此,因为一位顾客可能对折扣更敏感。或者,优先向一个群体而非另一个群体提供疫苗可能是有意义的,因为这些群体能从这种处理中获益更多。在这类情况下,个性化是关键。

实现这种个性化的一种方法是考虑效应异质性,这涉及到估计 \(CATE\)通过考虑每个单元的独特特征,你可以确定针对该特定情况的最有效处理方式:

\(E[Y_1 - Y_0|X]\)\(E[y^\prime(t)|X]\)

\(X\) 的条件限制意味着,现在你允许处理效应根据每个单元的协变量 \(X\) 所定义的特征而有所不同。同样,在这里,你认为并非所有实体对处理的反应都同样良好,并且你希望利用这种异质性。你希望只对合适的单元进行处理(在二元情况下),或者为每个单元找出最优的处理剂量(在连续情况下)

例如,如果你是一家银行,必须决定每位客户有资格获得的贷款,你可以肯定,把大量资金贷给所有人不是个好主意,尽管对一些人来说这可能是合理的。你必须明智地对待贷款金额这一 “处理”。或许,根据客户的信用评分,你可以算出合适的贷款金额。当然,你不必是大型机构才能利用个性化。适用个性化的例子不胜枚举。一年中哪些日子应该进行促销?一款产品应该定价多少?对每个人来说,多少运动量算是运动过量?

为什么预测并非答案

可以这样想:你有一群客户和一种处理方式(价格、折扣、贷款……),并且你想要对这种处理方式进行个性化设置 —— 例如,给不同的客户提供不同的折扣。假设你可以按照以下 “按结果划分的处理”图表来对客户进行分类:

你可以将个性化任务视为客户细分问题。你希望根据客户对处理的反应来创建客户群体。例如,假设你想找到对折扣反应良好的客户和对折扣反应不佳的客户。嗯,客户对处理的反应由条件处理效应 \(\frac{\delta Y}{\delta T}\) 给出。所以,如果你能以某种方式为每个客户估计出这个值,你就能将那些对处理反应很好(处理效应高)的客户和那些对处理反应不太好的客户分组。如果你这样做了,你对客户空间的划分会有点像下面的图片:

那会很棒,因为现在你将能够估计每个群体的不同处理效应。同样,由于效应就是处理响应函数的斜率,如果你能划分出斜率不同的群体,那么这些分组里的个体对处理会有不同的响应程度:

现在,将这与你用传统机器学习方法得到的结果进行对比。你可能会尝试预测 \(Y\),而不是每个单元的导数 \(\frac{\delta Y}{\delta T}\)。这本质上会在 \(Y\) 轴上对空间进行划分,假设你的预测模型能很好地近似目标。然而,这不一定会得到具有不同处理效应的群体。这就是为什么仅仅预测结果对于决策制定来说并不总是有用的:

好的,你可能会说,我明白我必须估计效应,而不只是预测结果,但这有点棘手。如果我看不到 \(\frac{\delta Sales}{\delta Discount}\) 这个斜率,我怎么能预测它呢?

这是个好问题。与原始结果 \(Y\) 不同,斜率或变化率在个体层面上本质上是无法直接观测到的。要看到个体的斜率,你必须观察每个单元在不同处理水平下的情况,并计算在这些不同处理下结果的变化情况:

\(\frac{\delta Y_i}{\delta T_i} \approx \frac{Y(T_i) - Y(T_i + \epsilon)}{T_i - (T_i + \epsilon)}\)

这又回到了因果推断的根本问题。你永远无法在不同的处理条件下观察同一个单元。那么,你能做些什么呢?


条件平均处理效应与个体处理效应

要记住,条件平均处理效应\((CATE)\)与个体处理效应\((ITE)\)是不同的。例如,假设你有两个组,\(x = 1\)\(x = 2\),每组各有 \(4\) 个单元,你想了解一种新药对一种通常会导致 \(50\%\) 患者死亡的疾病的疗效。对于组 \(x = 1\),这种药物对一名患者有害,导致其死亡,但救活了另一名患者。对于组 \(x = 2\),药物的效果为零,其中 \(1\) 人死亡(要记住,这种疾病的致死率为 50%)。在这两个组中,\(CATE\) 都是 \(0.5\),但没有单个单元的 \(ITE\)\(0.5\)


基于回归的条件平均处理效应

我想你可能已经猜到了:就像应用因果推断中的大多数内容一样,答案往往从线性回归开始。但在采用这种方法之前,让我们把事情弄得更具体一些。假设你在一家在全国范围内经营的连锁餐厅工作。这项业务的一个关键部分是了解 应该在什么时候向顾客提供折扣。出于这个原因,该公司开展了一项持续三年的全国性实验,在连锁旗下的六家不同餐厅中随机提供折扣。数据存储在以下数据框中:


import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
import seaborn as sns

import matplotlib

from cycler import cycler

color=['0.0', '0.4', '0.8']
default_cycler = (cycler(color=color))
linestyle=['-', '--', ':', '-.']
marker=['o', 'v', 'd', 'p']

plt.rc('axes', prop_cycle=default_cycler)

matplotlib.rcParams.update({'font.size': 18})
data = pd.read_csv("./data/daily_restaurant_sales.csv")
data.head()
##    rest_id         day  month  ...  competitors_price  discounts  sales
## 0        0  2016-01-01      1  ...               2.88          0   79.0
## 1        0  2016-01-02      1  ...               2.64          0   57.0
## 2        0  2016-01-03      1  ...               2.08          5  294.0
## 3        0  2016-01-04      1  ...               3.37         15  676.5
## 4        0  2016-01-05      1  ...               3.79          0   66.0
## 
## [5 rows x 11 columns]

你的目标是了解什么时候是提供折扣的最佳时机。在这些数据中,每个餐厅和日期的组合对应一行。这与本书中使用的大多数例子略有不同,那些例子中的单元是顾客。而现在,单元是 “日期 - 餐厅” 的组合。即便如此,你仍然可以应用之前的相同推理,只是不再是对顾客进行处理或者说提供折扣,而是对日期进行处理或者说提供折扣。你也可以在每个餐厅的每一天设置不同的价格,但为了简化问题,我们假设各餐厅的价格保持一致。

你可以将这个业务问题表述为一个条件平均处理效应估计问题。如果你能创建一个模型,输出每一天以及给定协变量下销售额对折扣的敏感度,即:

\(\frac{\partial}{\partial t} E\left[Sales(t) \mid X\right]\)

那么,你就可以使用这个模型来决定何时提供折扣以及提供多少折扣。


条件平均处理效应的识别

在本章中,你不必过于担心识别问题,因为在评估集中,处理是随机分配的。然而,估计条件平均处理效应的整个思路是基于使\(E\left[Sales(t) \mid X\right] = E\left[Sales \mid T = t, X\right]\)


现在你有了更具体的内容来着手处理,让我们看看回归如何能帮助你。回想一下,你之前处于一个复杂的情况中。你需要预测 \(\frac{\delta Y_i}{\delta T_i}\),但遗憾的是,它是不可观测的。所以,这并不像你可以简单地使用一个机器学习算法并将其作为目标变量那样。但也许你不需要观测 \(\frac{\delta Y_i}{\delta T_i}\) 就能对其进行预测。

例如,假设你用你的数据拟合了以下线性模型:

\(y_i = \beta_0 + \beta_1 t_i + \beta_2 X_i + e_i\)

如果你对处理变量求导,最终会得到:

\(\frac{\delta y_i}{\delta t_i} = \beta_1\)

在处理变量是随机分配的情况下,这就是 \(ATE\)

由于你可以估计前面的模型来得到 \(\widehat{\beta}_1\),你甚至可能进而认为,即使无法观测到斜率,你也能对其进行预测。在这个例子中,这是一个相当简单的预测。你在为每个人预测恒定值 \(\widehat{\beta}_1\)。这确实是有意义的,但并非你真正想要的。这是 \(ATE\),而不是 \(CATE\)。这对你弄清楚何时提供折扣的任务没有帮助,原因很简单,每个单元即日期与餐厅的组合得到的斜率预测都是相同的。

为了改进它,你可以做如下简单的修改:

\(y_i = \beta_0 + \beta_1 t_i + \beta_2 X_i + \beta_3 t_i X_i + e_i\)

相应地,这会给你如下的斜率预测:

\(\frac{\widehat{\delta y_i}}{\delta t_i} = \widehat{\beta_1} + \widehat{\beta_3} X_i\)

其中 \(\beta_3\) 是针对 \(X\) 中特征的向量系数。现在你取得进展了!由不同 \(X_i\) 定义的每个个体,都会有不同的斜率预测。换句话说,斜率预测会随着 \(X\) 的变化而变化。直观地说,纳入处理与协变量之间的交互项,能让模型了解效应如何随这些协变量而变化。这就是回归能为你提供一种估计 \(CATE\) 的方法的原因,即便你无法直接对其进行预测。

现在理论讲得够多了。让我们看看如何把这个用代码实现。首先,你需要定义协变量。在这个例子中,协变量基本上是特定日期的特征,比如月份、星期几,以及是否是节假日。我还纳入了竞争对手的平均价格,因为这可能会影响顾客对每家餐厅折扣的反应。

一旦你有了协变量,就需要将它们与处理变量进行交互。\(*\) 运算符正好能做到这一点。它会为左右两边创建一个加性项,再加上一个交互项。

例如,\(a*b\) 会在你的回归中包含项 \(a\)\(b\) 以及 \(a*b\)。在你的例子中,这会产生如下结果:

\(sales_i = \beta_0 + \beta_1 discount_i + \beta_2 X_i*discount_i + \beta_3 X_i + e_i\)

import statsmodels.formula.api as smf

X = ["C(month)", "C(weekday)", "is_holiday", "competitors_price"]
regr_cate = smf.ols(f"sales ~ discounts*({'+'.join(X)})",data=data).fit()
print(regr_cate.summary())
##                             OLS Regression Results                            
## ==============================================================================
## Dep. Variable:                  sales   R-squared:                       0.964
## Model:                            OLS   Adj. R-squared:                  0.964
## Method:                 Least Squares   F-statistic:                     5231.
## Date:                  周二, 30 9月 2025   Prob (F-statistic):               0.00
## Time:                        14:52:08   Log-Likelihood:                -42718.
## No. Observations:                7679   AIC:                         8.552e+04
## Df Residuals:                    7639   BIC:                         8.579e+04
## Df Model:                          39                                         
## Covariance Type:            nonrobust                                         
## ================================================================================================
##                                    coef    std err          t      P>|t|      [0.025      0.975]
## ------------------------------------------------------------------------------------------------
## Intercept                       41.9451      5.765      7.275      0.000      30.643      53.247
## C(month)[T.2]                   -5.4654      5.250     -1.041      0.298     -15.757       4.826
## C(month)[T.3]                  -11.7369      5.106     -2.299      0.022     -21.747      -1.727
## C(month)[T.4]                   -5.9218      5.174     -1.145      0.252     -16.064       4.221
## C(month)[T.5]                  -14.2763      5.097     -2.801      0.005     -24.268      -4.285
## C(month)[T.6]                  -15.6976      5.290     -2.967      0.003     -26.067      -5.328
## C(month)[T.7]                  -10.1173      5.178     -1.954      0.051     -20.267       0.033
## C(month)[T.8]                  -13.9666      5.189     -2.692      0.007     -24.138      -3.795
## C(month)[T.9]                   -9.2454      5.169     -1.789      0.074     -19.379       0.888
## C(month)[T.10]                   0.9170      5.103      0.180      0.857      -9.087      10.921
## C(month)[T.11]                   8.4826      5.152      1.646      0.100      -1.617      18.583
## C(month)[T.12]                  14.8578      5.148      2.886      0.004       4.767      24.949
## C(weekday)[T.1]                 -0.0569      4.030     -0.014      0.989      -7.957       7.843
## C(weekday)[T.2]                  0.7152      4.070      0.176      0.861      -7.264       8.694
## C(weekday)[T.3]                  6.6486      4.056      1.639      0.101      -1.302      14.599
## C(weekday)[T.4]                 11.9731      4.053      2.954      0.003       4.029      19.918
## C(weekday)[T.5]                 17.9391      4.023      4.459      0.000      10.052      25.826
## C(weekday)[T.6]                  7.9323      4.046      1.960      0.050       0.001      15.864
## is_holiday[T.True]              11.0338      6.852      1.610      0.107      -2.397      24.465
## discounts                       54.0127      0.413    130.846      0.000      53.203      54.822
## discounts:C(month)[T.2]         -6.2683      0.376    -16.686      0.000      -7.005      -5.532
## discounts:C(month)[T.3]         -6.6509      0.376    -17.695      0.000      -7.388      -5.914
## discounts:C(month)[T.4]         -7.9060      0.374    -21.163      0.000      -8.638      -7.174
## discounts:C(month)[T.5]        -10.2032      0.368    -27.718      0.000     -10.925      -9.482
## discounts:C(month)[T.6]        -10.8734      0.378    -28.778      0.000     -11.614     -10.133
## discounts:C(month)[T.7]         -6.5874      0.366    -18.007      0.000      -7.305      -5.870
## discounts:C(month)[T.8]         -9.3171      0.367    -25.365      0.000     -10.037      -8.597
## discounts:C(month)[T.9]         -6.9055      0.374    -18.464      0.000      -7.639      -6.172
## discounts:C(month)[T.10]        -2.1312      0.364     -5.860      0.000      -2.844      -1.418
## discounts:C(month)[T.11]        -0.0420      0.381     -0.110      0.912      -0.789       0.705
## discounts:C(month)[T.12]         1.6115      0.367      4.388      0.000       0.892       2.331
## discounts:C(weekday)[T.1]        0.2610      0.283      0.921      0.357      -0.295       0.817
## discounts:C(weekday)[T.2]        0.5658      0.285      1.988      0.047       0.008       1.124
## discounts:C(weekday)[T.3]        2.2163      0.287      7.735      0.000       1.655       2.778
## discounts:C(weekday)[T.4]        6.2810      0.283     22.196      0.000       5.726       6.836
## discounts:C(weekday)[T.5]        3.1135      0.284     10.963      0.000       2.557       3.670
## discounts:C(weekday)[T.6]        0.2112      0.285      0.740      0.459      -0.348       0.771
## discounts:is_holiday[T.True]     2.2529      0.492      4.577      0.000       1.288       3.218
## competitors_price                0.9388      0.595      1.577      0.115      -0.228       2.106
## discounts:competitors_price     -3.2097      0.042    -77.078      0.000      -3.291      -3.128
## ==============================================================================
## Omnibus:                     1320.497   Durbin-Watson:                   1.897
## Prob(Omnibus):                  0.000   Jarque-Bera (JB):            23624.005
## Skew:                           0.272   Prob(JB):                         0.00
## Kurtosis:                      11.575   Cond. No.                     1.61e+03
## ==============================================================================
## 
## Notes:
## [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
## [2] The condition number is large, 1.61e+03. This might indicate that there are
## strong multicollinearity or other numerical problems.

\(*\) \(:\) 运算符

如果你只想要乘法项,可以在公式中使用 \(:\) 运算符。


一旦你估计了模型,预测的斜率可以从参数估计值中提取:

\(\frac{\widehat{\delta sales_i}}{\delta discounts_i} = \widehat{\beta_1} + \widehat{\beta_3} X_i\)

其中 \(\beta_1\) 是折扣系数,\(\beta_3\) 是交互作用系数的向量。你可以直接从拟合的模型中提取这些参数,但获得斜率预测更简单的方法是使用导数的定义:

\(\frac{\delta y}{\delta t} = \frac{y(t + \epsilon) - y(t)}{(t + \epsilon) - t}\)

\(\epsilon\) 趋近于零时。你可以通过将 \(\epsilon\) 替换为 1 来近似这个定义:

\(\frac{\delta y}{\delta t} \approx \hat{y}(t + 1) - \hat{y}(t)\)

其中 \(\hat{y}\) 由你的模型预测值给出。由于这是一个线性模型,这种近似是精确的。

换句话说,你将用你的模型进行两次预测:一次使用原始数据,另一次使用原始数据,但处理变量增加一个单位。这些预测值之间的差异就是你的 \(CATE\) 预测值。下面是用一些代码展示的样子:


ols_cate_pred = regr_cate.predict(data.assign(discounts=data["discounts"]+1)) - regr_cate.predict(data)

好的,你已经有了 \(CATE\) 模型及其预测结果。但仍有一个潜在的问题:它有多好?换句话说,你如何评估这个模型?正如你可能意识到的,在这里比较实际值和预测值是行不通的,因为实际的处理效应在个体层面是无法观测到的。


实际案例:价格歧视

在微观经济文献中,本章所使用的例子就是所谓的价格歧视。尽管这个名字听起来不太好,但它仅仅意味着企业可以将消费者区分为那些愿意支付更多的人,并向他们收取更高的费用。一个非常著名的价格歧视例子是,航空公司根据机票的提前购买时间来调整机票价格:需要预订下周航班的顾客预计要比预订明年航班的顾客支付多得多的费用。这被称为跨期价格歧视,因为企业能够根据时间来区分顾客的价格敏感度。这与你在本章中看到的餐厅例子是非常相似的情况。

一个更声名狼藉的例子是,一家葡萄酒公司用两个不同的瓶子销售完全相同的葡萄酒,其中一个作为优质酒来营销,价格要高得多,另一个作为普通酒来营销,以更适中的价格出售。价格歧视的第三种方式是,向学生提供半价入场券。在这种情况下,公司知道学生平均收入较低,这意味着他们的消费能力较弱。


评估条件平均处理效应预测

如果你有传统数据科学的背景,你可能会发现这种 \(CATE\) 预测看起来很像常规的机器学习预测,但目标很狡猾,在个体层面是无法观测到的。这意味着,传统机器学习中使用的很多模型评估技术 —— 比如交叉验证 —— 在这里仍然适用,而其他一些技术则需要进行调整。

所以,为了遵循传统做法,我们把数据分成训练集和测试集。由于数据存在时间维度,我们就利用这一点。训练集将包含 \(2016\) 年和 \(2017\) 年的数据,测试集则包含 \(2018\) 年及以后的数据:


train = data.query("day<'2018-01-01'")
test = data.query("day>='2018-01-01'")

现在,让我们重新拟合之前用于 \(CATE\) 的回归模型,但只使用训练数据进行估计,并在测试集上进行预测:

X = ["C(month)", "C(weekday)", "is_holiday", "competitors_price"]
regr_model = smf.ols(f"sales ~ discounts*({'+'.join(X)})",data=train).fit()

cate_pred = regr_model.predict(test.assign(discounts=test["discounts"]+1)) - regr_model.predict(test)
print(cate_pred)
## 731     41.355802
## 732     44.743887
## 733     39.783798
## 734     40.770278
## 735     40.666949
##           ...    
## 7674    42.446960
## 7675    28.696779
## 7676    38.631153
## 7677    34.023585
## 7678    26.683185
## Length: 2562, dtype: float64

为了让事情更有意思,我们用一个纯粹的预测性机器学习模型来对这个回归模型进行基准测试。这个机器学习模型只是试图预测结果 \(Y\)

from sklearn.ensemble import GradientBoostingRegressor

X = ["month", "weekday", "is_holiday", "competitors_price", "discounts"]
y = "sales"

np.random.seed(1)
ml_model = GradientBoostingRegressor(n_estimators=50).fit(train[X],train[y])
ml_pred = ml_model.predict(test[X])
print(ml_pred)
## [236.3129596  470.21805035 429.18065202 ... 411.55117103 390.74557782
##  450.46119734]

最后,我们还要在比较中纳入一个非常糟糕的模型。这个模型只是输出介于 \(-1\)\(1\) 之间的随机数。

这显然是无意义的,但作为一个值得关注的基准很有意思。最终,你想知道通过 \(CATE\) 模型来分配处理措施是否会比完全随机分配更好,而最后这个模型做的就是完全随机分配。

为了方便起见,我把所有内容都存储在一个新的数据框 \(test\_pred\) 中:

np.random.seed(123)

test_pred = test.assign(
    ml_pred=ml_pred,
    cate_pred=cate_pred,
    rand_m_pred=np.random.uniform(-1, 1, len(test)),
)
test_pred.head()
##      rest_id         day  month  ...     ml_pred  cate_pred  rand_m_pred
## 731        0  2018-01-01      1  ...  236.312960  41.355802     0.392938
## 732        0  2018-01-02      1  ...  470.218050  44.743887    -0.427721
## 733        0  2018-01-03      1  ...  429.180652  39.783798    -0.546297
## 734        0  2018-01-04      1  ...  769.159322  40.770278     0.102630
## 735        0  2018-01-05      1  ...   83.426070  40.666949     0.438938
## 
## [5 rows x 14 columns]

一旦你有了模型,就该弄清楚如何评估和比较它们了。也就是说,你得应对真实情况不可观测这一事实。正如你很快会看到的,诀窍在于意识到,即使你无法测量单个个体的处理效应,你也可以对非常小的群体的处理效应进行估计。因此,如果你希望想出一种从 \(CATE\) 角度评估模型的方法,你将不得不依赖群体层面的指标。

按模型分位数的效应

构建 \(CATE\) 模型的想法源于这样的需求:找出哪些单元对处理更敏感,目的是更高效地分配处理。它源于对个性化的追求。如果这是目标,那么要是你能以某种方式将单元按从更敏感到较不敏感的顺序排列,会非常有用。而且,由于你有预测的 \(CATE\),你可以根据该预测对单元进行排序,并希望这也能按真实的 \(CATE\) 对它们进行排序。遗憾的是,你无法在单元层面评估这种排序。但要是你不需要在单元层面评估呢?要是转而评估由这种排序定义的群体呢

首先,回想一下,如果处理变量是随机分配的,你在这里就不必担心混杂偏差。估计一组单元的效应很容易。你所需要做的就是比较处理组和未处理组的结果。或者,更一般地说,在该组内对 \(Y\) 关于 \(T\) 进行简单回归:

\(y_i = \beta_0 + \beta_1 t_i + e_i \mid X = x\)

从简单线性回归的理论中,你知道:

\(\widehat{\beta_1} = \frac{\sum (t_i - \bar{t}) y_i}{\sum (t_i - \bar{t})^2}\)

其中 \(\bar{t}\) 是该组处理变量的样本平均值,\(\bar{y}\) 是该组结果的样本平均值。


柯里化

\(curry\) 装饰器是一种创建可部分应用函数的方式:

from toolz import curry

@curry
def addN(x, N):
    return x + N

ad5 = addN(N=5)
ad13 = addN(N=13)

print(ad5(5))
## 10
print(ad13(5))
## 18

要编写单变量回归的斜率参数估计代码,你可以使用 \(curry\)。当你需要创建仅接受数据框作为唯一参数的函数时,它非常有用:


from toolz import curry

@curry
def effect(data, y, t):
        return (np.sum((data[t] - data[t].mean())*data[y]) / np.sum((data[t] - data[t].mean())**2))

将这个函数应用于整个测试集,得到的是 \(ATE\)

effect(test, "sales", "discounts")
## 32.16196368039615

但那不是你想要的。相反,你想知道你刚刚拟合的模型是否能对数据进行划分,将对处理更敏感的单元和较不敏感的单元区分开来。为此,你可以根据模型预测的分位数对数据进行分段,并估计每个分位数的效应。如果每个分位数中的估计效应是有序的,你就知道该模型也能很好地对真实的 \(CATE\) 进行排序。


响应曲线形状

在这里,效应被定义为对 \(Y\) 关于 \(T\) 回归的估计斜率。如果你认为这不是一个好的效应指标,你可以使用其他指标。例如,如果你认为响应函数是凹的,你可以将效应定义为对 \(Y\) 关于 \(log(T)\)\(\sqrt{T}\) 回归的斜率。如果你的结果是二元的,使用逻辑回归的参数估计可能更有意义,而不是线性回归的参数估计。这里的关键是要明白,如果 \(T\) 是连续的,你必须将整个处理响应函数概括为一个单一的效应数值。


让我们编写一个函数来按分位数计算效应。它首先使用 \(pd.qcut\)\(q\) 个分位数(默认是 \(10\) 个)对数据进行分段。我用 \(pd.IntervalIndex\) 对其进行包装,以提取 \(pd.qcut\) 返回的每个组的中点。四舍五入只是为了让结果看起来更美观。

然后,我在数据中创建一个包含这些组的列,按这些组对数据进行划分,并估计每个划分中的效应。对于最后一步,我使用 \(pandas\)\(.apply(...)\) 方法。这个方法接受一个以数据框为输入并输出一个数值的函数:\(f(DataFrame) \rightarrow float\)。这就是你之前创建的 \(effect\) 函数发挥作用的地方。

你可以只传递结果和处理参数来调用它。这将返回一个部分应用的 \(effect\) 函数,该函数仅缺少数据框作为参数。这正是 \(.apply(...)\) 所期望的函数类型。

\(test\_pred\) 数据框中使用这个函数的结果是一个列,其中索引是模型预测的分位数,数值是该分位数中的处理效应

def effect_by_quantile(df, pred, y, t, q=10):
    
    # makes quantile partitions
    groups = np.round(pd.IntervalIndex(pd.qcut(df[pred], q=q)).mid, 2) 
    
    return (df
            .assign(**{f"{pred}_quantile": groups})
            .groupby(f"{pred}_quantile")
            # estimate the effect on each quantile
            .apply(effect(y=y, t=t))) 


effect_by_quantile(test_pred, "cate_pred", y="sales", t="discounts")
## <string>:7: FutureWarning: DataFrameGroupBy.apply operated on the grouping columns. This behavior is deprecated, and in a future version of pandas the grouping columns will be excluded from the operation. Either pass `include_groups=False` to exclude the groupings or explicitly select the grouping columns after groupby to silence this warning.
## cate_pred_quantile
## 17.50    20.494153
## 23.93    24.782101
## 26.85    27.494156
## 28.95    28.833993
## 30.81    29.604257
## 32.68    32.216500
## 34.65    35.889459
## 36.75    36.846889
## 39.40    39.125449
## 47.36    44.272549
## dtype: float64

注意,第 \(1\) 个分位数中的估计效应低于第 \(2\) 个分位数中的估计效应,第 \(2\) 个分位数中的估计效应又低于第 \(3\) 个分位数中的估计效应,依此类推。这证明你的 \(CATE\) 预测确实在对效应进行排序:预测值较低的日子对折扣的敏感度也较低,反之亦然。此外,每个分位数中的中点预测(前 \(1\) 列中的索引)与同一分位数的估计效应非常接近。这意味着你的 \(CATE\) 模型不仅能很好地对真实的 \(CATE\) 进行排序,还能正确地对其进行预测。换句话说,你有一个经过校准的 \(CATE\) 模型。

接下来,为了能与其他模型进行比较,你可以应用相同的函数,但传入预测性机器学习模型和随机模型。下面的图表展示了前面定义的 \(3\) 个模型按分位数的效应:

首先,看一下随机模型 \((rand\_m\_pred)\)。它在每个划分中的估计效应大致相同。仅通过看图表,你就能发现它对个性化没什么帮助,因为它无法区分对折扣敏感度高和低的日子。它所有划分中的效应都只是 \(ATE\)。接下来,考虑机器学习预测模型\((ml\_pred)\)。这个模型更有意思一些。看起来,销售预测高和销售预测低的群体都对折扣更敏感。不过,它并没有确切地生成一个排序分数,但你可以用它来进行个性化,比如在销售预测要么非常高要么非常低的时候提供更多折扣,因为这些情况表明处理折扣敏感度高。

最后,看一下你通过回归得到的 \(CATE\) 模型。\(CATE\)预测低的组,其实际 \(CATE\) 确实低于 \(CATE\) 预测高的组。看起来这个模型能很好地区分高效应和低效应。从分位数图中效应的阶梯形状就能看出来。一般来说,阶梯形状越陡,模型在对 \(CATE\) 进行排序方面就越好

在这个例子中,很明显哪个模型在对折扣敏感度进行排序方面更好。但如果有两个不错的模型,比较可能就没那么清晰了。而且,视觉验证虽然不错,但如果要进行模型选择(比如超参数调优或特征选择),就不是很理想了。理想情况下,你应该能用一个单一的数值来概括模型的质量。我们会讲到这一点,但要做到这一点,你首先需要了解累积效应曲线。

累积效应

如果你理解了按分位数的效应图,接下来这个就会很容易。同样,思路是用你的模型来定义组,并估计这些组内的效应。不过,不是按组估计效应,而是将一个组叠加在另一个组上进行累积。

首先,你需要用一个分数对数据进行排序 —— 通常是 \(CATE\) 模型,但实际上可以是任何东西。然后,根据这个排序,估计前 \(1\%\) 数据的效应。接下来,加入接下来的 \(1\%\),计算前 \(2\%\) 数据的效应,然后是前 \(3\%\),依此类推。结果会是一条按累积样本计算的效应曲线。下面是实现这一操作的简单代码:

def cumulative_effect_curve(dataset, prediction, y, t,
                            ascending=False, steps=100):
    size = len(dataset)
    ordered_df = (dataset
                  .sort_values(prediction, ascending=ascending)
                  .reset_index(drop=True))
    
    steps = np.linspace(size/steps, size, steps).round(0)
    
    return np.array([effect(ordered_df.query(f"index<={row}"), t=t, y=y) for row in steps])

cumulative_effect_curve(test_pred, "cate_pred", "sales", "discounts")
## array([49.65116279, 49.37712454, 46.20360341, 45.9770684 , 45.31711812,
##        45.23107566, 44.91075165, 44.74166167, 44.56309576, 44.27254859,
##        43.95937553, 43.66601716, 43.30710192, 42.71556602, 42.57717619,
##        42.14190207, 41.9120469 , 41.73012555, 41.59626585, 41.38073714,
##        41.11778111, 41.04864418, 40.95933655, 40.6949357 , 40.5585599 ,
##        40.32182723, 40.19431251, 40.0242349 , 39.92014581, 39.58537488,
##        39.44592589, 39.53382731, 39.41985355, 39.29287716, 39.27763551,
##        39.14056914, 39.04402217, 38.98421346, 38.9022713 , 38.71601845,
##        38.72423427, 38.53676256, 38.56709164, 38.36661748, 38.27685106,
##        38.06110588, 37.99479087, 37.89750017, 37.72699021, 37.55580363,
##        37.35825912, 37.28006651, 37.21139791, 37.06210857, 37.00753334,
##        36.87259686, 36.74299254, 36.52014459, 36.4873419 , 36.38463167,
##        36.21632077, 36.10908232, 35.90865968, 35.73915746, 35.56123101,
##        35.42248759, 35.37393053, 35.25994908, 35.24501879, 35.14156543,
##        35.09274947, 35.001377  , 34.96329233, 34.87674713, 34.79536682,
##        34.65873593, 34.59244929, 34.48071628, 34.38210781, 34.30316242,
##        34.08974584, 34.05306999, 33.92906318, 33.89927971, 33.73934419,
##        33.62832513, 33.47374484, 33.39942527, 33.25723224, 33.18300856,
##        33.08850259, 33.04965845, 33.00366002, 32.87840008, 32.81524103,
##        32.74069128, 32.6477962 , 32.46981935, 32.33428884, 32.16196368])

如果用于对数据排序的分数也能很好地对真实的条件平均处理效应(CATE)进行排序,那么得到的曲线会从很高的位置开始,然后逐渐下降到平均处理效应(ATE)。相反,一个糟糕的模型要么会快速收敛到 ATE,要么会一直围绕 ATE 波动。为了更好地理解这一点,下面是你创建的三个模型的累积效应曲线:

首先,注意回归条件平均处理效应(CATE)模型是如何从很高的位置开始,然后逐渐收敛到平均处理效应(ATE)的。例如,如果你用这个模型对数据进行排序,前 20% 数据的 ATE 约为 42,前 50% 数据的 ATE 大概是 37,而前 100% 数据的 ATE 就是处理的全局效应(ATE)。相反,一个只输出随机数的模型会围绕 ATE 波动,而一个对效应进行反向排序的模型会从低于 ATE 的位置开始。


排序不对称性

需要指出的是,这种排序是不对称的。也就是说,取一个分数并反转它,并不会简单地使曲线围绕 ATE 线翻转。


累积效应曲线比按分位数的效应曲线更好一些,因为它允许将其概括为一个单一的数值。例如,你可以计算曲线与 ATE 之间的面积,并利用这个面积来比较不同的模型。面积越大,模型越好。但这仍然存在一个缺点。如果这样做,曲线的起始部分会有最大的面积。但恰恰在这部分,由于样本量较小,不确定性是最大的。幸运的是,有一个非常简单的解决办法:累积增益曲线。

累积增益

如果采用与累积效应曲线完全相同的逻辑,但将每个点乘以累积样本数 \(N_{cum}/N\),你就得到了累积增益曲线。现在,即使曲线的起始部分(对于一个好的模型而言)会有最高的效应,它也会因相对规模小而被按比例缩小。

看一下代码,变化之处在于我现在在每次迭代时都将效应乘以 \((\text{row}/\text{size})\)。此外,我可以选择用 ATE 来对这条曲线进行归一化,这就是为什么我在每次迭代时也从效应中减去一个归一化因子:

def cumulative_gain_curve(df, prediction, y, t,
                          ascending=False, normalize=False, steps=100):
    
    effect_fn = effect(t=t, y=y)
    normalizer = effect_fn(df) if normalize else 0
    
    size = len(df)
    ordered_df = (df
                  .sort_values(prediction, ascending=ascending)
                  .reset_index(drop=True))
    
    steps = np.linspace(size/steps, size, steps).round(0)
    effects = [(effect_fn(ordered_df.query(f"index<={row}"))
                -normalizer)*(row/size) 
               for row in steps]

    return np.array([0] + effects)


cumulative_gain_curve(test_pred, "cate_pred", "sales", "discounts")
## array([ 0.        ,  0.50387597,  0.982917  ,  1.38863289,  1.83046877,
##         2.26408709,  2.71880783,  3.13779256,  3.58003148,  4.01798404,
##         4.42379877,  4.83861979,  5.23242282,  5.62890903,  5.98551452,
##         6.38159081,  6.74402024,  7.13257317,  7.5088165 ,  7.9068624 ,
##         8.26968674,  8.63441305,  9.03646968,  9.41649072,  9.76869065,
##        10.13172457, 10.48178647, 10.85654343, 11.20116176, 11.57715392,
##        11.88179285, 12.22484979, 12.65329367, 13.00147394, 13.35835129,
##        13.75177168, 14.08571614, 14.44720258, 14.8206963 , 15.16915263,
##        15.48942971, 15.87058782, 16.18483861, 16.5889676 , 16.87711862,
##        17.226077  , 17.51523959, 17.85547549, 18.1943502 , 18.48062947,
##        18.77790182, 19.05825319, 19.38214231, 19.7240743 , 20.00659491,
##        20.35269886, 20.65268403, 20.93862963, 21.18225404, 21.53351325,
##        21.82793867, 22.09450014, 22.38142964, 22.62161465, 22.87752468,
##        23.11063608, 23.3799479 , 23.70688475, 23.97456335, 24.32208947,
##        24.59360922, 24.91557817, 25.20590966, 25.51965521, 25.81042645,
##        26.10331578, 26.33901594, 26.6396965 , 26.89011363, 27.16213357,
##        27.44788562, 27.60976683, 27.92564404, 28.15503057, 28.47433643,
##        28.68239331, 28.91615935, 29.12294194, 29.39723028, 29.59660012,
##        29.8672981 , 30.10511301, 30.40516978, 30.69778369, 30.9021028 ,
##        31.17575982, 31.43719771, 31.66657828, 31.82346463, 32.00615008,
##        32.16196368])

另见

如果你不想费心实现所有这些函数,我和一些同事一直在开发一个 Python 库来为你处理这些。你可以直接从 fklearncausal 模块导入所有曲线及其曲线下面积(AUC):



from fklearn.causal.validation.auc import *
from fklearn.causal.validation.curves import *

三个模型的累积增益和归一化累积增益都显示在下图中。在这里,就对条件平均处理效应(CATE)的排序而言,更好的模型是在曲线与代表平均处理效应(ATE)的虚线之间具有最大面积的那个:

要将模型性能概括为一个单一的数值,你可以对归一化累积增益曲线的数值进行求和。在对条件平均处理效应(CATE)排序方面,数值最大的模型就是最好的。以下是到目前为止你评估的三个模型的曲线下面积(AUC)。注意,机器学习模型的面积为负,因为它对 CATE 进行了反向排序:

rand_m_pred 的 AUC:6.0745233598544495 ml_pred 的 AUC:-45.44063124684 cate_pred 的 AUC:181.74573239200615

再次说明,能够将模型的性能浓缩为一个单一的数值是很棒的,因为这允许进行自动化的模型选择。不过,尽管我很喜欢最后这个曲线,但在使用它们时,你需要注意一些注意事项。

首先,在你看到的所有曲线中,必须记住曲线中的每个点都是一个估计值,而不是真实值。它是对某个特定(有时规模非常小的)组进行回归斜率的估计。而且,由于它是一个回归估计,它依赖于 T 和 Y 之间的关系被正确设定。即使进行了随机化,如果处理和结果之间的关系(比如说)是一个对数函数,而把效应当作线性关系来估计,就会得到错误的结果。如果你知道处理响应函数的形状,你可以将效应函数调整为 \(y \sim log(t)\) 的斜率,而不是 \(y \sim t\) 的斜率。但要做到这一点,你需要知道正确的形状。

其次,这些曲线并不真正关心你是否正确得到了条件平均处理效应(CATE)。它们只关心排序是否正确。例如,如果你对任何一个模型的预测值减去 - 1000,它们的累积增益曲线将保持不变。因此,即使你对 CATE 有一个有偏的估计量,这种偏差也不会在这些曲线中体现出来。现在,如果你的全部关注点只是对处理进行优先排序,这可能不是问题。在这种情况下,排序就足够了。但如果你确实关心精确估计 CATE,这些曲线可能会产生误导。如果你有数据科学背景,你可以将累积增益曲线与 ROC 曲线进行类比。同样,一个 ROC - AUC 良好的模型不一定是经过校准的。

第三,或许也是最重要的一点,前面所有的方法都需要无混杂的数据。如果你有任何偏差,你所估计的效应 —— 无论是子组中的效应还是平均处理效应(ATE)—— 都会是错误的。如果处理不是随机分配的,理论上,只要你事先使用诸如正交化或逆概率加权(IPW)之类的方法对数据进行去偏,仍然可以使用这些评估技术。不过,我对此有点怀疑。相反,我强烈建议你投入一些实验数据,哪怕只是少量的,仅用于评估目的。这样,你就可以专注于效应异质性,而不必担心混杂因素的干扰。


另见

这里呈现的所有曲线,都是在处理为二元时,对传统上用于提升建模的曲线进行的一种泛化尝试。如果你想回顾相关文献,我推荐皮埃尔・古铁雷斯(Pierre Gutierrez)和让 - 伊夫・热拉迪(Jean - Yves Gérardy)的论文《因果推断与提升建模:文献综述》(“Causal Inference and Uplift Modeling a Review of the Literature”),以及迪维亚特・马哈詹(Divyat Mahajan)等人的论文《异质因果效应估计的模型选择实证分析》(“Empirical Analysis of Model Selection for Heterogeneous Causal Effect Estimation”)。


因果模型的评估是一个仍在发展的研究领域。因此,它仍然存在许多盲区。例如,到目前为止呈现的曲线只能告诉你一个模型在对条件平均处理效应(CATE)进行排序方面有多好。我还没有找到一个很好的方法来检验你的模型是否正确预测了 CATE。我喜欢做的一件事是将按分位数的效应图与累积增益曲线结合使用,因为前者能让我对模型的校准情况有一些了解,而后者能让我了解模型对 CATE 的排序效果如何。至于归一化累积增益,它只是一种放大手段,让可视化更容易。

但我得承认,这并不理想。如果你正在寻找像 \(R^2\) 或均方误差(MSE)这样的汇总指标 —— 这两者在预测模型中都很常用 —— 我很遗憾地说,在因果建模领域,我还没有找到任何能与它们相提并论的指标。不过,我确实找到了一个方法 —— 目标变换。

目标变换

事实证明,即使你无法观测到真实的处理效应 \(\tau(x_i)\),你也可以创建一个目标变量,使其在期望上近似于真实处理效应:

\(Y_i^* = \frac{(Y_i - \hat{\mu}_y(X_i))(T_i - \hat{\mu}_t(X_i))}{(T_i - \hat{\mu}_t(X_i))^2} = \frac{Y_i - \hat{\mu}_y(X_i)}{T_i - \hat{\mu}_t(X_i)}\)

其中 \(\mu_y\) 是结果的模型,\(\mu_t\) 是处理的模型。这个目标很有趣,因为 \(E[Y_i^*] = \tau_i\)。注意到它看起来很像回归系数的公式,分子是 Y 和 T 之间的协方差,分母是 T 的方差。不过,这里不是用期望来定义这些,而是在个体层面进行计算。

由于这个目标近似于真实的处理效应,你可以用它来计算偏差指标,比如均方误差(MSE)。如果你的条件平均处理效应(CATE)模型在预测个体层面的效应 \(\tau_i\) 方面表现良好,那么相对于这个目标,其预测的均方误差应该很小。

不过,这里有个问题。当接近处理平均值时,这个目标会非常嘈杂,因为分母会趋近于零。为了解决这个问题,你可以应用权重,对 \(T_i - \hat{\mu}_t(X_i)\) 较小的点赋予较低的重要性。例如,你可以用 \((T_i - \hat{\mu}_t(X_i))^2\) 对个体进行加权。


预示 R 损失

使用这些权重是有充分理论依据的。当我们在第 7 章讨论非参数双重 / 去偏机器学习时,你会了解到更多相关内容。目前,你先相信我说的话就好。


要编写这个目标(变量),你只需将结果模型和处理模型的残差相除即可:


X = ["C(month)", "C(weekday)", "is_holiday", "competitors_price"]

y_res = smf.ols(f"sales ~ {'+'.join(X)}", data=test).fit().resid
t_res = smf.ols(f"discounts ~ {'+'.join(X)}", data=test).fit().resid

tau_hat = y_res/t_res

接下来,你可以用它来计算所有模型的均方误差。请注意,我同样按照先前讨论的方式使用了权重:

from sklearn.metrics import mean_squared_error

for m in ["rand_m_pred", "ml_pred", "cate_pred"]:
    wmse = mean_squared_error(tau_hat, test_pred[m],sample_weight=t_res**2)
    print(f"MSE for {m}:", wmse)
## MSE for rand_m_pred: 1115.803515760459
## MSE for ml_pred: 576256.7425385384
## MSE for cate_pred: 42.90447405550294

根据这个加权均方误差(MSE),再次证明用于估计条件平均处理效应(CATE)的回归模型优于另外两个模型。此外,这里还有一个有趣的现象:机器学习模型的表现比随机模型更差。这并不令人意外,因为机器学习模型试图预测的是结果变量 \(Y\) 而非个体处理效应 \(\tau_i\)


另见

就像我之前说的,评估因果模型的文献仍处于起步阶段。这是一个相当令人兴奋的问题,新方法一直在被提出。例如,在论文《基于因果推断的消费贷款智能额度管理》中,蚂蚁金融服务集团的科学家们提出将个体划分为具有相似协变量的组(他们使用了超过 6000 个组!),假设结果是处理效应加上一些高斯随机噪声 \(\hat{y}_i = \hat{\tau}(x_i) + e_i\),计算每个组内的结果均方误差 \(N^{-1} \sum (y_i - \hat{y}_i)\),并使用每个组的样本量对结果进行平均。


只有当效应与结果相关时,对 Y 的预测在排序或预测 \(\tau_i\) 方面才会表现良好。一般情况下并非如此,但在某些情况下可能会出现这种情况。其中一些情况在商业中相当常见,所以值得研究一下。

预测模型何时有利于效应排序

就像我之前说的,对于一个预测 Y 的模型来说,要想也能很好地对条件平均处理效应(CATE) \(\tau(x_i)\) 进行排序,必须满足 Y 和 CATE \(\tau(x_i)\) 也是相关的。例如,在寻找餐厅顾客对折扣更有响应的日子的情境中,如果销售额高的日子与人们对折扣更敏感的日子相吻合,那么一个预测 Y 的模型也会很擅长对 T 对 Y 的效应进行排序。更一般地说,当处理响应函数是非线性的时候,这种情况可能会发生。

边际收益递减

当处理响应函数是凹函数时,额外的处理单位产生的效果会越来越小。这在商业中是非常常见的现象,因为事物往往存在饱和点。例如,即使你把折扣设为 100%,销售额也只能增长到一定程度,因为存在限制产量的因素。或者,你的营销预算的效果最终会趋于平稳,因为你能宣传触达的客户数量是有限的。

边际递减的处理响应函数大致是这样的:

很容易理解,在这种情况下,一个擅长预测结果 \(Y\) 的模型也能擅长对条件平均处理效应(CATE)进行排序:结果越高,效应越低。因此,如果你采用预测 Y 的模型,并根据这些预测的倒数对个体进行排序,你很可能会得到不错的 CATE 排序。

二元结果

另一种常见的情况是,当结果为二元时,预测 Y 的模型也能很好地对条件平均处理效应(CATE)进行排序。在这种情况下,\(E[Y|T]\) 呈 S 形,在 0 和 1 处趋于平缓:

在大多数商业应用中,数据会集中在这个 S 形函数的一端或另一端。例如,在银行业,只有一小部分客户会拖欠贷款,这意味着你大多处于这条曲线的左侧,那里看起来是指数型的。因此,如果你有一个预测客户拖欠情况的模型,预测值较高的客户很可能对处理也更敏感。直观地说,这些是处于不拖欠和拖欠之间临界点附近的客户。对他们来说,处理上的一个小变化可能会带来完全不同的结果。

相反,假设你在一家在线购物企业工作,大多数进入你网站的客户都会购买商品(转化)。在这种情况下,你更处于 S 形曲线的右侧。因此,如果你有一个预测转化的模型,这个模型很有可能也能对诸如折扣之类的因素的效应进行排序。转化的可能性越高,效应量就越低。这是因为在右侧,S 曲线看起来有点像你之前看到的边际收益递减的情况。

一般来说,当结果是二元的时候,你越接近中间 —— 也就是 \(E[Y|X] = 50\%\)—— 效应就会越高。


实际案例

疫苗优先级排序

你已经了解到二元结果如何在处理响应函数中产生非线性,这使得你可以利用结果的预测值来分配处理。在 2019 冠状病毒病(COVID - 19)大流行期间,就出现了对这一原理的一个非常有趣的应用。2021 年,全球首批获批的 COVID - 19 疫苗被交付给普通公众。当时,一个关键问题是,谁应该首先接种疫苗。不出所料,这是一个异质处理效应问题。政策制定者希望首先为那些受益最大的人接种疫苗。在这种情况下,处理效应是避免死亡或住院。那么,哪些人在接种疫苗后,死亡或住院的情况减少得最多呢?在大多数国家,是老年人和有既往健康问题(合并症)的人。现在,这些人是感染 COVID - 19 后更有可能死亡的群体。而且,幸好 COVID - 19 的死亡率远低于 50%,这使你处于逻辑函数的左侧。在这个区域,根据我们之前针对违约率的分析思路,优先为那些感染 COVID - 19 后死亡基线概率高的人接种疫苗是有道理的,而这些人正是前面提到的群体。这是巧合吗?也许吧。要记住,我不是健康专家,所以我在这里可能大错特错。但这个逻辑对我来说非常有道理。


当处理响应函数是非线性的,比如在二元结果的情况下,或者在结果呈边际递减的情况下,预测模型可能会对条件平均处理效应(CATE)产生良好的排序。不过,这并不意味着它会是最佳模型,也不意味着它不能被直接旨在预测 CATE 的模型超越。此外,尽管这样的模型可能会对处理效应进行排序,但它并不能预测该处理效应。只有当你只关心根据个体对处理的敏感性对其进行排序时,这才是可行的。但如果你的决策依赖于对 CATE 的正确估计,那么就需要进行额外的按组效应估计。


另见

有时,结果预测可能会优于 CATE 预测,因为 CATE 往往非常嘈杂。费尔南德斯 - 洛里亚(Fernández - Loria)和普罗沃斯特(Provost)在他们的论文《因果分类:处理效应估计与结果预测》(“Causal Classification: Treatment Effect Estimation vs. Outcome Prediction”)中进一步讨论了这一点。


说到这个,我认为有必要详细说明一下如何将条件平均处理效应(CATE)用于决策制定。你可能已经很清楚该怎么做了,但也许我能给出一些你没考虑过的建议。

用于决策制定的条件平均处理效应(CATE)

当处理是二元的时候,决策制定过程相当直接。本质上,你关心的是谁对处理有积极的响应。如果你的处理供应是无限的,那么你所要做的就是对所有 CATE 为正的人进行处理。如果您没有预测 CATE 的模型,但有一个对其进行排序的模型 —— 就像上一节中讨论的预测模型那样 —— 您可以使用按模型分位数的效应图。只需根据模型的分位数对数据进行划分,估计每个分位数中的处理效应,然后对效应仍为正的所有个体进行处理。

如果你的处理供应不是无限的,那么你需要增加第二条规则。你不仅要对那些有正效应的个体进行处理,还要对那些 CATE 最高的个体进行处理。例如,如果你只有 1000 个处理单位,根据某个 CATE 排序模型,你可能想要处理排名前 1000 的个体,前提是他们都有正效应。

如果处理是连续的或有序的,情况会变得稍微复杂一些。现在你不仅要决定对谁进行处理,还要决定处理的程度。这非常具有业务特异性。每个问题都会有自己的处理响应函数需要优化。这意味着我无法给你非常详细的操作指南,但我可以带你看一个典型的例子。

再考虑一下决定连锁餐厅每天给予多少折扣的问题。由于决定给予多少折扣只是决定收取什么价格的另一种方式(\(Price = Price_{base}*(1 - Discount)\)),我们把这个问题重新表述为一个价格优化问题。在所有商业问题中,都存在成本(即使不是货币成本)和收入函数。假设餐厅的收入由以下方程给出:

\(Demand_i = 50 - \tau(X_i)Price_i\)

\(Revenue_i = Demand_i Price_i\)

在第i天的收入就是价格乘以餐厅提供的餐食数量(需求)。然而,人们在某一天愿意购买的餐食数量与当天收取的价格成反比。也就是说,它有一个组成部分 \(-\tau(X_i)Price_i\),其中 \(\tau(X_i)\) 是顾客在那天对价格上涨的敏感程度(注意,这取决于特定日期的特征 \(X_i\))。换句话说,这种敏感性就是价格对需求的条件平均处理效应。

如果你针对不同的 \(\tau\) 值绘制需求曲线,你会发现 \(\tau\) 只不过是需求曲线的斜率。如果你将需求曲线乘以收入曲线,你会得到一个二次曲线形状。在这条曲线中,顾客对价格不太敏感(\(\tau = 1\))的那天,会在更高的价格值处达到峰值:

接下来,假设你花 3 美元来制作一份餐食。这意味着成本就是你生产的数量 q 乘以 3:

\(Costs(q_i) = 3q_i\)

要记住,成本方程并不直接依赖于处理效应,但如果你回想一下,生产的数量就是客户订单的数量 —— 也就是需求 —— 那么,价格越高,成本就越低,因为客户对餐食的需求会减少。

最后,一旦你有了收入和成本,就可以将它们结合起来,得到利润关于价格的函数:

\(Profit_i = Demand_i * Price_i - Cost(Demand_i)\)

如果你针对不同的 \(\tau_i\) 值绘制价格对应的利润,你会发现每个 \(\tau_i\) 都会产生不同的最优价格。\(\tau_i\) 越低,顾客对价格上涨的敏感度就越低,这使得餐厅可以提高价格以获取更多利润:

经济学家会很快意识到,这就是著名的企业问题。将边际成本等于边际收益并解出价格,就能得到利润最大化价格的数值解:

\((P(50 - \tau(X)P))' = (3(50 - \tau(X)P))'\)

\(P^* = \frac{3\tau(X) + 50}{2\tau(X)}\)

注意,在这种情况下,唯一的未知量是价格对需求的影响 \(\tau(X)\)。所以,如果你能用预测条件平均处理效应(CATE)的模型来估计它,就可以将该 CATE 预测转化为最优价格。

同样,这在很大程度上取决于收入和成本曲线的形式,而这又在很大程度上取决于你的业务。但总的来说,几乎任何你想要优化的处理都有好的一面 —— 在这个例子中是收入 —— 和坏的一面 —— 在这个例子中是成本。要使用条件平均处理效应(CATE)来决定连续处理的水平,你必须了解它如何影响这两个方面。


边界解

在一些罕见的情况下,能优化业务的处理水平要么是零,要么是允许的最高水平。例如,假设当地政府对你销售的产品设定了价格上限,而这个上限低于能使你利润最大化的价格。在这种情况下,最优价格就只是政府允许的最高价格。不过,这种情况很罕见。在大多数情况下,边界解的出现是因为没有考虑到隐性成本。例如,如果你试图优化交叉销售邮件,你可能会认为发送一封邮件的成本微不足道,所以应该直接发给所有人。但我会反驳说,你没有考虑到客户注意力方面的成本:如果你向客户发送垃圾邮件,最终他们会厌烦你并取消订阅你的邮件,这会让你失去未来通过邮件渠道获得的销售。这些隐性成本更难考虑进去,但这并不意味着它们不存在。事实上,找到这些成本的良好替代指标往往是一项极有价值的数据科学任务。


核心思想

本章介绍了处理异质性的概念。核心见解是,每个个体 i 可能有不同的处理效应 \(\tau_i\)。如果知道这种效应,就可以用它来更好地在个体间分配处理。遗憾的是,由于因果推断的基本问题,这种效应是无法观测的。不过,如果假设它依赖于个体的可观测特征 \(\tau(x_i)\),那么还是能取得一些进展的;也就是说,可以从估计平均处理效应转向估计条件平均处理效应(CATE):

\(CATE = \frac{\partial}{\partial t} E[Y(t)|X]\)

所以,即使在个体层面无法观测到处理效应,你仍然可以估计组效应。一种简单的方法是使用线性回归,纳入处理与协变量之间的交互项:

\(y_i = \beta_0 + \tau_0 T_i + \tau X_i T_i + \beta X_i + e_i\)

对这个模型进行估计,会得到以下条件平均处理效应(CATE)的估计值:

\(\widehat{CATE} = \tau_0 + \tau X_i\)

接下来,你了解了一些关于如何将交叉验证与条件平均处理效应(CATE)评估技术相结合,以评估你的 CATE 估计值的思路。由于 CATE 不是针对单个个体定义的,你必须依赖于特定组的指标,比如按分位数的效应曲线或累积增益曲线。如果这还不够,你也可以定义一个近似个体层面处理效应的目标,并使用它来计算偏差指标,比如均方误差(MSE)。

最后,值得强调的是,本章所讨论的一切都取决于这样一个事实:条件平均处理效应(CATE)作为一个因果量,可以从条件期望中识别出来,而条件期望是一个可从数据中恢复的统计量:

\(\frac{\partial}{\partial t} E[Y(t)|X] = \frac{\partial}{\partial t} E[Y|X, T = t]\)

如果没有这一点,将 CATE 作为可估计的组效应的想法就不再成立,这就是为什么随机数据对于 CATE 估计问题如此重要,即使只是为了评估你的处理异质性模型。