Rows: 1,209
Columns: 7
$ state <chr> "Rhode Island", "Tennessee", "Indiana", "Nevada", "Louisiana…
$ year <dbl> 1970, 1970, 1970, 1970, 1970, 1970, 1970, 1970, 1970, 1970, …
$ cigsale <dbl> 123.9, 99.8, 134.6, 189.5, 115.9, 108.4, 265.7, 93.8, 100.3,…
$ lnincome <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ beer <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ age15to24 <dbl> 0.1831579, 0.1780438, 0.1765159, 0.1615542, 0.1851852, 0.175…
$ retprice <dbl> 39.3, 39.9, 30.6, 38.9, 34.3, 38.4, 31.4, 37.3, 36.7, 28.8, …18 Synthetic Control in R and Stata
18.1 Synthetic control
Abadie et al. have a few papers on synthetic control (https://economics.mit.edu/files/11859, https://economics.mit.edu/files/17847). The key idea is for a case study situation. A single unit, say a state/firm/country is exposed to a shock, and we are interested in the effect of that shock. For example, Card et al studied the effect of minimum wage increase for the state of New Jersey. They used diff-in-diff, with Pennsylvania as the control state. The point of synthetic control is that maybe parallel trend assumption does not hold for the control unit. Is Pennsylvania a good control for New Jersey? Does it follow the same trend as New Jersey? Maybe not. Synthetic control is to generate a more comparable control with some combination of multiple control units, which can have a much closer parallel trend as the treated unit, than using any of the control units.
Now, how to construct this synthetic control unit? Basically it is a weighted average of multiple control units. There are two sets of weights to be determined. One set for the control units, one set to determine the importance of predictors. Abadie et al (2010) propose to choose a set of \(w_i\)’s so that the resulting synthetic control best resembles the pre-intervention values for the treated unit of predictors of the outcome variable, subject to the constraints that the weights are non-negative and sum to 1. Another set of weights \(v_h\) for predictor importance is determined such that it minimizes the mean squared prediction error (MSPE) of this synthetic control with respect to the outcome in the pre-treatment period. The intuition of these weights are to pick \(w_i\)’s to make the weighted sum of each predictor as close to the treated unit’s value of predictor as close as possible, given a set of \(v_h\)’s. Then \(v_h\)’s are determined by how important each predictor is to predict the outcome. We don’t use any information post-treatment to determine the weights. But we can use pre-treatment outcome to see the relative importance of each predictor in terms of predictor outcome. We can do this by out-of-sample validation. Basically divide the pre-treatment period into training and validation periods (that’s one reason that the pre-treatment period cannot be too short). Iterate the process of choosing \(w_i\) and \(v_h\) to achieve a minimization of MSPE in the validation period.
18.2 inference
Since synthetic control is for case study situation, we have essentially one treated unit and one synthetic control. The usual inference method would not work. Abadie et al suggested a permutation based approach. The basic idea is to permute the label of treatment; in other words, suppose one control unit is now a treated unit. Then follow the same procedure to get the synthetic control for that “treated unit”; look at the effect. Now if the treatment is truly for the treated unit only, then there would not be any effect for all the control unit. Therefore we can have a distribution of treatment effects after permutation of treatment label.
18.3 an example
The canonical synthetic control is implemented in both R and Stata (https://web.stanford.edu/~jhain/synthpage.html). We use a package called “tidysynth” here to study one of the original examples, which evaluates the impact of Proposition 99 on cigarette consumption in California.
smoking_out <-
smoking %>%
# initial the synthetic control object
synthetic_control(outcome = cigsale, # outcome
unit = state, # unit index in the panel data
time = year, # time index in the panel data
i_unit = "California", # unit where the intervention occurred
i_time = 1988, # time period when the intervention occurred
generate_placebos=T # generate placebo synthetic controls (for inference)
) %>%
# Generate the aggregate predictors used to fit the weights
# average log income, retail price of cigarettes, and proportion of the
# population between 15 and 24 years of age from 1980 - 1988
generate_predictor(time_window = 1980:1988,
ln_income = mean(lnincome, na.rm = T),
ret_price = mean(retprice, na.rm = T),
youth = mean(age15to24, na.rm = T)) %>%
# average beer consumption in the donor pool from 1984 - 1988
generate_predictor(time_window = 1984:1988,
beer_sales = mean(beer, na.rm = T)) %>%
# Lagged cigarette sales
generate_predictor(time_window = 1975,
cigsale_1975 = cigsale) %>%
generate_predictor(time_window = 1980,
cigsale_1980 = cigsale) %>%
generate_predictor(time_window = 1988,
cigsale_1988 = cigsale) %>%
# Generate the fitted weights for the synthetic control
generate_weights(optimization_window = 1970:1988, # time to use in the optimization task
margin_ipop = .02,sigf_ipop = 7,bound_ipop = 6 # optimizer options
) %>%
# Generate the synthetic control
generate_control()
smoking_out %>% plot_trends()
smoking_out %>% plot_differences()
smoking_out %>% plot_weights()
smoking_out %>% grab_balance_table()# A tibble: 7 × 4
variable California synthetic_California donor_sample
<chr> <dbl> <dbl> <dbl>
1 ln_income 10.1 9.85 9.83
2 ret_price 89.4 89.4 87.3
3 youth 0.174 0.174 0.173
4 beer_sales 24.3 24.2 23.7
5 cigsale_1975 127. 127. 137.
6 cigsale_1980 120. 120. 138.
7 cigsale_1988 90.1 91.4 114. smoking_out %>% plot_placebos()
smoking_out %>% plot_placebos(prune = FALSE)
smoking_out %>% plot_mspe_ratio()
smoking_out %>% grab_significance()# A tibble: 39 × 8
unit_name type pre_mspe post_mspe mspe_ratio rank fishers_exact_pvalue
<chr> <chr> <dbl> <dbl> <dbl> <int> <dbl>
1 California Trea… 3.17 392. 124. 1 0.0256
2 Georgia Donor 3.79 179. 47.2 2 0.0513
3 Indiana Donor 25.2 770. 30.6 3 0.0769
4 West Virginia Donor 9.52 284. 29.8 4 0.103
5 Wisconsin Donor 11.1 268. 24.1 5 0.128
6 Missouri Donor 3.03 67.8 22.4 6 0.154
7 Texas Donor 14.4 277. 19.3 7 0.179
8 South Carolina Donor 12.6 234. 18.6 8 0.205
9 Virginia Donor 9.81 96.4 9.83 9 0.231
10 Nebraska Donor 6.30 52.9 8.40 10 0.256
# ℹ 29 more rows
# ℹ 1 more variable: z_score <dbl>smoking_out %>% grab_synthetic_control()# A tibble: 31 × 3
time_unit real_y synth_y
<dbl> <dbl> <dbl>
1 1970 123 117.
2 1971 121 119.
3 1972 124. 124.
4 1973 124. 125.
5 1974 127. 127.
6 1975 127. 127.
7 1976 128 128.
8 1977 126. 126.
9 1978 126. 125.
10 1979 122. 123.
# ℹ 21 more rowssmoking_out %>% grab_synthetic_control(placebo = T)# A tibble: 1,209 × 5
.id .placebo time_unit real_y synth_y
<chr> <dbl> <dbl> <dbl> <dbl>
1 California 0 1970 123 117.
2 California 0 1971 121 119.
3 California 0 1972 124. 124.
4 California 0 1973 124. 125.
5 California 0 1974 127. 127.
6 California 0 1975 127. 127.
7 California 0 1976 128 128.
8 California 0 1977 126. 126.
9 California 0 1978 126. 125.
10 California 0 1979 122. 123.
# ℹ 1,199 more rows# A tibble: 1,482 × 50
.id .placebo .type time_unit California Alabama Arkansas Colorado
<chr> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
1 California 0 treated 1970 123 NA NA NA
2 California 0 treated 1971 121 NA NA NA
3 California 0 treated 1972 124. NA NA NA
4 California 0 treated 1973 124. NA NA NA
5 California 0 treated 1974 127. NA NA NA
6 California 0 treated 1975 127. NA NA NA
7 California 0 treated 1976 128 NA NA NA
8 California 0 treated 1977 126. NA NA NA
9 California 0 treated 1978 126. NA NA NA
10 California 0 treated 1979 122. NA NA NA
# ℹ 1,472 more rows
# ℹ 42 more variables: Connecticut <dbl>, Delaware <dbl>, Georgia <dbl>,
# Idaho <dbl>, Illinois <dbl>, Indiana <dbl>, Iowa <dbl>, Kansas <dbl>,
# Kentucky <dbl>, Louisiana <dbl>, Maine <dbl>, Minnesota <dbl>,
# Mississippi <dbl>, Missouri <dbl>, Montana <dbl>, Nebraska <dbl>,
# Nevada <dbl>, `New Hampshire` <dbl>, `New Mexico` <dbl>,
# `North Carolina` <dbl>, `North Dakota` <dbl>, Ohio <dbl>, Oklahoma <dbl>, …18.4 Generalized Synthetic Control
The idea of generalized SC is to combine both an interactive fixed effect model and synthetic control. It is basically a panel data synthetic control. It naturally takes multiple treated units and uses bootstrap to get standard errors. It generates ATT with standard errors.
Examples are here: https://yiqingxu.org/packages/gsynth/gsynth_examples.html
## for processing C++ code
require(Rcpp)
## for plotting
require(ggplot2)
require(GGally)
## for parallel computing
require(foreach)
require(future)
require(doParallel)
require(abind)
require(lfe)
library(gsynth)
## ## Syntax has been updated since v.1.2.0.
## ## Comments and suggestions -> yiqingxu@stanford.edu.
data(gsynth)
library(panelView)
## ## See bit.ly/panelview4r for more info.
## ## Report bugs -> yiqingxu@stanford.edu.
panelview(Y ~ D, data = simdata, index = c("id","time"), pre.post = TRUE)
system.time(
out <- gsynth(Y ~ D + X1 + X2, data = simdata, index = c("id","time"), force = "two-way", CV = TRUE, r = c(0, 5), se = TRUE, inference = "parametric", nboots = 100, parallel = FALSE)
) user system elapsed
1.375 0.019 1.402 print(out)Call:
gsynth(formula = Y ~ D + X1 + X2, data = simdata, index = c("id",
"time"), force = "two-way", r = c(0, 5), CV = TRUE, se = TRUE,
nboots = 100, inference = "parametric", parallel = FALSE,
vartype = "parametric")
Average Treatment Effect on the Treated:
ATT.avg S.E. CI.lower CI.upper p.value
[1,] 5.084 0.2918 4.512 5.656 0
~ by Period (including Pre-treatment Periods):
ATT S.E. CI.lower CI.upper p.value count
-19 0.05868 0.6322 -1.1804 1.29775 9.260e-01 5
-18 -0.11741 0.5355 -1.1670 0.93213 8.264e-01 5
-17 -0.89029 0.5042 -1.8784 0.09784 7.741e-02 5
-16 1.29462 0.4938 0.3267 2.26252 8.753e-03 5
-15 -1.20350 0.6416 -2.4610 0.05400 6.068e-02 5
-14 2.21992 0.6566 0.9330 3.50684 7.225e-04 5
-13 -0.18469 0.7023 -1.5612 1.19178 7.926e-01 5
-12 0.87413 0.5328 -0.1701 1.91835 1.009e-01 5
-11 -0.04824 0.5870 -1.1987 1.10220 9.345e-01 5
-10 0.84503 0.5443 -0.2218 1.91188 1.206e-01 5
-9 0.49754 0.5817 -0.6425 1.63762 3.924e-01 5
-8 -0.38001 0.6231 -1.6013 0.84130 5.420e-01 5
-7 -1.03493 0.5347 -2.0828 0.01297 5.290e-02 5
-6 -0.81532 0.6800 -2.1481 0.51743 2.305e-01 5
-5 -1.13329 0.5284 -2.1689 -0.09772 3.196e-02 5
-4 -0.17733 0.6114 -1.3757 1.02107 7.718e-01 5
-3 0.17671 0.8013 -1.3937 1.74714 8.255e-01 5
-2 -0.14376 0.5996 -1.3189 1.03139 8.105e-01 5
-1 -1.19697 0.4487 -2.0763 -0.31761 7.633e-03 5
0 1.35910 0.9498 -0.5024 3.22060 1.524e-01 5
1 0.34558 0.5937 -0.8180 1.50917 5.605e-01 5
2 0.15709 0.7348 -1.2831 1.59730 8.307e-01 5
3 3.64127 0.4977 2.6658 4.61673 2.549e-13 5
4 2.74975 0.6414 1.4925 4.00696 1.813e-05 5
5 4.45569 0.7330 3.0191 5.89233 1.211e-09 5
6 5.67231 0.8218 4.0616 7.28304 5.121e-12 5
7 8.17039 0.5637 7.0656 9.27519 0.000e+00 5
8 6.57085 0.7545 5.0921 8.04959 0.000e+00 5
9 9.11624 0.5965 7.9472 10.28534 0.000e+00 5
10 9.96205 0.4849 9.0117 10.91240 0.000e+00 5
Coefficients for the Covariates:
Coef S.E. CI.lower CI.upper p.value
X1 1.454 0.04187 1.372 1.537 0
X2 3.507 0.04288 3.423 3.591 0out$est.att ATT S.E. CI.lower CI.upper p.value count
-19 0.05868370 0.6321897 -1.1803854 1.29775282 9.260417e-01 5
-18 -0.11741144 0.5354905 -1.1669535 0.93213064 8.264478e-01 5
-17 -0.89029173 0.5041592 -1.8784256 0.09784213 7.741362e-02 5
-16 1.29461792 0.4938341 0.3267208 2.26251505 8.752727e-03 5
-15 -1.20349647 0.6415914 -2.4609926 0.05399966 6.068291e-02 5
-14 2.21991721 0.6566067 0.9329918 3.50684266 7.225056e-04 5
-13 -0.18468804 0.7022923 -1.5611556 1.19177952 7.925668e-01 5
-12 0.87413414 0.5327753 -0.1700862 1.91835448 1.008559e-01 5
-11 -0.04824208 0.5869693 -1.1986808 1.10219660 9.344969e-01 5
-10 0.84503247 0.5443195 -0.2218141 1.91187905 1.205530e-01 5
-9 0.49754000 0.5816863 -0.6425442 1.63762421 3.923625e-01 5
-8 -0.38001372 0.6231327 -1.6013314 0.84130399 5.419652e-01 5
-7 -1.03492639 0.5346505 -2.0828221 0.01296930 5.290366e-02 5
-6 -0.81531903 0.6799885 -2.1480720 0.51743399 2.305206e-01 5
-5 -1.13328625 0.5283589 -2.1688507 -0.09772183 3.195944e-02 5
-4 -0.17733004 0.6114381 -1.3757266 1.02106653 7.718000e-01 5
-3 0.17670638 0.8012576 -1.3937296 1.74714238 8.254535e-01 5
-2 -0.14376360 0.5995804 -1.3189195 1.03139231 8.105057e-01 5
-1 -1.19696676 0.4486581 -2.0763204 -0.31761309 7.633108e-03 5
0 1.35910371 0.9497582 -0.5023882 3.22059567 1.524303e-01 5
1 0.34558474 0.5936783 -0.8180034 1.50917289 5.604941e-01 5
2 0.15709192 0.7348156 -1.2831202 1.59730402 8.307154e-01 5
3 3.64127124 0.4976900 2.6658167 4.61672574 2.549072e-13 5
4 2.74975180 0.6414467 1.4925393 4.00696431 1.812680e-05 5
5 4.45569195 0.7329921 3.0190538 5.89233006 1.211063e-09 5
6 5.67231032 0.8218140 4.0615845 7.28303612 5.121015e-12 5
7 8.17038502 0.5636868 7.0655792 9.27519082 0.000000e+00 5
8 6.57084917 0.7544717 5.0921118 8.04958648 0.000000e+00 5
9 9.11624492 0.5964862 7.9471534 10.28533645 0.000000e+00 5
10 9.96204576 0.4848859 9.0116869 10.91240458 0.000000e+00 5out$est.avg ATT.avg S.E. CI.lower CI.upper p.value
[1,] 5.084123 0.2918341 4.512138 5.656107 0out$est.beta Coef S.E. CI.lower CI.upper p.value
X1 1.454447 0.04186583 1.372392 1.536503 0
X2 3.507197 0.04287653 3.423161 3.591234 0system.time(
out <- gsynth(Y ~ D + X1 + X2, data = simdata, index = c("id","time"),
force = "two-way", CV = TRUE, r = c(0, 5), se = TRUE,
inference = "parametric", nboots = 100, parallel = TRUE, cores = 4)
) user system elapsed
7.173 0.340 22.550 # cumuEff was removed in newer gsynth versions (now a wrapper for fect)
# cumu1 <- cumuEff(out, cumu = TRUE, id = NULL, period = c(0,5))
# cumu1$est.catt
plot(out)
Systematic treatment: Julia.