Rastros del polvo blanco
Published:
Fecha de última actualización: 15 de marzo de 2024
TL;DR: La evidencia indica que tras el cierre de la Base de Manta en 2009, las provincias claves en rutas de narcotráfico, especialmente en la costa, experimentaron un aumento del 75,48% en las tasas de muertes violentas de intención indeterminada en 2010. Este efecto se invierte en el tiempo, lo que sugiere la apertura y/o consolidación de rutas logísticas en provincias amazónicas como Sucumbíos y Napo. La metodología se basa en un estudio de eventos aplicado a un modelo de Diferencia-en-Diferencias (DiD), utilizando 252 observaciones anuales a nivel provincial entre 2007 y 2018.
Un fantasma recorre las filas de la Revolución Ciudadana: el fantasma del narcotráfico. En un [post anterior], señalé inconsistencias en el principal argumento defensivo del expresidente Rafael Correa sobre las bajas tasas de homicidio durante su gobierno. Específicamente, hubo un aumento simultáneo de muertes violentas de intención indeterminada—como estrangulamientos y envenenamientos en los que "la información disponible no permite distinguir entre accidente, lesión autoinfligida y agresión" [Clasificación Internacional de Enfermedades]—, una categoría que podría ocultar asesinatos relacionados con el narcotráfico, dado que su tasa alcanzó cifras inéditas en la historia del Ecuador.
Desde el cierre de la base militar estadounidense en Manta en septiembre de 2009, se observó un marcado incremento en la tasa de muertes violentas indeterminadas, alcanzando la cifra récord en 2014. Para intentar comprender este fenómeno, mi [anterior análisis] se centró en explorar una posible correlación entre esta métrica y las incautaciones de cocaína, que pueden servir como un indicador del nivel de permisividad o indulgencia del gobierno hacia el narcotráfico. Comparé los períodos 1999-2007 (pre-Correa) y 2010-2018 (post-Manta hasta la ruptura Moreno-Correa). Mi intuición teórica sugería que estas dos variables debían moverse en direcciones opuestas: si las incautaciones bajan, entonces esta tasa sube. Para ponerla a prueba, realicé análisis estadísticos que efectivamente indicaron la existencia de una relación significativa, pero únicamente durante el segundo período.
Dado que las incautaciones de cocaína cayeron un 80% en 2010 tras el cese de operaciones de la base, esto podría constituir un factor explicativo importante. De hecho, el análisis de regresión con quiebres estructurales de Manuel González concluyó que: "[El] análisis contrafactual indicaría que si la Base de Manta se hubiera mantenido, la incautación de drogas sería al menos 75% más alta que lo que actualmente observamos" [Link].
La improbable hipótesis de que la caída permanente haya sido producto de falta de droga para confiscar, debido a una disminución en la producción de coca en Colombia, también fue refutada por González. Su análisis, controlando por el cierre de la base, arrojó una asociación nula en donde "solo la constante y el coeficiente de la variable dummy son estadísticamente significativos (el valor p del coeficiente de producción de coca es 0.85)" [Link]. En otras palabras, es probable que las rutas narcotraficantes ecuatorianas fueran para la cocaína lo que el balón era para los pies de Maradona (no pun intended).
Aún si excluimos factores del lado de la oferta, la relación previamente detectada entre las incautaciones y la tasa de muertes violentas indeterminadas podría ser espuria debido a la escasez de datos nacionales. Para remediarlo, decidí llevar a cabo el análisis a escala provincial en lugar de nacional. Mi premisa era que las estimaciones previas deberían ser más contundentes en zonas donde el narcotráfico es más intenso, especialmente en provincias fronterizas con Colombia como Esmeraldas.
Ahora, si la cocaína transita activamente por rutas logísticas, entonces la visualización espacial de las nuevas estimaciones debería permitirnos rastrear sus huellas. Aplicando la [misma estimación] a nivel provincial, y utilizando el logaritmo natural para asegurar comparabilidad, empleé los valores p como indicadores de la evidencia contra la hipótesis nula (H0) de no asociación entre estas dos variables a lo largo del período de 9 años, desde 2010 hasta 2018, en la provincia en cuestión. Con el "termómetro" en mano, es posible ilustrar los resultados mediante el siguiente mapa de calor:
Para mi sorpresa, este mapa coincide con notable precisión con las rutas de la cocaína recientemente identificadas por la Policía Nacional, tal como se muestra en la siguiente imagen del diario digital [Primicias]. Al contrastar los mapas, la evidencia estadística parece capturar detalles muy sutiles. Un ejemplo de esto es la hipótesis que sugiere que la cocaína se transporta desde Pastaza hacia Bolívar, evitando detenerse en Tungurahua; o el caso de Santo Domingo de los Tsáchilas, que permanece relativamente intacta a pesar de estar circundada por las rutas tanto del Pacífico como de la Sierra.
Los indicios en mi [post anterior] no solo que parecen no ser espurios al identificar las correspondencias a nivel provincial, sino que, figurativamente, el mapa resultante se "triangula" a través de las coincidencias con las rutas de cocaína delineadas por la inteligencia policial. En este sentido, podríamos utilizar el mapa de calor para elaborar un diseño metodológico cuasi-experimental que identifique de manera potencialmente causal el impacto del cierre de la base militar sobre este tipo de muerte violenta. Asignaríamos cada provincia a un grupo de Tratamiento o Control según los patrones territoriales. Esto nos permite una clasificación inicial donde la pertenencia a cada grupo segmenta las provincias en términos de alta (Tratamiento) versus baja (Control) exposición a rutas narcotraficantes:
| Grupo | Provincias | Cantidad |
|---|---|---|
| Tratamiento | Esmeraldas, Manabí, Guayas, Santa Elena, Carchi, Imbabura, Pichincha, Cotopaxi, Los Ríos, Bolívar, El Oro, Pastaza | 12 |
| Control | Santo Domingo, Sucumbíos, Napo, Orellana, Tungurahua, Chimborazo, Cañar, Morona Santiago, Azuay, Loja, Zamora Chinchipe | 11 |
Conceptualmente, la clasificación de provincias nos permite utilizar una metodología que combine el enfoque de Diferencia-en-Diferencias (DiD) con un estudio de eventos. De ese modo podríamos aislar los efectos dinámicos derivados del cierre de la Base de Manta en 2009; es decir, no solo detectar su impacto inmediato en 2010 (que es el prioritario en este blog post), pero también intentar monitorear su evolución en el tiempo. Mediante la base de datos, que la comparto en el siguiente [link], y aplicando el código en R que se encuentra al final de este documento, obtuve estimaciones a partir de 252 observaciones a nivel provincial a lo largo de un período de 12 años, desde 2007 hasta 2018, cuyas estadísticas descriptivas son las siguientes:
| Grupo | Mean | Std. Dev | Min | Max | Obs |
|---|---|---|---|---|---|
| Tratamiento | 1.36 | 0.800 | 0 | 3.12 | 132 |
| Control | 1.10 | 0.795 | 0 | 2.87 | 120 |
La similitud en las desviaciones estándar entre los grupos (Std. Dev: 0.800 versus 0.795) indica que son comparables en términos de la consistencia de sus distribuciones. A su vez, esto asegura comparabilidad en la variabilidad de las tasas de muertes violentas de intención no determinada; es decir, es una buena indicación de que cualquier diferencia en promedios (Mean) no se debe a que el grupo de Tratamiento sea simplemente más variable que el grupo de Control.
Al final de este post especifico mayores detalles metodológicos. Allí también abordo el supuesto de "tendencias paralelas" sobre el que descansa la validez del estimador DiD, y realizo pruebas de falsificación asignando aleatoriamente un tratamiento placebo. Utilizando errores estándar robustos y agrupados a nivel de provincia, se obtienen los siguientes resultados:
Al final de este post especifico mayores detalles metodológicos. Allí también abordo el supuesto de "tendencias paralelas" sobre el que descansa la validez del estimador DiD, y realizo pruebas de falsificación asignando aleatoriamente un tratamiento placebo. Utilizando errores estándar robustos y agrupados a nivel de provincia, se obtienen los siguientes resultados:
Las estimaciones muestran que el cierre de la Base de Manta en 2009 ocasionó un aumento del 75,48% en la tasa de muertes violentas durante 2010 en el grupo de Tratamiento (provincias fuertemente vinculadas a rutas del narcotráfico), en comparación con las del grupo de Control. Este efecto se estima manteniendo constantes otras características que pudieran incidir, puesto que el diseño metodológico neutraliza diferencias preexistentes, cambios provinciales específicos y tendencias comunes previas entre grupos.
Adicionalmente, el análisis gráfico indica que el impacto empieza a disminuir hacia el final del segundo mandato de Correa, entre 2015 y 2017, y se invierte con la ruptura política entre Correa y su sucesor y exaliado, Lenin Moreno, ocurrida en 2017-2018. En ese período, Correa acusó a Moreno de traición y corrupción [Link]. La reversión del efecto coincidió con un aumento de las medidas gubernamentales contra el narcotráfico, particularmente en provincias costeras, destacando un incremento récord en las incautaciones de cocaína [Link]. Al respecto, el expresidente Moreno arguyó que si anteriormente "vivíamos una paz" era porque "se permitía que la droga transite por nuestro territorio" [Link].
Reflexionando retrospectivamente, se pueden teorizar mecanismos posibles. Por ejemplo, mayores controles en Esmeraldas pudieron fomentar nuevas rutas por zonas amazónicas controladas. Ello explicaría la aparición ahora de Napo en informes policiales y el aumento en las tasas, a pesar de no estar antes señalada en el Otro Mapa del Ecuador (2010-2018). Otra posibilidad, más directa, es que una parte importante de estas muertes violentas, especialmente en la costa, empezaran a ser finalmente reconocidas por lo que eran: homicidios. O quizá el efecto observado sea producto de una combinación de ambos factores.
Con cautela pero firmeza, el incremento en muertes violentas en 2010 se presenta como un posible efecto causal del cierre de la Base de Manta en 2009. La coincidencia entre grupos y rutas narcotraficantes refuerza la plausibilidad de una conexión. Si bien los hallazgos no prueban de forma definitiva, sí brindan un respaldo estadístico a ciertas conclusiones de la [Comisión de Transparencia y Verdad] que investigó el caso Angostura en 2008. En su informe oficial de diciembre 2009, la Comisión señaló que "las FARC y el narcotráfico han logrado penetrar varias instancias de los organismos políticos, judiciales, culturales y sociales" [Link].
Investigaciones futuras con más datos, especificaciones detalladas y controles adicionales, podrían robustecer la validez de estos resultados. Ahora bien, alcanzar certeza estadística es imposible por definición. Es más, quizá no exista certidumbre absoluta en ningún ámbito, como lo exhibe el principio legal de duda razonable para determinar culpabilidad—un principio que algunos actores políticos intentan desplazar hacia los márgenes, como quien barre bajo la alfombra ciertos rastros de un polvo blanco.
Detalles metodológicos.
1. Tendencias paralelas.
Específicamente, la ecuación:
Específicamente, la ecuación:
$$ \begin{align*} \ln(\text{Tasa}_t)&=\beta_0+\beta_1\cdot\text{Grupo}+\beta_2\cdot D_t+{\color{red}{\beta_{3}}}\cdot\text{Grupo}\cdot D_{t}+\text{Prov}+\varepsilon_t\notag\\\\ &\text{donde } D_t:=[t-2009] \end{align*} $$
$\qquad$ se enfoca en explorar el efecto diferencial que el cierre de la base ejerció sobre el grupo tratado a lo largo del tiempo, mediante el coeficiente ${\color{red}{\beta_{3}}}.$ Este coeficiente de interacción grupo-tiempo captura fluctuaciones específicas del impacto más allá del promedio. Las demás variables controlan por factores fijos de provincia y tendencias comunes entre grupos.
Antes de ponernos a ejecutar regresiones cual monos con cuchillo, debemos abordar el supuesto de "tendencias paralelas" sobre el que descansa la validez del estimador DiD. Este supuesto implica que, antes del evento de 2009, las tasas de los grupos de Tratamiento y Control evolucionaron de manera similar. La figura siguiente, que muestra las tasas entre 2007 y 2010, revela un comportamiento casi idéntico entre los grupos durante 2007-2009, seguido de una divergencia marcada precisamente en 2010. La observación de tendencias paralelas entre 2007 y 2008 suele ser suficiente para satisfacer el supuesto. La continuidad de esta similitud hasta el año 2009, justo antes del cierre de la Base de Manta en septiembre, proporciona un respaldo adicional al minimizar la posibilidad de que otros factores hayan alterado estas tendencias.
Antes de ponernos a ejecutar regresiones cual monos con cuchillo, debemos abordar el supuesto de "tendencias paralelas" sobre el que descansa la validez del estimador DiD. Este supuesto implica que, antes del evento de 2009, las tasas de los grupos de Tratamiento y Control evolucionaron de manera similar. La figura siguiente, que muestra las tasas entre 2007 y 2010, revela un comportamiento casi idéntico entre los grupos durante 2007-2009, seguido de una divergencia marcada precisamente en 2010. La observación de tendencias paralelas entre 2007 y 2008 suele ser suficiente para satisfacer el supuesto. La continuidad de esta similitud hasta el año 2009, justo antes del cierre de la Base de Manta en septiembre, proporciona un respaldo adicional al minimizar la posibilidad de que otros factores hayan alterado estas tendencias.
Una preocupación válida radica en la creación de las provincias de Santa Elena y Santo Domingo de los Tsáchilas en noviembre de 2007. Omitirlas del análisis adopta una postura más prudente, especialmente dada la relevancia de Santa Elena en las rutas de narcotráfico. Exclusión mediante [ver sección 2.2], el supuesto de tendencias paralalelas no solo continúa superando la verificación visual, pero produce una notable mejora:
También se superan pruebas estadísticas formales. En concreto, los coeficientes estimados para el período previo a 2010, en relación con la pertenencia a los grupos, resultan no significativos. Utilizando el comando en R mostrado en la siguiente sección, obtenemos:
Tabla 1. Test formal: Tendencias paralelas
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Grupo:[2008] | 0.122216 | 0.341932 | 0.3574 | 0.72275 |
| Grupo:[2009] | 0.139178 | 0.395636 | 0.3518 | 0.72695 |
Residual standard error: 0.5608 on 38 degrees of freedom
Multiple R-squared: 0.6995, Adjusted R-squared: 0.5097
F-statistic: 3.686 on 24 and 38 DF, p-value: < 0.0001664
2. Modelo central y pruebas de falsificación.
Respecto al gráfico previo, los resultados se obtuvieron con errores estándar robustos HC1 agrupados por provincia. Las estimaciones se verifican consistentes también tras aplicar el estimador de Newey-West con cuatro rezagos, que ajusta heterocedasticidad y autocorrelaciones. Las estimaciones del gráfico son númericamente las siguientes:
Tabla 2. Resultados númericos del modelo central
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Grupo:[-2] | -0.0971708 | 0.4026634 | -0.2413 | 0.8095448 |
| Grupo:[-1] | -0.0169621 | 0.2767917 | -0.0613 | 0.9511943 |
| Grupo:[1] | 0.5623433 | 0.2665823 | 2.1095 | 0.0361013* |
| Grupo:[2] | -0.7166077 | 0.4007359 | -1.7882 | 0.0751946 |
| Grupo:[3] | -0.1591828 | 0.3385036 | -0.4703 | 0.6386661 |
| Grupo:[4] | -0.5492535 | 0.3159830 | -1.7382 | 0.0836488 |
| Grupo:[5] | -0.6257075 | 0.3430680 | -1.8239 | 0.0696085 |
| Grupo:[6] | -0.7187493 | 0.3282678 | -2.1895 | 0.0296713* |
| Grupo:[7] | -0.6404073 | 0.3433707 | -1.8651 | 0.0635811 |
| Grupo:[8] | -0.9249716 | 0.4061550 | -2.2774 | 0.0237800* |
| Grupo:[9] | -1.1923565 | 0.2906964 | -4.1017 | 5.844e-5 *** |
Residual standard error: 0.568 on 208 degrees of freedom
Multiple R-squared: 0.5887, Adjusted R-squared: 0.5036
F-statistic: 6.922 on 43 and 208 DF, p-value: < 2.2e-16
2.1. Prueba de falsificación
Realizar un test de placebo mediante la asignación aleatoria de un tratamiento ficticio, sin efectos reales, ayuda a verificar la robustez de las estimaciones ante la presencia de sesgos o variables omitidas en el modelo central. Si bajo esta asignación se identificaran efectos significativos, esto indicaría la posible presencia de sesgos. Su ausencia, en cambio, reforzaría la credibilidad de que los efectos se deben al cierre de la Base de Manta en 2009. Aplicando el comando en R que se presenta en la siguiente sección, observamos que el modelo central pasa satisfactoriamente la prueba de placebo:
Realizar un test de placebo mediante la asignación aleatoria de un tratamiento ficticio, sin efectos reales, ayuda a verificar la robustez de las estimaciones ante la presencia de sesgos o variables omitidas en el modelo central. Si bajo esta asignación se identificaran efectos significativos, esto indicaría la posible presencia de sesgos. Su ausencia, en cambio, reforzaría la credibilidad de que los efectos se deben al cierre de la Base de Manta en 2009. Aplicando el comando en R que se presenta en la siguiente sección, observamos que el modelo central pasa satisfactoriamente la prueba de placebo:
Tabla 3. Resultados de la prueba de falsificación
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Placebo:[-2] | -0.470525 | 0.408685 | -1.151 | 0.25093 |
| Placebo:[-1] | 0.059447 | 0.394333 | 0.151 | 0.88032 |
| Placebo:[1] | -0.035625 | 0.396637 | -0.090 | 0.92852 |
| Placebo:[2] | 0.348469 | 0.390852 | 0.892 | 0.37366 |
| Placebo:[3] | 0.177191 | 0.391193 | 0.453 | 0.65106 |
| Placebo:[4] | 0.292472 | 0.390889 | 0.748 | 0.45517 |
| Placebo:[5] | 0.337489 | 0.388983 | 0.868 | 0.38660 |
| Placebo:[6] | -0.054863 | 0.404951 | -0.135 | 0.89236 |
| Placebo:[7] | 0.403286 | 0.412249 | 0.978 | 0.32908 |
| Placebo:[8] | 0.514318 | 0.39744 | 1.294 | 0.19707 |
| Placebo:[9] | 0.071283 | 0.394338 | 0.181 | 0.85673 |
Residual standard error: 0.6057 on 208 degrees of freedom
Multiple R-squared: 0.5322, Adjusted R-squared: 0.4355
F-statistic:5.504 on 43 and 208 DF, p-value: < 2.2e-16
2.2. Robustez: Santa Elena y Santo Domingo de los Tsáchilas
Para resolver cualquier inquietud sobre la exclusión de estas dos provincias, esta sección repite el análisis incorporándolas. En primer lugar, los siguientes resultados corroboran que el supuesto de tendencias paralelas previas al tratamiento se sigue manteniendo:
Para resolver cualquier inquietud sobre la exclusión de estas dos provincias, esta sección repite el análisis incorporándolas. En primer lugar, los siguientes resultados corroboran que el supuesto de tendencias paralelas previas al tratamiento se sigue manteniendo:
Tabla 4. Test formal: Tendencias paralelas
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Grupo:[2008] | 0.2109490 | 0.3329251 | 0.6336 | 0.5297625 |
| Grupo:[2009] | 0.2202526 | 0.3725037 | 0.5913 | 0.5575044 |
Residual standard error: 0.5607 on 42 degrees of freedom
Multiple R-squared: 0.6959, Adjusted R-squared: 0.5077
F-statistic: 3.697 on 26 and 42 DF, p-value: < 8.251e-05
También se mantienen los resultados del modelo central:
Tabla 5. Resultados númericos del modelo alterno
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Grupo:[-2] | -0.2202526 | 0.3707599 | -0.5941 | 0.5530552 |
| Grupo:[-1] | -0.0093036 | 0.2599900 | -0.0358 | 0.9714852 |
| Grupo:[1] | 0.4849058 | 0.2460748 | 1.9706 | 0.0499682* |
| Grupo:[2] | -0.6476735 | 0.3780252 | -1.7133 | 0.0879976. |
| Grupo:[3] | -0.0651889 | 0.3266296 | -0.1996 | 0.8419842 |
| Grupo:[4] | -0.5973396 | 0.2888526 | -2.0680 | 0.0397561* |
| Grupo:[5] | -0.7654834 | 0.3207095 | -2.3868 | 0.0177996* |
| Grupo:[6] | -0.8096031 | 0.3044824 | -2.6589 | 0.0083868** |
| Grupo:[7] | -0.7674772 | 0.3184151 | -2.4103 | 0.0167219* |
| Grupo:[8] | -0.9823449 | 0.3696537 | -2.6575 | 0.0084225** |
| Grupo:[9] | -1.2553163 | 0.2663352 | -4.7133 | 4.213e-06*** |
Residual standard error: 0.5601 on 231 degrees of freedom
Multiple R-squared: 0.5911, Adjusted R-squared: 0.5133
F-statistic: 7.591 on 44 and 231 DF, p-value: < 2.2e-16
Y también se continúa superando las pruebas de falsificación:
Tabla 6. Resultados de la prueba de falsificación (modelo alterno)
| Coefficient | Estimate | Std. Error | t value | Pr(>|t|) |
|---|---|---|---|---|
| Placebo:[-2] | -0.4672519 | 0.4800489 | -0.9733 | 0.331406 |
| Placebo:[-1] | 0.2482002 | 0.3151349 | 0.7876 | 0.431742 |
| Placebo:[1] | 0.0017733 | 0.3665619 | 0.0048 | 0.996144 |
| Placebo:[2] | 0.4635328 | 0.3916779 | 1.1835 | 0.237851 |
| Placebo:[3] | 0.4215750 | 0.3662935 | 1.1509 | 0.250960 |
| Placebo:[4] | 0.3365370 | 0.3079688 | 1.0928 | 0.275641 |
| Placebo:[5] | 0.5124930 | 0.3422538 | 1.4974 | 0.135659 |
| Placebo:[6] | 0.1438626 | 0.3549316 | 0.4053 | 0.685615 |
| Placebo:[7] | 0.5838254 | 0.3231690 | 1.8066 | 0.072137 |
| Placebo:[8] | 0.5588859 | 0.3997009 | 1.3983 | 0.163382 |
| Placebo:[9] | 0.1167071 | 0.3407552 | 0.3425 | 0.732291 |
Residual standard error: 0.6013 on 230 degrees of freedom
Multiple R-squared: 0.5309, Adjusted R-squared: 0.4391
F-statistic: 5.784 on 45 and 230 DF, p-value: < 2.2e-16
Replicación: Comandos en R
A. La base de datos la pueden encontrar haciendo click [aquí].
B. Lo más tedioso al construir esta base de datos fue el año 2007, ya que la categoría específica de "Eventos de intención no determinada" no está incluida (o al menos yo no la encontré) en el anuario del Registro Estadístico de Defunciones Generales del INEC. Sin embargo, se puede construir esta variable utilizando la base de datos sin procesar [aquí], cuyos códigos pertinentes van desde van desde Y10 hasta Y34 conforme la [CIE-10].
B. Lo más tedioso al construir esta base de datos fue el año 2007, ya que la categoría específica de "Eventos de intención no determinada" no está incluida (o al menos yo no la encontré) en el anuario del Registro Estadístico de Defunciones Generales del INEC. Sin embargo, se puede construir esta variable utilizando la base de datos sin procesar [aquí], cuyos códigos pertinentes van desde van desde Y10 hasta Y34 conforme la [CIE-10].
####### 1. Chequeo formal de tendencias paralelas #######
# Filtrar datos para el periodo previo al tratamiento (2007-2009)
pre_treatment <- subset(data, Year <= 2009)
# Ajustar modelo con interacción entre TreatmentGroup y el año, con efectos fijos de Provincia
model_pre_trend <- lm(LnTasa ~ TreatmentGroup * factor(Year) + as.factor(Provincia), data = pre_treatment)
# Calcular errores estándar robustos (HC1) agrupados por Provincia
robust_se_pre_trend <- vcovHC(model_pre_trend, type = "HC1", cluster = ~Provincia)
# Realizar prueba de coeficientes con errores estándar robustos
coeftest(model_pre_trend, vcov = robust_se_pre_trend)
####### 2. Ajuste del modelo principal #######
# Excluir provincias específicas de los datos
data <- subset(data, !Provincia %in% c("StDomingo", "StElena", "Nacional"))
# Ajustar 'eventtime' como factor, tratando '0' por separado
data$eventtime <- factor(data$eventtime, levels = c("0", unique(as.character(data$eventtime[data$eventtime != 0]))))
# Ajustar modelo con LnTasa, interacción entre TreatmentGroup y eventtime, con efectos fijos de Provincia
model <- lm(LnTasa ~ TreatmentGroup * eventtime + as.factor(Provincia), data = data)
# Calcular errores estándar robustos (HC1) agrupados por Provincia
robust_se <- vcovHC(model, type = "HC1", cluster = ~Provincia)
# Prueba de coeficientes con errores estándar robustos
coeftest(model, vcov = robust_se)
####### 2.1. Tests de falsificación #######
# Crear variable de tratamiento placebo
set.seed(123) # Asegurar reproducibilidad
data$Placebo <- ifelse(runif(nrow(data)) > 0.5, 1, 0)
# Ajustar modelo usando la variable de tratamiento placebo, con efectos fijos de Provincia
model_placebo <- lm(LnTasa ~ Placebo * factor(eventtime) + as.factor(Provincia), data = data)
# Calcular errores estándar robustos (HC1) agrupados por Provincia
robust_se_placebo <- vcovHC(model_placebo, type = "HC1", cluster = ~Provincia)
# Prueba de coeficientes con errores estándar robustos
coeftest(model_placebo, vcov = robust_se_placebo)