-
Notifications
You must be signed in to change notification settings - Fork 45
/
Copy pathREADME.Rmd
261 lines (172 loc) · 15.4 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
---
title: "HonestDiD"
output: github_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning = FALSE)
```
The HonestDiD R package implements the tools for robust inference and sensitivity analysis for differences-in-differences and event study designs developed in [Rambachan and Roth (2022)](https://asheshrambachan.github.io/assets/files/hpt-draft.pdf). There is also an [HonestDiD Stata package](https://github.com/mcaceresb/stata-honestdid#honestdid), and a [Shiny app](https://ccfang2.shinyapps.io/HonestDiDSenAnlys/) developed by Chengcheng Fang.
## Background
The robust inference approach in Rambachan and Roth formalizes the intuition that pre-trends are informative about violations of parallel trends. They provide a few different ways of formalizing what this means.
**Bounds on relative magnitudes.** One way of formalizing this idea is to say that the violations of parallel trends in the post-treatment period cannot be much bigger than those in the pre-treatment period. This can be formalized by imposing that the post-treatment violation of parallel trends is no more than some constant <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> larger than the maximum violation of parallel trends in the pre-treatment period. The value of <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> = 1, for instance, imposes that the post-treatment violation of parallel trends is no longer than the worst pre-treatment violation of parallel trends (between consecutive periods). Likewise, setting <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> = 2 implies that the post-treatment violation of parallel trends is no more than twice that in the pre-treatment period.
**Smoothness restrictions.** A second way of formalizing this is to say that the post-treatment violations of parallel trends cannot deviate too much from a linear extrapolation of the pre-trend. In particular, we can impose that the slope of the pre-trend can change by no more than *M* across consecutive periods, as shown in the figure below for an example with three periods.
![diagram-smoothness-restriction](deltaSD.png)
Thus, imposing a smoothness restriction with *M* = 0 implies that the counterfactual difference in trends is exactly linear, whereas larger values of *M* allow for more non-linearity.
**Other restrictions**. The Rambachan and Roth framework allows for a variety of other restrictions on the differences in trends as well. For example, the smoothness restrictions and relative magnitudes ideas can be combined to impose that the non-linearity in the post-treatment period is no more than <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> times larger than that in the pre-treatment periods. The researcher can also impose monotonicity or sign restrictions on the differences in trends as well.
**Robust confidence intervals**. Given restrictions of the type described above, Rambachan and Roth provide methods for creating robust confidence intervals that are guaranteed to include the true parameter at least 95% of the time when the imposed restrictions on satisfied. These confidence intervals account for the fact that there is estimation error both in the treatment effects estimates and our estimates of the pre-trends.
**Sensitivity analysis**. The approach described above naturally lends itself to sensitivity analysis. That is, the researcher can report confidence intervals under different assumptions about how bad the post-treatment violation of parallel trends can be (e.g., different values of <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> or *M*.) They can also report the "breakdown value" of <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> (or *M*) for a particular conclusion -- e.g. the largest value of <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> for which the effect is still significant.
## Package installation
The package may be installed by using the function `install_github()` from the `remotes` package:
```{r,eval = FALSE}
## Installation
# Install remotes package if not installed
install.packages("remotes")
# Turn off warning-error-conversion, because the tiniest warning stops installation
Sys.setenv("R_REMOTES_NO_ERRORS_FROM_WARNINGS" = "true")
# install from github
remotes::install_github("asheshrambachan/HonestDiD")
```
## Example usage -- Medicaid expansions
As an illustration of the package, we will examine the effects of Medicaid expansions on insurance coverage using publicly-available data derived from the ACS. We first load the data and packages relevant for the analysis.
```{r, message = F, warning=F}
#Install here, dplyr, did, haven, ggplot2, fixest packages from CRAN if not yet installed
#install.packages(c("here", "dplyr", "did", "haven", "ggplot2", "fixest"))
library(here)
library(dplyr)
library(did)
library(haven)
library(ggplot2)
library(fixest)
library(HonestDiD)
df <- read_dta("https://raw.githubusercontent.com/Mixtape-Sessions/Advanced-DID/main/Exercises/Data/ehec_data.dta")
head(df,5)
```
The data is a state-level panel with information on health insurance coverage and Medicaid expansion. The variable `dins` shows the share of low-income childless adults with health insurance in the state. The variable `yexp2` gives the year that a state expanded Medicaid coverage under the Affordable Care Act, and is missing if the state never expanded.
### Estimate the baseline DiD
For simplicity, we will first focus on assessing sensitivity to violations of parallel trends in a non-staggered DiD (see below regarding methods for staggered timing). We therefore restrict the sample to the years 2015 and earlier, and drop the small number of states who are first treated in 2015. We are now left with a panel dataset where some units are first treated in 2014 and the remaining units are not treated during the sample period. We can then estimate the effects of Medicaid expansion using a canonical two-way fixed effects event-study specification,
<img src= "https://latex.codecogs.com/svg.image?Y_%7Bit%7D%20=%20%5Calpha_i%20+%20%5Clambda_t%20+%20%5Csum_%7Bs%20%5Cneq%202013%7D%201%5Bs=t%5D%20%5Ctimes%20D_i%20%5Ctimes%20%5Cbeta_s%20+%20u_%7Bit%7D%20" title = "TWFE" />
where D is 1 if a unit is first treated in 2014 and 0 otherwise.
```{r, message = F, warning=F}
df <- read_dta("https://raw.githubusercontent.com/Mixtape-Sessions/Advanced-DID/main/Exercises/Data/ehec_data.dta")
#Keep years before 2016. Drop the 2016 cohort
df_nonstaggered <- df %>% filter(year < 2016 &
(is.na(yexp2)| yexp2 != 2015) )
#Create a treatment dummy
df_nonstaggered <- df_nonstaggered %>% mutate(D = case_when( yexp2 == 2014 ~ 1,
T ~ 0))
#Run the TWFE spec
twfe_results <- fixest::feols(dins ~ i(year, D, ref = 2013) | stfips + year,
cluster = "stfips",
data = df_nonstaggered)
betahat <- summary(twfe_results)$coefficients #save the coefficients
sigma <- summary(twfe_results)$cov.scaled #save the covariance matrix
fixest::iplot(twfe_results)
```
## Sensitivity analysis using relative magnitudes restrictions
We are now ready to apply the HonestDiD package to do sensitivity analysis. Suppose we're interested in assessing the sensitivity of the estimate for 2014, the first year after treatment.
```{r, cache = TRUE, warning=FALSE}
delta_rm_results <-
HonestDiD::createSensitivityResults_relativeMagnitudes(
betahat = betahat, #coefficients
sigma = sigma, #covariance matrix
numPrePeriods = 5, #num. of pre-treatment coefs
numPostPeriods = 2, #num. of post-treatment coefs
Mbarvec = seq(0.5,2,by=0.5) #values of Mbar
)
delta_rm_results
```
The output of the previous command shows a robust confidence interval for different values of <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" />. We see that the "breakdown value" for a significant effect is <img src="https://latex.codecogs.com/svg.image?%5Cbar%7BM%7D" title="Mbar" /> = 2, meaning that the significant result is robust to allowing for violations of parallel trends up to twice as big as the max violation in the pre-treatment period.
We can also visualize the sensitivity analysis using the `createSensitivityPlot_relativeMagnitudes`. To do this, we first have to calculate the CI for the original OLS estimates using the `constructOriginalCS` command. We then pass our sensitivity analysis and the original results to the `createSensitivityPlot_relativeMagnitudes` command.
```{r}
originalResults <- HonestDiD::constructOriginalCS(betahat = betahat,
sigma = sigma,
numPrePeriods = 5,
numPostPeriods = 2)
HonestDiD::createSensitivityPlot_relativeMagnitudes(delta_rm_results, originalResults)
```
## Sensitivity Analysis Using Smoothness Restrictions
We can also do a sensitivity analysis based on smoothness restrictions -- i.e. imposing that the slope of the difference in trends changes by no more than *M* between periods.
```{r, cache=TRUE}
delta_sd_results <-
HonestDiD::createSensitivityResults(betahat = betahat,
sigma = sigma,
numPrePeriods = 5,
numPostPeriods = 2,
Mvec = seq(from = 0, to = 0.05, by =0.01))
delta_sd_results
createSensitivityPlot(delta_sd_results, originalResults)
```
We see that the breakdown value for a significant effect is *M* ≈ 0.03, meaning that we can reject a null effect unless we are willing to allow for the linear extrapolation across consecutive periods to be off by more than 0.03 percentage points.
## Sensitivity Analysis for Average Effects or Other Periods
So far we have focused on the effect for the first post-treatment period, which is the default in HonestDiD. If we are instead interested in the average over the two post-treatment periods, we can use the option `l_vec = c(0.5,0.5)`. More generally, the package accommodates inference on any scalar parameter of the form *θ* = *l*<sub>*v**e**c*</sub>′*τ*<sub>*p**o**s**t*</sub>, where *τ*<sub>*p**o**s**t*</sub> = (*τ*<sub>1</sub>,...,*τ*<sub>*T̄*</sub>)′ is the vector of dynamic treatment effects. Thus, for example, setting `l_vec = basisVector(2,numPostPeriods)` allows us to do inference on the effect for the second period after treatment.
```{r, cache=TRUE}
delta_rm_results_avg <-
HonestDiD::createSensitivityResults_relativeMagnitudes(betahat = betahat,
sigma = sigma,
numPrePeriods = 5,
numPostPeriods = 2, Mbarvec = seq(0,2,by=0.5),
l_vec = c(0.5,0.5))
originalResults_avg <- HonestDiD::constructOriginalCS(betahat = betahat,
sigma = sigma,
numPrePeriods = 5,
numPostPeriods = 2,
l_vec = c(0.5,0.5))
HonestDiD::createSensitivityPlot_relativeMagnitudes(delta_rm_results_avg, originalResults_avg)
```
## Regressions with controls
The parameters `betahat` and `sigma` should be the coefficients and variance-covariance matrix for the **event-study** plot only. Sometimes we may have regression results that contain both the event-study coefficients and coefficients on other auxilliary variables (e.g. controls). In this case, we should subset `betahat` and `sigma` to the relevant coefficients corresponding to the event-study. An example is given below:
```{r, cache=T}
# CREATE A FAKE CONTROL FOR ILLUSTRATION
set.seed(0)
df_nonstaggered$control <- rnorm(NROW(df_nonstaggered), 0, 1)
#Run the TWFE spec with the control added
# (Note that TWFEs with controls may not yield ATT under het effects; see Abadie 2005)
twfe_results <- fixest::feols(dins ~ i(year, D, ref = 2013) + control | stfips + year,
cluster = "stfips",
data = df_nonstaggered)
betahat <- summary(twfe_results)$coefficients #save the coefficients
sigma <- summary(twfe_results)$cov.scaled #save the covariance matrix
#Subset the coefficients to exclude the control, which here is in position 8
betahat <- betahat[-8]
sigma <-sigma[-8,-8]
HonestDiD::createSensitivityResults(betahat = betahat,
sigma = sigma,
numPrePeriods = 5,
numPostPeriods = 2,
Mvec = seq(from = 0, to = 0.05, by =0.01))
```
## Staggered timing
So far we have focused on a simple case without staggered timing. Fortunately, the HonestDiD approach works well with recently-introduced methods for DiD under staggered treatment timing. Below, we show how the package can be used with the [did package](https://github.com/bcallaway11/did#difference-in-differences-) implementing Callaway and Sant'Anna. (See, also, the example on the did package [website](https://github.com/pedrohcgs/CS_RR)). We are hoping to more formally integrate the did and HonestDiD packages in the future -- stay tuned!
First, we import the function Pedro Sant'Anna created for formatting did output for HonestDiD:
```{r, code = readLines("R/honest_did.R")}
```
```{r, cache = TRUE}
###
# Run the CS event-study with 'universal' base-period option
## Note that universal base period normalizes the event-time minus 1 coef to 0
cs_results <- did::att_gt(yname = "dins",
tname = "year",
idname = "stfips",
gname = "yexp2",
data = df %>% mutate(yexp2 = ifelse(is.na(yexp2), 3000, yexp2)),
control_group = "notyettreated",
base_period = "universal")
es <- did::aggte(cs_results, type = "dynamic",
min_e = -5, max_e = 5)
#Run sensitivity analysis for relative magnitudes
sensitivity_results <-
honest_did(es,
e=0,
type="relative_magnitude",
Mbarvec=seq(from = 0.5, to = 2, by = 0.5))
HonestDiD::createSensitivityPlot_relativeMagnitudes(sensitivity_results$robust_ci,
sensitivity_results$orig_ci)
```
## Additional options and resources
See the previous package [vignette] for additional examples and package options, including incorporating sign and monotonicity restrictions, and combining relative magnitudes and smoothness restrictions.
You can also view a video presentation about this paper [here](https://www.youtube.com/watch?v=6-NkiA2jN7U).
## Authors
- [Ashesh Rambachan](https://asheshrambachan.github.io/)
- [Jonathan Roth](https://www.jonathandroth.com/)
## Acknowledgements
This software package is based upon work supported by the National Science Foundation Graduate Research Fellowship under Grant DGE1745303 (Rambachan) and Grant DGE1144152 (Roth). We thank [Mauricio Cáceres Bravo](https://mcaceresb.github.io/) for his help in developing the package.