Stata:双层优化匹配估计量-blopmatch

发布时间:2022-12-19 阅读 954

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

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

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc, ihelp, rdbalance, gitee, installpkg

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

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

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者:陈广宇 (中山大学)
邮箱cgyyy_sysu@qq.com


目录


1. 简介

在经济学中,我们常常希望评估某项公共政策实施后的效应。为了解决样本自选择问题,并得到政策的平均处理效应 (average treatment effect),匹配估计量在评价文献中被广泛使用。

常用的匹配方法有最近邻特征匹配 (k-Nearest neighbors matching)、马氏距离匹配 (Mahalanobis matching)、半径匹配 (Radius matching)、倾向得分匹配 (Propensity score match, PSM) 等。其中,倾向得分匹配又需要基于前面提到的任意一种匹配的方法。

最流行匹配方法的莫过于基于最近邻匹配的倾向得分匹配。然而,在进行最邻近匹配的过程中,匹配的结果和效果会受到研究人员给定的近邻数量影响。双层优化匹配 (BLOP Matching) 通过将匹配和权重的选择与匹配后协变量平衡的改进联系起来,能够内生决定最佳近邻数量,且在协变量平衡方面表现出色。

接下来,本文主要介绍双层优化匹配 (BLOP Matching) 的理论内容和 Stata 操作。

2. 理论背景

2.1 基本定义

通常一个匹配问题包含三个关键变量:

  • 结果变量 Yi:即我们关心的实验结果,或者理解为回归中的被解释变量;
  • 分组变量 Wi:即样本是否被处理,通常用哑变量表示;
  • 协变量 Xi:即匹配时需要被控制,使得处理组和控制组尽可能接近的变量。当需要控制的协变量有 K 个时,就构成一个 K 维向量 Xi

以下是关于双层优化匹配问题的一些重要定义:

2.1.1 反事实估计

Y^ib:表示 Yi 的反事实估计,即在控制协变量、未被处理 (被处理) 的情况下,处理组 (对照组) 的估计结果

2.1.2 平均处理效应

τ^b:全样本的平均处理效应 (average treatment effect,ATE)。τ^tre b:处理组的平均处理效应 (average treatment effect on the treated,ATET)。其中 N1 是处理组样本的个数, N0 是控制组样本的个数,N=N1+N0 是总样本的个数。

2.1.3 凸组合和凸壳

凸组合 (Convex combination):λ1X1+...+λNXN,其中 λi 非负且 λ1+...+λN=1。凸组合是多个协变量的一种特殊的线性组合。

凸壳 (Convex hull):co{X1,,XN}={j=1NλjXj,(λ1,,λN)ΔN},其中ΔN={(λ1,,λN)R+N,j=1Nλj=1}。凸壳是所有凸组合的集合,形成了 K 维向量空间的一个壳。

2.1.4 距离和投影

距离:用 ||·|| 表示,可以是欧氏距离、马氏距离等。

投影:用 Proj(X) 表示,是凸壳 (Convex hull) 上到向量 X 距离最短的向量。

2.2 双层优化匹配

双层优化匹配将匹配问题分解为两个问题:第一个问题是投影匹配问题,即在凸壳上求解待匹配协变量 X 的投影 Proj(X),也就是在所有控制组向量组成的凸壳里寻找凸组合,使得这个凸组合距离向量 X 最短。用数学公式表示为:

第二个问题是优化问题。由于第一个问题给出的投影 Proj(X) 多数情况下不唯一,需要用一个描述协变量平衡好坏程度的指标 j=1N0λjXiXj 来确定最优的匹配估计量。用数学公式表示为:

解出上述两个线性优化问题能够对每个处理组样本求解出一组权重 λi=(λ1,...,λN),使得凸组合 λ1X1+...+λNXN 到待匹配协变量 Xi 距离最短,且使得协变量能够最好地进行平衡。这相当于对于处理组里面的每一个样本,都找到了最优的邻居数,并进行了匹配。

类似地,对于控制组的个体,将凸壳从 ΔN0 换成 ΔN1 即可。最终,反事实估计量可以由如下式子表示:

从而通过下式计算出 ATE 和 ATT:

值得一提的是,BLOP 估计量具有一致性和渐近正态性。

3. 命令介绍

blopmatch 命令是 Díaz et al. (2021) 根据以上估计方法所编写的 Stata 新命令。

命令安装:

net install st0632.pkg, replace

命令语法:

blopmatch (ovar omvarlist) (tvar) [if] [in] [weight] [, stat options]

其中,必填项包括:

  • ovar:结果变量 (outcome) ,可以是二进制、计数、连续、分数或非负结果;
  • omvarlist:指定结果模型中的协变量 (covariate);
  • tvar:表示处理级别,必须是整数,只允许二值变量 (通常是哑变量)。

可选项 [stat] 包括:

  • ate:默认,估计全样本的平均处理效应 (average treatment effect);
  • atet:估计处理组的平均处理效应 (average treatment effect on the treated)。

可选项 [options] 包括:

  • ematch(varlist):对于某些变量采用精确匹配;
  • level:设置置信度;
  • dmvariables:输出匹配变量的名字;
  • controltvar 控制组的值,默认是 tvar 较小的,比如哑变量中的 0;
  • tleveltvar 实验组的值,默认是 tvar 较大的,比如哑变量中的 1;
  • metric:指定多维度下距离的度量方式,默认为马氏距离。

4. Stata 实例

4.1 实例演示

接下来,我们使用 LaLonde (1986) 文章中美国国家劳动培训计划 (National Supported Work Demonstration) 的数据来演示双层优化匹配 (BLOP Matching) 的 Stata 操作。

. * 数据描述
. lxhuse nsw.dta, clear
. describe

 Observations:   722   Sample from the National Supported Work Demonstration
    Variables:    10   15 Nov 2019 12:24
                        (_dta has notes)
-----------------------------------------------------------------------------
Variable      Storage   Display    Value
    name         type    format    label      Variable label
-----------------------------------------------------------------------------
data_id         str14   %14s                  
treat           byte    %9.0g      treat      treatment indicator
age             byte    %8.0g                 age (in years)
education       byte    %8.0g                 schooling (in years)
black           byte    %9.0g      black      1(black)
hispanic        byte    %9.0g      hispanic   1(hispanic)
married         byte    %9.0g      married    1(married)
nodegree        byte    %21.0g     nodegree   1(no high school degree)
re75            float   %9.0g                 earnings in 1975 (in 1982 USD)
re78            float   %9.0g                 earnings in 1978 (in 1982 USD)
-----------------------------------------------------------------------------

其中,re78 是工人 1978 年的收入,也是我们感兴趣的结果。协变量包括年龄、教育、人种哑变量 (以白种人为基准组)、是否接受过中学教育、工人 1975 年的收入 re75。对照组与实验组用 treat 区分,1 为实验组,0 为对照组。

. bysort treat: fsum

-> treat = otherwise
  Variable |        N     Mean       SD      Min      Max                                                   
-----------+---------------------------------------------
     treat |      425     0.00     0.00     0.00     0.00  
       age |      425    24.45     6.59    17.00    55.00  
 education |      425    10.19     1.62     3.00    14.00  
     black |      425     0.80     0.40     0.00     1.00  
  hispanic |      425     0.11     0.32     0.00     1.00  
   married |      425     0.16     0.36     0.00     1.00  
  nodegree |      425     0.81     0.39     0.00     1.00  
      re75 |      425  3026.68  5201.25     0.00 36941.27  
      re78 |      425  5090.05  5718.09     0.00 39483.53  
---------------------------------------------------------

-> treat = treated
  Variable |        N     Mean       SD      Min      Max                                                     
-----------+---------------------------------------------
     treat |      297     1.00     0.00     1.00     1.00  
       age |      297    24.63     6.69    17.00    49.00  
 education |      297    10.38     1.82     4.00    16.00  
     black |      297     0.80     0.40     0.00     1.00  
  hispanic |      297     0.09     0.29     0.00     1.00  
   married |      297     0.17     0.37     0.00     1.00  
  nodegree |      297     0.73     0.44     0.00     1.00  
      re75 |      297  3066.10  4874.89     0.00 37431.66  
      re78 |      297  5976.35  6923.80     0.00 60307.93  

可以看出,样本中属于控制组的有 425 人,属于处理组的有 297 人。控制组和处理组在年龄、教育、人种、是否结婚、是否是否接受过中学教育和 1975 年的工资均有较为相似的分布。

. egen stdre75 = std(re75) // 对 re75 进行标准化

. * 匹配并计算平均处理效应 (ATE)
. global Y = "re78"
. global W = "treat"
. global X = "age education black hispanic married nodegree stdre75"
. blopmatch ($Y $X) ($W), ate
                                                                    
Treatment-effects estimation                                        
Estimator      : blop matching                 Number of obs        =      722
Outcome model  : matching                      Control group size   =      425
Distance Metric: Mahalanobis                   Treatment group size =      297
------------------------------------------------------------------------------
        re78 | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         ate |    975.625    503.853     1.94   0.053      -11.910    1963.159
------------------------------------------------------------------------------

根据 BLOP matching,所有样本都得到了匹配。ATE 显示,美国国家劳动培训计划 (National Supported Work Demonstration) 平均给工人 1978 年的工资带来了 975.625 美元的提升,该效应在 10% 的统计水平上显著。

4.2 结果对比

接下来我们对比 BLOP matching 和最近邻匹配的结果:


. * 首先是 k=1,即进行 1:1 匹配
. teffects nnmatch ($Y $X) ($W), ate nneighbor(1)  // ATE = 744.979

Treatment-effects estimation                   Number of obs      =        722
Estimator      : nearest-neighbor matching     Matches: requested =          1
Outcome model  : matching                                     min =          1
Distance metric: Mahalanobis                                  max =          8
-----------------------------------------------------------------------------------------
                        |              AI robust
                   re78 | Coefficient  std. err.      z    P>|z|     [95% conf. interval]
------------------------+----------------------------------------------------------------
ATE                     |
                  treat |
(treated vs otherwise)  |    744.979    565.032     1.32   0.187     -362.464    1852.422
-----------------------------------------------------------------------------------------

. * 更改 k=16
. teffects nnmatch ($Y $X) ($W), ate nneighbor(16) // ATE = 968.917

Treatment-effects estimation                   Number of obs      =        722
Estimator      : nearest-neighbor matching     Matches: requested =         16
Outcome model  : matching                                     min =         16
Distance metric: Mahalanobis                                  max =         23
-----------------------------------------------------------------------------------------
                        |              AI robust
                   re78 | Coefficient  std. err.      z    P>|z|     [95% conf. interval]
------------------------+----------------------------------------------------------------
ATE                     |
                  treat |
(treated vs otherwise)  |    968.917    517.388     1.87   0.061      -45.145    1982.978
-----------------------------------------------------------------------------------------

可以看出,随着 k 的改变,ATE 的计算结果也会改变。而怎样才能在确定最佳的 k 近邻数量,做好协变量平衡,最近邻匹配没有给出答案。

5. 参考资料

  • Díaz, J. D., Gutiérrez, I., & Rivera, J. (2021). Implementing blopmatching in Stata. The Stata Journal, 21(1), 180–194. -Link- -PDF-
  • Díaz, J. D., Rau, T. & Rivera, J. (2015). A Matching Estimator Based on a Bilevel Optimization Problem. The Review of Economics and Statistics, 97(4), 803-812. -Link- -PDF-

6. 相关推文

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

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,700+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxhsongbl 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh