---
title: "Ayudantía 9"
date: "Semestre otoño 2021"
output: html_document
---

## Evaluación de impacto: Regresión Discontinua

Siguiendo con el ejemplo de las ayudantías anteriores, el HISP tiene como propósito central el reducir el gasto en salud de las familias pobres (`health_expenditures`). Ahora bien, pongámonos en el escenario de que la asignación de los hogares al programa depende del índice de pobreza (`poverty_index`) , el cual es un indicador que oscila entre 0 y 100, en donde los hogares más pobres son los hogares con menor puntaje en el índice.

Por lo anterior, se decide que:

* Participa (P=1) si índice_pobreza<58
* Participa (P=0) si índice_pobreza>=58

Dado que se utiliza un índice para seleccionar a aquellos que participan y quienes no, se utiliza como método de evaluación el diseño de regresión discontinua (RD). Para realizar la estimación del impacto mediante este método, cargaremos los paquetes estadísticos y nuestra base de datos (instalar paquete `rddensity` y `rdrobust`):
```{r}
library(tidyverse)
library(haven)
library(rddensity)
library(rdrobust)
rm(list=ls()) #vacía todas las variables del entorno de R
HISP<-read_dta("D:/Dropbox/Cursos/Uchile/Evaluación de Políticas/Rodrigo/Material Apoyo Evaluación  de Impacto en la Práctica/Technical Companion/evaluation.dta")
```

En primer lugar, es importante destacar que para que una RD produzca una estimación no sesgada de la puntuación límite, es importante que el índice de elegibilidad no sea manipulado en la cercanía de la puntuación límite de modo que un individuo pueda cambiar su condición de tratamiento o control. Para ello,verificaremos si existe o no manipulación del índice de pobreza:
```{r}
ronda0<-HISP %>% 
  filter(round==0 & treatment_locality==1)

#Opción 1
ggplot(ronda0, aes(x = poverty_index)) + 
  geom_histogram(aes(y = ..density..),
                 colour = 1, fill = "light green") +
  geom_density() +
  geom_vline (xintercept = 58, color="red")+
  labs(x = "Poverty Index", y = "Density")

#Opción 2
ggplot(ronda0, aes(x = poverty_index)) +
  geom_histogram(color = "yellow", boundary = 58) + 
  geom_vline(xintercept = 58, color="red") + 
  labs(x = "Poverty Index", y = "Hogares")
```

El gráfico resultante muestra la distribución de los hogares según su índice de línea de base cuando no hay manipulación. En este caso, la densidad de los hogares en torno al límite (58) es continua (o fluida). Sin perjuicio de ello, hay una pequeña diferencia visible entre la altura de las barras justo antes y justo después del umbral de 58 puntos, pero parece seguir la forma general de la distribución general. Podemos verificar si ese salto es estadísticamente significativo con una prueba de densidad de McCrary.

Análisis Gráfico de discontinuidad con Kernel e intervalos de confianza

```{r}
rdplotdensity(rdd = rddensity(ronda0$poverty_index , c = 58), X = ronda0$poverty_index)
```


Esta prueba estadística establece que:

* H0= No hay diferencia significativa entre las dos líneas de normalidad.
* H1= Hay diferencia significativa entre las dos líneas de normalidad.
```{r}
test_density <- rddensity(ronda0$poverty_index, c = 58)
summary(test_density)
```
El valor p para el tamaño de esa superposición es 0.6454, que es mucho mayor que 0.05, por lo que se puede afirmar que no existe evidencia de que haya una diferencia significativa entre las dos líneas, es decir, probablemente no existe evidencia de manipulación en el puntaje de asignación a tratamiento (pobreza).

En segundo lugar, aun cuando verificamos que no existe evidencia de manipulación en el índice de pobreza, puede que todavía persista un problema si los hogares no respetan su asignación al grupo de tratamiento o de comparación. Cuando todas las unidades cumplen con la asignación que les corresponde sobre la base de su índice de elegibilidad, se dice que el DRD es “nítido”. Verifiquemos aquello:
```{r}
ggplot(ronda0, aes(x = poverty_index, y = enrolled, color = enrolled)) +
  geom_point(size = 0.5, alpha = 0.5, 
             position = position_jitter(width = 0, height = 0.25, seed = 1234))+
  geom_vline(xintercept = 58)+ 
  labs(x = "Poverty Index", y = "Enrolled household")+
  guides(color = FALSE)
```

A partir del gráfico, se observa que todos los hogares elegibles participaron del programa y, por contra parte, todos los hogares no elegibles no participaron de él.

En tercer lugar, ahora que sabemos que este es un diseño nítido y que no hay un montón de puntajes de pruebas alrededor del umbral de 58 puntos, debemos verificar la discontinuidad en el resultado en la variable de ejecución, es decir, en el gasto en salud. Lo anterior según los siguientes códigos:
```{r}
ggplot(ronda0, aes(x = poverty_index, y = health_expenditures, color=enrolled)) +
  geom_point(size = 0.5, alpha = 0.5) + 
  geom_smooth(data = filter(ronda0, poverty_index <=58), method = "lm") +
  geom_smooth(data = filter(ronda0, poverty_index >58), method = "lm") +
  geom_vline(xintercept = 58) +
  labs(x = "Poverty Index", y = "Health Expenditures")
```

Según el gráfico, hay una clara discontinuidad en el gasto en salud, lo que nos permite intuir que la participación en el programa HISP logró dismunuir el gasto en salud de los hogares.
**no olvidar testear igualdad en la línea base de ambos grupos en variables observables disponibles en la base de datos para sujetos con y sin tratamiento incluidos en la banda creada (test de hipótesis de diferencia de medias / proporciones)**.

# ESTIMACION DE EFECTO REGRESION DISCONTINUA MANUAL
Finalmente, solo nos queda la estimación de efecto en Regresión Discontinua. Para ello, escojeremos un “ancho de banda” alrededor del punto de corte establecido por el programa. Por regla general la mejor idea es generar puntos de corte que permitan incluir en las estimaciones un número suficiente de observaciones (100 por grupo aprox.). En nuestro ejemplo escogeremos un ancho de banda de 10 unidades a ambos lados del punto de corte, y crearemos una variable de "participación" en el programa a partir del poverty_index

Crea Base filtrada (ronda = 0, sujetos treatment_locality=1, con poverty_index en la banda 58 +-10)
```{r}
en_banda <- HISP %>% 
  filter(round==1 & treatment_locality==1 & poverty_index>=48 & poverty_index<=68)
en_banda$participa[en_banda$poverty_index>=48 & en_banda$poverty_index<=58] <-1
en_banda$participa[en_banda$poverty_index> 58 & en_banda$poverty_index<=68] <-0
```

Cantidad de integrantes grupo control y tratamiento
```{r}
table(en_banda$participa)
```

Efectúa estimación de efecto mediante test de hipótesis y mediante regresión
```{r}
t.test(health_expenditures ~ participa , data = en_banda, var.equal=TRUE)

summary(lm(en_banda$health_expenditures ~ en_banda$participa))
```
t.test : el t de 36,882 y el pvalue 0 permiten rechazar la hipótesis nula de igualdad de gasto en salud entre los dos grupos

regresión : el t de la variable explicativa en_banda$participa (-36,88) y el pvalue 0 permiten rechazar la hipótesis nula de igualdad de gasto en salud entre los dos grupos

CONCLUSION: tanto el test de hipótesis de diferencia arrojan exactamente el mismo resultado: EL PROGRAMA REDUCE EL GASTO EN SALUD


# ESTIMACION EFECTO REGRESION DISCONTINUA CATTANEO
```{r}
ronda1<- HISP %>% 
  filter(round==1 & treatment_locality==1)

library(rdrobust)
rdrobust(y = ronda1$health_expenditures, x = ronda1$poverty_index, c = 58, h = 10) %>% 
  summary()
```