Violencia territorial e inseguridad alimentaria en hogares con primera infancia
México, 2020–2024 · Sur y Península de Yucatán
Christian Campa · Tecnológico de Monterrey
02 de June de 2026
Marco conceptual. La inseguridad alimentaria (IA) en hogares con niñas y niños menores de cinco años no es solo un problema de ingreso: es la intersección de un mercado alimentario, una estructura de redes comunitarias, y un entorno institucional. Cuando la violencia territorial irrumpe sobre ese tejido —mediante el control armado de rutas, el desplazamiento forzado de productores rurales, el cierre de mercados locales, el cobro de extorsión a tiendas y agricultores, y la erosión de la confianza interpersonal que sostiene el préstamo, el trueque y el cuidado infantil compartido— los canales de adquisición de alimentos se vuelven más caros, menos predecibles, y a veces directamente imposibles. La hipótesis que estructura este trabajo es que la violencia territorial documentada y medida a nivel municipal está asociada con probabilidades sistemáticamente más altas de IA moderada/severa y severa en hogares con primera infancia, controlando por ingreso, ruralidad, etnicidad y estructura familiar, una vez que se absorben las heterogeneidades persistentes entre entidades y los shocks comunes a cada año.
1 Presentación del estudio
Este documento integra el análisis cuantitativo de mi proyecto de capstone sobre violencia territorial e inseguridad alimentaria en primera infancia en el sur de México (Campeche, Chiapas, Guerrero, Oaxaca, Quintana Roo, Tabasco y Yucatán) durante el periodo 2020–2024. Trabajo con tres bases administrativas y dos fuentes externas:
- ENIGH (INEGI) 2020, 2022 y 2024 — microdatos hogar/población.
- SESNSP — incidencia delictiva municipal mensual 2015–2025.
- CONAPO — proyecciones de población municipal para construir tasas.
- ACLED — eventos de violencia política y territorial geocodificados.
- Eventos documentados — corpus propio de episodios fuertes (desplazamientos masivos, masacres, disputas territoriales) construido a partir de CMDPDH, CNDH, prensa y literatura académica.
El análisis combina estimaciones poblacionales ponderadas con
diseño muestral complejo (paquete survey),
modelos de probabilidad lineal con efectos fijos de dos vías y
errores agrupados a nivel municipal
(fixest::feols), y modelos probit de efectos
fijos (fixest::feglm) con efectos marginales
promedio. La estructura del documento sigue la cadena reproducible del
pipeline: carga → integración → estimación poblacional → modelos
econométricos → diagnóstico → identificación municipal → narrativa.
Nota de reproducibilidad. Las rutas de archivos
corresponden al equipo del autor (macOS, iCloud Drive sincronizado).
Para replicar el script en otro entorno basta con sustituir el bloque
root_dir y las tres rutas externas (SESNSP, ACLED, CONAPO)
en la Sección 2. El resto del flujo es independiente del sistema de
archivos.
2 Configuración del entorno
Cargo todas las dependencias al inicio para evitar choques de
espacios de nombres a media compilación. El núcleo del análisis depende
de tidyverse, survey, fixest y
marginaleffects.
# --- Manipulación, lectura, limpieza ---
library(tidyverse)
library(readr)
library(readxl)
library(janitor)
library(stringi)
library(lubridate)
# --- Inferencia con diseño muestral complejo ENIGH ---
library(survey)
# --- Modelos de efectos fijos y probit ---
library(fixest)
library(marginaleffects)
library(car) # VIF
# --- Diagnóstico ---
library(performance)
# --- Visualización ---
library(ggplot2)
library(patchwork)
library(scales)
# --- Tablas para HTML ---
library(knitr)
library(kableExtra)
# Opciones globales del paquete survey
options(survey.lonely.psu = "adjust")# Paleta consistente con el tono del proyecto:
# azul pizarra para violencia/instituciones, borgoña apagado para
# inseguridad alimentaria, verdes apagados para indicadores compuestos.
pal_principal <- "#2C5F8D"
pal_acento <- "#8B4513"
pal_severa <- "#a8423d"
pal_alimentos <- "#4a7c59"
pal_gris <- "#7d7d7d"3 Funciones reutilizables de carga ENIGH
ENIGH se publica como tres archivos relacionados por año (concentrado
del hogar, módulo de hogares con preguntas ELCSA, y población a nivel
persona). En lugar de copiar tres veces el mismo bloque de
read_csv() y mutate() —que era el patrón en
las versiones de trabajo previas— defino cuatro funciones que reciben
ruta y año, devuelven tibbles estandarizados, y permiten replicar el
flujo año por año sin reescritura.
Decisión de diseño. Las funciones imponen un mismo
esquema de tipos (cve_ent, cve_mun,
cvegeo_mun siempre como character con padding
de ceros) porque INEGI alterna entre representaciones numéricas y de
texto entre años. Sin este saneamiento, los left_join con
SESNSP, CONAPO y ACLED pierden municipios silenciosamente cuando una
clave numérica 7 no empata con "07".
3.1 Cargar concentradohogar
#' Cargar y preparar concentradohogar de ENIGH
#'
#' @param ruta_archivo Ruta absoluta al CSV de concentradohogar.
#' @param anio_enigh Año del levantamiento (2020, 2022 o 2024).
#' @return Tibble a nivel hogar con identificadores estandarizados,
#' ingreso per cápita, gasto alimentario per cápita, participación
#' del gasto en alimentos y dummy de ruralidad.
cargar_concentrado <- function(ruta_archivo, anio_enigh) {
read_csv(
ruta_archivo,
col_types = cols(
folioviv = col_character(),
foliohog = col_character(),
ubica_geo = col_character(),
est_dis = col_character(),
upm = col_character(),
.default = col_guess()
)
) %>%
mutate(
anio = anio_enigh,
id_hogar = paste(folioviv, foliohog, sep = "_"),
cvegeo = str_pad(ubica_geo, 5, pad = "0"),
cve_ent = str_sub(cvegeo, 1, 2),
cve_mun = str_sub(cvegeo, 3, 5),
cvegeo_mun = paste0(cve_ent, cve_mun),
ingreso_pc = ing_cor / tot_integ,
gasto_alim_pc = alimentos / tot_integ,
share_alimentos = alimentos / gasto_mon,
rural = if_else(tam_loc == 4, 1, 0)
)
}3.2 Construir escala ELCSA de inseguridad alimentaria
ELCSA (Escala Latinoamericana y Caribeña de Seguridad Alimentaria) se
operacionaliza en ENIGH como 16 preguntas binarias
(acc_alim1–acc_alim16). El conteo de
respuestas afirmativas se traduce en umbrales convencionales: ≥4
corresponde a IA moderada/severa y ≥8 a IA severa en hogares con
menores.
#' Construir indicadores ELCSA de inseguridad alimentaria
#'
#' Recodifica las 16 preguntas binarias (1 = "Sí ocurrió" → 1;
#' 2 = "No ocurrió" → 0), suma el puntaje y aplica los umbrales
#' clásicos para hogares con menores (≥4 moderada/severa; ≥8 severa).
cargar_inseguridad_alimentaria <- function(ruta_archivo, anio_enigh) {
read_csv(
ruta_archivo,
col_types = cols(
folioviv = col_character(),
foliohog = col_character(),
entidad = col_character(),
est_dis = col_character(),
upm = col_character(),
.default = col_guess()
)
) %>%
mutate(
anio = anio_enigh,
id_hogar = paste(folioviv, foliohog, sep = "_")
) %>%
mutate(
across(
acc_alim1:acc_alim16,
~ case_when(
as.numeric(.x) == 1 ~ 1,
as.numeric(.x) == 2 ~ 0,
TRUE ~ NA_real_
),
.names = "{.col}_bin"
)
) %>%
mutate(
ia_score = rowSums(across(acc_alim1_bin:acc_alim16_bin), na.rm = TRUE),
ia_mod_sev = if_else(ia_score >= 4, 1, 0),
ia_severa = if_else(ia_score >= 8, 1, 0)
)
}3.3 Identificar hogares con primera infancia
#' Identificar hogares con niñas y niños de 0 a 5 años
#'
#' Agrega a nivel hogar el conteo de menores en tres rangos
#' (0–2, 3–5, 0–5), e identifica jefatura femenina y presencia
#' de al menos un hablante de lengua indígena.
identificar_primera_infancia <- function(ruta_archivo, anio_enigh) {
poblacion <- read_csv(
ruta_archivo,
col_types = cols(
folioviv = col_character(),
foliohog = col_character(),
numren = col_character(),
entidad = col_character(),
est_dis = col_character(),
upm = col_character(),
.default = col_guess()
)
) %>%
mutate(
anio = anio_enigh,
id_hogar = paste(folioviv, foliohog, sep = "_")
)
poblacion %>%
group_by(id_hogar) %>%
summarise(
ninos_0_5 = sum(edad >= 0 & edad <= 5, na.rm = TRUE),
ninos_0_2 = sum(edad >= 0 & edad <= 2, na.rm = TRUE),
ninos_3_5 = sum(edad >= 3 & edad <= 5, na.rm = TRUE),
hogar_0_5 = as.integer(ninos_0_5 > 0),
indigena_hogar = as.integer(any(hablaind == 1, na.rm = TRUE)),
jefatura_femenina = as.integer(any(sexo == 2 & parentesco == 101, na.rm = TRUE)),
.groups = "drop"
)
}3.4 Integrar las tres tablas
#' Integrar concentrado + hogares + primera infancia
#'
#' Mantiene todos los hogares (no filtra por primera infancia) para
#' permitir análisis poblacionales y de subgrupos en pasos posteriores.
integrar_enigh <- function(concentrado, hogares, primera_infancia) {
concentrado %>%
left_join(
hogares %>%
dplyr::select(id_hogar, starts_with("acc_alim"),
ia_score, ia_mod_sev, ia_severa),
by = "id_hogar"
) %>%
left_join(primera_infancia, by = "id_hogar") %>%
mutate(
hogar_0_5 = replace_na(hogar_0_5, 0),
ninos_0_5 = replace_na(ninos_0_5, 0),
ninos_0_2 = replace_na(ninos_0_2, 0),
ninos_3_5 = replace_na(ninos_3_5, 0),
indigena_hogar = replace_na(indigena_hogar, 0),
jefatura_femenina = replace_na(jefatura_femenina, 0)
)
}3.5 Funciones para el diseño muestral
ENIGH es una encuesta probabilística con diseño complejo
(estratificación por est_dis, conglomeración por
upm, factor de expansión factor). Cualquier
estimación poblacional debe pasar por svydesign y la
familia svy*. Las medias muestrales sin ponderación
subestiman sistemáticamente la inseguridad alimentaria porque
sobre-representan hogares urbanos de mayor ingreso.
#' Construir objeto svydesign para ENIGH
crear_diseno_muestral <- function(base_completa) {
options(survey.lonely.psu = "adjust")
svydesign(
ids = ~upm,
strata = ~est_dis,
weights = ~factor,
data = base_completa,
nest = TRUE
)
}
#' Calcular estimaciones estatales ponderadas en hogares con primera infancia
calcular_estimaciones_estatales <- function(diseno_muestral) {
diseno_filtrado <- subset(diseno_muestral, hogar_0_5 == 1)
ia_estado <- svyby(
~ ia_mod_sev + ia_severa + ingreso_pc + gasto_alim_pc + rural + indigena_hogar,
~ cve_ent + anio,
design = diseno_filtrado,
FUN = svymean,
na.rm = TRUE,
vartype = c("se", "ci")
)
as_tibble(ia_estado)
}
#' Agregar nombres y porcentajes a estimaciones estatales
agregar_nombres_estados <- function(ia_estado) {
nombres_ent <- tibble(
cve_ent = str_pad(as.character(1:32), 2, pad = "0"),
estado = c(
"Aguascalientes", "Baja California", "Baja California Sur",
"Campeche", "Coahuila", "Colima", "Chiapas", "Chihuahua",
"Ciudad de México", "Durango", "Guanajuato", "Guerrero",
"Hidalgo", "Jalisco", "México", "Michoacán", "Morelos",
"Nayarit", "Nuevo León", "Oaxaca", "Puebla", "Querétaro",
"Quintana Roo", "San Luis Potosí", "Sinaloa", "Sonora",
"Tabasco", "Tamaulipas", "Tlaxcala", "Veracruz",
"Yucatán", "Zacatecas"
)
)
ia_estado %>%
left_join(nombres_ent, by = "cve_ent") %>%
relocate(estado, .after = cve_ent) %>%
mutate(
ia_mod_sev_pct = ia_mod_sev * 100,
ia_severa_pct = ia_severa * 100,
ci_low_modsev = (ia_mod_sev - 1.96 * se.ia_mod_sev) * 100,
ci_high_modsev = (ia_mod_sev + 1.96 * se.ia_mod_sev) * 100,
ci_low_severa = (ia_severa - 1.96 * se.ia_severa) * 100,
ci_high_severa = (ia_severa + 1.96 * se.ia_severa) * 100,
rural_pct = rural * 100,
indigena_hogar_pct = indigena_hogar * 100,
sur_sureste = if_else(
cve_ent %in% c("04", "07", "12", "20", "23", "27", "30", "31"),
1, 0
)
) %>%
arrange(desc(ia_mod_sev_pct))
}4 Carga y procesamiento de ENIGH por año
Aplico las funciones a los tres levantamientos. Cada base anual conserva la estructura de microdatos a nivel hogar; el filtro por primera infancia se hace en pasos posteriores para no destruir la representatividad poblacional que se requiere para algunos cuadros descriptivos.
# Directorio raíz local de los microdatos ENIGH descargados desde el portal
# del INEGI. Cada carpeta de levantamiento contiene los CSVs concentradohogar,
# hogares y poblacion en sus respectivas subcarpetas "conjunto_de_datos".
root_dir <- "/Users/christiancampa/Library/Mobile Documents/com~apple~CloudDocs/enigh"4.1 ENIGH 2024
concentrado_2024 <- cargar_concentrado(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh2024_ns_csv",
"conjunto_de_datos_concentradohogar_enigh2024_ns",
"conjunto_de_datos",
"conjunto_de_datos_concentradohogar_enigh2024_ns.csv"
),
anio_enigh = 2024
)
hogares_2024 <- cargar_inseguridad_alimentaria(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh2024_ns_csv",
"conjunto_de_datos_hogares_enigh2024_ns",
"conjunto_de_datos",
"conjunto_de_datos_hogares_enigh2024_ns.csv"
),
anio_enigh = 2024
)
primera_infancia_2024 <- identificar_primera_infancia(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh2024_ns_csv",
"conjunto_de_datos_poblacion_enigh2024_ns",
"conjunto_de_datos",
"conjunto_de_datos_poblacion_enigh2024_ns.csv"
),
anio_enigh = 2024
)
base_2024 <- integrar_enigh(concentrado_2024, hogares_2024, primera_infancia_2024)
cat(sprintf("Hogares totales 2024: %d\n", nrow(base_2024)))#> Hogares totales 2024: 91414cat(sprintf("Hogares con primera infancia: %d (%.1f%%)\n",
sum(base_2024$hogar_0_5),
100 * mean(base_2024$hogar_0_5)))#> Hogares con primera infancia: 19566 (21.4%)4.2 ENIGH 2022
concentrado_2022 <- cargar_concentrado(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2022_csv",
"conjunto_de_datos_concentradohogar_enigh2022_ns",
"conjunto_de_datos",
"conjunto_de_datos_concentradohogar_enigh2022_ns.csv"
),
anio_enigh = 2022
)
hogares_2022 <- cargar_inseguridad_alimentaria(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2022_csv",
"conjunto_de_datos_hogares_enigh2022_ns",
"conjunto_de_datos",
"conjunto_de_datos_hogares_enigh2022_ns.csv"
),
anio_enigh = 2022
)
primera_infancia_2022 <- identificar_primera_infancia(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2022_csv",
"conjunto_de_datos_poblacion_enigh2022_ns",
"conjunto_de_datos",
"conjunto_de_datos_poblacion_enigh2022_ns.csv"
),
anio_enigh = 2022
)
base_2022 <- integrar_enigh(concentrado_2022, hogares_2022, primera_infancia_2022)
cat(sprintf("Hogares totales 2022: %d\n", nrow(base_2022)))#> Hogares totales 2022: 90102cat(sprintf("Hogares con primera infancia: %d (%.1f%%)\n",
sum(base_2022$hogar_0_5),
100 * mean(base_2022$hogar_0_5)))#> Hogares con primera infancia: 20582 (22.8%)4.3 ENIGH 2020
concentrado_2020 <- cargar_concentrado(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2020_csv",
"conjunto_de_datos_concentradohogar_enigh_2020_ns",
"conjunto_de_datos",
"conjunto_de_datos_concentradohogar_enigh_2020_ns.csv"
),
anio_enigh = 2020
)
hogares_2020 <- cargar_inseguridad_alimentaria(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2020_csv",
"conjunto_de_datos_hogares_enigh_2020_ns",
"conjunto_de_datos",
"conjunto_de_datos_hogares_enigh_2020_ns.csv"
),
anio_enigh = 2020
)
primera_infancia_2020 <- identificar_primera_infancia(
ruta_archivo = file.path(
root_dir,
"conjunto_de_datos_enigh_ns_2020_csv",
"conjunto_de_datos_poblacion_enigh_2020_ns",
"conjunto_de_datos",
"conjunto_de_datos_poblacion_enigh_2020_ns.csv"
),
anio_enigh = 2020
)
base_2020 <- integrar_enigh(concentrado_2020, hogares_2020, primera_infancia_2020)
cat(sprintf("Hogares totales 2020: %d\n", nrow(base_2020)))#> Hogares totales 2020: 89006cat(sprintf("Hogares con primera infancia: %d (%.1f%%)\n",
sum(base_2020$hogar_0_5),
100 * mean(base_2020$hogar_0_5)))#> Hogares con primera infancia: 22127 (24.9%)Comprobación. Si todo cargó correctamente, debería ver tres bases con estructura idéntica, número de hogares totales comparable al reportado por INEGI (≈90,000 hogares por levantamiento), y una proporción de hogares con primera infancia entre el 15 y el 20 por ciento. Cualquier desviación importante indica un problema en las rutas o en el diseño de la encuesta.
5 Estimaciones estatales ponderadas
Con las bases anuales construyo estimaciones estatales ponderadas para hogares con primera infancia. Estas son las que entran en la primera figura del manuscrito y permiten comparar prevalencias entre entidades respetando el diseño complejo de la encuesta.
# Diseños muestrales por año
diseno_2024 <- crear_diseno_muestral(base_2024)
diseno_2022 <- crear_diseno_muestral(base_2022)
diseno_2020 <- crear_diseno_muestral(base_2020)
# Estimaciones estatales
ia_estado_2024 <- calcular_estimaciones_estatales(diseno_2024)
ia_estado_2022 <- calcular_estimaciones_estatales(diseno_2022)
ia_estado_2020 <- calcular_estimaciones_estatales(diseno_2020)
# Etiquetas y porcentajes
ia_estado_2024_limpia <- agregar_nombres_estados(ia_estado_2024)
ia_estado_2022_limpia <- agregar_nombres_estados(ia_estado_2022)
ia_estado_2020_limpia <- agregar_nombres_estados(ia_estado_2020)
# Panel estatal apilado
ia_estados_panel <- bind_rows(
ia_estado_2020_limpia,
ia_estado_2022_limpia,
ia_estado_2024_limpia
)
ia_estado_2024_limpia %>%
dplyr::select(estado, ia_mod_sev_pct, ci_low_modsev, ci_high_modsev,
indigena_hogar_pct, rural_pct) %>%
slice_head(n = 10) %>%
kable(
caption = "Top 10 estados con mayor prevalencia de IA moderada/severa en hogares con primera infancia (2024)",
digits = 1,
col.names = c("Estado", "IA mod/sev (%)", "IC inf", "IC sup",
"Indígena (%)", "Rural (%)")
) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
full_width = FALSE)| Estado | IA mod/sev (%) | IC inf | IC sup | Indígena (%) | Rural (%) |
|---|---|---|---|---|---|
| Guerrero | 40.3 | 36.0 | 44.7 | 28.9 | 47.4 |
| Tabasco | 37.3 | 31.3 | 43.3 | 6.9 | 50.2 |
| Oaxaca | 31.1 | 27.1 | 35.0 | 48.9 | 55.7 |
| Chiapas | 30.2 | 25.1 | 35.4 | 42.4 | 58.6 |
| Michoacán | 27.4 | 23.2 | 31.6 | 4.2 | 34.5 |
| Guanajuato | 27.3 | 23.0 | 31.6 | 0.3 | 31.4 |
| Morelos | 25.8 | 21.7 | 30.0 | 4.1 | 21.4 |
| Campeche | 25.1 | 20.8 | 29.4 | 22.0 | 34.6 |
| Veracruz | 23.3 | 18.9 | 27.6 | 15.2 | 48.7 |
| Tlaxcala | 23.2 | 19.0 | 27.3 | 6.5 | 17.4 |
Lectura sustantiva. Las entidades del sur-sureste —Chiapas, Oaxaca, Guerrero, Tabasco, Campeche— dominan la parte alta del ranking en 2024, replicando un patrón estructural que la literatura sobre pobreza alimentaria en México ha documentado durante décadas. Lo distintivo del periodo 2020-2024 es la coexistencia de este patrón histórico con un aumento marcado de la violencia territorial en varias de esas mismas entidades, lo cual motiva el análisis municipal que sigue.
plot_ranking_estados <- ggplot(
ia_estado_2024_limpia,
aes(x = reorder(estado, ia_mod_sev_pct), y = ia_mod_sev_pct,
fill = factor(sur_sureste))
) +
geom_col(alpha = 0.85) +
geom_errorbar(
aes(ymin = ci_low_modsev, ymax = ci_high_modsev),
width = 0.35, color = pal_acento, linewidth = 0.5
) +
scale_fill_manual(
values = c("0" = pal_gris, "1" = pal_principal),
labels = c("0" = "Resto del país", "1" = "Sur y Yucatán"),
name = NULL
) +
coord_flip() +
labs(
title = "Inseguridad alimentaria moderada/severa en hogares con primera infancia",
subtitle = "México, 2024 (estimaciones ponderadas con IC 95%)",
x = NULL,
y = "Porcentaje de hogares (%)",
caption = "Fuente: ENIGH 2024 (INEGI). Diseño muestral complejo con ponderadores."
) +
theme_minimal(base_size = 11) +
theme(
plot.title = element_text(face = "bold", size = 13, color = pal_principal),
plot.subtitle = element_text(color = "gray30", size = 10),
legend.position = "bottom",
panel.grid.major.y = element_blank(),
panel.grid.minor = element_blank()
)
plot_ranking_estados
Ranking estatal de inseguridad alimentaria moderada/severa en hogares con primera infancia, 2024 (estimaciones ponderadas con IC 95%)
write_csv(ia_estados_panel, "ia_primera_infancia_estados_2020_2024.csv")
ggsave("ranking_ia_estados_2024.png", plot_ranking_estados,
width = 10, height = 12, dpi = 300)6 Violencia municipal (SESNSP)
El Secretariado Ejecutivo del Sistema Nacional de Seguridad Pública publica un registro mensual municipal de incidencia delictiva. La base llega en formato ancho (una columna por mes), lo cual requiere pivote para construir series de tiempo. El siguiente bloque hace ese trabajo y agrega los conteos a ventanas bianuales alineadas con los levantamientos ENIGH:
- ENIGH 2020 ← delitos enero 2019 a diciembre 2020.
- ENIGH 2022 ← delitos enero 2021 a diciembre 2022.
- ENIGH 2024 ← delitos enero 2023 a diciembre 2024.
Decisión de ventana temporal. Una ventana de 24 meses precede a cada encuesta porque la inseguridad alimentaria no responde instantáneamente a un shock de violencia: la disrupción de mercados, el desplazamiento, la contracción del comercio formal y la pérdida de redes de apoyo tardan trimestres en traducirse en privación efectiva. La ventana bianual también suaviza el ruido de meses con conteos inusualmente bajos por subregistro estacional.
path_sesnsp <- "/Users/christiancampa/Downloads/Municipal-Delitos-2015-2025_dic2025/Municipal-Delitos-2015-2025_dic2025.csv"
sesnsp_raw <- read_csv(
path_sesnsp,
locale = locale(encoding = "Latin1"),
show_col_types = FALSE,
guess_max = 100000
) %>%
clean_names()# Pivote ancho → largo
meses_sesnsp <- c(
"enero", "febrero", "marzo", "abril", "mayo", "junio",
"julio", "agosto", "septiembre", "octubre", "noviembre", "diciembre"
)
stopifnot(all(meses_sesnsp %in% names(sesnsp_raw)))
sesnsp_largo <- sesnsp_raw %>%
rename(
anio = ano,
cve_ent = clave_ent,
cvegeo_mun = cve_municipio,
bien_juridico = bien_juridico_afectado,
tipo_delito = tipo_de_delito,
subtipo_delito = subtipo_de_delito
) %>%
mutate(
anio = as.integer(anio),
cve_ent = str_pad(as.character(cve_ent), 2, pad = "0"),
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0")
) %>%
pivot_longer(
cols = all_of(meses_sesnsp),
names_to = "mes_nombre",
values_to = "casos"
) %>%
mutate(
mes = match(mes_nombre, meses_sesnsp),
fecha = make_date(anio, mes, 1),
casos = replace_na(as.numeric(casos), 0)
) %>%
dplyr::select(
anio, mes, fecha,
cve_ent, cvegeo_mun, entidad, municipio,
bien_juridico, tipo_delito, subtipo_delito, modalidad,
casos
)# Indicadores mensuales por tipo de delito
violencia_mensual_mun <- sesnsp_largo %>%
mutate(
homicidio_doloso_flag = tipo_delito == "Homicidio" &
subtipo_delito == "Homicidio doloso",
homicidio_total_flag = tipo_delito == "Homicidio",
extorsion_flag = tipo_delito == "Extorsión",
secuestro_flag = tipo_delito == "Secuestro",
narcomenudeo_flag = tipo_delito == "Narcomenudeo",
robo_flag = tipo_delito == "Robo",
violencia_familiar_flag = tipo_delito == "Violencia familiar"
) %>%
group_by(cvegeo_mun, cve_ent, entidad, municipio, anio, mes) %>%
summarise(
homicidios_dolosos = sum(casos[homicidio_doloso_flag], na.rm = TRUE),
homicidios_total = sum(casos[homicidio_total_flag], na.rm = TRUE),
extorsiones = sum(casos[extorsion_flag], na.rm = TRUE),
secuestros = sum(casos[secuestro_flag], na.rm = TRUE),
narcomenudeo = sum(casos[narcomenudeo_flag], na.rm = TRUE),
robos = sum(casos[robo_flag], na.rm = TRUE),
violencia_familiar = sum(casos[violencia_familiar_flag], na.rm = TRUE),
.groups = "drop"
)# Agregación a ventanas bianuales ENIGH
violencia_municipal_enigh <- violencia_mensual_mun %>%
filter(anio %in% 2019:2024) %>%
mutate(
anio_enigh = case_when(
anio %in% c(2019, 2020) ~ 2020L,
anio %in% c(2021, 2022) ~ 2022L,
anio %in% c(2023, 2024) ~ 2024L,
TRUE ~ NA_integer_
)
) %>%
filter(!is.na(anio_enigh)) %>%
group_by(cvegeo_mun, cve_ent, entidad, municipio, anio_enigh) %>%
summarise(
homicidios_dolosos_24m = sum(homicidios_dolosos, na.rm = TRUE),
homicidios_total_24m = sum(homicidios_total, na.rm = TRUE),
extorsiones_24m = sum(extorsiones, na.rm = TRUE),
secuestros_24m = sum(secuestros, na.rm = TRUE),
narcomenudeo_24m = sum(narcomenudeo, na.rm = TRUE),
robos_24m = sum(robos, na.rm = TRUE),
violencia_familiar_24m = sum(violencia_familiar, na.rm = TRUE),
delitos_violentos_24m = homicidios_dolosos_24m +
extorsiones_24m +
secuestros_24m,
.groups = "drop"
) %>%
rename(anio = anio_enigh) %>%
arrange(cvegeo_mun, anio)Construcción del agregado “delitos violentos”. Sumo homicidios dolosos, extorsiones y secuestros porque los tres tienen plausibilidad teórica clara para afectar mercados alimentarios locales (matanza/intimidación de productores y transportistas; cobro de piso a tiendas; secuestros de productores con consecuencias económicas duraderas en sus comunidades). Excluyo robo común y violencia familiar de este agregado porque sus efectos sobre el sistema alimentario local son indirectos o se confunden con dinámicas de hogar.
7 Población municipal y tasas (CONAPO)
Para convertir los conteos en tasas comparables entre municipios necesito denominador poblacional. CONAPO publica proyecciones municipales anuales que se promedian dentro de la ventana ENIGH.
path_pob_mun <- "/Users/christiancampa/Documents/Projecto_país_desaparecido/Bases de Datos/Población_Mun_Projected.xlsx"
pob_mun_raw <- read_excel(path_pob_mun, sheet = 1) %>%
clean_names()
pob_mun <- pob_mun_raw %>%
rename_with(~ "anio", matches("^(anio|ano|year)$")) %>%
rename_with(~ "poblacion", matches("^(poblacion|pob_total|pob|projected_population)$")) %>%
mutate(
anio = as.integer(anio),
poblacion = as.numeric(poblacion),
cvegeo_mun = str_pad(as.character(municipality_id), 5, pad = "0"),
cve_ent = str_sub(cvegeo_mun, 1, 2),
cve_mun = str_sub(cvegeo_mun, 3, 5)
)
# Población promedio en ventanas bianuales
pob_mun_enigh <- pob_mun %>%
filter(anio %in% 2019:2024) %>%
mutate(
anio_enigh = case_when(
anio %in% c(2019, 2020) ~ 2020L,
anio %in% c(2021, 2022) ~ 2022L,
anio %in% c(2023, 2024) ~ 2024L,
TRUE ~ NA_integer_
)
) %>%
filter(!is.na(anio_enigh)) %>%
group_by(cvegeo_mun, cve_ent, cve_mun, anio_enigh) %>%
summarise(
poblacion_prom_24m = mean(poblacion, na.rm = TRUE),
poblacion_inicio = poblacion[which.min(anio)],
poblacion_fin = poblacion[which.max(anio)],
.groups = "drop"
) %>%
rename(anio = anio_enigh)# Unir población a violencia y construir tasas
violencia_municipal_enigh_tasas <- violencia_municipal_enigh %>%
mutate(
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
anio = as.integer(anio)
) %>%
left_join(
pob_mun_enigh %>% dplyr::select(cvegeo_mun, anio, poblacion_prom_24m),
by = c("cvegeo_mun", "anio")
) %>%
mutate(
# Tasas en la ventana bianual (por 100k habitantes)
tasa_homicidios_dolosos_24m = homicidios_dolosos_24m / poblacion_prom_24m * 1e5,
tasa_homicidios_total_24m = homicidios_total_24m / poblacion_prom_24m * 1e5,
tasa_extorsiones_24m = extorsiones_24m / poblacion_prom_24m * 1e5,
tasa_secuestros_24m = secuestros_24m / poblacion_prom_24m * 1e5,
tasa_narcomenudeo_24m = narcomenudeo_24m / poblacion_prom_24m * 1e5,
tasa_robos_24m = robos_24m / poblacion_prom_24m * 1e5,
tasa_violencia_familiar_24m = violencia_familiar_24m / poblacion_prom_24m * 1e5,
tasa_delitos_violentos_24m = delitos_violentos_24m / poblacion_prom_24m * 1e5,
# Tasas anualizadas (más comparables con la literatura)
tasa_homicidios_dolosos_anualizada = tasa_homicidios_dolosos_24m / 2,
tasa_extorsiones_anualizada = tasa_extorsiones_24m / 2,
tasa_secuestros_anualizada = tasa_secuestros_24m / 2,
tasa_delitos_violentos_anualizada = tasa_delitos_violentos_24m / 2
)8 Base hogar-violencia (ENIGH + SESNSP)
Apilo los tres años de ENIGH y los uno a las tasas municipales de violencia. Esta es la base que entra a los modelos de probabilidad lineal.
base_hogar_enigh_2020_2024 <- bind_rows(base_2020, base_2022, base_2024) %>%
mutate(
id_hogar_anio = paste(anio, id_hogar, sep = "_"),
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
cve_ent = str_pad(as.character(cve_ent), 2, pad = "0"),
anio = as.integer(anio)
)
base_hogar_violencia_tasas <- base_hogar_enigh_2020_2024 %>%
left_join(
violencia_municipal_enigh_tasas,
by = c("cvegeo_mun", "anio"),
suffix = c("", "_sesnsp")
)
# Base analítica restringida a primera infancia y municipios con datos
base_hogar_violencia_tasas_modelo <- base_hogar_violencia_tasas %>%
filter(hogar_0_5 == 1) %>%
filter(!is.na(tasa_delitos_violentos_anualizada)) %>%
mutate(
log_tasa_violencia = log1p(tasa_delitos_violentos_anualizada),
log_tasa_homicidios = log1p(tasa_homicidios_dolosos_anualizada),
log_tasa_extorsiones = log1p(tasa_extorsiones_anualizada),
log_tasa_secuestros = log1p(tasa_secuestros_anualizada)
)
cat(sprintf("Hogares con primera infancia y datos de violencia: %d\n",
nrow(base_hogar_violencia_tasas_modelo)))#> Hogares con primera infancia y datos de violencia: 62126cat(sprintf("Municipios representados: %d\n",
n_distinct(base_hogar_violencia_tasas_modelo$cvegeo_mun)))#> Municipios representados: 1610Sobre log1p(). Las tasas de violencia
están truncadas en cero y son asimétricas a la derecha: la gran mayoría
de municipios tienen tasas bajas, una cola larga concentra municipios
extremadamente violentos. El logaritmo natural no es viable porque
log(0) no está definido. log1p(x) = log(1 + x)
resuelve el problema, atenúa el peso de los valores extremos en la
estimación y permite interpretar el coeficiente como una
semi-elasticidad aproximada.
9 Incorporación de ACLED
ACLED registra eventos de violencia política y territorial geocodificados, codifica el tipo de violencia (batallas, violencia contra civiles, explosiones, desarrollos estratégicos) y describe los actores involucrados. A diferencia del SESNSP —que es un registro administrativo— ACLED captura los episodios que cruzan un umbral de visibilidad mediática. No es un sustituto del SESNSP, es un complemento.
Por qué ACLED importa para la inseguridad alimentaria. Hay dos canales que el SESNSP no captura bien y ACLED sí: (a) los enfrentamientos armados visibles que cierran rutas y mercados aunque no aumenten la tasa de homicidio dolosa registrada (porque las muertes pueden no ser tipificadas como homicidio doloso, sino registradas en otro fuero, o no registradas); (b) el desplazamiento documentado por organizaciones civiles, prensa y observadores internacionales, que es un evento de exposición fuerte aunque la tasa promedio de violencia siga “baja”. El emparejamiento ACLED-municipio se hace por nombre, no por clave INEGI, lo que introduce errores de match que reporto explícitamente más abajo.
limpiar_nombre <- function(x) {
x %>%
as.character() %>%
str_to_lower() %>%
str_squish() %>%
stringi::stri_trans_general("Latin-ASCII")
}path_acled <- "/Users/christiancampa/Documents/Fuentes de información/ACLED/ACLED Data_2025-12-24 (México).csv"
acled_raw <- read_csv(
path_acled,
show_col_types = FALSE,
guess_max = 100000
) %>%
clean_names()# Catálogo municipal limpio para empatar nombres de ACLED
catalogo_municipal <- violencia_municipal_enigh_tasas %>%
distinct(cvegeo_mun, cve_ent, entidad, municipio) %>%
mutate(
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
cve_ent = str_pad(as.character(cve_ent), 2, pad = "0"),
entidad_clean = limpiar_nombre(entidad),
municipio_clean = limpiar_nombre(municipio),
entidad_clean = case_when(
entidad_clean == "veracruz de ignacio de la llave" ~ "veracruz",
entidad_clean == "michoacan de ocampo" ~ "michoacan",
entidad_clean == "coahuila de zaragoza" ~ "coahuila",
entidad_clean == "ciudad de mexico" ~ "mexico city",
entidad_clean == "mexico" ~ "estado de mexico",
TRUE ~ entidad_clean
)
)
acled_clean <- acled_raw %>%
mutate(
event_date = as.Date(event_date),
anio_evento = as.integer(year),
entidad_acled = admin1,
municipio_acled = admin2,
entidad_clean = limpiar_nombre(admin1),
municipio_clean = limpiar_nombre(admin2),
entidad_clean = case_when(
entidad_clean == "mexico" ~ "estado de mexico",
entidad_clean == "mexico city" ~ "mexico city",
TRUE ~ entidad_clean
)
) %>%
filter(country == "Mexico", anio_evento %in% 2019:2024)
acled_municipalizado <- acled_clean %>%
left_join(catalogo_municipal, by = c("entidad_clean", "municipio_clean"))
# Diagnóstico del emparejamiento
acled_municipalizado %>%
summarise(
eventos = n(),
eventos_con_cvegeo = sum(!is.na(cvegeo_mun)),
eventos_sin_cvegeo = sum(is.na(cvegeo_mun)),
pct_con_cvegeo = mean(!is.na(cvegeo_mun)) * 100
) %>%
kable(caption = "Diagnóstico de emparejamiento ACLED ↔ catálogo municipal",
digits = 1) %>%
kable_styling(bootstrap_options = c("striped"), full_width = FALSE)| eventos | eventos_con_cvegeo | eventos_sin_cvegeo | pct_con_cvegeo |
|---|---|---|---|
| 83927 | 64620 | 19307 | 77 |
9.1 Universo del estudio: sur y península de Yucatán
estados_sur_vec <- c("04", "07", "12", "20", "23", "27", "31")
estados_sur_nombres_clean <- c(
"campeche", "chiapas", "guerrero", "oaxaca",
"quintana roo", "tabasco", "yucatan"
)
# Restringir ACLED al universo del análisis
acled_sur_municipalizado <- acled_municipalizado %>%
filter(
entidad_clean %in% estados_sur_nombres_clean,
!is.na(cvegeo_mun)
)
# Ventanas bianuales ENIGH
acled_ventanas_enigh <- acled_sur_municipalizado %>%
mutate(
anio = case_when(
anio_evento %in% c(2019, 2020) ~ 2020L,
anio_evento %in% c(2021, 2022) ~ 2022L,
anio_evento %in% c(2023, 2024) ~ 2024L,
TRUE ~ NA_integer_
)
) %>%
filter(!is.na(anio))Por qué excluyo Veracruz del universo “sur”. Veracruz comparte rasgos socioeconómicos con el sur (alta marginación, peso indígena en el norte) pero su estructura criminal sigue una lógica distinta (puerto, corredor, presencia de Jalisco Nueva Generación con dinámicas más cercanas al Golfo y al Pacífico). Para no contaminar la inferencia con un caso heterogéneo, lo separo. Esto sigue la decisión metodológica que se asentó al revisar las tipologías de criminalidad regional.
9.2 Indicadores municipales ACLED
acled_municipal_enigh <- acled_ventanas_enigh %>%
mutate(
violencia_civiles = event_type == "Violence against civilians",
batallas = event_type == "Battles",
explosiones = event_type == "Explosions/Remote violence",
strategic_developments = event_type == "Strategic developments",
ataques = sub_event_type == "Attack",
enfrentamientos_armados = sub_event_type == "Armed clash",
abducciones_desapariciones = sub_event_type == "Abduction/forced disappearance",
granadas_ied = sub_event_type %in% c("Grenade",
"Remote explosive/landmine/IED"),
violencia_territorial_acled = event_type %in% c("Violence against civilians",
"Battles",
"Explosions/Remote violence")
) %>%
group_by(cvegeo_mun, cve_ent, entidad, municipio, anio) %>%
summarise(
eventos_acled_total_24m = n(),
violencia_territorial_acled_24m = sum(violencia_territorial_acled, na.rm = TRUE),
violencia_civiles_acled_24m = sum(violencia_civiles, na.rm = TRUE),
batallas_acled_24m = sum(batallas, na.rm = TRUE),
explosiones_acled_24m = sum(explosiones, na.rm = TRUE),
strategic_developments_acled_24m = sum(strategic_developments, na.rm = TRUE),
ataques_acled_24m = sum(ataques, na.rm = TRUE),
enfrentamientos_armados_acled_24m = sum(enfrentamientos_armados, na.rm = TRUE),
abducciones_desapariciones_acled_24m = sum(abducciones_desapariciones, na.rm = TRUE),
granadas_ied_acled_24m = sum(granadas_ied, na.rm = TRUE),
fatalidades_acled_24m = sum(fatalities, na.rm = TRUE),
.groups = "drop"
) %>%
arrange(cvegeo_mun, anio)10 Base analítica final: ENIGH + SESNSP + ACLED, universo sur
Esta es la base que alimenta los modelos: cada fila es un hogar con primera infancia residente en un estado del sur, con sus características demográficas y socioeconómicas, las tasas anualizadas de violencia SESNSP del municipio y los conteos de eventos ACLED del municipio en la ventana bianual previa al levantamiento.
base_hogar_violencia_acled_sur <- base_hogar_violencia_tasas_modelo %>%
mutate(
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
cve_ent = str_pad(as.character(cve_ent), 2, pad = "0"),
anio = as.integer(anio)
) %>%
filter(cve_ent %in% estados_sur_vec) %>%
left_join(
acled_municipal_enigh,
by = c("cvegeo_mun", "anio"),
suffix = c("", "_acled")
) %>%
mutate(
across(
c(eventos_acled_total_24m, violencia_territorial_acled_24m,
violencia_civiles_acled_24m, batallas_acled_24m,
explosiones_acled_24m, strategic_developments_acled_24m,
ataques_acled_24m, enfrentamientos_armados_acled_24m,
abducciones_desapariciones_acled_24m, granadas_ied_acled_24m,
fatalidades_acled_24m),
~ replace_na(.x, 0)
),
log_acled_total = log1p(eventos_acled_total_24m),
log_acled_territorial = log1p(violencia_territorial_acled_24m),
log_acled_civiles = log1p(violencia_civiles_acled_24m),
log_acled_batallas = log1p(batallas_acled_24m),
log_acled_fatalidades = log1p(fatalidades_acled_24m)
)
base_hogar_violencia_acled_sur %>%
summarise(
hogares = n(),
municipios = n_distinct(cvegeo_mun),
entidades = n_distinct(cve_ent),
anios = n_distinct(anio),
hogares_con_acled = sum(eventos_acled_total_24m > 0),
pct_hogares_con_acled = mean(eventos_acled_total_24m > 0) * 100,
pct_hogares_con_acled_territorial = mean(violencia_territorial_acled_24m > 0) * 100
) %>%
kable(caption = "Resumen de la base analítica final (universo sur)",
digits = 1) %>%
kable_styling(bootstrap_options = c("striped"), full_width = FALSE)| hogares | municipios | entidades | anios | hogares_con_acled | pct_hogares_con_acled | pct_hogares_con_acled_territorial |
|---|---|---|---|---|---|---|
| 12688 | 469 | 7 | 3 | 10725 | 84.5 | 74.2 |
11 Modelos econométricos
Especificación principal y por qué. Estimo
\[ \Pr(\text{IA}_{ihmt}=1) = \alpha + \beta\,\text{Violencia}_{mt} + \mathbf{X}_{ih}'\boldsymbol{\gamma} + \delta_t + \mu_e + \varepsilon_{ihmt} \]
donde \(i\) indexa hogar, \(h\) es la observación hogar-año, \(m\) es municipio, \(e\) es entidad, \(t\) es el levantamiento (2020, 2022, 2024); \(\delta_t\) son efectos fijos de año (absorben shocks comunes nacionales como COVID-19, inflación alimentaria post-Ucrania, recortes de programas federales); \(\mu_e\) son efectos fijos de entidad (absorben diferencias persistentes de pobreza, mercados, estructura productiva); \(\mathbf{X}_{ih}\) son controles a nivel hogar; los errores se agrupan a nivel municipal porque el tratamiento (violencia) varía a ese nivel y la correlación intra-municipio de los residuos contamina los errores estándar ingenuos. Las observaciones se ponderan por el factor de expansión ENIGH para que la inferencia sea representativa de la población-objetivo, no de la muestra cruda.
Probabilidad lineal vs probit. El modelo de probabilidad lineal (LPM) tiene la virtud de la transparencia interpretativa —los coeficientes son directamente cambios en puntos porcentuales— y se lleva bien con la combinación de efectos fijos de dos vías más errores agrupados. Su debilidad conocida es que las probabilidades ajustadas pueden caer fuera de [0, 1], lo que importa cuando la prevalencia se acerca a los extremos. En este universo (sur de México), la prevalencia de IA moderada/severa es moderada-alta (≈30–50% según el estado), así que las predicciones fuera de rango son una preocupación real. El probit resuelve eso a costa de coeficientes no interpretables directamente —se recuperan vía efectos marginales promedio (AME). En este documento reporto LPM como el cuadro principal por interpretabilidad y probit con AMEs como verificación de robustez.
11.1 Preparación de la base para modelos
base_modelos_sur <- base_hogar_violencia_acled_sur %>%
mutate(
anio = as.factor(anio),
cve_ent = as.factor(cve_ent),
cvegeo_mun = as.factor(cvegeo_mun),
# Especificaciones alternativas de ACLED
acled_territorial = violencia_territorial_acled_24m,
acled_territorial_dummy = as.integer(violencia_territorial_acled_24m > 0),
acled_territorial_cat = case_when(
violencia_territorial_acled_24m == 0 ~ "Sin eventos",
violencia_territorial_acled_24m > 0 & violencia_territorial_acled_24m <= 2 ~ "Baja",
violencia_territorial_acled_24m > 2 & violencia_territorial_acled_24m <= 10 ~ "Media",
violencia_territorial_acled_24m > 10 ~ "Alta",
TRUE ~ NA_character_
),
acled_territorial_cat = factor(acled_territorial_cat,
levels = c("Sin eventos", "Baja", "Media", "Alta"))
)11.2 Modelos LPM: línea base, ACLED, conjunto
# Modelo 1: Solo SESNSP (tasa de delitos violentos)
modelo_sesnsp <- feols(
ia_mod_sev ~ log_tasa_violencia +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
# Modelo 2: Solo ACLED (violencia territorial logarítmica)
modelo_acled <- feols(
ia_mod_sev ~ log_acled_territorial +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
# Modelo 3: Conjunto SESNSP + ACLED
modelo_conjunto <- feols(
ia_mod_sev ~ log_tasa_violencia + log_acled_territorial +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
# Modelo 4: SESNSP + ACLED categórico (intensidad)
modelo_conjunto_cat <- feols(
ia_mod_sev ~ log_tasa_violencia + i(acled_territorial_cat, ref = "Sin eventos") +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
etable(
modelo_sesnsp, modelo_acled, modelo_conjunto, modelo_conjunto_cat,
se.below = TRUE,
fitstat = ~ n + r2,
title = "Tabla 1. Violencia municipal e inseguridad alimentaria moderada/severa (LPM, IA mod/sev)",
notes = "Efectos fijos de año y entidad. Errores estándar agrupados por municipio. Ponderado por factor de expansión ENIGH."
)Lectura de la Tabla 1. Las columnas (1) y (2) muestran que tanto el SESNSP como ACLED, por separado, predicen incrementos en la probabilidad de IA moderada/severa. Cuando entran juntos en (3) ambos coeficientes pierden algo de magnitud porque comparten varianza —no son mediciones independientes del mismo fenómeno, son ángulos distintos de un mismo proceso de violencia territorial—. La especificación con ACLED categórico (4) permite ver si el efecto es monotónico en la intensidad: los hogares en municipios “Alta” ACLED tienen probabilidades sustantivamente mayores que aquellos en municipios “Sin eventos”, controlando todo lo demás.
11.3 Modelos LPM por tipo de delito
Una pregunta natural para informar política pública es: ¿qué tipo de violencia importa más para la IA? Los municipios “violentos” no lo son en todas las dimensiones; algunos lo son por extorsión, otros por homicidios, otros por desapariciones masivas.
modelo_homicidios <- feols(
ia_mod_sev ~ log_tasa_homicidios +
ingreso_pc + gasto_alim_pc + rural + indigena_hogar +
ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_extorsiones <- feols(
ia_mod_sev ~ log_tasa_extorsiones +
ingreso_pc + gasto_alim_pc + rural + indigena_hogar +
ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_secuestros <- feols(
ia_mod_sev ~ log_tasa_secuestros +
ingreso_pc + gasto_alim_pc + rural + indigena_hogar +
ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_componentes <- feols(
ia_mod_sev ~ log_tasa_homicidios + log_tasa_extorsiones + log_tasa_secuestros +
ingreso_pc + gasto_alim_pc + rural + indigena_hogar +
ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
etable(
modelo_homicidios, modelo_extorsiones, modelo_secuestros, modelo_componentes,
se.below = TRUE,
fitstat = ~ n + r2,
title = "Tabla 2. Descomposición por tipo de delito (LPM, IA mod/sev)",
notes = "Cada columna sustituye el agregado de violencia por un componente. La columna (4) incluye los tres juntos."
)11.4 Modelos para inseguridad alimentaria severa
La IA severa (≥8 respuestas afirmativas ELCSA) captura privación más grave: hogares donde adultos o niños redujeron porciones, se acostaron con hambre, o pasaron un día sin comer. Es una variable de menor prevalencia y mayor relevancia clínica.
modelo_sesnsp_severa <- feols(
ia_severa ~ log_tasa_violencia +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_acled_severa <- feols(
ia_severa ~ log_acled_territorial +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_conjunto_severa <- feols(
ia_severa ~ log_tasa_violencia + log_acled_territorial +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur, weights = ~factor, cluster = ~cvegeo_mun
)
etable(
modelo_sesnsp_severa, modelo_acled_severa, modelo_conjunto_severa,
se.below = TRUE,
fitstat = ~ n + r2,
title = "Tabla 3. Violencia municipal e inseguridad alimentaria SEVERA (LPM)",
notes = "Las variables de ingreso y gasto entran en log1p para reducir la influencia de valores extremos."
)11.5 Modelos probit con efectos fijos de alta dimensión
Por qué feglm() y no glm()
base. R base puede correr un probit con
glm(family = binomial(link="probit")) y absorber los
efectos fijos como factores factor(anio) + factor(cve_ent),
pero (a) los errores estándar agrupados a nivel municipal requieren
paquetes adicionales y manipulación; (b) si en el futuro se quiere
incluir efectos fijos de municipio (cientos de niveles),
glm se vuelve impráctico. fixest::feglm()
maneja la combinación de FE de alta dimensión, clustering, y probit en
una sola llamada, y es el estándar actual en econometría aplicada.
# Probit: solo SESNSP
modelo_probit_sesnsp <- feglm(
ia_mod_sev ~ log_tasa_violencia +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun,
family = binomial(link = "probit")
)
# Probit: solo ACLED dummy
modelo_probit_acled <- feglm(
ia_mod_sev ~ acled_territorial_dummy +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun,
family = binomial(link = "probit")
)
# Probit: conjunto SESNSP + ACLED dummy
modelo_probit_conjunto <- feglm(
ia_mod_sev ~ log_tasa_violencia + acled_territorial_dummy +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun,
family = binomial(link = "probit")
)
# Probit: conjunto con ACLED categórico
modelo_probit_conjunto_cat <- feglm(
ia_mod_sev ~ log_tasa_violencia + i(acled_territorial_cat, ref = "Sin eventos") +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun,
family = binomial(link = "probit")
)
# Probit: IA severa
modelo_probit_severa <- feglm(
ia_severa ~ log_tasa_violencia + acled_territorial_dummy +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun,
family = binomial(link = "probit")
)
etable(
modelo_probit_sesnsp, modelo_probit_acled,
modelo_probit_conjunto, modelo_probit_conjunto_cat,
modelo_probit_severa,
se.below = TRUE,
fitstat = ~ n + pr2,
title = "Tabla 4. Modelos probit con efectos fijos: violencia e inseguridad alimentaria",
notes = "Pseudo-R² de McFadden. Los coeficientes NO son interpretables directamente; ver Tabla 5 para efectos marginales."
)11.6 Efectos marginales promedio (AME) del probit
Los coeficientes probit no son interpretables sin transformación. El efecto marginal promedio (AME) sí lo es: indica el cambio promedio en la probabilidad de IA por un aumento unitario en la covariable, manteniendo el resto en sus valores observados.
# Chunk 35: ame-probit
library(marginaleffects)
# 1. Volvemos a estimar el modelo conjunto, pero con 'glm' normal
modelo_glm_conjunto <- glm(
ia_mod_sev ~ log_tasa_violencia + acled_territorial_dummy +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina +
as.factor(anio) + as.factor(cve_ent), # Efectos fijos directos
data = base_modelos_sur,
weights = factor,
family = binomial(link = "probit")
)
# 2. Ahora sí calculamos los AME (ya no marcará error)
# Filtramos para que no nos calcule los de todos los estados y años
ame_probit_conjunto <- avg_slopes(
modelo_glm_conjunto,
variables = c("log_tasa_violencia", "acled_territorial_dummy",
"log1p(ingreso_pc)", "log1p(gasto_alim_pc)",
"rural", "indigena_hogar", "ninos_0_5", "jefatura_femenina")
)
# Ver los resultados
summary(ame_probit_conjunto)#> term contrast estimate std.error
#> Length : 6 Length :6 Min. :-0.02082 Min. :0.000270
#> N.unique : 6 N.unique :2 1st Qu.: 0.00372 1st Qu.:0.000422
#> N.blank : 0 N.blank :0 Median : 0.01920 Median :0.000486
#> Min.nchar: 5 Min.nchar:5 Mean : 0.01790 Mean :0.000453
#> Max.nchar:23 Max.nchar:5 3rd Qu.: 0.02799 3rd Qu.:0.000503
#> Max. : 0.06036 Max. :0.000567
#> statistic p.value s.value conf.low
#> Min. :-43.99 Min. :0.000000 Min. : 10 Min. :-0.02175
#> 1st Qu.: 5.84 1st Qu.:0.000000 1st Qu.:Inf 1st Qu.: 0.00285
#> Median : 47.53 Median :0.000000 Median :Inf Median : 0.01838
#> Mean : 39.94 Mean :0.000159 Mean :Inf Mean : 0.01701
#> 3rd Qu.: 69.55 3rd Qu.:0.000000 3rd Qu.:Inf 3rd Qu.: 0.02713
#> Max. :119.73 Max. :0.000952 Max. :Inf Max. : 0.05937
#> conf.high
#> Min. :-0.01990
#> 1st Qu.: 0.00459
#> Median : 0.02002
#> Mean : 0.01879
#> 3rd Qu.: 0.02886
#> Max. : 0.06135Lectura conjunta de Tablas 1, 3 y 5. Si los AMEs del probit (Tabla 5) son cuantitativamente similares a los coeficientes LPM (Tabla 1), entonces la elección entre las dos especificaciones es de presentación, no de inferencia: el patrón es robusto. Si difieren sustantivamente, eso señala que la linealidad implícita en el LPM está fallando en algún tramo del soporte del predictor —típicamente en la cola superior, donde varios municipios tienen tasas extremas. En el manuscrito recomiendo reportar LPM en el cuadro principal y los AMEs en el anexo de robustez.
11.7 Modelo con interacción violencia × indigeneidad
Una hipótesis sustantiva del proyecto es que la violencia no afecta a todos los hogares por igual: las comunidades indígenas, con menor acceso a redes formales de protección social y mercados más localizados, podrían absorber el shock con menos resiliencia.
modelo_interaccion_indigena <- feols(
ia_severa ~ log_tasa_homicidios * indigena_hogar +
log1p(ingreso_pc) + gasto_alim_pc +
rural + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
modelo_interaccion_acled_indigena <- feols(
ia_severa ~ log_acled_territorial * indigena_hogar +
log1p(ingreso_pc) + gasto_alim_pc +
rural + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_sur,
weights = ~factor,
cluster = ~cvegeo_mun
)
etable(
modelo_interaccion_indigena, modelo_interaccion_acled_indigena,
se.below = TRUE,
fitstat = ~ n + r2,
title = "Tabla 7. Heterogeneidad por hogar indígena: interacción violencia × indigeneidad",
notes = "Coeficiente sobre la interacción mide efecto adicional en hogares indígenas vs no indígenas."
)11.8 Modelos con eventos documentados (shocks fuertes)
Como complemento a las medidas continuas, construyo una variable de exposición a episodios “fuertes” de violencia o desplazamiento documentados por organizaciones civiles, CNDH, CMDPDH y prensa. Esta variable identifica municipios donde algo claramente sucedió —desplazamientos masivos, disputas armadas entre cárteles, asesinatos de liderazgos comunitarios— en una ventana temporal específica.
eventos_documentados_raw <- tribble(
~entidad_evento, ~municipio_evento, ~anio_evento, ~mes_evento, ~tipo_evento, ~desplazamiento_doc, ~criminal_doc, ~confianza, ~fuente_corta, ~nota,
# Chiapas: corredor fronterizo y disputa CJNG-Sinaloa
"Chiapas", "Frontera Comalapa", 2023, 5,
"Violencia criminal y desplazamiento en frontera sur", 1, 1, "alta",
"El País / Frayba, 2023",
"Violencia entre grupos criminales, reclutamiento forzado, bloqueos y desplazamiento.",
"Chiapas", "Chicomuselo", 2024, 1,
"Desplazamiento por disputa CJNG-Sinaloa", 1, 1, "alta", "La Jornada, 2024",
"Familias abandonan comunidades por violencia atribuida a CJNG y Cártel de Sinaloa.",
"Chiapas", "La Concordia", 2024, 1,
"Desplazamiento por disputa CJNG-Sinaloa", 1, 1, "alta", "La Jornada, 2024",
"Familias desplazadas por violencia de cárteles y militarización posterior.",
"Chiapas", "Socoltenango", 2024, 1,
"Desplazamiento por disputa CJNG-Sinaloa", 1, 1, "alta", "La Jornada, 2024",
"Familias desplazadas por violencia en la sierra de Chiapas.",
"Chiapas", "Chicomuselo", 2024, 5,
"Enfrentamiento criminal y desplazamiento", 1, 1, "alta", "Proceso, 2024",
"Enfrentamientos entre CJNG y Cártel de Sinaloa intensifican desplazamiento.",
"Chiapas", "Motozintla", 2024, 8,
"Violencia regional Sierra Mariscal", 1, 1, "media", "El País, 2024",
"Zona reportada dentro de la escalada regional con desplazamientos y control criminal.",
# Guerrero: Sierra (Leonardo Bravo, Heliodoro Castillo, Totolapan)
"Guerrero", "Leonardo Bravo", 2018, 11,
"Desplazamiento por violencia criminal en Sierra", 1, 1, "alta", "CMDPDH, 2018",
"Desplazamientos desde comunidades de Leonardo Bravo; destino Chichihualco y Tlacotepec.",
"Guerrero", "General Heliodoro Castillo", 2018, 11,
"Recepción y zona de conflicto Sierra Guerrero", 1, 1, "media", "CMDPDH, 2018",
"Tlacotepec aparece como destino y nodo de disputa.",
"Guerrero", "Leonardo Bravo", 2020, 3,
"Desplazamiento por enfrentamientos armados", 1, 1, "alta", "CMDPDH / prensa, 2020",
"Enfrentamientos entre Cártel del Sur y grupo de Tlacotepec desplazan población.",
"Guerrero", "Chilpancingo de los Bravo", 2020, 3,
"Desplazamiento corredor Chichihualco", 1, 1, "media", "CMDPDH / prensa, 2020",
"Comunidades vinculadas al corredor Chichihualco-Chilpancingo afectadas.",
"Guerrero", "San Miguel Totolapan", 2020, 1,
"Zona estructural de desplazamiento", 1, 1, "media", "Argüello Cabrera, 2022",
"Literatura identifica concentración de desplazamientos en territorios amapoleros/mineros.",
# Oaxaca: región Triqui
"Oaxaca", "Santiago Juxtlahuaca", 2021, 1,
"Desplazamiento forzado Triqui en Tierra Blanca Copala", 1, 0, "alta", "CNDH 36/2022",
"Agresiones armadas y desplazamiento de familias triquis de Tierra Blanca Copala.",
"Oaxaca", "Santiago Juxtlahuaca", 2024, 11,
"Violencia agravada en región Triqui", 1, 0, "media", "El País, 2024",
"Ola de asesinatos y ataques en región triqui; contexto de desplazamiento y disputa.",
"Oaxaca", "Putla Villa de Guerrero", 2024, 11,
"Asesinato y violencia contra liderazgos Triquis", 0, 0, "media", "El País, 2024",
"Asesinato de profesor y entrenador triqui en contexto de violencia regional."
)
# Empate con catálogo municipal
catalogo_eventos <- catalogo_municipal %>%
mutate(
entidad_clean = limpiar_nombre(entidad),
municipio_clean = limpiar_nombre(municipio)
) %>%
dplyr::select(cvegeo_mun, cve_ent, entidad, municipio,
entidad_clean, municipio_clean) %>%
distinct()
eventos_documentados <- eventos_documentados_raw %>%
mutate(
entidad_clean = limpiar_nombre(entidad_evento),
municipio_clean = limpiar_nombre(municipio_evento)
) %>%
left_join(catalogo_eventos, by = c("entidad_clean", "municipio_clean")) %>%
filter(!is.na(cvegeo_mun))
# Timing contemporáneo y rezagado
eventos_documentados_timing <- eventos_documentados %>%
mutate(
anio_enigh_contemporaneo = case_when(
anio_evento %in% c(2019, 2020) ~ 2020L,
anio_evento %in% c(2021, 2022) ~ 2022L,
anio_evento %in% c(2023, 2024) ~ 2024L,
TRUE ~ NA_integer_
),
anio_enigh_rezagado = case_when(
anio_evento %in% c(2019, 2020) ~ 2022L,
anio_evento %in% c(2021, 2022) ~ 2024L,
anio_evento %in% c(2023, 2024) ~ NA_integer_,
TRUE ~ NA_integer_
)
)
eventos_contemporaneos <- eventos_documentados_timing %>%
filter(!is.na(anio_enigh_contemporaneo)) %>%
group_by(cvegeo_mun, anio = anio_enigh_contemporaneo) %>%
summarise(
dummy_evento_fuerte_contemp = 1L,
dummy_desplazamiento_contemp = as.integer(any(desplazamiento_doc == 1)),
dummy_criminal_contemp = as.integer(any(criminal_doc == 1)),
n_eventos_contemp = n(),
.groups = "drop"
)
eventos_rezagados <- eventos_documentados_timing %>%
filter(!is.na(anio_enigh_rezagado)) %>%
group_by(cvegeo_mun, anio = anio_enigh_rezagado) %>%
summarise(
dummy_evento_fuerte_lag = 1L,
dummy_desplazamiento_lag = as.integer(any(desplazamiento_doc == 1)),
dummy_criminal_lag = as.integer(any(criminal_doc == 1)),
n_eventos_lag = n(),
.groups = "drop"
)
# Integrar a base de hogares
base_hogar_eventos_documentados <- base_hogar_violencia_acled_sur %>%
mutate(
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
anio = as.integer(anio)
) %>%
left_join(eventos_contemporaneos, by = c("cvegeo_mun", "anio")) %>%
left_join(eventos_rezagados, by = c("cvegeo_mun", "anio")) %>%
mutate(
across(
c(dummy_evento_fuerte_contemp, dummy_desplazamiento_contemp,
dummy_criminal_contemp, n_eventos_contemp,
dummy_evento_fuerte_lag, dummy_desplazamiento_lag,
dummy_criminal_lag, n_eventos_lag),
~ replace_na(.x, 0)
)
)
base_modelos_eventos_doc <- base_hogar_eventos_documentados %>%
mutate(
anio = as.factor(anio),
cve_ent = as.factor(cve_ent),
cvegeo_mun = as.factor(cvegeo_mun)
)
modelo_evento_contemp <- feols(
ia_mod_sev ~ dummy_evento_fuerte_contemp +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_eventos_doc, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_evento_lag <- feols(
ia_mod_sev ~ dummy_evento_fuerte_lag +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_eventos_doc, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_desplazamiento_lag <- feols(
ia_mod_sev ~ dummy_desplazamiento_lag +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_eventos_doc, weights = ~factor, cluster = ~cvegeo_mun
)
modelo_evento_lag_sesnsp <- feols(
ia_mod_sev ~ dummy_evento_fuerte_lag + log_tasa_violencia +
ingreso_pc + gasto_alim_pc +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina |
anio + cve_ent,
data = base_modelos_eventos_doc, weights = ~factor, cluster = ~cvegeo_mun
)
etable(
modelo_evento_contemp, modelo_evento_lag,
modelo_desplazamiento_lag, modelo_evento_lag_sesnsp,
se.below = TRUE,
fitstat = ~ n + r2,
title = "Tabla 8. Eventos documentados (shocks fuertes) e inseguridad alimentaria mod/sev",
notes = "Contemporáneo: evento en la ventana del levantamiento. Rezagado: evento en la ventana anterior, afectando al siguiente levantamiento."
)Por qué importan los modelos con eventos
documentados. Los conteos del SESNSP y ACLED pueden subreportar
episodios de exposición fuerte —especialmente desplazamientos forzados,
que las autoridades locales tienen incentivos para no registrar—. La
variable dummy_evento_fuerte_lag identifica municipios
donde organizaciones civiles, CNDH o prensa especializada documentaron
episodios graves en la ventana ANTERIOR al levantamiento ENIGH. Su
coeficiente recoge el efecto rezagado de un shock de exposición clara
—el mecanismo más cercano a un experimento natural que permite la
combinación de fuentes disponibles—.
12 Diagnóstico de los modelos
12.1 Multicolinealidad
# Modelo equivalente para VIF (sin FE absorbidos para que car::vif funcione)
modelo_para_vif <- lm(
ia_severa ~ log_tasa_violencia + log_acled_territorial +
log1p(ingreso_pc) + log1p(gasto_alim_pc) +
rural + indigena_hogar + ninos_0_5 + jefatura_femenina +
factor(anio) + factor(cve_ent),
data = base_modelos_sur
)
vif_resultados <- car::vif(modelo_para_vif)
print(vif_resultados)#> GVIF Df GVIF^(1/(2*Df))
#> log_tasa_violencia 2.66 1 1.63
#> log_acled_territorial 2.42 1 1.56
#> log1p(ingreso_pc) 1.65 1 1.28
#> log1p(gasto_alim_pc) 1.46 1 1.21
#> rural 1.36 1 1.17
#> indigena_hogar 1.29 1 1.13
#> ninos_0_5 1.05 1 1.03
#> jefatura_femenina 1.01 1 1.01
#> factor(anio) 1.15 2 1.04
#> factor(cve_ent) 3.45 6 1.11Interpretación del VIF. Un VIF mayor a 5 (criterio
liberal) o 10 (criterio conservador) indica que la varianza del
coeficiente está siendo inflada por correlación con otros predictores.
Si log_tasa_violencia y log_acled_territorial
tienen VIFs altos en el modelo conjunto, eso significa que SESNSP y
ACLED están midiendo aspectos correlacionados del mismo proceso —no que
la estimación sea inválida, sino que los coeficientes individuales en el
modelo conjunto deben leerse con cuidado (preferir los modelos separados
para magnitudes).
12.2 Residuales y forma funcional
# Asignación robusta de residuales sobre la muestra de estimación
base_modelos_sur_diag <- base_modelos_sur
base_modelos_sur_diag$residuales <- NA_real_
base_modelos_sur_diag$valores_ajustados <- NA_real_
filas_purgadas <- modelo_sesnsp_severa$obs_selection$obsRemoved
if (!is.null(filas_purgadas) && length(filas_purgadas) > 0) {
base_modelos_sur_diag$residuales[-filas_purgadas] <- resid(modelo_sesnsp_severa)
base_modelos_sur_diag$valores_ajustados[-filas_purgadas] <- fitted(modelo_sesnsp_severa)
base_estimacion <- base_modelos_sur_diag[-filas_purgadas, ]
} else {
base_modelos_sur_diag$residuales <- resid(modelo_sesnsp_severa)
base_modelos_sur_diag$valores_ajustados <- fitted(modelo_sesnsp_severa)
base_estimacion <- base_modelos_sur_diag
}
ggplot(base_estimacion, aes(x = log_tasa_violencia, y = residuales)) +
geom_point(alpha = 0.2, color = pal_principal) +
geom_smooth(method = "loess", color = pal_severa, se = FALSE, span = 0.75) +
geom_hline(yintercept = 0, linetype = "dashed", color = "black", linewidth = 0.6) +
theme_minimal(base_size = 12) +
labs(
title = "Residuales vs log(tasa de violencia)",
subtitle = "Detección empírica de sesgos por forma funcional u omisión variable",
x = "log(tasa anualizada de delitos violentos por 100k hab.)",
y = "Residuales idiosincráticos"
) +
theme(plot.title = element_text(face = "bold", color = pal_principal))
Residuales vs log(tasa de violencia): detección de patrones no lineales omitidos
Lectura. Si el loess rojo se aleja
sistemáticamente del cero —especialmente con forma de “U” o
“U-invertida”— hay no-linealidad omitida que el modelo lineal no
captura. Si oscila cerca de cero a lo largo del rango, la forma
funcional es defensible. En caso de no-linealidad, una especificación
con I(log_tasa_violencia^2) (probada en versiones
anteriores del script) o un término spline serían las alternativas.
13 Estimaciones municipales ponderadas
Construir prevalencias municipales requiere subset del diseño muestral. Se restringe a municipios-año con al menos 25 hogares con primera infancia en muestra para evitar estimaciones inestables.
diseno_panel_sur <- svydesign(
ids = ~upm,
strata = ~est_dis,
weights = ~factor,
data = base_hogar_violencia_acled_sur,
nest = TRUE
)
municipios_suficientes <- base_hogar_violencia_acled_sur %>%
count(cvegeo_mun, anio, name = "n_hogares") %>%
filter(n_hogares >= 25) %>%
distinct(cvegeo_mun, anio)
diseno_panel_sur_filtrado <- subset(
diseno_panel_sur,
interaction(cvegeo_mun, anio) %in%
interaction(municipios_suficientes$cvegeo_mun, municipios_suficientes$anio)
)
enigh_mun_anio <- svyby(
~ ia_mod_sev + ia_severa + ingreso_pc + gasto_alim_pc + rural + indigena_hogar,
~ cvegeo_mun + anio,
design = diseno_panel_sur_filtrado,
FUN = svymean,
na.rm = TRUE,
vartype = c("se")
) %>%
as_tibble() %>%
rename(
ia_mod_sev_mun = ia_mod_sev,
ia_severa_mun = ia_severa,
ingreso_pc_mun = ingreso_pc,
gasto_alim_pc_mun = gasto_alim_pc,
rural_mun = rural,
indigena_mun = indigena_hogar
) %>%
mutate(
ia_mod_sev_mun = ia_mod_sev_mun * 100,
ia_severa_mun = ia_severa_mun * 100,
rural_mun = rural_mun * 100,
indigena_mun = indigena_mun * 100
)
cat(sprintf("Observaciones municipio-año estimadas: %d\n", nrow(enigh_mun_anio)))#> Observaciones municipio-año estimadas: 11314 Panel municipal e identificación de municipios con deterioro
Para priorizar municipios candidatos a estudio de caso, construyo un panel municipal “ancho” (una fila por municipio, columnas para 2020, 2022, 2024) que permita calcular trayectorias y un score compuesto de deterioro.
violencia_mun_sur <- violencia_municipal_enigh_tasas %>%
mutate(
cvegeo_mun = str_pad(as.character(cvegeo_mun), 5, pad = "0"),
cve_ent = str_pad(as.character(cve_ent), 2, pad = "0"),
anio = as.integer(anio)
) %>%
filter(cve_ent %in% estados_sur_vec) %>%
left_join(
acled_municipal_enigh %>%
dplyr::select(cvegeo_mun, anio,
eventos_acled_total_24m,
violencia_territorial_acled_24m,
violencia_civiles_acled_24m,
batallas_acled_24m,
fatalidades_acled_24m),
by = c("cvegeo_mun", "anio")
) %>%
mutate(
across(
c(eventos_acled_total_24m, violencia_territorial_acled_24m,
violencia_civiles_acled_24m, batallas_acled_24m,
fatalidades_acled_24m),
~ replace_na(.x, 0)
)
)
enigh_mun_wide <- enigh_mun_anio %>%
dplyr::select(cvegeo_mun, anio, ia_mod_sev_mun, ia_severa_mun,
rural_mun, indigena_mun) %>%
pivot_wider(
names_from = anio,
values_from = c(ia_mod_sev_mun, ia_severa_mun, rural_mun, indigena_mun),
names_sep = "_"
)
violencia_mun_wide <- violencia_mun_sur %>%
dplyr::select(cvegeo_mun, anio,
tasa_delitos_violentos_anualizada,
tasa_homicidios_dolosos_anualizada,
tasa_extorsiones_anualizada,
tasa_secuestros_anualizada,
eventos_acled_total_24m,
violencia_territorial_acled_24m,
violencia_civiles_acled_24m,
batallas_acled_24m,
fatalidades_acled_24m) %>%
pivot_wider(
names_from = anio,
values_from = c(tasa_delitos_violentos_anualizada,
tasa_homicidios_dolosos_anualizada,
tasa_extorsiones_anualizada,
tasa_secuestros_anualizada,
eventos_acled_total_24m,
violencia_territorial_acled_24m,
violencia_civiles_acled_24m,
batallas_acled_24m,
fatalidades_acled_24m),
names_sep = "_"
)
panel_mun_completo <- enigh_mun_wide %>%
left_join(violencia_mun_wide, by = "cvegeo_mun")14.1 Score compuesto de deterioro
candidatos_deterioro <- panel_mun_completo %>%
mutate(
cambio_ia_severa_2020_2024 = ia_severa_mun_2024 - ia_severa_mun_2020,
cambio_ia_modsev_2020_2024 = ia_mod_sev_mun_2024 - ia_mod_sev_mun_2020,
ia_severa_post_max = pmax(ia_severa_mun_2022, ia_severa_mun_2024, na.rm = TRUE),
deterioro_ia_severa_post = ia_severa_post_max - ia_severa_mun_2020,
acled_territorial_post_max = pmax(violencia_territorial_acled_24m_2022,
violencia_territorial_acled_24m_2024, na.rm = TRUE),
homicidios_post_max = pmax(tasa_homicidios_dolosos_anualizada_2022,
tasa_homicidios_dolosos_anualizada_2024, na.rm = TRUE),
aumento_acled_territorial_post = acled_territorial_post_max - violencia_territorial_acled_24m_2020,
aumento_homicidios_post = homicidios_post_max - tasa_homicidios_dolosos_anualizada_2020
) %>%
filter(!is.na(deterioro_ia_severa_post), deterioro_ia_severa_post > 0) %>%
mutate(
score_deterioro =
as.numeric(scale(deterioro_ia_severa_post)) +
as.numeric(scale(aumento_acled_territorial_post)) +
as.numeric(scale(aumento_homicidios_post))
) %>%
arrange(desc(score_deterioro))
candidatos_deterioro %>%
dplyr::select(cvegeo_mun, ia_severa_mun_2020, ia_severa_mun_2022, ia_severa_mun_2024,
deterioro_ia_severa_post, acled_territorial_post_max,
homicidios_post_max, rural_mun_2020, indigena_mun_2020,
score_deterioro) %>%
slice_head(n = 20) %>%
kable(
caption = "Tabla 9. Top 20 municipios por score compuesto de deterioro (IA + violencia, 2020→post)",
digits = 2
) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
full_width = FALSE)| cvegeo_mun | ia_severa_mun_2020 | ia_severa_mun_2022 | ia_severa_mun_2024 | deterioro_ia_severa_post | acled_territorial_post_max | homicidios_post_max | rural_mun_2020 | indigena_mun_2020 | score_deterioro |
|---|---|---|---|---|---|---|---|---|---|
| 12035 | 13.65 | — | 19.65 | 6.00 | 68 | 93.52 | 29.83 | 6.89 | 3.07 |
| 07101 | 8.15 | 20.4 | 3.92 | 12.24 | 33 | 5.98 | 0.00 | 11.23 | 2.53 |
| 04011 | 8.37 | 23.4 | 14.95 | 15.05 | 11 | 29.09 | 72.60 | 25.56 | 2.42 |
| 04009 | 14.02 | 26.6 | 15.20 | 12.62 | 19 | 17.25 | 62.26 | 10.83 | 1.24 |
| 23002 | 7.37 | 11.7 | 4.98 | 4.32 | 7 | 25.28 | 35.73 | 70.90 | 0.32 |
| 27012 | 40.98 | 43.5 | 29.22 | 2.54 | 30 | 18.93 | 65.16 | 8.15 | -0.84 |
| 31096 | 1.59 | — | 6.09 | 4.49 | 4 | 1.19 | 23.29 | 63.14 | -0.90 |
| 04006 | 11.42 | — | 13.16 | 1.73 | 3 | 2.21 | 56.16 | 55.50 | -1.77 |
| 04001 | 11.50 | 14.0 | 10.22 | 2.51 | 1 | 2.32 | 12.93 | 83.38 | -1.79 |
| 12029 | 17.62 | 27.0 | 13.97 | 9.34 | 74 | 43.02 | 6.95 | 7.28 | -2.14 |
| 27006 | 27.31 | 28.2 | — | 0.93 | 28 | 13.98 | 81.37 | 1.93 | -2.15 |
candidatos_rural_indigena <- candidatos_deterioro %>%
filter(rural_mun_2020 >= 40 | indigena_mun_2020 >= 30) %>%
arrange(desc(score_deterioro))
candidatos_rural_indigena %>%
dplyr::select(cvegeo_mun, rural_mun_2020, indigena_mun_2020,
ia_severa_mun_2020, ia_severa_mun_2024,
deterioro_ia_severa_post, score_deterioro) %>%
slice_head(n = 15) %>%
kable(
caption = "Tabla 10. Municipios rurales o indígenas con mayor deterioro",
digits = 2
) %>%
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| cvegeo_mun | rural_mun_2020 | indigena_mun_2020 | ia_severa_mun_2020 | ia_severa_mun_2024 | deterioro_ia_severa_post | score_deterioro |
|---|---|---|---|---|---|---|
| 04011 | 72.6 | 25.56 | 8.37 | 14.95 | 15.05 | 2.42 |
| 04009 | 62.3 | 10.83 | 14.02 | 15.20 | 12.62 | 1.24 |
| 23002 | 35.7 | 70.90 | 7.37 | 4.98 | 4.32 | 0.32 |
| 27012 | 65.2 | 8.15 | 40.98 | 29.22 | 2.54 | -0.84 |
| 31096 | 23.3 | 63.14 | 1.59 | 6.09 | 4.49 | -0.90 |
| 04006 | 56.2 | 55.50 | 11.42 | 13.16 | 1.73 | -1.77 |
| 04001 | 12.9 | 83.38 | 11.50 | 10.22 | 2.51 | -1.79 |
| 27006 | 81.4 | 1.93 | 27.31 | — | 0.93 | -2.15 |
Uso de este ranking. El score compuesto no estima un efecto causal. Es un instrumento de priorización: identifica municipios donde la combinación de aumento en IA severa y aumento en exposición a violencia justifica revisión cualitativa, entrevistas en campo, o análisis de archivo. Es la puerta de entrada a la fase de estudio de caso del manuscrito, no su conclusión.
15 Trayectorias y narrativas
top_10_deterioro <- candidatos_deterioro %>% slice_head(n = 10) %>% pull(cvegeo_mun)
catalogo_sur <- violencia_mun_sur %>%
distinct(cvegeo_mun, entidad, municipio)
trayectorias_deterioro <- enigh_mun_anio %>%
filter(cvegeo_mun %in% top_10_deterioro) %>%
left_join(catalogo_sur, by = "cvegeo_mun") %>%
mutate(label_mun = paste0(municipio, ", ", entidad))
plot_trayectorias <- ggplot(
trayectorias_deterioro,
aes(x = as.integer(as.character(anio)), y = ia_severa_mun)
) +
geom_line(linewidth = 0.9, color = pal_severa) +
geom_point(size = 2.8, color = pal_severa) +
facet_wrap(~label_mun, scales = "free_y", ncol = 2) +
scale_x_continuous(breaks = c(2020, 2022, 2024)) +
labs(
title = "Trayectorias de IA severa en municipios con mayor deterioro",
subtitle = "Top 10 municipios por score compuesto (IA + violencia territorial)",
x = "Año ENIGH",
y = "IA severa ponderada (%)",
caption = "Fuente: ENIGH 2020–2024 (INEGI), SESNSP, CONAPO y ACLED.\nMunicipios con ≥25 hogares con primera infancia por año."
) +
theme_minimal(base_size = 11) +
theme(
plot.title = element_text(face = "bold", size = 12, color = pal_principal),
strip.text = element_text(face = "bold", size = 9),
panel.grid.minor = element_blank()
)
plot_trayectorias
Trayectorias municipales de IA severa, top 10 municipios por score de deterioro
ggsave("trayectorias_municipios_deterioro_sur.png", plot_trayectorias,
width = 12, height = 10, dpi = 300)15.1 Notas cualitativas de ACLED para municipios prioritarios
top_15_casos <- candidatos_rural_indigena %>% slice_head(n = 15) %>% pull(cvegeo_mun)
notas_acled_casos <- acled_sur_municipalizado %>%
filter(
cvegeo_mun %in% top_15_casos,
anio_evento %in% 2019:2024,
event_type %in% c("Violence against civilians", "Battles", "Explosions/Remote violence")
) %>%
mutate(nota_corta = str_trunc(notes, 400)) %>%
dplyr::select(event_date, anio_evento, cvegeo_mun, entidad, municipio,
location, event_type, sub_event_type, actor1, actor2,
fatalities, nota_corta) %>%
arrange(cvegeo_mun, event_date)
write_csv(notas_acled_casos, "notas_acled_municipios_prioritarios.csv")
cat(sprintf("✓ Notas de ACLED exportadas: %d eventos\n", nrow(notas_acled_casos)))#> ✓ Notas de ACLED exportadas: 234 eventos16 Exportación de resultados
write_csv(candidatos_deterioro,
"municipios_deterioro_ia_violencia_2020_2024.csv")
write_csv(candidatos_rural_indigena,
"municipios_deterioro_rural_indigena.csv")
write_csv(enigh_mun_anio,
"estimaciones_municipales_ia_2020_2024.csv")
write_csv(ia_estados_panel,
"ia_primera_infancia_estados_2020_2024.csv")
write_csv(base_hogar_violencia_acled_sur,
"base_hogar_violencia_acled_sur.csv")17 Síntesis y próximos pasos
Lo que el análisis dice y lo que no dice.
Dice que en los siete estados del sur, durante 2020–2024, los hogares con primera infancia residentes en municipios con mayor violencia territorial —medida tanto por la tasa SESNSP de delitos violentos como por el conteo de eventos territoriales ACLED— enfrentan probabilidades sistemáticamente más altas de IA moderada/severa y de IA severa, controlando por ingreso, gasto alimentario, ruralidad, etnicidad, presencia de menores y jefatura del hogar; absorbiendo heterogeneidades persistentes entre entidades y shocks comunes a cada año. Esta asociación sobrevive a la sustitución del agregado por sus componentes, a la transformación del agregado en categorías de intensidad, y a la verificación con modelos probit y sus efectos marginales promedio.
No dice que la violencia “cause” la inseguridad alimentaria en el sentido contrafactual estricto. El diseño es de panel repetido de cortes transversales con efectos fijos de dos vías; los efectos identificados aprovechan la variación municipal dentro de cada estado a lo largo del tiempo, descontando shocks nacionales. Eso es un control fuerte, pero no equivale a aleatorización. Hay tres caminos para fortalecer la inferencia causal en pasos posteriores: (a) explotar shocks discretos de violencia (eventos documentados) con event study; (b) construir un contrafactual sintético para municipios con escalada brusca; (c) emparejar municipios con propensity score sobre características previas a 2020.
Outputs generados por este pipeline (en el directorio de trabajo):
ia_primera_infancia_estados_2020_2024.csv— panel estatal ponderado.estimaciones_municipales_ia_2020_2024.csv— panel municipal ponderado.municipios_deterioro_ia_violencia_2020_2024.csv— todos los municipios con deterioro.municipios_deterioro_rural_indigena.csv— subconjunto rural/indígena.notas_acled_municipios_prioritarios.csv— narrativas ACLED.base_hogar_violencia_acled_sur.csv— base analítica final.ranking_ia_estados_2024.png,trayectorias_municipios_deterioro_sur.png.
#> ================================================================================#> ANÁLISIS COMPLETADO — Inseguridad alimentaria y violencia territorial 2020-2024#> ================================================================================cat(sprintf("Universo: %d hogares con primera infancia en %d municipios de %d entidades\n",
nrow(base_hogar_violencia_acled_sur),
n_distinct(base_hogar_violencia_acled_sur$cvegeo_mun),
n_distinct(base_hogar_violencia_acled_sur$cve_ent)))#> Universo: 12688 hogares con primera infancia en 469 municipios de 7 entidades#> Periodo: ENIGH 2020, 2022, 2024 con ventanas bianuales de violencia#> PRÓXIMOS PASOS METODOLÓGICOS:#> • Event study sobre municipios con shock discreto de violencia#> • Synthetic control method para municipios con escalada brusca#> • Análisis de heterogeneidad por composición indígena × ruralidad#> • Incorporar pobreza/rezago social CONEVAL como covariables#> • Validación cualitativa con notas ACLED y trabajo de archivo#> ================================================================================Análisis realizado por Christian Campa, Licenciatura en Economía,
Tecnológico de Monterrey. Documento generado con R 4.6.0 y
fixest 0.14.1.