ebalance：基于熵平衡法的协变量平衡性检验

发布时间：2021-05-05 阅读 15126

Stata连享会主页 || 视频 || 推文 || 知乎

温馨提示： 定期清理浏览器缓存，可以获得最佳浏览体验。

New！ lianxh 命令发布了：
随时搜索推文、Stata 资源。安装命令如下：
. ssc install lianxh
详情参见帮助文件 (有惊喜)：
. help lianxh

✌ 课程详情： https://gitee.com/lianxh/Course

⛳ 课程主页： https://gitee.com/lianxh/Course

⛳ Stata 系列推文：

全部 | Stata入门 | Stata教程 | Stata资源 | Stata命令
计量专题 | 论文写作 | 数据分享 | 专题课程
结果输出 | Stata绘图 | 数据处理 | Stata程序
回归分析 | 面板数据 | 交乘项-调节 | IV-GMM
内生性-因果推断 | 倍分法DID | 断点回归RDD | PSM-Matching | 合成控制法
Probit-Logit | 时间序列 | 空间计量 | 分位数回归 | 生存分析 | SFA-DEA
文本分析-爬虫 | Python-R-Matlab | 机器学习
Markdown | 工具软件 | 其它

☝ PDF下载 - 推文合集

作者：王文韬 (山东大学)
邮箱：190138445@qq.com

1. 引言

在进行随机实验时，我们需要对处理组和对照组的协变量进行平衡性检验，以保证实验的随机性或外生性，增强回归结果的有效性。Hainmueller (2012) 构建了「熵平衡方法」 (Entropy Balancing)，使研究者能够同时控制处理组与对照组样本协变量多维平衡性，如同时考虑协变量的一阶矩、二阶矩、三阶矩和交叉矩等，进而最大程度上使两组样本在实现精确匹配。

需要强调的是，熵平衡方法并不能完全克服内生性问题。在实际工作中，研究者仍需搭配其他计量工具进行系统检验，以保证回归结果的稳健性和合理性。

本文的主要目的是介绍熵平衡方法理论、及 Stata 中的实现命令 ebalance。

2. 熵平衡

2.1 基本含义

在高中物理课堂上，我们曾学习过 “熵” 的概念。它是热力学中表征物质状态的参数之一，能够度量一个体系的混乱程度。在物理热力学中，熵越大表示该体系中的热量越高。在计量经济学中，“熵平衡” 则是一个基于最大化熵的权重计算方法。

所谓最大化熵，指的是研究者进行随机性实验研究时，在对数据进行预处理的过程中，最大化处理组与对照组的协变量匹配效率。

2.2 总体平均处理效应与熵平衡

在进行随机性实验时，研究者主要关注的结果是总体平均处理效应 (the Population Average Treament Effect on the Treated)，可以由以下公式表示：

τ = E [Y (1) | D = 1] - E [Y (0) | D = 1]

其中， $E [Y (1) | D = 1]$ 表示处理组实验发生后的效应大小， $E [Y (0) | D = 1]$ 表示处理组在没有发生实验后的效应大小。

显然， $E [Y (0) | D = 1]$ 是一个反事实指标，在现实中不可能存在对应的数据。为此，学者尽可能寻找与处理组协变量分布相似的对照组样本，来近似估计这一项。

目前，我们熟知的协变量匹配方法包括最近邻匹配 (NNM)、粗化精确匹配 (CEM)、倾向得分匹配 (PSM) 等。以倾向得分匹配为例， $E [Y (0) | D = 1]$ 可以通过以下表达式计算得到：

E [Y (0) | D = 1] = \frac{\sum_{i | D = 0} Y_{i} d_{i}}{\sum_{i | D = 0} d_{i}}

其中， $d_{i} = \frac{p (x_{i})}{1 - p (x_{i})}$ ， $p_{i}$ 为倾向得分值。

不过，这些匹配方法面临的一个共同问题是：研究者无法保证联合并平衡所有协变量，难以避免倾向得分模型被错误指定的可能。为了解决这一问题，学者通常利 Logit 回归挑选可以通过平衡性检验的协变量组合。但是，这样做费时费力，同时无法完全保证最终得到高水平的协变量组合。

在上述基础上，Hainmueller (2012) 推导出熵平衡方法，以熵权重 $ω_{i}$ 代替 $d_{i}$ ，得到反事实指标 $E [Y (0) | D = 1]$ 的表达式如下所示：

E [Y (0) | D = 1] = \frac{\sum_{i | D = 0} Y_{i} ω_{i}}{\sum_{i | D = 0} ω_{i}}

权重 $ω_{i}$ 通过以下公式决定：

m i n H (ω_{i}) = \sum_{i | D = 0} ω_{i} l o g (ω_{i} / q_{i})

其中， $q_{i} = 1 / n$ ， $n$ 表示对照组样本数量。 $ω_{i}$ 具有一定约束条件，其核心思想是，给予对照组协变量增加一组矩约束，使其与处理组的协变量相平衡。矩约束包括均值 (一阶矩)、方差 (二阶矩) 和偏度 (三阶矩)。有关最小化熵距离和 $ω_{i}$ 的相关数理验证这里就不展开描述了，对此感兴趣的读者可以自行阅读 Hainmueller 的相关文献。

综上所述，Hainmueller (2012) 创造的熵平衡方法为研究者在平衡处理组与对照组的协变量分布上提供了新的科学方法。接下来，我们使用 Hainmueller 和 Xu (2013) 提供的案例及其数据来说明 ebalance 命令在 Stata 中的具体应用。

3. Stata 实现

3.1 命令安装和数据介绍

首先，安装 ebalance 命令，以及下载所需要的数据，具体命令如下：

. net describe ebalance, from(http://fmwww.bc.edu/RePEc/bocode/e)

. ssc install ebalance, all replace

. net get ebalance.pkg, replace // 下载数据到当前工作路径

我们使用的案例数据是 Dehejia 和 Wahba (1999) 提供的「cps1re74.dta」。它包括 185 名参与 NSW (the National Supported Work Demonstration，一个关于补贴资助工作的随机评估) 项目的成员和 15992 名未参与项目的成员 (来自美国社会保障管理档案 CPS-1)。其主要变量如下所示：


. use "cps1re74.dta", clear     // 载入数据
. des

Contains data from cps1re74.dta
  obs:        16,177             DW Subset of LaLonde Data
 vars:            12             6 Oct 2011 14:57
-------------------------------------------------------------
              storage   display                              
variable name   type    format   variable label              
-------------------------------------------------------------
re78            double  %9.0g    real earnings 78
treat           long    %9.0g    1 if treated, 0 control
age             long    %9.0g    age in years
educ            long    %9.0g    years of schooling
black           long    %9.0g    indicator for black
hispan          long    %9.0g    indicator for hispanic
married         long    %9.0g    indicator for married
nodegree        long    %9.0g    indicator for no HS degree
re74            double  %9.0g    real earnings 74
re75            double  %9.0g    real earnings 75
u74             double  %9.0g    indicator for unemployed 74
u75             double  %9.0g    indicator for unemployed 75
-------------------------------------------------------------

本文衡量政策效果的变量为 re78，根据 Dehejia 和 Wahba (1999) 的研究分析，NSW 这一补贴资助计划平均提高了 1794 美元的收入水平，95% 的置信区间为 [551, 3038]。如果我们直接使用 OLS 回归进行检验，结果如下所示：

. use cps1re74.dta, clear
(DW Subset of LaLonde Data)

. reg re78 treat age-u75

   Source |       SS           df       MS      Number of obs   =    16,177
----------+----------------------------------   F(11, 16165)    =   1343.88
    Model |  7.2418e+11        11  6.5835e+10   Prob > F        =    0.0000
 Residual |  7.9190e+11    16,165  48988567.3   R-squared       =    0.4777
----------+----------------------------------   Adj R-squared   =    0.4773
    Total |  1.5161e+12    16,176  93724175.2   Root MSE        =    6999.2
---------------------------------------------------------------------------
     re78 |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
----------+----------------------------------------------------------------
    treat |   1067.546    554.060     1.93   0.054      -18.472    2153.564
      age |    -94.541      6.000   -15.76   0.000     -106.302     -82.780
     educ |    175.225     28.697     6.11   0.000      118.977     231.474
    black |   -811.089    212.849    -3.81   0.000    -1228.296    -393.881
   hispan |   -230.535    218.610    -1.05   0.292     -659.034     197.965
  married |    153.228    142.775     1.07   0.283     -126.626     433.083
 nodegree |    342.927    177.878     1.93   0.054       -5.734     691.587
     re74 |      0.291      0.013    22.89   0.000        0.266       0.316
     re75 |      0.443      0.013    34.35   0.000        0.417       0.468
      u74 |    355.556    231.600     1.54   0.125      -98.406     809.519
      u75 |  -1612.758    239.803    -6.73   0.000    -2082.798   -1142.717
    _cons |   5762.180    445.614    12.93   0.000     4888.726    6635.634
---------------------------------------------------------------------------

可以发现，OLS 估计的处理效应为 1067.546，这大大低于真正的政策实施效果。于是我们考虑利用熵平衡对数据进行预处理，看看是否能使我们更准确地评估实际情况。

3.2 ebalance 的基本语法

ebalance的基本语法如下所示：

ebalance [treat] [covar], targets(numlist)

其中，treat 代表是否为处理组变量，covar 代表需要在熵平衡中进行调整的协变量组合，targets(numlist) 是最为关键的选择指标，numlist 对应协变量组合需要调整的矩阶数。具体来看：

. ebalance treat age black educ, targets(3) //表示对 age、black 和 educ 三个协变量的一阶、二阶和三阶矩进行调整

Data Setup
Treatment variable:   treat
Covariate adjustment: age black educ (1st order). age black educ (2nd order). age black educ (3rd order).

Optimizing...
Iteration 1: Max Difference = 580799.347
Iteration 2: Max Difference = 213665.688
Iteration 3: Max Difference = 78604.7628
...
Iteration 15: Max Difference = .420310244
Iteration 16: Max Difference = .037151116
Iteration 17: Max Difference = .008791339
maximum difference smaller than the tolerance level; convergence achieved


Treated units: 185     total of weights: 185
Control units: 15992   total of weights: 185


Before: without weighting

      |         Treat             |         Control            
      |  mean  variance  skewness |   mean  variance  skewness 
------+---------------------------+----------------------------
  age | 25.82     51.19     1.115 |  33.23       122     .3478 
black | .8432     .1329    -1.888 | .07354    .06813     3.268 
 educ | 10.35     4.043    -.7212 |  12.03     8.242    -.4233 


After:  _webal as the weighting variable

      |         Treat             |        Control            
      |  mean  variance  skewness |  mean  variance  skewness 
------+---------------------------+---------------------------
  age | 25.82     51.19     1.115 | 25.75     51.22      1.14 
black | .8432     .1329    -1.888 | .8423     .1328    -1.879 
 educ | 10.35     4.043    -.7212 | 10.35     4.039    -.7224

可以看出，在调整后，处理组和加权调整后的对照组的年龄均值、方差和偏度都非常相似。需要指出的是，ebalance treat age, targets(2) 与 ebalance treat age age2, targets(1) 表示的执行指令相同，这里 age2 代表年龄变量 age 的平方项。

3.3 平衡性检验

为了验证某一个协变量匹配后是否平衡，可以使用 tabstat 进行验证，以年龄 age 为例：

. tabstat age [aweight=_webal], by(treat) s(N me v) nototal //检验是否平衡

Summary for variables: age
     by categories of: treat (1 if treated, 0 control)

   treat |         N      mean  variance
---------+------------------------------
       0 |     15992  25.75072  51.22003
       1 |       185  25.81622   51.1943
----------------------------------------

默认情况下，ebalance 命令执行后，Stata 会自动将匹配权重储存为变量 _webal。从上面的检验结果看，经过熵平衡调整之后，age 变量在对照组的分布几乎与处理组一致。

3.4 变量交互项的匹配

在实际研究中，我们还需要考虑协变量之间的交互项。例如，为了调整协变量组，我们应保证黑人和非黑人的年龄均值相似，这里便要引入 black 和 age 的交互项。

. gen ageXblack = age*black  //设置 age 和 black 的交互项 ageXblack

. ebalance treat age educ black ageXblack, targets(3 2 1 1)  //进行熵平衡

Data Setup
Treatment variable:   treat
Covariate adjustment: age educ black ageXblack (1st order). age educ (2nd order). age (3rd order).

Optimizing...
Iteration 1: Max Difference = 573647.601
Iteration 2: Max Difference = 211032.079
Iteration 3: Max Difference = 77633.2843
...
Iteration 15: Max Difference = .277671087
Iteration 16: Max Difference = .015938697
Iteration 17: Max Difference = .000055717
maximum difference smaller than the tolerance level; convergence achieved

Treated units: 185     total of weights: 185
Control units: 15992   total of weights: 185

Before: without weighting

          |         Treat             |         Control            
          |  mean  variance  skewness |   mean  variance  skewness 
----------+---------------------------+----------------------------
      age | 25.82     51.19     1.115 |  33.23       122     .3478 
     educ | 10.35     4.043    -.7212 |  12.03     8.242    -.4233 
    black | .8432     .1329    -1.888 | .07354    .06813     3.268 
ageXblack | 21.91     134.6    -.4435 |  2.402     81.55     3.893 


After:  _webal as the weighting variable

          |         Treat             |        Control            
          |  mean  variance  skewness |  mean  variance  skewness 
----------+---------------------------+---------------------------
      age | 25.82     51.19     1.115 | 25.82     51.19     1.115 
     educ | 10.35     4.043    -.7212 | 10.35     4.043    -.7231 
    black | .8432     .1329    -1.888 | .8432     .1322    -1.888 
ageXblack | 21.91     134.6    -.4435 | 21.91     133.2    -.4572

. bysort black: tabstat age [aweight=_webal], by(treat) s(N me v) nototal //根据 _webal 验证年龄均值在黑人和非黑人之间是平衡的

--------------------------------------------------------
-> black = 0

Summary for variables: age
     by categories of: treat (1 if treated, 0 control)

   treat |         N      mean  variance
---------+------------------------------
       0 |     14816  24.93109  45.28087
       1 |        29  24.93103  40.49507
----------------------------------------

--------------------------------------------------------
-> black = 1

Summary for variables: age
     by categories of: treat (1 if treated, 0 control)

   treat |         N      mean  variance
---------+------------------------------
       0 |      1176  25.98077  52.16224
       1 |       156  25.98077   53.2835
----------------------------------------

3.5 NSW 实验效果的无偏估计

为最大程度的平衡样本，我们使用所有协变量的一阶、二阶、三阶以及所有一阶交互项进行 ebalance 匹配。

. use cps1re74.dta, clear
(DW Subset of LaLonde Data)

. *协变量组合生成二阶矩及一阶交互项
. foreach v in age educ black hispan married nodegree re74 re75 u74 u75 {
    foreach m in age educ black hispan married nodegree re74 re75 u74 u75 {
      gen `v'X`m'=`v'*`m'
    }
  }

. *age、educ、re74 和 re75 作为连续变量，还需设置他们的三阶矩
. foreach v in age educ re74 re75 {
      gen `v'X`v'X`v' = `v'^3
  }

. *进行熵平衡，并将结果保存到 baltable.dta 
. ebalance treat-re75Xre75Xre75, keep(baltable) replace

. reg re78 treat [pweight=_webal] 

Linear regression               Number of obs     =     16,177
                                F(1, 16175)       =       5.58
                                Prob > F          =     0.0182
                                R-squared         =     0.0161
                                Root MSE          =     6889.6

---------------------------------------------------------------
       |            Robust                                     
  re78 |    Coef.  Std. Err.   t   P>|t|   [95% Conf. Interval]
-------+-------------------------------------------------------
 treat | 1761.344   745.533  2.36  0.018    300.017    3222.671
 _cons | 4587.800   472.243  9.71  0.000   3662.151    5513.449
---------------------------------------------------------------

可以发现，上述回归结果显示的处理效果为 1761 美元，这表明熵平衡预处理步骤能够让我们非常接近随机实验的实际评估结果，95% 置信区间为 [300, 3223]。

4. 结论

熵平衡可以与回归方法相结合。研究者首先通过调整协变量对协变量进行重新加权，然后将加权后的数据应用于回归模型。这一过程可应用于回归前对处理组和对照组数据的预处理工作以及后续的稳健性检验。

在实证回归中，ebalance 命令的关键点是如何设置阶数 (numlist)：

对于虚拟变量及包含虚拟变量的交互项，只需进行一阶矩加权调整即可实现平衡；
对于年龄等连续变量，通常选择其一阶、二阶、三阶矩以及一阶与其他协变量一阶矩的交互项进行加权调整；
在进行熵平衡之后，权重值自动储存为变量 _webal。

熵平衡不是万能的。在实际研究中，为保证检验结果的稳健性，研究者最好搭配 NNM、PSM、CEM 等其他匹配方法一起进行实证分析。Hainmueller 和 Xu (2013) 也提出，未来会考虑如何将熵平衡与 Stata 中其他匹配方法相结合，进一步提高匹配效率。

5. 参考文献

Hainmueller J. Entropy balancing for causal effects: A multivariate reweighting method to produce balanced samples in observational studies[J]. Political analysis, 2012: 25-46. -PDF-
Hainmueller J, Xu Y. Ebalance: A Stata package for entropy balancing[J]. Journal of Statistical Software, 2013, 54(7). -PDF-
Dehejia R H, Wahba S. Causal effects in nonexperimental studies: Reevaluating the evaluation of training programs[J]. Journal of the American statistical Association, 1999, 94(448): 1053-1062. -PDF-
张海峰, 梁若冰, 林细细. 子女数量对农村家庭经济决策的影响——兼谈对 “二孩政策” 的启示[J]. 中国经济问题, 2019 (3): 68-80. -Link-

6. 相关推文

Note：产生如下推文列表的 Stata 命令为：
lianxh 匹配, m
安装最新版 lianxh 命令：
ssc install lianxh, replace

专题	嘉宾	直播/回看视频
⭐ 最新专题		文本分析、机器学习、效率专题、生存分析等
研究设计	连玉君	我的特斯拉-实证研究设计，-幻灯片-
面板模型	连玉君	动态面板模型，-幻灯片-
面板模型	连玉君	直击面板数据模型 [免费公开课，2小时]