Taller 8: Diff-in-Diff
2023-11-01
Resumen de contenido
Nuestra meta
Entender que són las interacciones y su en las estimaciones de MCO.
Entender Diff-in-Diff de forma intuitiva.
Aprender a estimar Diff-in-Diff con función lm.
Librerías
Usaremos las siguientes librerías:
Interacciones en regresión lineal
Botones y Deslizantes
Esperanza de vida e ingreso por habitante
\[ LifExp_i=\alpha+ \beta_1 log(GPDpc_i)+\beta_2Europa+\epsilon_i \]
data<-gapminder %>%
filter(year==2007) %>%
mutate(europa=ifelse(continent=="Europe",1,0))
model1<-lm(lifeExp~log(gdpPercap)+europa, data=data)
tidy(model1)# A tibble: 3 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 6.51 4.35 1.50 1.36e- 1
2 log(gdpPercap) 6.99 0.520 13.4 6.97e-27
3 europa 1.35 1.72 0.784 4.34e- 1GPD per capita es continua (y log): Cada 1% adicional de inrgeso per capita se asocia a \(\frac{\beta_1}{100}\)% años extra de esperanza de vida.
Europa es binaria: Estar en Europa se asocia con un aumento de \(\beta_2\) años en la esperanza de vida.
Sin control binario
\[ LifExp_i=\alpha+ \beta_1 log(GPDpc_i)+\epsilon_i \]
data %>%
ggplot(aes(x=gdpPercap,y=lifeExp))+
geom_point(aes(color=(continent=="Europe")))+
scale_color_manual(values = c("gray60","purple"))+
# Grafico modelo sencillo (lm of lifexp~log_gdppc)
geom_smooth(method="lm",se = F,color="black")+
scale_x_log10()+
theme_minimal()+
labs(y="Esperanza de vida al nacer",
x="PIB per capita (log scale)",
subtitle = "Sin 'Europa' en regresion",
color="¿Europa?")Sin control binario
Control binario: Cambio de intercepto
broom::augment_columns(model1,newdata = data) %>%
ggplot(aes(x=gdpPercap,y=lifeExp))+
geom_point(aes(color=(continent=="Europe")))+
scale_color_manual(values = c("gray60","purple"))+
# our model
geom_line(aes(y=.fitted,group=continent))+
scale_x_log10()+
theme_minimal()+
labs(y="Esperanza de vida al nacer",
x="PIB per capita (log scale)",
subtitle = "Con 'Europa' en regresion",
color="¿Europa?")Control binario: Cambio de intercepto
Término de interacción
\[ LifExp_i=\alpha+ \beta_1 log(GPDpc_i)+\beta_2Europa+\beta_3 (log(GPDpc_i)\times Europa)+\epsilon_i \]
# A tibble: 4 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 5.22 4.45 1.17 2.43e- 1
2 log(gdpPercap) 7.15 0.534 13.4 1.03e-26
3 europa 30.2 22.9 1.32 1.89e- 1
4 log(gdpPercap):europa -2.92 2.31 -1.27 2.08e- 1gdpPercap x europa es una interacción: En Europa, cada incremento de 1% en ingreso per capita se asocia a \(\frac{\beta_1+\beta_3}{100}\)% años extra de esperanza de vida.
Sin logs, la asociacion fuera de \(\beta_1+\beta_3\) por cada dolar adicional.
Imágen de cambio de pendiente para grupo específico
broom::augment_columns(model2,newdata = data) %>%
ggplot(aes(x=gdpPercap,y=lifeExp))+
geom_point(aes(color=(continent=="Europe")))+
scale_color_manual(values = c("gray60","purple"))+
# General lm of lifexp~gdppc
geom_line(aes(y=.fitted,group=continent))+
scale_x_log10()+
theme_minimal()+
labs(y="Esperanza de vida al nacer",
x="PIB per capita (log scale)",
subtitle = "Con interacción 'log(PIBpc)xEuropa' en regresion",
color="¿Europa?")Imágen de cambio de pendiente para grupo específico
¿Qué pasaría si corremos esto?
Imaginemos que olvidé poner las variables por separado en la regresión.
# A tibble: 4 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 5.22 4.45 1.17 2.43e- 1
2 log(gdpPercap) 7.15 0.534 13.4 1.03e-26
3 europa 30.2 22.9 1.32 1.89e- 1
4 log(gdpPercap):europa -2.92 2.31 -1.27 2.08e- 1Nada. Todo pelota, siga!
¿Qué pasaría si la variable tiene varios niveles?
# A tibble: 10 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 22.9 6.22 3.68 0.000337
2 log(gdpPercap) 4.26 0.824 5.17 0.000000836
3 continentAmericas 10.2 15.7 0.648 0.518
4 continentAsia 2.74 9.81 0.280 0.780
5 continentEurope 12.5 19.9 0.629 0.530
6 continentOceania 23.9 279. 0.0855 0.932
7 log(gdpPercap):continentAmericas 0.233 1.79 0.130 0.897
8 log(gdpPercap):continentAsia 0.896 1.19 0.753 0.453
9 log(gdpPercap):continentEurope -0.0340 2.06 -0.0165 0.987
10 log(gdpPercap):continentOceania -0.965 27.2 -0.0355 0.972 Tenemos una pendiente única para cada continente.
¿Qué pasaría si la variable tiene varios niveles?
broom::augment_columns(model4,newdata = data) %>%
ggplot(aes(x=gdpPercap,y=lifeExp))+
geom_point(aes(color=continent))+
scale_color_manual(values = c("gray60","purple","darkgreen","blue","red"))+
geom_line(aes(y=.fitted,group=continent, color=continent))+
scale_x_log10()+
theme_minimal()+
labs(y="Esperanza de vida al nacer",
x="PIB per capita (log scale)",
subtitle = "Con interacción 'log(PIBpc)xContinente' en regresion",
color="Continente")Tenemos una pendiente única para cada continente.
¿Qué pasaría si la variable tiene varios niveles?
Idea general de las interacciones
Los términos de interacción muestran el cambio adicional (promedio) que observamos al cumplirse dos condiciones simultáneamente.
Vimos:
Efecto del GDP per cápita.
Efecto de Europa.
Efecto del GDP per cápita si el país está en Europa (o en cualquier otro continente en la interacción).
Meyer, Viscusi y Durbin (1995)
“Workers’ Compensation and Injury Duration: Evidence from a Natural Experiment” Una aplicación del método de Diferencias en Diferencias (Diff-in-Diff).
Meyer, Viscusi y Durbin (1995)
En 1980 las autoridades de Kentuky aumentaron en 60% la compensación por reposo percibido por los trabajadores de ingresos altos.
Los autores encuentran datos a nivel de trabajador. Pueden ver:
quienes eran de altos ingresos (beneficiarios de la politica).
quienes quienes estaban de reposo cuando el cambio ocurrió.
El paper busca ver si cuando aumenta la compensación, los trabajadores de altos ingresos demoran mas en volver a trabajar. Oh sorpresa: Si demoran más.
Lectura y limpieza de datos
Leo los datos, almaceno la versión cruda en un objeto. Luego restrinjo las observaciones a Kentucky y cambio el nombre de algunas variables
Exploración de datos (Control vs. Tratamiento)
Veamos diferencias en duración de desempleo entre grupos de control y tratamiento:
La mayoría tiene duración entre 0 y 16, con algunas exepciones destacables. Mejor usemos escala log.
Exploración de datos (Control vs. Tratamiento)
Exploración de datos (Control vs. Tratamiento)
Usando logs, y resaltando el promedio de cada grupo:
Veo una diferencia, pero ¿cómo saber si no viene de antes?
Exploración de datos (Control vs. Tratamiento)
Exploración de datos (Before & After)
Hacemos lo mismo para ver la duración del desempleo antes y despues de la medida.
ggplot(data = injury, mapping = aes(x = log_duration)) +
geom_histogram(binwidth = 0.5, color = "white", boundary = 0) +
geom_vline(data=injury %>%
group_by(after_1980) %>%
summarise(mean=mean(log_duration)),
aes(xintercept = mean), color="red", linetype="dashed")+
facet_wrap(vars(ifelse(after_1980==1,"b. After 1980","a. Before 1980")))Exploración de datos (Before & After)
Exploración de datos (Tratamiento vs Control, Before vs After)
No pude encontrar indicios de nada. ¿Que pasa si busco la distribución en todos los grupos?
(plot<-ggplot(injury, aes(x = ifelse(highearn==0,"a. Ingresos bajos","b. Ingresos Altos"),
y = log_duration)) +
# calculate the mean and 95% confidence interval:
stat_summary(geom = "pointrange", size = 1, color = "red",
fun.data = "mean_se", fun.args = list(mult = 1.96)) +
facet_wrap(vars(ifelse(after_1980==1,"b. After 1980","a. Before 1980")))+
labs(x=NULL))¿Ustedes ven lo mismo que yo?
Exploración de datos (Tratamiento vs Control, Before vs After)
Diseño Diff in Diff
Diff-in-Diff con regresión 00
Estamos buscando (D − C) − (B − A). En esta regresión implica encontrar \(\delta\).
\[ log(duration)=\alpha+\beta highearn +\gamma after_1980 + \delta(highearn\times after_1980)+epsilon \]
model_small <- lm(log_duration ~ highearn + after_1980 + highearn * after_1980,data = injury)
tidy(model_small)# A tibble: 4 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 1.13 0.0307 36.6 1.62e-263
2 highearn 0.256 0.0474 5.41 6.72e- 8
3 after_1980 0.00766 0.0447 0.171 8.64e- 1
4 highearn:after_1980 0.191 0.0685 2.78 5.42e- 3Recibir el tratamiento causaría un incremento promedio del 19% en la duración del desempleo.
Comprobemos que \(\delta\) nos da lo que buscamos
Calculo el promedio de la duración para los cuatro grupos.
diffs <- injury %>%
group_by(after_1980, highearn) %>%
summarize(mean_duration = mean(log_duration),
# Calculate average with regular duration too, just for fun
mean_duration_for_humans = mean(duration))
diffs# A tibble: 4 × 4
# Groups: after_1980 [2]
after_1980 highearn mean_duration mean_duration_for_humans
<dbl> <dbl> <dbl> <dbl>
1 0 0 1.13 6.27
2 0 1 1.38 11.2
3 1 0 1.13 7.04
4 1 1 1.58 12.9 Comprobemos que \(\delta\) nos da lo que buscamos
Aplico la formula (C-D) - (B-A)
# promedios de d
A_before_control <- diffs %>%
filter(after_1980 == 0, highearn == 0) %>%
pull(mean_duration)
C_before_treatment <- diffs %>%
filter(after_1980 == 0, highearn == 1) %>%
pull(mean_duration)
B_after_control <- diffs %>%
filter(after_1980 == 1, highearn == 0) %>%
pull(mean_duration)
D_after_treatment <- diffs %>%
filter(after_1980 == 1, highearn == 1) %>%
pull(mean_duration)
# (D-C) - (B-A)
diff_in_treatment <- D_after_treatment - C_before_treatment
diff_in_control <- B_after_control - A_before_control
print(diff_in_treatment-diff_in_control)[1] 0.1906Diff-in-Diff con regresión, añadiendo controles
\(\epsilon\) pudiera correlacionarse con el tratamiento de contener variables como género, industria, y edad del trabajador, violando así el supuesto de exogeneidad.
Lo bueno de la regresión es que nos permite remover estas variables del término de error al añadirlas como controles. Con ello haremos nuestra estimación mas creíble.
Preparamos las variables para la regresión:
Diff-in-Diff con regresión, añadiendo controles
Añadimos género, estado civil, edad, hospitalización, código de industria, típo de lesión, y salario antes de estar desempleado.
model_big <- lm(log_duration ~ highearn + after_1980 + highearn * after_1980 + male + married + age + hosp + indust + injtype + lprewage,
data = injury_fixed)
tidy(model_big) # A tibble: 18 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) -1.53 0.422 -3.62 2.98e- 4
2 highearn -0.152 0.0891 -1.70 8.86e- 2
3 after_1980 0.0495 0.0413 1.20 2.31e- 1
4 male -0.0843 0.0423 -1.99 4.64e- 2
5 married 0.0567 0.0373 1.52 1.29e- 1
6 age 0.00651 0.00134 4.86 1.19e- 6
7 hosp 1.13 0.0370 30.5 5.20e-189
8 indust2 0.184 0.0541 3.40 6.87e- 4
9 indust3 0.163 0.0379 4.32 1.60e- 5
10 injtype2 0.935 0.144 6.51 8.29e- 11
11 injtype3 0.635 0.0854 7.44 1.19e- 13
12 injtype4 0.555 0.0928 5.97 2.49e- 9
13 injtype5 0.641 0.0854 7.51 7.15e- 14
14 injtype6 0.615 0.0863 7.13 1.17e- 12
15 injtype7 0.991 0.191 5.20 2.03e- 7
16 injtype8 0.434 0.119 3.65 2.64e- 4
17 lprewage 0.284 0.0801 3.55 3.83e- 4
18 highearn:after_1980 0.169 0.0640 2.64 8.38e- 3Controlando por las variables mencionadas, recibir el tratamiento causaría un incremento promedio del 17% en la duración del desempleo.