Stata+R：合成DID原理及实现-sdid

发布时间：2022-07-13 阅读 6644

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 系列推文：

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

☝ PDF下载 - 推文合集

作者：陈卓然 (中山大学)
邮箱：chenzhr25@mail2.sysu.edu.cn

编者按：本文主要摘译自下文，特此致谢！
Source：Arkhangelsky D, Athey S, Hirshberg D A, et al. Synthetic difference-in-differences[J]. American Economic Review, 2021, 111(12): 4088-4118. -PDF- -Replicate- -sdid-

DID 和 SC 在现在的实证研究中占据着主导地位，但是二者的适用性不同：DID 要求处理组的个体要有很多，并且施加平行趋势假设；而在 SC 当中只有一个或者极少数个体被处理，通过对控制组个体进行加权平均来解决不满足平行趋势时的问题。SDID 试图将二者结合起来：

类似于 SC，SDID 通过重新赋权的办法匹配处理前趋势，从而弱化对平行趋势的依赖性；
类似于 DID，SDID 对于加和型个体的水平移动不敏感 (invariant to additive unit-level shifts)。

1. 模型简介

考虑一个平衡面板： $N$ 个个体， $T$ 个时期，个体 $i$ 在时期 $t$ 的结果变量为 $Y_{i t}$ ，二元处理变量为 $W_{i t} \in {0, 1}$ 。同时假定前 $N_{c o}$ 个个体不被处理，后 $N_{t r} = N - N_{c o}$ 个个体在 $T_{pre}$ 会被处理，类似 SC。

首先，我们寻找权重 ${\hat{ω}}^{sdid}$ ，以使得处理组的事前趋势和控制组的事前趋势尽可能相似，即 $\sum_{i = 1}^{N_{c o}} {\hat{ω}}_{i}^{sdid} Y_{i t} \approx N_{t r}^{- 1} \sum_{i = N_{c o} + 1}^{N} Y_{i t}$ ，对于所有的 $t = 1, \dots, T_{pre}$ 都成立。然后，寻找权重 ${\hat{λ}}_{t}^{s d i d}$ ，以平衡处理前的时间趋势和处理后的时间趋势。最后，在二维固定效应模型下估计平均处理效应 ( $τ$ )：

({\hat{τ}}^{s d i d}, \hat{μ}, \hat{α}, \hat{β}) = \underset{τ, μ, α, β}{\arg min} {\sum_{i = 1}^{N} \sum_{t = 1}^{T} {(Y_{i t} - μ - α_{i} - β_{t} - W_{i t} τ)}^{2} {\hat{ω}}_{i}^{sdid} {\hat{λ}}_{t}^{sddid}} (1)

相比之下，DID 是通过下述不含任何时间和个体权重的二维固定效应模型来估计处理效应的：

({\hat{τ}}^{d i d}, \hat{μ}, \hat{α}, \hat{β}) = \underset{α, β, μ, τ}{\arg min} {\sum_{i = 1}^{N} \sum_{t = 1}^{T} {(Y_{i t} - μ - α_{i} - β_{t} - W_{i t} τ)}^{2}} (2)

因此 SDID 估计量事实上使得二维固定效应更加 local，也就是说它将更多的权重放在与处理组个体更相似的个体身上，以及将更多的权重放在与处理期更相似的时期上。与 SDID 不同，SC 估计量事实上求解如下的优化问题：

({\hat{τ}}^{s c}, \hat{μ}, \hat{β}) = \underset{μ, β, τ}{\arg min} {\sum_{i = 1}^{N} \sum_{t = 1}^{T} {(Y_{i t} - μ - β_{t} - W_{i t} τ)}^{2} {\hat{ω}}_{i}^{s c}} (3)

可以看到，SC 仅考虑了截面维度上的权重，忽略了时间维度上的权重。SDID 通过加入时间权重可以将哪些和处理期相差甚远的时期排除掉，以减小估计的偏误、提高估计的精度。

2. 加州禁烟

为了更好地比较 ${\hat{τ}}^{d i d}$ ， ${\hat{τ}}^{s c}$ 和 ${\hat{τ}}^{s d i d}$ ，作者回顾了Abadie 的成名之作：加州禁烟，也就是估计美国 California 州增加香烟税之后对其吸烟量的影响。Abadie 考虑了包括 California 在内的 39 个州在 1970 年到 2000 年的观测值。其中，California 州在 1989 年通过了 Proposition 99 禁烟法案，因此在我们的设定中 $T_{pre} = 19$ ， $T_{post} = T - T_{pre} = 12$ ， $N_{c o} = 38$ ， $N_{t r} = 1$ 。

上述 SDID 优化问题的核心在于估计参数： ${\hat{ω}}^{sdid}$ 和 ${\hat{λ}}^{sdid}$ 。首先对于个体权重 (unit weights) ${\hat{ω}}^{sdid}$ ，我们通过求解如下优化问题：

({\hat{ω}}_{0}, {\hat{ω}}^{sdid}) = \underset{ω_{0} \in R, ω \in Ω}{\arg min} ℓ_{unit} (ω_{0}, ω) (4)

其中，

ℓ_{unit} (ω_{0}, ω) = \sum_{t = 1}^{T_{pre}} {(ω_{0} + \sum_{i = 1}^{N_{c o}} ω_{i} Y_{i t} - \frac{1}{N_{t r}} \sum_{i = N_{c o} + 1}^{N} Y_{i t})}^{2} + ζ^{2} T_{p r e} ∥ ω ∥_{2}^{2}

Ω = {ω \in R_{+}^{N} : \sum_{i = 1}^{N_{c o}} ω_{i} = 1, ω_{i} = N_{t r}^{- 1} for all i = N_{c o} + 1, \dots, N}

正则化参数 $ζ$ 等于

ζ = {(N_{t r} T_{p o s t})}^{1 / 4} \hat{σ} (5)

{\hat{σ}}^{2} = \frac{1}{N_{c o} (T_{pre} - 1)} \sum_{i = 1}^{N_{c o}} \sum_{t = 1}^{T_{pre} - 1} {(Δ_{i t} - \bar{Δ})}^{2}

Δ_{i t} = Y_{i (t + 1)} - Y_{i t}

\bar{Δ} = \frac{1}{N_{c o} (T_{p r e} - 1)} \sum_{i = 1}^{N_{c o}} \sum_{t = 1}^{T_{p r e} - 1} Δ_{i t}

SDID 的权重与 Abadie 在合成 DID 中使用的权重主要有以下两点不同：

SDID 中考虑了截距项 $ω_{0}$ ，从而使得权重 ${\hat{ω}}^{sdid}$ 不再需要控制组的事前趋势和处理组的事前趋势完美重合，或者说这一截距项的加入使得在事前只需要控制组和处理组的趋势保持平行即可。之所以可以这样灵活地选择权重，是因为 SDID 在模型设定中加入了个体固定效应 $α_{i}$ ，从而可以吸收个体之间不变的差异。
SDID 中加入了正则化惩罚项来增加权重的分散性 (dispersion)，并确保其唯一性 (uniqueness)。

其次，对于时间权重 ${\hat{λ}}^{s d i d}$ ，我们通过求解如下最优化问题：

({\hat{λ}}_{0}, {\hat{λ}}^{sdid}) = \underset{λ_{0} \in R, λ \in Λ}{\arg min} ℓ_{time} (λ_{0}, λ) (6)

ℓ_{time} (λ_{0}, λ) = \sum_{i = 1}^{N_{c o}} {(λ_{0} + \sum_{t = 1}^{T_{pre}} λ_{t} Y_{i t} - \frac{1}{T_{post}} \sum_{t = T_{pre} + 1}^{T} Y_{i t})}^{2}

Λ = {λ \in R_{+}^{T} : \sum_{t = 1}^{T_{pre}} λ_{t} = 1, λ_{t} = T_{post}^{- 1} for all t = T_{pre} + 1, \dots, T}

SDID 在求解时间权重的时并没有加入任何正则化惩罚项。这反映了在 SDID 中允许同一个个体在不同时间内的观测值是相关的，但是不允许不同个体在同一时期内是相关的。SDID 的具体算法如下：

最终的估计结果如下表所示：

可以看出，DID 要明显高估加州禁烟的政策 (-27.3)。相比之下，SC 的估计结果明显要更可信 (-19.6)，而 SDID 估计出来的结果更加保守。但是作者指出当 ${\hat{τ}}^{s c}$ 和 ${\hat{τ}}^{s d i d}$ 有差异时，后者要更为可信。此外我们也可以看到 SDID 估计出来的方差也是更小的，这是由于 SDID 是局部拟合导致的。

为比较 DID、SC、SDID 三者得到的估计量，我们可以将其分别写成如下的形式：

\hat{τ} = {\hat{δ}}_{t r} - \sum_{i = 1}^{N_{c o}} {\hat{ω}}_{i} {\hat{δ}}_{i} (7)

其中，

{\hat{δ}}_{t r} = \frac{1}{N_{t r}} \sum_{i = N_{c o} + 1}^{N} {\hat{δ}}_{i}

具体而言，

{\hat{δ}}_{i}^{s c} = \frac{1}{T_{post}} \sum_{t = T_{pre} + 1}^{T} Y_{i t} (8)

{\hat{δ}}_{i}^{did} = \frac{1}{T_{post}} \sum_{t = T_{pre} + 1}^{T} Y_{i t} - \frac{1}{T_{pre}} \sum_{t = 1}^{T_{pre}} Y_{i t},

{\hat{δ}}_{i}^{s d i d} = \frac{1}{T_{post}} \sum_{t = T_{p r e} + 1}^{T} Y_{i t} - \sum_{t = 1}^{T_{pre}} {\hat{λ}}_{t}^{s d i d} Y_{i t}

可以看出，传统的标准 DID 事前平行趋势假设并不满足。下面我们比较每一种方法下的 ${\hat{δ}}_{t r} - {\hat{δ}}_{i}$ 。对于每一个控制组的州，点的大小对应于其权重，当权重为 0 时，用 $\times$ 表示。

可以看出，合成控制法下的权重是非常稀疏的。并且，DID 和 SC 都使得某几个州有极高的影响力，也就是说这些州有极高的 ${\hat{ω}}_{i} ({\hat{δ}}_{t r} - {\hat{δ}}_{i})$ 。但是 SDID 并没有给某个州极高的影响力，这意味着在赋权之后，SDID 基本上保证了平行趋势假设的成立。

3. Stata 实操

3.1 命令介绍

命令安装：

ssc install sdid, replace // https://github.com/Daniel-Pailanir/sdid

命令语法：

sdid Y S T D [if] [in], vce(method) seed(#) reps(#) covariates(varlist [, method]) 
    graph g1_opt(string) g2_opt(string) unstandardized graph_export([stub] , type)

其中，

Y：产出变量，只能是数值型；
S：个体变量，可以是数值型或字符串；
T：时间变量，只能是数值型；
D：处理变量，当个体被处理时取值为 1，否则取值为 0；
vce(method)：有三种计算标准误的方法，即 bootstrap、jackknife、placebo；
seed()：设定随机数的种子；
reps：设定 bootstrap 和 placebo 的抽样次数；
covariates(varlist [, method])：用来调整 $Y$ 的控制变量，调整方法有两种。一种是 Arkhangelsky 等提出的 optimized (默认)，另一种是 Kranz (2021) 提出的 projected，后者运算速度要更快；
graph：指定这一选项将会绘制出第 2 部分图形；
g1_opt() 和 g2_opt()：一些画图的选项 (也就是 twoway_options 中的一些选项)；
unstandardized：如果指定这一选项，控制变量会被标准化，避免了在最优化的过程中控制变量过度分散。如果不指定这一选项，则控制变量将以原始形态进入回归当中；
graph_export([stub] , type)：输出图片，命名格式为 weightsYYYY 和 trendsYYYY。其中，YYYY 指的是处理时期，如果处理时期有多起，它将会对每一个处理时期输出上述两张图。在这一选项中 type 是必须指定的，其类型可以是 Stata 支持的任何一种格式 (eps、pdf、png 等)。当然也可以指定图片名字的前缀 stub。

3.2 具体应用

3.2.1 单期应用

. webuse set www.damianclarke.net/stata/
. webuse prop99_example.dta, clear
. sdid packspercapita state year treated, vce(placebo) seed(1213) g1_opt(xtitle("") ///
>     ylabel(-35(5)10) scheme(white_tableau)) g2_opt( ytitle("Packs per capita")    ///
>     xtitle("") scheme(white_tableau)) graph graph_export(lianxh, .png)


Synthetic Difference-in-Differences Estimator
-----------------------------------------------------------------------------
packsperca~a |     ATT     Std. Err.     t      P>|t|    [95% Conf. Interval]
-------------+---------------------------------------------------------------
   treatment | -15.60383    9.00803    -1.73    0.083   -33.25924     2.05158
-----------------------------------------------------------------------------
95% CIs and p-values are based on Large-Sample approximations.
Refer to Arkhangelsky et al., (2020) for theoretical derivations.

3.2.2 多期应用

本部分以 Bhalotra 等 (2020) 为例，在这篇文章中，Bhalotra 等认为增加议会中女性比例可以显著减少女性分娩时的死亡率。作者利用的事件是在 90 年代中期，席卷整个发展中国家的议会性别配额的立法。下图展现了议会中女性比例随时间变化的趋势，以及女性生育致死率的时间趋势。

我们节选了其中的部分数据，在这份数据中一共有 6 个变量，分别为：

womparl：议会中女性的占比；
lnmmrt：女性的生育致死率的自然对数；
country：国家；
year：年份；
quota：处理变量，当一个国家通过了一个性别配额的立法，取值为 1，否则取值为 0。

我们将议会中女性占比作为结果变量，采用 sdid 的方法探究是否通过性别配额立法可以显著增加议会中女性的占比。

. webuse set www.damianclarke.net/stata/
. webuse quota_example.dta, clear
. sdid womparl country year quota, vce(bootstrap) seed(1213) graph graph_export(lianxh,.png)

Synthetic Difference-in-Differences Estimator
-----------------------------------------------------------------------------
     womparl |     ATT     Std. Err.     t      P>|t|    [95% Conf. Interval]
-------------+---------------------------------------------------------------
   treatment |   8.03410    3.74040     2.15    0.032     0.70305    15.36516
-----------------------------------------------------------------------------
95% CIs and p-values are based on Large-Sample approximations.
Refer to Arkhangelsky et al., (2020) for theoretical derivations.

可以看出，系数在 5% 的水平上显著为正。下面将 GDP 的自然对数值作为控制变量加入回归中。需要注意的是，由于 SDID 只能处理平衡面板数据，我们需要先将 lngdp 缺失的观测值删掉。

. drop if lngdp==.
. sdid womparl country year quota, vce(bootstrap) seed(1213) covariates(lngdp, projected)

Synthetic Difference-in-Differences Estimator
-----------------------------------------------------------------------------
     womparl |     ATT     Std. Err.     t      P>|t|    [95% Conf. Interval]
-------------+---------------------------------------------------------------
   treatment |   8.05927    3.11913     2.58    0.010     1.94589    14.17264
-----------------------------------------------------------------------------
95% CIs and p-values are based on Large-Sample approximations.
Refer to Arkhangelsky et al., (2020) for theoretical derivations.

可以看出，系数依然显著，表明通过性别配额法案确实可以有效地增加议会中女性占比。当然，我们也可以使用 jackknife 对标准误进行调整。在使用 jackknife 时，需要保证每一个处理时期内至少要有两个处理个体。

. drop if country=="Algeria" | country=="Kenya" | country=="Samoa"  | ///
>     country=="Swaziland" | country=="Tanzania"
. sdid womparl country year quota, vce(jackknife)

Synthetic Difference-in-Differences Estimator
-----------------------------------------------------------------------------
     womparl |     ATT     Std. Err.     t      P>|t|    [95% Conf. Interval]
-------------+---------------------------------------------------------------
   treatment |  10.41200    6.05456     1.72    0.085    -1.45472    22.27872
-----------------------------------------------------------------------------
95% CIs and p-values are based on Large-Sample approximations.
Refer to Arkhangelsky et al., (2020) for theoretical derivations.

4. R 实操

4.1 基础篇

# 导入 R 包
library(synthdid)
library(ggplot2)
set.seed(12345)

# 导入数据并进行估计
data('california_prop99')
setup = panel.matrices(california_prop99)
tau.hat = synthdid_estimate(setup$Y, setup$N0, setup$T0)
print(summary(tau.hat))

se = sqrt(vcov(tau.hat, method='placebo'))
sprintf('point estimate: %1.2f', tau.hat)
sprintf('95%% CI (%1.2f, %1.2f)', tau.hat - 1.96 * se, tau.hat + 1.96 * se)

# 绘制平行趋势图
plot(tau.hat, se.method='placebo')

其中，深色箭头表示处理效应的下边界，浅色的箭头表示处理效应的上边界。如果在 plot 中不加入 se.method 的选项，则不会绘制浅色箭头。

4.1.1 绘制控制组个体贡献图

正如前文所言，实际处理效应是处理组的结果变量，与加权平均得到的合成控制组的结果变量间的差异，即 $\hat{τ} = {\hat{δ}}_{t r} - \sum_{i = 1}^{N_{c o}} {\hat{ω}}_{i} {\hat{δ}}_{i}$ 。下图展示的是每一个控制个体的贡献 ${\hat{δ}}_{t r} - {\hat{δ}}_{i}$ ，相应的权重反映在点的大小上面。中间深色的实现是估计出来的处理效应，上下两条浅色的部分是处理效应的 95% 置信区间。

sdid_units_plot(tau.hat,se.method="placebo")

4.1.2 检查处理前平行趋势

我们可以将平行趋势图中的两条线上下平移，使之接近或者重合，从而更好地比较处理前是否满足平行趋势。

plot(tau.hat,overlay=1,se.method = 'placebo')

plot(tau.hat,overlay=.8,se.method = 'placebo')

4.1.3 与其他估计量的比较

我们将合成 DID、合成控制法 SC、双重差分 DID 进行比较，对于后两者我们分别采用 sc_estimate 和 did_estimate 进行估计。

# 采用 did 和 sc 的方法去进行估计
tau.sc = sc_estimate(setup$Y,setup$N0,setup$T0)
tau.did = did_estimate(setup$Y,setup$N0,setup$T0)
estimates = list(tau.sc,tau.did,tau.hat)
names(estimates) = c('Synthetic Control',
                     'Diff-in-Diff',
                     'Synthetic-DID')
print(unlist(estimates))
Synthetic Control      Diff-in-Diff     Synthetic-DID 
        -19.61966         -27.34911         -15.60383 

# 三种平行趋势图的比较
synthdid_plot(estimates,se.method = 'placebo')

# 三种单位个体贡献图的比较
synthdid_units_plot(estimates,se.method = 'placebo')

4.1.4 美化图形

由于 synthdid 的绘图是基于 ggplot2 来进行的，我们可以加入一些 ggplot2 的选项使其更加美观。

synthdid_plot(estimates, facet.vertical=FALSE,
              control.name='control', treated.name='california',
              lambda.comparable=TRUE, se.method = 'none',
              trajectory.linetype = 1, line.width=.75, effect.curvature=-.4,
              trajectory.alpha=.7, effect.alpha=.7,
              diagram.alpha=1, onset.alpha=.7) +
              theme(legend.position=c(.21,.07), legend.direction='horizontal',
                    legend.key=element_blank(), legend.background=element_blank(),
                    strip.background=element_blank(), strip.text.x = element_blank())

synthdid_units_plot(estimates, se.method='none') +
    theme(legend.background=element_blank(), legend.title = element_blank(),
          legend.direction='horizontal', legend.position=c(.17,.07),
          strip.background=element_blank(), strip.text.x = element_blank())

4.2 进阶篇

有时我们希望将平行趋势图的横坐标设定为具体的时间而非仅仅是某一年，这时可以使用 as.Date 函数将一个日期的字符串转化为日期格式，作为一个变量储存在数据框中。

data(california_prop99)
california_prop99$date = as.Date(sprintf('%04d/%02d/%02d', california_prop99$Year, 1, 1))
setup = panel.matrices(california_prop99[! california_prop99$Year %in%  c(1974:1977, 1989:1992),], time='date') 
estimate = synthdid_estimate(setup$Y, setup$N0, setup$T0)
plot(estimate, se.method='placebo')

4.2.1 绘制意大利面条图

setup = panel.matrices(california_prop99)
estimate = synthdid_estimate(setup$Y, setup$N0, setup$T0)
top.controls = synthdid_controls(estimate)[1:10, , drop=FALSE]
plot(estimate, spaghetti.units=rownames(top.controls))

我们挑选了控制组中权重最大的 10 个州，然后绘制了他们各自的轨迹。除了将单独的每个州的轨迹绘制出来以外，我们还可以对某些州的轨迹进行加总，或者对其中几个重要的州的轨迹取平均值。

northern.new.england = c('New Hampshire', 'Vermont', 'Maine')
spaghetti.matrices = rbind(colMeans(setup$Y[rownames(setup$Y) %in% northern.new.england, ]),
                           colMeans(setup$Y[rownames(setup$Y) %in% rownames(top.controls), ]))
rownames(spaghetti.matrices) = c('Northern New England States Average', 'Top-10 Control Average')
plot(estimate, spaghetti.matrices=list(spaghetti.matrices), spaghetti.line.alpha=.4)
ggsave("Output/Spaghettiaverage.png")

4.2.2 控制组个体贡献图的简化

当有很多控制组个体时，在个体贡献图中全部将其列出来显得有些拥挤，我们可以采用以下的方式来挑选一部分的个体列示。

synthdid_units_plot(estimate,units = rownames(top.control))

synthdid_units_plot(estimates,units = rownames(top.control))

4.2.3 将图中的一部分放大

estimate.sc = sc_estimate(setup$Y, setup$N0, setup$T0)
with.overlay = function(est, s) { attr(est,'overlay') = s; est }
estimators = function(s) {
  estimator.list = list(with.overlay(estimate, s), estimate.sc, estimate)
  names(estimator.list)=c('synth. diff-in-diff', 'synth. control', '')
  estimator.list
}

plot.estimators = function(ests, alpha.multiplier) {
  p = synthdid_plot(ests, se.method='none',
                alpha.multiplier=alpha.multiplier, facet=rep(1,length(ests)),
            trajectory.linetype = 1, effect.curvature=-.4,
            trajectory.alpha=.5, effect.alpha=.5, diagram.alpha=1)
  suppressMessages(p + scale_alpha(range=c(0,1), guide='none'))
  # scale alpha so alpha=0 means totally invisible, which is unusual but useful
  # for hiding our invisible estimate. We have to suppress a warning that 
  # we're overriding an extant alpha scale that's added in synthdid_plot 
}

# set up the box we zoom in on in plot 5
lambda = attr(estimate, 'weights')$lambda
time = as.integer(timesteps(setup$Y))
xbox.ind = c(which(lambda > .01)[1], setup$T0+4)
xbox = time[xbox.ind] + c(-.5,.5)
ybox = range(setup$Y[setup$N0+1, min(xbox.ind):(max(xbox.ind))]) + c(-4,4)


p1 = plot.estimators(estimators(0),   alpha.multiplier=c(1,.1,0))
p2 = plot.estimators(estimators(.75), alpha.multiplier=c(1,.1,0))
p3 = plot.estimators(estimators(1),   alpha.multiplier=c(1,.1,0))
p4 = plot.estimators(estimators(1),   alpha.multiplier=c(1, 1,0))
p4.zoom = p4 + coord_cartesian(xlim=xbox, ylim=ybox) + xlab('') + ylab('') +
    theme(axis.ticks.x= element_blank(), axis.text.x = element_blank(),
          axis.ticks.y=element_blank(), axis.text.y=element_blank(),
          legend.position='off')
p5 = p4 + annotation_custom(ggplotGrob(p4.zoom), xmin = 1968,   # location by
                            xmax = 1984.7,   ymin=2, ymax=95) + # trial and error
            geom_rect(aes(xmin=min(xbox), xmax=max(xbox),
                          ymin=min(ybox), ymax=max(ybox)),
                      color=alpha('black', .25), size=.3, fill=NA)


plot.theme = theme(legend.position=c(.9,.85), legend.direction='vertical',
                   legend.key=element_blank(), legend.background=element_blank())

p5 + plot.theme

5. 总结

本文简要介绍了最近十分受欢迎的合成 DID 理论，以及在 Stata 和 R 中实现这一方法。但是这篇推文并没有对 SDID 更深层的原理进行介绍，感兴趣的读者可以自行查阅原文以及 Athey Susan 的讲座视频：Synthetic Difference in Difference。

6. 相关推文

Note：产生如下推文列表的 Stata 命令为：
lianxh did 合成控制法, m
安装最新版 lianxh 命令：
ssc install lianxh, replace

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