高维情况下基于倾向性评分的因果推断方法*

山西医科大学卫生统计学教研室(030001) 高 倩 王 彤

【提 要】 目的 观察性研究中使用倾向性评分法估计因果效应依赖于不存在未测量混杂假设。随着大数据时代的到来，越来越多的处理前协变量可被测量并纳入分析中，使得这一假设更合理，但同时也导致标准的倾向性评分法不再适用，这就促使了高维情况下因果推断方法的发展。本研究对现有的高维情况下基于倾向性评分的因果推断方法进行概述，为实际应用提供参考依据。方法 从变量选择和协变量组间均衡性两个角度出发，对现有方法做归类和总结，系统回顾各方法的基本原理、方法步骤及优缺点等。结果 基于变量选择的方法可分为Screening、惩罚和贝叶斯三类；基于均衡性的方法有明确的机制确保协变量的组间均衡性。结论 实际应用中研究者应根据自身数据特点结合各方法的优缺点选择合适的方法进行分析，从而得到一致有效的因果效应估计量。

因果推断的金标准是随机对照试验，然而完全随机化试验费时费力，且在很多场景下不可行[1]。相对而言，观察性研究容易实现且更接近现实情况，因此有必要发展统计学方法以从观察性研究中获得因果效应估计值。Rosenbaum和Rubin于1983年提出的倾向性评分(propensity score，PS)方法广泛应用于观察性研究的因果推断[2]。常用的基于PS的方法主要包括匹配、分层、回归校正及逆概率加权(inverse propensity weighting，IPW)，其中IPW方法应用最广，从属于因果模型家族中的边际结构模型[3-8]。将IPW与回归方法结合发展出双稳健估计方法(doubly robust，DR)，这类方法的优点是只要PS模型和结局模型之一正确设定，那么所得估计量就是一致且渐近正态的[9]。

使用基于PS的方法估计因果效应值需要满足不存在未测量混杂假设，即校正充分混杂变量集后，处理组和对照组个体是可比的。这一假设无法检验，一般认为分析中包含的协变量越多那么这一假设越合理，然而因果效应估计值对PS模型中纳入的协变量敏感，如遗漏重要的混杂变量会导致偏倚；模型中包含仅预测处理因素的协变量或无关变量会导致因果效应估计值的方差膨胀而不能降低偏倚[10-11]。从这个意义上说，选择一组适当的混杂因素进行控制对于可靠的因果推断是至关重要的。另一方面，观察性研究中因果效应估计值的偏倚主要是由协变量在处理组和对照组间不均衡引起的，从这一角度讲，确保协变量实现充分的组间均衡性是获得因果效应一致估计的另一关键因素[12-13]。随着大数据时代的到来，成百上千甚至数以万计的处理前协变量可被收集测量，在这种情况下如何使用基于PS的方法得到因果效应的一致有效估计是一个亟待解决的问题。目前已针对这一问题发展出诸多方法，本文将从变量选择及协变量均衡性两个方面介绍高维情况下基于PS的因果推断方法。

基于变量选择的方法

遗漏重要的混杂变量会导致因果效应估计值有偏，因此，研究者们一直提倡将所有观测到的变量均加入到PS模型中[14]。但在高维情况下这一策略不再适用甚至不可行，这是因为将大量无关的变量加入到PS模型中会导致有效性的损失，同时可能由于维度过高而无法拟合PS模型。因此，当存在大量可能的混杂因素时，需要某种类型的变量选择以实现无偏的有效的估计。Brookhart等的模拟研究表明，校正所有的混杂变量足以去除混杂偏倚，且额外的校正仅预测结局的协变量可改善方差估计[15]。这意味着，一种有效的因果推断变量选择方法应同时考虑协变量与结局及处理因素之间的关系。基于此，Schneeweiss等于2009年提出高维倾向性评分方法(high-dimensional propensity score，HDPS)，该方法首先根据协变量可能造成的偏倚对协变量进行排序，随后取前k个协变量构建PS模型并估计因果效应[16]。这一方法存在两个缺陷：一是没有明确的标准以确定进入PS模型中的最优协变量数，为了解决这一问题，Wyss等提出将HDPS方法与超级算法(super learner)或与协作目标最大似然估计(collaborative targeted maximum likelihood estimation，C-TMLE)联合使用[17]；二是计算协变量可能造成的偏倚时依赖于协变量与结局及协变量与处理因素间的边际相关，这可能把与结局条件独立的协变量排在前面。在实际中，HDPS类方法常用于分析电子医疗数据。

无论是screening方法还是惩罚类的方法，最终都是利用一个统计模型估计因果效应，这一过程并没有考虑“校正不确定性”，即模型中包含哪些协变量可恰当地校正混杂偏倚所带来的不确定性。贝叶斯模型平均(Bayesian model averaging，BMA)类方法试图通过对整个模型空间的效应估计进行平均来考虑这一不确定性[28]。简单来说，这类方法将协变量是否包含在模型中转化为一个指示变量α，并将其视为未知的冗余参数。随后计算各模型对应估计值的加权平均作为最终估计值，权重依赖于每个模型的后验分布概率[29-30]。传统的BMA方法在效应估计上存在局限性。Dominici和Parmigiani指出，BMA方法对包含所有混杂变量的模型子空间和没有包含所有混杂变量的模型子空间进行平均(两类模型估计的效应值解释不同)，这导致了估计偏倚[31]。BAC(Bayesian adjustment for confounding)方法通过侧重考虑校正了所有必要混杂变量的模型来避免这一问题[32]。这一方法的核心是定义了一个依赖先验(dependence prior)参数ω，先验地假定如果一个协变量对处理因素具有高度的预测性，那么相同的协变量将有很大的先验概率被纳入结局模型。这使得在计算效应估计值的加权平均时，对可能包含所有必要混杂因素的模型赋予较高的权重。Zigler和Dominici提出的贝叶斯模型平均方法借用PS回归校正的思想，通过将PS以线性预测变量的形式加入结局模型中来控制混杂[33]。此外，还将包含在PS模型中的变量也以线性形式加入到结局模型中以解决“feedback”问题。该方法通过使用联合贝叶斯PS估计(joint Bayesian PS estimation)将模型选择的不确定性整合入因果效应估计的后验分布中。无论是BAC还是Zigler和Dominici提出的贝叶斯模型平均方法都依赖于参数假设，如假设协变量以线性组合形式加入结局模型中，这导致二者对模型误设不稳健。为解决这一问题，Cefalu等人提出了MA-DR(model averaged double robust)方法[34]。与BAC方法类似，MA-DR也定义了一个依赖先验，将PS模型中包含的协变量强制为结局模型中包含的协变量的子集。不同之处在于MA-DR使用DR法估计因果效应值。Wilson和Reich于2014年提出Bayesian credible region方法，主要思想是在所有的可行模型中识别最简单的模型[28]。可行模型定义为包含在PS模型和结局模型(PS模型和结局模型均为贝叶斯回归全模型)中协变量系数的(1-α)×100%的后验可信区域内的模型。通过惩罚不包含混杂因素的可行模型识别最简单的模型。Bayesian credible region方法同时适用于处理因素为二分类和连续型的情况，但值得注意的是该方法是保守的，它能选出更多的协变量[11]。

基于均衡性的方法

基于PS的方法旨在通过均衡处理组和对照组间的协变量实现因果效应的无偏估计，然而传统的PS估计方法，如最大似然估计，并没有正式的机制可以保证协变量最终得到平衡[35]。针对这一问题，有学者提出使用calibrated estimation(CAL)替代最大似然估计拟合PS模型。CAL的主要思想是构建一个新的损失函数lCAL使其对应的估计方程为协变量均衡条件，即处理组或对照组子样本协变量的加权均数等于样本协变量均数[36]。Tan通过在lCAL上增加一个LASSO惩罚项提出RCAL(regularized calibrated estimation)，将CAL推广到高维情况[35]。显然，CAL方法除了可以得到令人满意的协变量组间均衡性外，对模型误设也更加稳健。但RCAL中由于引入了惩罚项导致PS估计值不能严格满足协变量组间均衡性，这可能引发估计偏倚。Zhao提出了一个分析框架(covariate balancing scoring rules，CBSR)用于统一和推广基于协变量均衡性拟合PS的方法[13]。在实际应用中，研究者可根据感兴趣的因果参数及PS模型的链接函数推导对应的得分方程(score function)。与RCAL及CBSR不同，高维协变量均衡倾向性评分(high-dimensional covariate balancing propensity score，hdCBPS)通过最优化协变量的均衡性对PS模型进行修正，包含四步：第一步构建初始的PS模型，该过程使用最大化惩罚广义伪似然函数估计参数；第二步采用加权最小二乘法拟合结局模型；前两步的权重函数至关重要，影响模型错误设定时hdCBPS的表现，作者给出一组权重可供使用；第三步是通过均衡结局模型的预测变量校正PS模型，这一步保证PS的估计值满足弱的协变量均衡性；最后采用IPW估计平均处理效应[37]。与RCAL类似，hdCBPS对模型误设也相对稳健。使用惩罚回归方法校正高维混杂会导致大的偏倚，一方面是由于遗漏弱混杂变量引起的，另一方面是由于惩罚估计是有偏的。Athey等发现可通过实现协变量组间近似平衡去除这一偏倚，并基于此提出approximate residual balancing算法，包括两步[12]，第一步是使用弹性网或LASSO拟合结局模型；第二步是对第一步的残差进行加权，该权重使得所有协变量在处理组与对照组间近似均衡。这一算法将回归与加权的方法结合起来，克服单独使用其中一种方法去除混杂偏倚的局限性。然而，尽管该算法不要求PS模型是可估计的，但其一致估计依赖于结局模型的线性稀疏假设，这一假设在高维情况下很难被满足。

总 结

在观察性研究中，为了估计因果效应，研究者有时需要考虑大量的处理前协变量以使得不存在未测量混杂假设合理。随着可使用数据的增多，如电子病历的出现等，高维情况下如何使用基于PS的方法估计因果效应成为如今研究的热点。本研究从变量选择和协变量均衡性两方面对现有的方法进行了阐述。两类方法有各自的优势和不足。基于变量选择的方法大多从现有的方法发展而来，易于理解且有夯实的理论基础，但因果效应的一致估计往往依赖于混杂变量的正确识别及其进入模型的函数形式。基于均衡性的方法有明确的机制确保协变量在处理组和对照组间实现均衡，且对模型误设相对稳健。然而，在不同高维数据结构下，这两类方法的表现如何仍需进一步研究，这对于实际应用中如何选择最优的方法估计因果效应值有指导意义。此外，在实际应用尤其是公共卫生研究中，我们所关注的处理因素往往不是二分类的，因此，将上述方法扩展到处理因素为连续型的情况、处理因素为随时间变化的情况等将是未来的一个研究方向。