双重稳健学习

这是什么?

双重稳健学习(Doubly Robust Learning)类似于双重机器学习(Double Machine Learning),是一种用于估计(异质)处理效应的方法。它适用于处理是分类变量,且所有潜在的混杂因素/控制变量(即同时对收集数据中的处理决策和观测结果有直接影响的因素)都被观测到,但这些变量数量过多(高维)以至于经典统计方法无法适用,或者它们对处理和结果的影响无法通过参数函数(非参数)令人满意地建模的情况。后者的两个问题都可以通过机器学习技术解决(例如参见[Chernozhukov2016], [Foster2019])。该方法可以追溯到 [Robins1994], [Bang] 的早期工作(更多细节请参见 [Tsiatis]),这些工作主要将该方法应用于平均处理效应的估计。在此库中,我们实现了对双重稳健方法的最新修改,这些修改允许估计异质处理效应(例如参见 [Foster2019])。该方法最近也在策略学习领域被大量使用(例如参见 [Dudik2014], [Athey2017])。

它将问题简化为首先估计两个预测任务

  1. 根据处理和控制变量预测结果,

  2. 根据控制变量预测处理;

因此,与双重机器学习不同,第一个模型是根据处理和控制变量预测结果,而不是仅根据控制变量。然后,该方法在最后阶段的估计中结合这两个预测模型,从而创建异质处理效应的模型。该方法允许对这两个预测任务使用任意机器学习算法,同时保持与最终模型相关的许多有利的统计属性(例如,小的均方误差、渐近正态性、置信区间的构建)。当这两个预测任务中的任意一个达到小的均方误差时,后者有利的统计属性就成立(因此得名双重稳健)。

我们的软件包为最终模型估计提供了几种变体。其中许多变体还提供有效推断(置信区间构建),用于衡量学习模型的 不确定性。

有哪些相关的估计器类别?

本节描述了以下类中实现的方法:DRLearner, LinearDRLearner, SparseLinearDRLearner, ForestDRLearner。点击这些链接可以查看每个类的详细模块文档和输入参数。

何时使用它?

假设您有观测数据(或来自 A/B 测试的实验数据),其中从一组有限的处理中选择了某个处理/干预/行动 \(T\),并观测到某个结果 \(Y\),并且数据集中记录了所有可能影响 \(T\) 选择,同时又可能直接影响结果 \(Y\) 的变量 \(W\)(也称为控制变量或混杂因素)。

如果您的目标是理解每种处理对结果的影响,以及这种影响如何作为处理样本的一组可观测特征 \(X\) 的函数而变化,那么可以使用此方法。例如调用

from econml.dr import LinearDRLearner
est = LinearDRLearner()
est.fit(y, T, X=X, W=W)
est.effect(X, T0=t0, T1=t1)

通过这种方式,可以通过简单地检查对于哪些 \(X\) 效应是正的来学习最优处理策略。

正式方法论概述

该模型的假设最好用潜在结果的语言来解释。如果我们用 \(Y^{(t)}\) 表示如果对样本应用处理 \(T=t\) 我们可能观测到的潜在结果,则该方法假设

\[\begin{split}Y^{(t)} =~& g_t(X, W) + \epsilon_t ~~~&~~~ \E[\epsilon | X, W] = 0 \\ \Pr[T = t | X, W] =~& p_t(X, W) & \\ \{Y^{(t)}\}_{t=1}^{n_t} \perp T | X, W\end{split}\]

它不对 \(g_t\)\(p_t\) 做进一步的结构性假设,并使用任意非参数机器学习方法对其进行非参数估计。我们的目标是估计与每个可能的处理 \(t \in \{1, \ldots, n_t\}\) 相关的 CATE,与某个基线处理 \(t=0\) 相比,即

\[\theta_t(X) = \E[Y^{(t)} - Y^{(0)} | X] = \E[g_t(X, W) - g_0(X, W) | X]\]

估计 \(\theta_t(X)\) 的一种方法是直接方法 (DM) 方法,我们简单地通过回归 \(Y\)\(T, X, W\) 来估计一个回归模型,以学习 \(g_T(X, W) = \E[Y | T, X, W]\) 的模型,然后通过回归来评估 \(\theta_t(X)\)

\[Y_{i, t}^{DM} = g_t(X_i, W_i) - g_0(X_i, W_i)\]

基于 \(X\)。这种方法的主要问题在于它严重依赖于通过回归拟合模型隐含完成的基于模型的推断。本质上,当我们在特征为 \(X, W\) 的样本上评估 \(g_t(X, W)\),而该样本接受了其他处理 \(T=t'\) 时,我们是在从接受了处理 \(T=t\) 的具有相似 \(X, W\) 的其他样本进行推断。然而,“相似性”的定义非常依赖于模型,在某些情况下,我们甚至可能从非常远的样本点进行推断(例如,如果我们拟合线性回归模型)。

另一种不受上述问题困扰的方法是逆倾向得分 (IPS) 方法。这种方法源于这样一个事实:由于无混杂假设,我们可以通过将每个样本按其观测到的处理分配概率的倒数进行重新加权(即,对具有“令人惊讶”处理分配的样本给予更高的权重),从而创建每个潜在结果的无偏估计。更具体地说,如果我们令

\[Y_{i, t}^{IPS} = \frac{Y_i 1\{T_i=t\}}{\Pr[T_i=t | X_i, W_i]} = \frac{Y_i 1\{T_i=t\}}{p_t(X_i, W_i)}\]

那么成立

\[\begin{split}\E[Y_{i, t}^{IPS} | X, W] =~& \E\left[\frac{Y_i 1\{T_i=t\}}{p_t(X_i, W_i)} | X_i, W_i\right] = \E\left[\frac{Y_i^{(t)} 1\{T_i=t\}}{p_t(X_i, W_i)} | X_i, W_i\right]\\ =~& \E\left[\frac{Y_i^{(t)} \E[1\{T_i=t\} | X_i, W_i]}{p_t(X_i, W_i)} | X_i, W_i\right] = \E\left[Y_i^{(t)} | X_i, W_i\right]\end{split}\]

因此,我们可以通过对 \(Y_{i, t}^{IPS} - Y_{i, 0}^{IPS}\) 基于 \(X\) 进行回归来估计 \(\theta_t(X)\)。这种方法有两个缺点:1) 首先,即使我们知道处理的概率 \(p_t(X, W)\),这种方法的方差很高,因为我们将观测值除以一个相对较小的数(尤其是在 \(X, W\) 的某些区域,某些处理的可能性很小),2) 其次,在观测数据中,我们通常不知道处理的概率,因此也需要估计处理概率的模型。这对应于一个多类分类任务,当 \(X, W\) 是高维的或者我们使用像随机森林这样的非线性模型时,估计速率可能会很慢。这种方法将继承这些速率。此外,如果我们使用机器学习来拟合这些倾向模型,则很难描述我们估计量的极限分布,从而提供有效的置信区间。

双重稳健方法通过结合上述两种方法来避免其缺点。特别是,它首先拟合一个直接回归模型,然后通过对该模型的残差应用逆倾向方法来消除其偏差,即它构建了潜在结果的以下估计

\[Y_{i, t}^{DR} = g_t(X_i, W_i) + \frac{Y_i -g_t(X_i, W_i)}{p_t(X_i, W_i)} \cdot 1\{T_i=t\}\]

然后,我们可以通过对 \(Y_{i, t}^{DR} - Y_{i, 0}^{DR}\) 基于 \(X_i\) 进行回归来学习 \(\theta_t(X)\)

这就得到了整体算法:首先学习一个回归模型 \(\hat{g}_t(X, W)\),通过对 \(Y\) 基于 \(T, X, W\) 进行回归;以及一个倾向模型 \(\hat{p}_t(X, W)\),通过对 \(T\) 基于 \(X, W\) 进行分类预测。然后构建如上所述的双重稳健随机变量,并将其基于 \(X\) 进行回归。

双重稳健方法的主要优点是,最终估计量 \(\theta_t(X)\) 的均方误差仅受回归估计量 \(\hat{g}_t(X, W)\) 和倾向估计量 \(\hat{p}_t(X, W)\) 的均方误差乘积影响。因此,只要其中一个准确,最终模型就是正确的。例如,只要两者都没有以慢于 \(n^{-1/4}\) 的速率收敛,那么最终模型就能达到参数速率 \(n^{-1/2}\)。此外,在对最终阶段使用的估计算法做一些进一步假设的情况下,最终估计量是渐近正态的,并且可以构建有效的置信区间。为了使该定理成立,辅助估计量需要以交叉拟合的方式进行拟合(参见_OrthoLearner)。后者的稳健性属性源于这样一个事实:对应于最终最小二乘估计(即平方损失的梯度)的矩方程,对于辅助参数 \(q, f\) 满足 Neyman 正交条件。关于 Neyman 正交性如何导致稳健性的更详细阐述,我们建议读者参考[Chernozhukov2016], [Mackey2017], [Nie2017], [Chernozhukov2017], [Chernozhukov2018], [Foster2019]。实际上,双重稳健估计量满足比 Neyman 正交性稍强的属性,这就是为什么它具有更强的稳健性保证,即只有第一阶段模型的两个均方误差的乘积会影响最终估计量的误差和分布性质。

双重稳健方法相对于 DML 方法的另一个优点是,即使我们在最终回归中最小化损失的函数空间不包含真实的 CATE 函数,最终回归结果仍然有意义。在这种情况下,该方法将估计 CATE 函数在我们最终回归中优化的模型空间上的投影。例如,这允许对 CATE 函数的最佳线性投影进行推断,或者对潜在可能产生异质性的一组特征上的最佳 CATE 函数进行推断。例如,可以使用 DR 方法配合非参数最终模型(如 Honest Forest),并对相对于单个特征的边际处理效应异质性进行推断,而无需对处理效应异质性的形式做出任何进一步的假设。

DR 方法相对于 DML 方法的缺点是,它通常具有更高的方差,尤其是在控制变量空间 \(X, W\) 中,某些处理被分配的概率很小(在文献中通常称为“小重叠”)的区域。在这种情况下,DML 方法可能具有更好的外推能力,因为它只需要“平均”具有良好的重叠就能达到较低的均方误差。

类层次结构

在此库中,我们实现了双重稳健方法的几种变体,具体取决于最终阶段选择的估计算法类型。在所有这些变体中,用户都可以为第一阶段模型选择任何回归/分类方法。实现的 CATE 估计器的层次结构如下。

Inheritance diagram of econml.dr.DRLearner, econml.dr.LinearDRLearner, econml.dr.SparseLinearDRLearner, econml.dr.ForestDRLearner

下面我们简要描述这些类

  • DRLearner.DRLearner 对每个结果 \(i\) 和处理 \(t\) 的效应模型不作任何假设。任何 scikit-learn 回归器都可用于最终阶段的估计。类似地,任何 scikit-learn 回归器都可用于回归模型,任何 scikit-learn 分类器都可用于倾向模型

    from econml.dr import DRLearner
    from sklearn.ensemble import GradientBoostingRegressor, GradientBoostingClassifier
    est = DRLearner(model_regression=GradientBoostingRegressor(),
                    model_propensity=GradientBoostingClassifier(),
                    model_final=GradientBoostingRegressor())
    est.fit(y, T, X=X, W=W)
    point = est.effect(X, T0=T0, T1=T1)
    

    模型的例子包括随机森林 (RandomForestRegressor)、梯度提升森林 (GradientBoostingRegressor) 和支持向量机 (SVC)。此外,甚至可以使用交叉验证的估计器,对这些模型中的每一个执行自动模型选择。

    from econml.dr import DRLearner
    from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
    from sklearn.model_selection import GridSearchCV
    model_reg = lambda: GridSearchCV(
                    estimator=RandomForestRegressor(),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=10, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    model_clf = lambda: GridSearchCV(
                    estimator=RandomForestClassifier(min_samples_leaf=10),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=10, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    est = DRLearner(model_regression=model_reg(), model_propensity=model_clf(),
                    model_final=model_reg(), cv=5)
    est.fit(y, T, X=X, W=W)
    point = est.effect(X, T0=T0, T1=T1)
    

    从这个角度来看,这个估计器也是一个元学习器,因为估计的所有步骤都使用了现成的机器学习算法。更多信息,请查阅元学习器用户指南。这种通用方法由 [Foster2019] 提出。

    • LinearDRLearner. 子类 LinearDRLearner 使用未正则化的最终线性模型,基本上只在特征向量 \(\phi(X)\) 是低维时有效。鉴于它是一个未正则化的低维最终模型,该类还通过渐近正态性论证提供置信区间。这主要通过使用 StatsModelsLinearRegression(它是 scikit-learn LinearRegression 估计器的扩展,也支持推断功能)作为最终模型来实现。该类的理论基础主要遵循 [Chernozhukov2016] 中的论证。例如,要获取从基线处理(假定为处理 0)到任何其他处理 T1 的效应的置信区间,只需调用

      from econml.dr import LinearDRLearner
      est = LinearDRLearner()
      est.fit(y, T, X=X, W=W)
      point = est.effect(X, T1=t1)
      lb, ub = est.effect_interval(X, T1=t1, alpha=0.05)
      # Get CATE for all treatments
      point = est.const_marginal_effect(X)
      lb, ub = est.const_marginal_effect_interval(X, alpha=0.05)
      

      还可以通过设置 inference=’bootstrap’ 构建基于自助法的置信区间。

    • SparseLinearDRLearner. 子类 SparseLinearDRLearner 使用 \(\ell_1\) 正则化的最终模型。特别是,它使用了我们实现的 DebiasedLasso 算法 [Buhlmann2011](参见 DebiasedLasso)。利用去偏 Lasso 的渐近正态性性质,该类也提供了基于渐近正态性的置信区间。该类的理论基础主要遵循 [Chernozhukov2017], [Chernozhukov2018] 中的论证。例如,要获取从任何处理 T0 到任何其他处理 T1 的效应的置信区间,只需调用

      from econml.dr import SparseLinearDRLearner
      est = SparseLinearDRLearner()
      est.fit(y, T, X=X, W=W)
      point = est.effect(X, T1=T1)
      lb, ub = est.effect_interval(X, T1=T1, alpha=0.05)
      # Get CATE for all treatments
      point = est.const_marginal_effect(X)
      lb, ub = est.const_marginal_effect_interval(X, alpha=0.05)
      
    • ForestDRLearner. 子类 ForestDRLearner 使用子采样诚实森林回归器作为最终模型(参见 [Wager2018][Athey2019])。子采样诚实森林作为 RandomForestRegressor 的 scikit-learn 扩展,实现在我们的库中的 RegressionForest 类中。该估计器通过 [Athey2019] 中描述的 Bootstrap-of-Little-Bags 方法提供置信区间。利用此功能,我们也可以构建 CATE 的置信区间。

      from econml.dr import ForestDRLearner
      from sklearn.ensemble import GradientBoostingRegressor, GradientBoostingClassifier
      est = ForestDRLearner(model_regression=GradientBoostingRegressor(),
                            model_propensity=GradientBoostingClassifier())
      est.fit(y, T, X=X, W=W)
      point = est.effect(X, T0=T0, T1=T1)
      lb, ub = est.effect_interval(X, T0=T0, T1=T1, alpha=0.05)
      

      此方法与 DROrthoForest 相关,您可以查看 [Oprescu2019] 以获取更多技术细节;主要区别在于如何在目标 \(X=x\) 处构建用于 CATE 估计的辅助模型。有关基于森林的 CATE 模型和 CausalForestDML 的其他替代方案的更多信息,请查阅森林估计器用户指南

使用常见问题解答

  • 如果我想要置信区间怎么办?

    对于有效的置信区间,如果用于分析异质性的特征 \(X\) 数量相对于样本数量较少,请使用 LinearDRLearner,例如:

    from econml.dr import LinearDRLearner
    est = LinearDRLearner()
    est.fit(y, T, X=X, W=W)
    lb, ub = est.const_marginal_effect_interval(X, alpha=.05)
    lb, ub = est.coef__interval(T=1, alpha=.05)
    lb, ub = est.effect_interval(X, T0=T0, T1=T1, alpha=.05)
    

    如果特征数量与样本数量相当甚至更多,那么请使用 SparseLinearDRLearner,并设置 inference='debiasedlasso'。如果您想要非线性模型,则使用 ForestDRLearner 并设置 inference='blb'

  • 如果我对异质性的形式一无所知怎么办?

    要么使用灵活的特征化器,例如具有多个自由度的多项式特征化器,并使用 SparseLinearDRLearner

    from econml.dr import SparseLinearDRLearner
    from sklearn.preprocessing import PolynomialFeatures
    est = SparseLinearDRLearner(featurizer=PolynomialFeatures(degree=3, include_bias=False))
    est.fit(y, T, X=X, W=W)
    lb, ub = est.const_marginal_effect_interval(X, alpha=.05)
    lb, ub = est.coef__interval(T=1, alpha=.05)
    lb, ub = est.effect_interval(X, T0=T0, T1=T1, alpha=.05)
    

    或者,您也可以使用基于森林的估计器,例如 ForestDRLearner。这个估计器也可以处理许多特征,尽管通常比稀疏线性 DRLearner 处理的特征数量少。此外,这个估计器本质上执行自动特征化,并且可以拟合非线性模型。

    from econml.dr import ForestDRLearner
    from sklearn.ensemble import GradientBoostingRegressor, GradientBoostingClassifier
    est = ForestDRLearner(model_regression=GradientBoostingRegressor(),
                          model_propensity=GradientBoostingClassifier())
    est.fit(y, T, X=X, W=W)
    point = est.effect(X, T0=T0, T1=T1)
    lb, ub = est.effect_interval(X, T0=T0, T1=T1, alpha=0.05)
    lb, ub = est.const_marginal_effect_interval(X, alpha=0.05)
    

    如果您更关心均方误差而不是置信区间和假设检验,那么请使用 DRLearner 类,并选择交叉验证的最终模型(请查阅森林学习器 Jupyter notebook 获取此类示例)。另外,请查阅正交随机森林用户指南元学习器用户指南

  • 如果我有很多可能产生异质性的特征怎么办?

    使用 SparseLinearDRLearnerForestDRLearnerDRLearner。(见上文)。

  • 如果我有很多需要控制的特征怎么办?

    使用适用于高维特征的第一阶段模型。例如,Lasso 或 ElasticNet 或梯度提升森林都是不错的选择(后者允许模型中存在非线性,但通常处理的特征数量少于前者),例如:

    from econml.dr import SparseLinearDRLearner
    from sklearn.linear_model import LassoCV, LogisticRegressionCV, ElasticNetCV
    from sklearn.ensemble import GradientBoostingRegressor
    est = SparseLinearDRLearner(model_regression=LassoCV(),
                                model_propensity=LogisticRegressionCV())
    est = SparseLinearDRLearner(model_regression=ElasticNetCV(),
                                model_propensity=LogisticRegressionCV())
    est = SparseLinearDRLearner(model_regression=GradientBoostingRegressor(),
                                model_propensity=GradientBoostingClassifier())
    

    只要这些第一阶段模型达到小的均方误差,置信区间仍然有效。

  • 我应该用什么进行第一阶段估计?

    见上文。第一阶段问题是纯预测任务,因此任何与您的预测问题相关的机器学习方法都是好的。

  • 如何选择第一阶段模型或最终模型的超参数?

    您可以使用自动选择超参数的交叉验证模型,例如使用 LassoCV 代替 Lasso。类似地,对于基于森林的估计器,您可以用网格搜索交叉验证 GridSearchCV 包裹它们,例如:

    from econml.dr import DRLearner
    from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
    from sklearn.model_selection import GridSearchCV
    model_reg = lambda: GridSearchCV(
                    estimator=RandomForestRegressor(),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=5, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    model_clf = lambda: GridSearchCV(
                    estimator=RandomForestClassifier(min_samples_leaf=10),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=5, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    est = DRLearner(model_regression=model_reg(), model_propensity=model_clf(),
                    model_final=model_reg(), cv=5)
    est.fit(y, T, X=X, W=W)
    point = est.effect(X, T0=T0, T1=T1)
    

    另一种方法是,您可以在 EconML 框架之外选择最佳的第一阶段模型,并将选定的模型传递给 EconML。这可以节省运行时和计算资源。此外,由于所有数据都用于超参数调优,而不是 DML 算法内部的单个折叠(只要您选择的超参数值数量不是样本数量的指数,这种方法在统计上是有效的),因此它在统计上更稳定。例如:

    from econml.dr import DRLearner
    from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
    from sklearn.model_selection import GridSearchCV
    model_reg = lambda: GridSearchCV(
                    estimator=RandomForestRegressor(),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=5, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    model_clf = lambda: GridSearchCV(
                    estimator=RandomForestClassifier(min_samples_leaf=10),
                    param_grid={
                            'max_depth': [3, None],
                            'n_estimators': (10, 50, 100)
                        }, cv=5, n_jobs=-1, scoring='neg_mean_squared_error'
                    )
    XW = np.hstack([X, W])
    model_regression = model_reg().fit(XW, Y).best_estimator_
    model_propensity = model_clf().fit(XW, T).best_estimator_
    est = DRLearner(model_regression=model_regression,
                    model_propensity=model_propensity,
                    model_final=model_regression, cv=5)
    est.fit(y, T, X=X, W=W)
    point = est.effect(X, T0=T0, T1=T1)
    
  • 例如:

    如果我有很多处理怎么办?

  • 该方法允许存在多个离散(分类)处理,并将为每个处理估计一个 CATE 模型。

    如何评估 CATE 模型的性能?

    from econml.dr import DRLearner
    from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
    est = DRLearner(model_regression=RandomForestRegressor(oob_score=True),
                    model_propensity=RandomForestClassifier(min_samples_leaf=10, oob_score=True),
                    model_final=RandomForestRegressor())
    est.fit(y, T, X=X, W=W)
    est.score_
    

    每个 DRLearner 类在拟合后都有一个属性 score_。因此,可以通过访问该属性并比较不同建模参数下的性能(分数越低越好)。

    est.score(Y_val, T_val, X_val, W_val)
    

    这主要衡量基于最终阶段损失的得分。此外,可以通过在未用于训练的单独验证样本上调用 score 方法来评估样本外得分。

    [mdl.oob_score_ for mdls in est.models_regression for mdl in mdls]
    

    此外,可以通过检查已拟合的模型来独立检查已拟合的第一阶段模型的拟合优度。您可以通过方法 models_tmodels_y 访问已拟合的第一阶段模型的嵌套列表(每个交叉拟合折叠对应一个)。然后,如果这些模型也有关联的得分属性,则可以用作第一阶段性能的指标。例如,在上述示例中的随机森林第一阶段情况下,如果 oob_score 设置为 True,则估计器具有拟合后的性能衡量指标。

  • 如果使用交叉验证估计器作为第一阶段,则第一阶段模型的模型选择是自动执行的。

    我应该如何设置参数 cv

此参数定义了为以交叉拟合方式拟合第一阶段而创建的数据分区的数量(参见 _OrthoLearner)。默认值为 2,这是最小值。然而,较大的值(如 5 或 6)可以提高方法的统计稳定性,尤其是在样本数量较少的情况下。因此,我们建议对于小型数据集应提高此值。这会增加计算成本,因为需要拟合更多的第一阶段模型。

使用示例