library(tidyverse) # ggplot(), %>%, mutate(), and friends
library(broom) # Convert models to data frames
library(rdrobust) # For robust nonparametric regression discontinuity
library(rddensity) # For nonparametric regression discontinuity density tests
library(modelsummary) # Create side-by-side regression tablesDiscontinuidad en Regresión (Ejemplo)
La regresión discontinua se utiliza generalmente cuando la asignación al tratamiento se basa en un punto de corte de una variable continua. Al igual que Diff-in-Diff es un método cuasi-experimental.
Para este ejercicio utilizaremos el archivo denominado tutoring_program.csv
Antecedentes del programa
En este ejemplo hipotético, los estudiantes presentan un examen de admisión al inicio del año escolar. Quienes obtienen una puntuación de 70 o menos se inscriben automáticamente en un programa de tutoría gratuito y reciben apoyo durante todo el año. Al final del año escolar, los estudiantes presentan un examen final (con un máximo de 100 puntos) para evaluar su aprendizaje general. Recuerden que este es un ejemplo hipotético y que este tipo de exámenes no existen, pero simplemente tómenlo como referencia.
Tienes un conjunto de datos con cuatro columnas:
id:El ID del estudianteentrance_exam:La puntuación del examen de ingreso del estudiante (sobre 100)exit_exam:Puntuación del examen de salida del estudiante (sobre 100)tutoring:Una variable indicadora que muestra si el estudiante estuvo inscrito en el programa de tutoría.
Cargar y limpiar datos
Primero, descarguemos el conjunto de datos (si aún no lo has hecho), póngalo en una carpeta llamada datay carguelo:
# Load the data.
# It'd be a good idea to click on the "tutoring" object in the Environment
# panel in RStudio to see what the data looks like after you load it
tutoring <- read_csv("data/tutoring_program.csv")Paso 1: Determinar si el proceso de asignación del tratamiento se basa en reglas
Para unirse al programa de tutorías, los estudiantes deben obtener una puntuación de 70 puntos o menos en el examen de admisión. Quienes obtengan una puntuación superior a 70 no podrán participar. Dado que tenemos una regla clara de 70 puntos, podemos asumir que el proceso de participación en el programa de tutorías se rige por normas.
Paso 2: Determinar si el diseño es difuso o claro
Dado que sabemos que el programa se aplicó con base en una regla, ahora queremos determinar cuán estrictamente se aplicó. El umbral era de 70 puntos en la prueba: ¿quienes obtuvieron 68 puntos se evadieron de los trámites burocráticos y no participaron, o quienes obtuvieron 73 puntos se colaron en el programa? La forma más fácil de comprobarlo es con un gráfico, y podemos obtener cifras exactas con una tabla.
ggplot(tutoring, aes(x = entrance_exam, y = tutoring, color = tutoring)) +
# Make points small and semi-transparent since there are lots of them
geom_point(size = 0.5, alpha = 0.5,
position = position_jitter(width = 0, height = 0.25, seed = 1234)) +
# Add vertical line
geom_vline(xintercept = 70) +
# Add labels
labs(x = "Entrance exam score", y = "Participated in tutoring program") +
# Turn off the color legend, since it's redundant
guides(color = "none")
Esto parece bastante claro: no parece que las personas con una puntuación inferior a 70 hayan participado en el programa. Podemos comprobarlo con una tabla. No hay personas donde ” entrance_exames mayor que 70” y tutoringsea falso, ni personas donde ” entrance_exames menor que 70” y tutoringsea verdadero.
tutoring %>%
group_by(tutoring, entrance_exam <= 70) %>%
summarize(count = n())`summarise()` has grouped output by 'tutoring'. You can override using the
`.groups` argument.# A tibble: 2 × 3
# Groups: tutoring [2]
tutoring `entrance_exam <= 70` count
<lgl> <lgl> <int>
1 FALSE FALSE 759
2 TRUE TRUE 241Paso 3: Verificar si hay discontinuidad en la variable de ejecución alrededor del punto de corte
A continuación, necesitamos ver si hubo alguna manipulación en la variable móvil; quizá mucha gente se agrupó alrededor de 70 debido a la calificación del examen (por ejemplo, los estudiantes querían ingresar al programa, así que obtuvieron malos resultados a propósito para obtener menos de 70). Podemos hacer esto de dos maneras diferentes. Primero, crearemos un histograma de la variable móvil (puntuaciones del examen) y veremos si hay saltos significativos alrededor del umbral:
ggplot(tutoring, aes(x = entrance_exam, fill = tutoring)) +
geom_histogram(binwidth = 2, color = "white", boundary = 70) +
geom_vline(xintercept = 70) +
labs(x = "Entrance exam score", y = "Count", fill = "In program")
Aquí no parece que haya un salto alrededor del corte—hay una pequeña diferencia visible entre la altura de las barras justo antes y justo después del umbral de 70 puntos, pero parece seguir la forma general de la distribución total. Podemos verificar si ese salto es estadísticamente significativo con una prueba de densidad de McCrary. Esta pone los datos en contenedores como un histograma y luego traza los promedios e intervalos de confianza de esos contenedores. Si los intervalos de confianza de las líneas de densidad no se superponen, entonces es probable que haya algo sistemáticamente incorrecto con cómo se calificó la prueba (es decir, demasiadas personas obtuvieron 69 frente a 71). Si los intervalos de confianza se superponen, no hay ninguna diferencia significativa alrededor del umbral y estamos bien.
Primero, realizamos rddensity()la prueba estadística. Debes introducir dos datos: la variable de ejecución y el valor de corte.
# Notice this tutoring$entrance_exam syntax. This is one way for R to access
# columns in data frames---it means "use the entrance_exam column in the
# tutoring data frame". The general syntax for it is data_frame$column_name
test_density <- rddensity(tutoring$entrance_exam, c = 70)
summary(test_density)Manipulation testing using local polynomial density estimation.
Number of obs = 1000
Model = unrestricted
Kernel = triangular
BW method = estimated
VCE method = jackknife
c = 70 Left of c Right of c
Number of obs 237 763
Eff. Number of obs 208 577
Order est. (p) 2 2
Order bias (q) 3 3
BW est. (h) 22.444 19.966
Method T P > |T|
Robust -0.5521 0.5809 Luego podemos graficar esa prueba de densidad:
# The syntax for rdplotdensity is kinda wonky here. You have to feed it the
# rddensity() test, and then you have to specify x, which is your running
# variable (again!). The type argument tells the plot to show both points and
# lines---without it, it'll only show lines.
#
# Finally, notice how I assigned the output of rdplotdensity to a variable named
# plot_density_test. In theory, this should make it show nothing---all the
# output should go to that object. Because of a bug in rdplotdensity, though, it
# will show a plot automatically even if assigning it to a variable. If we don't
# assign it to a variable you'll see two identical plots when knitting, which is
# annoying. So we save it as a variable to hide the output, but get the output
# for a single plot anyway. Ugh.
plot_density_test <- rdplotdensity(rdd = test_density,
X = tutoring$entrance_exam,
type = "both") # This adds both points and lines
Hay muchos resultados, pero lo que nos interesa es la línea que empieza con “Robusto”, que muestra la prueba t para la diferencia entre los dos puntos a ambos lados del punto de corte en el gráfico. Observe que los intervalos de confianza se superponen considerablemente. El valor p para la magnitud de dicha superposición es 0,5809, mucho mayor que 0,05, por lo que no tenemos evidencia sólida de que exista una diferencia significativa entre las dos líneas. Con base en este gráfico y el estadístico t, probablemente podamos afirmar con seguridad que no hay manipulación ni agrupamiento.
Paso 4: Verificar si hay discontinuidad en los resultados a lo largo de la variable en ejecución
Ahora que sabemos que este es un diseño preciso y que no hay acumulación de puntajes en torno al umbral de 70 puntos, finalmente podemos ver si existe una discontinuidad en los puntajes finales según la participación en el programa de tutoría. Grafica la variable móvil en el eje x, la variable de resultado en el eje y, y colorea los puntos según si participaron en el programa.
ggplot(tutoring, aes(x = entrance_exam, y = exit_exam, color = tutoring)) +
geom_point(size = 0.5, alpha = 0.5) +
# Add a line based on a linear model for the people scoring 70 or less
geom_smooth(data = filter(tutoring, entrance_exam <= 70), method = "lm") +
# Add a line based on a linear model for the people scoring more than 70
geom_smooth(data = filter(tutoring, entrance_exam > 70), method = "lm") +
geom_vline(xintercept = 70) +
labs(x = "Entrance exam score", y = "Exit exam score", color = "Used tutoring")
Según este gráfico, ¡hay una clara discontinuidad! Parece que la participación en el programa de tutorías mejoró las calificaciones de los exámenes finales.
Paso 5: Mide el tamaño del efecto
Hay una discontinuidad, pero ¿de qué magnitud es? ¿Es estadísticamente significativa?
Podemos comprobar el tamaño de dos maneras: paramétricamente (es decir, utilizando lm()parámetros y coeficientes específicos) y no paramétricamente (es decir, sin usar lm()ningún tipo de línea recta, sino dibujando líneas que se ajusten a los datos con mayor precisión). Lo haremos de ambas maneras.
Estimación paramétrica
Primero, lo haremos de forma paramétrica mediante regresión lineal. Aquí queremos explicar la variación en las puntuaciones finales según la puntuación del examen de admisión y la participación en el programa de tutorías:
\[ \text{Exit exam} = \beta_0 + \beta_1 \text{Entrance exam score}_\text{centered} + \beta_2 \text{Tutoring program} + \epsilon \]
Para facilitar la interpretación de los coeficientes, podemos centrar la columna del examen de admisión para que, en lugar de mostrar la puntuación real, muestre cuántos puntos por encima o por debajo de 70 obtuvo el estudiante. De esta forma, podemos usar el coeficiente para tutoringel efecto causal.
tutoring_centered <- tutoring %>%
mutate(entrance_centered = entrance_exam - 70)
model_simple <- lm(exit_exam ~ entrance_centered + tutoring,
data = tutoring_centered)
tidy(model_simple)# A tibble: 3 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 59.4 0.442 134. 0
2 entrance_centered 0.510 0.0269 18.9 1.40e-68
3 tutoringTRUE 10.8 0.800 13.5 3.12e-38Esto es lo que significan estos coeficientes:
\(\beta_0\) es la intersección. Dado que centramos las puntuaciones del examen de admisión, muestra la puntuación promedio del examen de egreso en el umbral de 70 puntos. Las personas que obtuvieron 70.001 puntos en el examen de admisión obtuvieron un promedio de 59.4 puntos en el examen de egreso. (Otra forma de verlo es que muestra la puntuación prevista del examen de egreso cuando
entrance_centeredes 0 (es decir, 70) y cuandotutoringesFALSE).\(\beta_1\) es el coeficiente de [número
entrance_centered]. Por cada punto por encima de 70 que las personas obtienen en el examen de admisión, obtienen 0,51 puntos más en el examen de egreso. Realmente no nos importa mucho esta cifra.\(\beta_2\) es el coeficiente del programa de tutoría, y es el que más nos importa . Es el desplazamiento en el intercepto cuando
tutoringes verdadero, o la diferencia entre las puntuaciones en el umbral. Participar en el programa de tutoría aumenta la puntuación del examen de egreso en 10,8 puntos.
Una ventaja de usar un enfoque paramétrico es que permite incluir otras covariables, como las demográficas. También se puede usar la regresión polinómica e incluir términos como ” entrance_centered²o entrance_centered³” o “par” entrance_centered⁴para que la línea se ajuste lo más posible a los datos.
Aquí ajustamos el modelo a todos los datos, pero en la práctica, nos importan más las observaciones cercanas al umbral, lo que significa que este modelo es incorrecto . Las puntuaciones muy altas o muy bajas no deberían influir en el tamaño del efecto, ya que solo nos importan las personas con puntuaciones ligeramente por debajo y ligeramente por encima de 70.
Podemos ajustar el mismo modelo pero restringirlo a personas dentro de una ventana o ancho de banda más pequeño, como ±10 puntos o ±5 puntos:
model_bw_10 <- lm(exit_exam ~ entrance_centered + tutoring,
data = filter(tutoring_centered,
entrance_centered >= -10 &
entrance_centered <= 10))
tidy(model_bw_10)# A tibble: 3 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 60.4 0.752 80.3 2.99e-249
2 entrance_centered 0.388 0.114 3.40 7.45e- 4
3 tutoringTRUE 9.27 1.31 7.09 6.27e- 12model_bw_5 <- lm(exit_exam ~ entrance_centered + tutoring,
data = filter(tutoring_centered,
entrance_centered >= -5 &
entrance_centered <= 5))
tidy(model_bw_5)# A tibble: 3 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 60.6 1.12 54.3 4.78e-118
2 entrance_centered 0.380 0.331 1.15 2.53e- 1
3 tutoringTRUE 9.12 1.91 4.77 3.66e- 6Podemos comparar todos estos modelos simultáneamente con modelsummary():
modelsummary(list("Full data" = model_simple,
"Bandwidth = 10" = model_bw_10,
"Bandwidth = 5" = model_bw_5))| Datos completos | Banda = 10 | Banda = 5 | |
|---|---|---|---|
| + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001 | |||
| (Intercept) | 59.411*** | 60.377*** | 60.631*** |
| (0.442) | (0.752) | (1.117) | |
| entrance_centered | 0.510*** | 0.388*** | 0.380 |
| (0.027) | (0.114) | (0.331) | |
| tutoringTRUE | 10.800*** | 9.273*** | 9.122*** |
| (0.800) | (1.309) | (1.912) | |
| Num.Obs. | 1000 | 404 | 194 |
| R2 | 0.268 | 0.162 | 0.222 |
| R2 Adj. | 0.267 | 0.158 | 0.214 |
El efecto tutoringvaría mucho entre estos distintos modelos, de 9.1 a 10.8. ¿Cuál es el correcto? No lo sé. Definitivamente no es el de datos completos.
Observe también cuánto varía el tamaño de la muestra (N) en los modelos. A medida que se reduce el ancho de banda, se observan cada vez menos observaciones.
ggplot(tutoring, aes(x = entrance_exam, y = exit_exam, color = tutoring)) +
geom_point(size = 0.5, alpha = 0.2) +
# Add lines for the full model (model_simple)
geom_smooth(data = filter(tutoring, entrance_exam <= 70),
method = "lm", se = FALSE, linetype = "dotted", linewidth = 1) +
geom_smooth(data = filter(tutoring, entrance_exam > 70),
method = "lm", se = FALSE, linetype = "dotted", linewidth = 1) +
# Add lines for bandwidth = 10
geom_smooth(data = filter(tutoring, entrance_exam <= 70, entrance_exam >= 60),
method = "lm", se = FALSE, linetype = "dashed", linewidth = 1) +
geom_smooth(data = filter(tutoring, entrance_exam > 70, entrance_exam <= 80),
method = "lm", se = FALSE, linetype = "dashed", linewidth = 1) +
# Add lines for bandwidth = 5
geom_smooth(data = filter(tutoring, entrance_exam <= 70, entrance_exam >= 65),
method = "lm", se = FALSE, linewidth = 2) +
geom_smooth(data = filter(tutoring, entrance_exam > 70, entrance_exam <= 75),
method = "lm", se = FALSE, linewidth = 2) +
geom_vline(xintercept = 70) +
# Zoom in
coord_cartesian(xlim = c(50, 90), ylim = c(55, 75)) +
labs(x = "Entrance exam score", y = "Exit exam score", color = "Used tutoring")
Estimación no paramétrica
En lugar de usar regresión lineal para medir la magnitud de la discontinuidad, podemos usar métodos no paramétricos. En esencia, esto significa que R no intentará ajustar una línea recta a los datos, sino que se curvará alrededor de los puntos e intentará ajustar todo de la forma más uniforme posible.
La rdrobust()función facilita enormemente la medición de la brecha en el punto de corte mediante estimación no paramétrica. Aquí está la versión más sencilla:
# Notice how we have to use the tutoring$exit_exam syntax here. Also make sure
# you set the cutoff with c
rdrobust(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70) %>%
summary()Warning in rdrobust(y = tutoring$exit_exam, x = tutoring$entrance_exam, : Mass
points detected in the running variable.Sharp RD estimates using local polynomial regression.
Number of Obs. 1000
BW type mserd
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Eff. Number of Obs. 144 256
Order est. (p) 1 1
Order bias (q) 2 2
BW est. (h) 9.969 9.969
BW bias (b) 14.661 14.661
rho (h/b) 0.680 0.680
Unique Obs. 155 262
=====================================================================
Point Robust Inference
Estimate z P>|z| [ 95% C.I. ]
---------------------------------------------------------------------
RD Effect -8.578 -4.352 0.000 [-12.101 , -4.587]
=====================================================================Hay algunos datos importantes que se deben tener en cuenta en este resultado:
Lo que más le importa es el tamaño real del efecto. Este es el coeficiente en la tabla inferior, indicado con el método “Convencional”. Aquí es -8,578, lo que significa que el programa de tutoría provoca un cambio de 8 puntos en las calificaciones del examen de egreso. La tabla inferior también incluye errores estándar, valores p e intervalos de confianza para el coeficiente, tanto estimaciones normales (convencionales) como estimaciones robustas (robustas). Según ambos tipos de estimaciones, este aumento de 8 puntos es estadísticamente significativo (p < 0,001; el intervalo de confianza del 95 % definitivamente no incluye el 0).
Es importante destacar que el coeficiente aquí es negativo (-8,578), mientras que nuestras estimaciones paramétricas previas fueron todas positivas. Esto no significa que el programa provoque una caída en las puntuaciones de las pruebas. Ese valor negativo es simplemente un efecto secundario de cómo
rdrobust()se mide la brecha. Analiza el valor del grupo de tratamiento justo antes del umbral y luego muestra que las puntuaciones de las pruebas disminuyen al pasar al grupo de control (mientras que antes hablamos de la situación opuesta: las puntuaciones de las pruebas aumentan al pasar del grupo de control al grupo de tratamiento). No hay forma de invertir el signo dentro de [número]rdrobust(), por lo que es necesario observar el gráfico y ver qué está haciendo realmente la brecha.El modelo utilizó un ancho de banda de 9,969 (
BW est. (h)en la salida), lo que significa que solo consideró personas con puntajes de prueba de 70 ± 10. Decidió este ancho de banda automáticamente, pero puedes cambiarlo a lo que desees.El modelo utilizó un kernel triangular. Esta es la parte más compleja del modelo: el kernel decide la ponderación de las observaciones cercanas al punto de corte. Puntuaciones como 69,99 o 70,01 están muy cerca de 70, por lo que reciben la mayor ponderación. Puntuaciones como 67 o 73 están un poco más alejadas, por lo que su importancia es menor. Puntuaciones como 64 o 76 son aún menos importantes, por lo que su importancia es aún menor. También se pueden usar diferentes kernels, y [Wikipedia tiene un gráfico que muestra las formas de estos diferentes kernels](https://en.wikipedia.org/wiki/Kernel_(statistics)#Kernel_functions_in_common_use) y cómo asignan diferentes ponderaciones a las observaciones a diferentes distancias.
Podemos graficar este modelo no paramétrico con rdplot().
rdplot(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70)[1] "Mass points detected in the running variable."
¡Hay un salto de 8,5 puntos en el 70!
Debemos tener en cuenta que los puntos aquí no son realmente las observaciones del conjunto de datos. La rdplot()función crea intervalos de puntos (como un histograma) y luego muestra el resultado promedio dentro de cada intervalo. Puedes controlar cuántos intervalos se utilizan en el eje x con los argumentos nbins``` o ` ` (ejecútalo en tu consola para ver todas las opciones posibles para graficar discontinuidades no paramétricas).binselectrdplot()?rdplot .
De forma predeterminada, rdrobust()se elige automáticamente el tamaño del ancho de banda basándose en algoritmos sofisticados desarrollados por economistas. Puedes usar [nombre del método] rdbwselect()para ver cuál es ese ancho de banda y, si incluyes el all = TRUEargumento, puedes ver varios anchos de banda potenciales diferentes basados en diversos algoritmos.
# This says to use the mserd version, which according to the help file for
# rdbwselect means "the mean squared error-optimal bandwidth selector for RD
# treatment effects". Sounds great.
rdbwselect(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70) %>%
summary()Call: rdbwselect
Number of Obs. 1000
BW type mserd
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Order est. (p) 1 1
Order bias (q) 2 2
Unique Obs. 155 262
=======================================================
BW est. (h) BW bias (b)
Left of c Right of c Left of c Right of c
=======================================================
mserd 9.969 9.969 14.661 14.661
=======================================================¿Cuáles son los diferentes anchos de banda posibles que podríamos utilizar?
rdbwselect(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70, all = TRUE) %>%
summary()Call: rdbwselect
Number of Obs. 1000
BW type All
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Order est. (p) 1 1
Order bias (q) 2 2
Unique Obs. 155 262
=======================================================
BW est. (h) BW bias (b)
Left of c Right of c Left of c Right of c
=======================================================
mserd 9.969 9.969 14.661 14.661
msetwo 11.521 10.054 17.067 14.907
msesum 12.044 12.044 17.631 17.631
msecomb1 9.969 9.969 14.661 14.661
msecomb2 11.521 10.054 17.067 14.907
cerrd 7.058 7.058 14.661 14.661
certwo 8.156 7.118 17.067 14.907
cersum 8.526 8.526 17.631 17.631
cercomb1 7.058 7.058 14.661 14.661
cercomb2 8.156 7.118 17.067 14.907
=======================================================Hay muchos aquí. El mejor fue mserd±9,69. Algunos dicen ±7, otros ±12, y otros son asimétricos y dicen -11,5 antes y +10 después. Prueba con varios anchos de banda diferentes como parte de tu análisis de sensibilidad y comprueba si el tamaño del efecto cambia sustancialmente.
Otro enfoque común para el análisis de sensibilidad consiste en utilizar el ancho de banda ideal, el doble y la mitad (en nuestro caso, 10, 20 y 5) y comprobar si la estimación cambia sustancialmente. Utilice el hargumento para especificar su propio ancho de banda.
rdrobust(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70, h = 9.969) %>%
summary()Sharp RD estimates using local polynomial regression.
Number of Obs. 1000
BW type Manual
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Eff. Number of Obs. 144 256
Order est. (p) 1 1
Order bias (q) 2 2
BW est. (h) 9.969 9.969
BW bias (b) 9.969 9.969
rho (h/b) 1.000 1.000
Unique Obs. 237 763
=====================================================================
Point Robust Inference
Estimate z P>|z| [ 95% C.I. ]
---------------------------------------------------------------------
RD Effect -8.578 -3.276 0.001 [-12.483 , -3.138]
=====================================================================rdrobust(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70, h = 9.969 * 2) %>%
summary()Sharp RD estimates using local polynomial regression.
Number of Obs. 1000
BW type Manual
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Eff. Number of Obs. 206 577
Order est. (p) 1 1
Order bias (q) 2 2
BW est. (h) 19.938 19.938
BW bias (b) 19.938 19.938
rho (h/b) 1.000 1.000
Unique Obs. 237 763
=====================================================================
Point Robust Inference
Estimate z P>|z| [ 95% C.I. ]
---------------------------------------------------------------------
RD Effect -9.151 -4.980 0.000 [-11.670 , -5.078]
=====================================================================rdrobust(y = tutoring$exit_exam, x = tutoring$entrance_exam, c = 70, h = 9.969 / 2) %>%
summary()Sharp RD estimates using local polynomial regression.
Number of Obs. 1000
BW type Manual
Kernel Triangular
VCE method NN
Number of Obs. 237 763
Eff. Number of Obs. 82 109
Order est. (p) 1 1
Order bias (q) 2 2
BW est. (h) 4.984 4.984
BW bias (b) 4.984 4.984
rho (h/b) 1.000 1.000
Unique Obs. 237 763
=====================================================================
Point Robust Inference
Estimate z P>|z| [ 95% C.I. ]
---------------------------------------------------------------------
RD Effect -8.201 -2.032 0.042 [-13.618 , -0.246]
=====================================================================Aquí los coeficientes cambian un poco :
| Bandwidth | Effect size |
|---|---|
| 9.969 (ideal) | 8.578 |
| 19.938 (twice) | 9.151 |
| 4.984 (half) | 8.201 |
También puede ajustar el kernel. Por defecto, rd_robustse usa un kernel triangular (las observaciones más distantes tienen menor ponderación linealmente), pero puede cambiarlo a Epanechnikov (las observaciones más distantes tienen menor ponderación siguiendo una curva) o uniforme (las observaciones más distantes tienen la misma ponderación que las observaciones más cercanas; esto no está ponderado).
Nuevamente los coeficientes cambian ligeramente:
| Kernel | Effect size |
|---|---|
| Triangular | 8.578 |
| Epanechnikov | 8.389 |
| Uniform | 8.175 |
¿Cuál es mejor? ¯\_(ツ)_/¯. Lo importante es que la magnitud y la dirección del efecto no cambien. Sigue siendo positivo y se mantiene en el rango de 8 a 9 puntos.
Paso 6: Comparar todos los efectos
Acabamos de estimar una gran cantidad de tamaños de efecto. En la práctica, generalmente se reportaría uno de estos como efecto final, pero se ejecutarían los diferentes modelos paramétricos y no paramétricos para comprobar la fiabilidad y solidez de los hallazgos. Aquí está todo lo que encontramos:
| Method | Bandwidth | Kernel | Estimate |
|---|---|---|---|
| Parametric | Full data | Unweighted | 10.8 |
| Parametric | 10 | Unweighted | 9.273 |
| Parametric | 5 | Unweighted | 9.122 |
| Nonparametric | 9.969 | Triangular | 8.578 |
| Nonparametric | 19.938 | Triangular | 9.151 |
| Nonparametric | 4.984 | Triangular | 8.201 |
| Nonparametric | 8.201 | Epanechnikov | 8.389 |
| Nonparametric | 7.346 | Uniform (unweighted) | 8.175 |
En la vida real, probablemente informaría el más simple (fila 4: no paramétrico, ancho de banda automático, kernel triangular), pero es útil saber cuánto varía el efecto según las especificaciones del modelo.