Stata：广义精确匹配-Coarsened Exact Matching (CEM)

「匹配」是一种文献中流行的方法，但是，匹配成功与否取决于是否通过平衡性检验。常用的 PSM 方法平衡性检验十分费力，而且严格来说，即使发表的论文，协变量在匹配后也不一定更加平衡，往往是某些变量的平衡性得到提升，而另一部分变量的平衡性有所下降。因此，我们介绍一种无需检查协变量平衡性、模型依赖度更低的匹配方法 ——Coarsened Exact Matching (广义精确匹配 \ 粗粒度精确匹配)。

0. 命令汇总

可以先安装外部命令cem：

ssc install cem

cem 的基本语法格式是：

cem varname1 [(cutpoints1)] [varname2 [(cutpoints2)]] ... [if] [in] [, options]

varname# 是一组协变量(匹配时关注的特征，比如，进行个人特征匹配时，可能是 高矮、胖瘦、年龄、性别等)，cutpoints 是事先 coarsening 分组的标准，options 比较多，介绍一个好用的 ，保证匹配后的处理组和控制组有相等的观测值。

1. CEM 原理简介

匹配本质上是一种控制处理组和控制组处理前差异的非参数方法，关键在于筛选样本，提高平衡性，换言之，处理组和控制组的协变量经验分布更加接近。常用的 PSM 方法往往无法确保在匹配后提升平衡性，而 Coarsened Exact Matching(CEM) 可以通过控制观测数据中混杂因素对政策结果影响使处理组与控制组的协变量的分布尽可能保持平衡，从而增强两组数据之间的可比性。在 CEM 中，处理组和控制组之间的最大不平衡可以由用户事先选择，而不是通过那些费力的事后检查和重复重新评估的过程来发现，并且调整一个变量上的不平衡不会影响任何其他变量的不平衡。

CEM 通过预先的用户选择，将变量根据用户的标准分组，在保留基本信息的同时，提升了匹配的效率，限制了模型依赖程度和平均处理效应估计误差。它能够将数据限制在共同支持范围内，符合一致性原则，对测量误差的敏感度较低，能够很好地用多重插补方法处理缺失值数据。此外，它可以完全自动化地进行这系列操作，即使对于非常大的数据集，计算处理速度也很快。

CEM 有一个非常重要的性质——单调不平衡边界。它通过变量事先分组的标准，确保匹配提升样本的平衡性，而且限制了模型依赖程度和平均处理效应估计误差，因此无需匹配后检查协变量的平衡性。

2. CEM 操作步骤

首先根据用户事先的截断点，将协变量进行分层。其次运用精确匹配的算法根据样本的经验分布将每层中的研究对象进行精确匹配，保证每层中至少有一个处理组和一个对照组的研究对象，否则将该研究对象删除。最后保留匹配成功的研究对象，用匹配后的数据，研究政策变量对结果变量的影响。

匹配中有一个重要的变量L1，L1的取值范围是[0,1]。若L1=0，则说明两组数据完全平衡，若L1=1，则说明两组数据完全不平衡，越接近1则说明不平衡程度越大。一般来说，匹配后 L1 较匹配前的 L1 有所下降，CEM 的匹配效果较好。CEM 匹配后两组数据的样本量可能会不相等，因此 CEM 匹配过程中会产生权重变量(Weight)，以此来平衡每层中参与组和控制组的人数。

3. CEM 演示数据

我们使用 Lalonde (1986) 的数据进行演示：这个数据是关于一个美国 12-18 月的就业培训项目，研究的目的是探究接受培训是否提升了 1978 年的工资。首先我们从 NBER 中下载 nsw 这个数据集（Source: https://users.nber.org/~rdehejia/data/nswdata2.html）, 然后导入Stata，我们看一下样本的分布，

. use "https://users.nber.org/~rdehejia/data/nsw_dw.dta", clear
. tab treat

      treat |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |        260       58.43       58.43
          1 |        185       41.57      100.00
------------+-----------------------------------
      Total |        445      100.00

可以看到样本中有 185 个处理组观测值和 260 个控制组观测值。接着，我们可以通过一个简单线性回归来找到两组之间均值的差异，

. reg re78 treat

      Source |       SS           df       MS      Number of obs   =       445
-------------+----------------------------------   F(1, 443)       =      8.04
       Model |   348013183         1   348013183   Prob > F        =    0.0048
    Residual |  1.9178e+10       443  43290369.3   R-squared       =    0.0178
-------------+----------------------------------   Adj R-squared   =    0.0156
       Total |  1.9526e+10       444  43976681.9   Root MSE        =    6579.5

------------------------------------------------------------------------------
        re78 |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
       treat |   1794.342    632.853     2.84   0.005      550.574    3038.110
       _cons |   4554.801    408.046    11.16   0.000     3752.855    5356.747
------------------------------------------------------------------------------

根据回归结果，我们得到样本的平均处理效应是 1794.342。显然，这个结果是有偏的，受到两组之间很多潜在因素的影响，因此我们考虑协变量年龄 age，受教育年限 education，种族 black，无高中文凭 nodegree，1975 年的工资 re75。我们用 L1 统计量来判断两组之间协变量的不平衡程度：

. imb age education black nodegree re74, treatment(treat)

Multivariate L1 distance: .48284823

Univariate imbalance:

                L1     mean      min      25%      50%      75%      max
      age   .10915   .76237        0        1        1        1       -7
education   .16933   .25748        1        0        1        1        2
    black   .01632   .01632        0        0        0        0        0
 nodegree   .12651  -.12651        0       -1        0        0        0
     re74   .05894  -11.453        0        0        0   1291.5  -4530.6

从表中可以看出，两组协变量的 L1 统计量是 0.4828，这为我们的匹配提供了一个基准组，匹配后的不平衡程度可以与该值进行比对。下面的表从第二列开始分别报告了两组在最小值、25 分位数、50 分位数、75 分位数和最大值上的差异。数据显示，1975 年工资的不平衡性最高，年龄的不平衡性最低但不体现在分位数上。所以，均值上的差异并不能完整地反映不平衡性。这从侧面体现出，L1 统计量从整体度量不平衡性的重要性。然后，我们依据协变量进行 CEM 匹配，

. cem age education black nodegree re75, tr(treat)
(using the scott break method for imbalance)

Matching Summary:
-----------------
Number of strata: 121
Number of matched strata: 38

             0    1
      All  260  185
  Matched  204  133
Unmatched   56   52


Multivariate L1 distance: .41924802

Univariate imbalance:

                 L1      mean       min       25%       50%       75%       max
      age    .06116     -.066         0         1         0        -1        -1
education     .0662   -.01983        -1        -1         0         0         0
    black   1.7e-16   2.2e-16         0         0         0         0         0
 nodegree   5.6e-17  -1.1e-16         0         0         0         0         0
     re75    .07594    112.24         0         0         0    225.39   -505.18

我们可以看到，控制组中的 185 个观测值中有 133 个得到匹配 ，处理组中对应的数字则是 260:204。L1 统计量下降为 0.419，其它分位数上的不平衡性也得到改善，这说明我们的匹配是成功的。但是如前面所说，CEM 可以通过变量事先分组，确保匹配提升样本的平衡性，限制模型依赖程度和平均处理效应估计误差，因此无需费力地进行协变量平衡性检验。在美国的教育体系中，有这样的一个年限的划分：

Grade school		0-6
Middle school		7-8
High school		9-12
College			13-16
Graduate school		>16

我们把这个作为 coarsening 的依据，给受教育年限这个协变量进行分组，然后进行改进的 CEM 的匹配

. cem age education (0 6 8 12 16) black nodegree re75, tr(treat)

Matching Summary:
-----------------
Number of strata: 99
Number of matched strata: 39

             0    1
      All  260  185
  Matched  223  150
Unmatched   37   35


Multivariate L1 distance: .49501281

Univariate imbalance:

                L1     mean      min      25%      50%      75%      max
      age   .06147  -.07964        0        0        0       -1       -1
education   .07539  -.05256       -1       -1        0        0        1
    black  7.8e-16  1.2e-15        0        0        0        0        0
 nodegree  1.1e-15  1.6e-15        0        0        0        0        0
     re75   .09364   71.046        0        0        0   355.65  -505.18

新的 L1 统计量为 0.495，似乎变得更加糟糕了。对比单个变量的匹配情况，不平衡性反而加大了。为此，我们仍然选用此前的设定方法进行匹配，进而使用其返回的权重矩阵来执行估计，以便得到样本平均处理效应的估计值。

. qui cem age education black nodegree re75, tr(treat)
. reg re78 treat [iweight=cem_weights]

      Source |       SS           df       MS      Number of obs   =       337
-------------+----------------------------------   F(1, 335)       =      5.26
       Model |   239267421         1   239267421   Prob > F        =    0.0225
    Residual |  1.5249e+10       335    45519065   R-squared       =    0.0154
-------------+----------------------------------   Adj R-squared   =    0.0125
       Total |  1.5488e+10       336  46095697.1   Root MSE        =    6746.8

------------------------------------------------------------------------------
        re78 |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
       treat |   1723.915    751.918     2.29   0.022      244.838    3202.992
       _cons |   4469.369    472.369     9.46   0.000     3540.185    5398.552
------------------------------------------------------------------------------

非常简单，只需要把匹配的权重代入我们的简单回归就可以了，任何支持权重 (aweight or iweight) 的回归命令都可以和 CEM 匹配兼容。我们的例子中，新的样本平均处理效应是 1723.915，未进行匹配的 OLS 估计结果 (1794.342) 高估了就业培训对工资的作用。