Introducción

El conflicto armado en Colombia constituye una de las expresiones más persistentes y territorializadas de violencia en América Latina. Desde mediados del siglo XX, ha emergido como resultado de una compleja interacción entre la desigualdad en el acceso a la tierra, la debilidad del Estado en regiones, la expansión de economías ilegales y la presencia intermitente de grupos armados. Si bien el Acuerdo de Paz firmado en 2016 entre el Gobierno Nacional y las FARC-EP representó un hito histórico para la desmovilización de uno de los principales actores armados, aun así, la violencia no desapareció de manera homogénea en el territorio colombiano.

La firma del Acuerdo Final de Paz en 2016 entre el Gobierno Nacional y las FARC-EP representó un hito histórico al permitir la desmovilización del grupo armado más antiguo del hemisferio. Sin embargo, el fin del conflicto con las FARC no significó el cese de la violencia. Por el contrario, muchos territorios han experimentado una recomposición del conflicto, con la reconfiguración de actores armados y la persistencia de amenazas a la seguridad humana.

En este escenario, el Estado colombiano creó los Programas de Desarrollo con Enfoque Territorial (PDET) en 2017, como estrategia integral para transformar 170 municipios priorizados por su histórica afectación por el conflicto armado y la pobreza estructural. Según cifras oficiales, estos municipios concentran el 36% de las víctimas del conflicto armado (UARIV, 2023) y el 40% de los cultivos ilícitos del país (UNODC, 2023), lo que evidencia su centralidad en la política de paz territorial. Entre 2018 y 2021 se destinaron más de 4,2 billones de pesos a su implementación (Contraloría, 2022), buscando cerrar brechas en infraestructura, salud, educación y desarrollo rural. Sin embargo, cinco años después del inicio de su implementación, diversos estudios y evaluaciones independientes han evidenciado resultados dispares. Mientras municipios como Briceño (Antioquia) presentan reducciones significativas del conflicto —con una disminución del 56% en el Índice de Incidencia del Conflicto Armado (IICA)— otros como Tumaco (Nariño) mantienen niveles críticos de violencia, lo que evidencia que los resultados han tenido un impacto desigual en el territorio.

Esta variabilidad territorial plantea interrogantes sustanciales: • ¿Por qué municipios con niveles similares de inversión presentan trayectorias opuestas en términos de violencia? • ¿Hasta qué punto las dinámicas del conflicto están condicionadas por factores espaciales, más allá de las decisiones locales o de la asignación de recursos?

En este escenario cobra relevancia el uso del Índice de Incidencia del Conflicto Armado (IICA), desarrollado metodológicamente por el Departamento Nacional de Planeación (DNP) y adoptado en instrumentos como el CONPES 3932 de 2018. El IICA permite condensar la magnitud del conflicto a nivel municipal mediante unas variables clave: acciones armadas, homicidios, secuestros, desplazamiento forzado y presencia de minas antipersonal, integradas a través de técnicas de normalización y ponderación. No obstante, si bien el IICA captura la intensidad del conflicto de forma puntual en cada municipio, poco se ha analizado sobre su comportamiento espacial y su relación con los municipios vecinos. ¿Es el conflicto un fenómeno territorialmente aislado o responde a lógicas regionales de propagación? Para abordar esta pregunta, este trabajo propone un enfoque de estadística espacial, que permite evaluar patrones de autocorrelación y dependencia geográfica del conflicto armado entre los años 2017 a 2021. A través de herramientas como el índice local de Moran (LISA) y el modelo espacial autorregresivo (SLM), se busca evidenciar si el conflicto persiste y se propaga con base en relaciones espaciales significativas, o si su comportamiento es aleatorio a nivel municipal.

Este estudio busca entonces, combinando el uso del IICA con técnicas de análisis espacial, el objetivo de identificar cómo se distribuye el conflicto armado en el territorio nacional, qué patrones de concentración o dispersión se observan, y en qué medida la pertenencia a los PDET incide en estas configuraciones. A través de un enfoque espaciotemporal y utilizando datos oficiales entre 2017 a 2021, se pretende aportar evidencia sobre la estructura territorial del conflicto armado y su evolución reciente.

Datos y preparación

library(sf)
## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
library(readxl)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(stringr)
library(spdep)
## Loading required package: spData
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
library(tmap)
library(leaflet)
library(plotly)
## Loading required package: ggplot2
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(tibble)
library(spatialreg)
## Loading required package: Matrix
## 
## Attaching package: 'spatialreg'
## The following objects are masked from 'package:spdep':
## 
##     get.ClusterOption, get.coresOption, get.mcOption,
##     get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
##     set.coresOption, set.mcOption, set.VerboseOption,
##     set.ZeroPolicyOption
# rutas
path_shape <- "/Users/daniel/Desktop/MGN_2021_COLOMBIA/ADMINISTRATIVO/MGN_MPIO_POLITICO.shp"
path_excel <- "/Users/daniel/Desktop/Base IICA.xlsx"

# shape municipal DANE
Shape <- st_read(path_shape, quiet = TRUE) |>
  mutate(cod_dane = str_pad(as.character(MPIO_CDPMP), 5, pad = "0"))

# base IICA (2017–2021)
Base_IICA <- read_excel(path_excel)

# llave común
if ("Divipola #" %in% names(Base_IICA)) {
  Base_IICA <- Base_IICA |>
    mutate(cod_dane = sprintf("%05.0f", as.numeric(`Divipola #`)))
} else if ("Divipola" %in% names(Base_IICA)) {
  Base_IICA <- Base_IICA |>
    mutate(cod_dane = str_pad(str_replace_all(Divipola, "\\D", ""), 5, pad = "0"))
}

# exclusión de islas (San Andrés y Providencia)
drop_divipola <- c("88001","88564")
Shape <- Shape |> filter(!cod_dane %in% drop_divipola)
Base_IICA <- Base_IICA |> filter(!cod_dane %in% drop_divipola)

# unión shape + datos IICA
panel_wide <- Shape |> left_join(Base_IICA, by = "cod_dane")

# chequeo general
panel_wide |> st_drop_geometry() |>
  summarise(across(starts_with("IICA"), ~sum(!is.na(.x))))
##   IICA 2017 IICA 2018 IICA 2019 IICA 2020 IICA 2021 IICA 2017-2021
## 1      1119      1119      1119      1119      1119           1119
##   IICA 2017 Categoría IICA 2018\r\nCategoría IICA 2019\r\nCategoría
## 1                1119                   1119                   1119
##   IICA 2020\r\nCategoría IICA 2021\r\nCategoría IICA 2017-2021 Categoría
## 1                   1119                   1119                     1119

En el proceso de construcción de la base, se decidió excluir los municipios correspondientes a San Andrés (código 88001) y Providencia y Santa Catalina (código 88564) de los análisis estadísticos espaciales. Esta decisión responde a una razón técnica: ambos territorios insulares se encuentran aislados geográficamente del territorio continental colombiano, lo cual impide la conformación de una estructura de vecindad adecuada bajo el criterio de los k vecinos más cercanos (KNN). Dado que los métodos de autocorrelación espacial —como el índice de Moran o el modelo SAR— requieren una matriz de pesos espaciales conectada, mantener municipios sin vecinos genera “islas” en la matriz W que distorsionan los cálculos de dependencia espacial y afectan la estabilidad numérica del modelo.

Por otra parte, San Andrés y Providencia presentan una dinámica del conflicto armado sustancialmente diferente al resto del país. Históricamente, el archipiélago ha mostrado una baja incidencia del conflicto y no ha sido escenario de confrontaciones armadas de los actores tradicionales, por lo que su inclusión no aporta información relevante sobre los patrones de concentración, persistencia o propagación de la violencia en el territorio continental. Por el contrario, su permanencia en la base podría sesgar la interpretación de los índices globales de autocorrelación al incorporar valores atípicos cercanos a cero en un contexto de alta dispersión continental.

Finalmente, la exclusión de estos dos municipios no altera la validez del análisis ni la representatividad territorial del estudio, dado que los 1.119 municipios restantes abarcan la totalidad del territorio continental y concentran la dinámica real del conflicto armado colombiano. Así, la decisión metodológica permite garantizar una matriz espacial conectada, consistente y analíticamente robusta, preservando la comparabilidad de los resultados tanto en el cálculo de los índices globales (Moran I) como en los modelos econométricos espaciales (SAR y SEM).

Matriz espacial

muni_proj <- st_transform(Shape, 3116)
centroids <- st_centroid(st_geometry(muni_proj))
coords    <- st_coordinates(centroids)# vecinos más cercanos
nb_full <- knn2nb(knearneigh(coords, k = 4))
W_full  <- nb2listw(nb_full, style = "W")

Para el análisis de autocorrelación espacial fue necesario definir una estructura que representara las relaciones de proximidad e influencia entre los municipios. Esta estructura se formalizó mediante una matriz de pesos espaciales (W), la cual constituye la base sobre la que se calculan tanto los índices de Moran (global y local) como los modelos econométricos espaciales (SAR y SEM). Se optó por construir dicha matriz utilizando el método de los k vecinos más cercanos (K-nearest neighbors, KNN), con un valor de k = 4, lo que significa que cada municipio se conecta con sus cuatro municipios más próximos en términos de distancia euclidiana.

Autores como Cliff y Ord (1981) y Anselin (1988) han destacado que la estructura de pesos espaciales constituye una hipótesis explícita sobre el espacio: define la forma en que los fenómenos locales influyen entre sí. En este caso, asumir k = 4 implica que el conflicto en un municipio puede estar parcialmente determinado por lo que ocurre en sus cuatro vecinos más próximos, lo cual resulta especialmente pertinente para fenómenos que se propagan territorialmente, como el conflicto armado.

Finalmente, la matriz W obtenida es simétrica y completamente conectada, cumpliendo los requisitos para su uso en pruebas de autocorrelación (Moran global y local) y en modelos espaciales autorregresivos. De este modo, se garantiza que cada municipio tenga influencia recíproca con su entorno y que los resultados derivados de los análisis sean estadísticamente estables y territorialmente consistentes.

Mapa de conexiones

vec_lines <- st_as_sf(spdep::nb2lines(nb_full, coords = coords,
                                      proj4string = st_crs(muni_proj)$proj4string))
pts <- st_as_sf(data.frame(x = coords[,1], y = coords[,2]),
                coords = c("x","y"), crs = st_crs(muni_proj))

tmap_mode("plot")
## ℹ tmap modes "plot" - "view"
## ℹ toggle with `tmap::ttm()`
map_w_k4 <- tm_shape(muni_proj) +
  tm_borders(lwd = 0.2, col = "grey80") +
  tm_shape(vec_lines) + tm_lines(lwd = 0.4, col = "steelblue", alpha = 0.55) +
  tm_shape(pts) + tm_dots(size = 0.05, col = "firebrick3") +
  tm_layout(title = "matriz de vecindad knn (k = 4)", legend.show = FALSE, frame = FALSE)
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_lines()`: use `col_alpha` instead of `alpha`.[v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
map_w_k4

Este mapa permite observar patrones regionales de interacción territorial, donde ciertas zonas presentan una estructura espacial más compacta y otras, en cambio, más dispersa. Por ejemplo, en regiones como el Eje Cafetero, el Valle del Cauca o la zona andina central, la red de conexiones es más densa debido a la alta proximidad entre municipios, lo que sugiere una mayor capacidad de propagación espacial de los fenómenos sociales, como el conflicto armado. En contraste, zonas como la Orinoquía y la Amazonía presentan una red más laxa, con distancias mayores entre centroides municipales, lo cual puede limitar la transmisión espacial de la violencia o hacerla más localizada.

La visualización de la matriz K=4 ofrece una lectura territorial del conflicto armado desde la perspectiva de las interdependencias espaciales. Al representar gráficamente las conexiones entre municipios, se revela el potencial de contagio territorial del conflicto, es decir, la posibilidad de que los altos niveles de violencia en un municipio se asocien con comportamientos similares en sus vecinos inmediatos. Esta visión espacial constituye el punto de partida para los análisis de autocorrelación global y local desarrollados en los apartados siguientes.

Mapa de vecindades y su interpretación

panel_wide <- panel_wide |>
  mutate(IICA_mean = rowMeans(across(c("IICA 2017","IICA 2018","IICA 2019",
                                       "IICA 2020","IICA 2021")), na.rm = TRUE))

cod_top10 <- panel_wide |> arrange(desc(IICA_mean)) |> slice(1:10) |> pull(cod_dane)

plot_ego_municipio_zoom_lbl <- function(target_cod,
                                        buffer_km = 42,
                                        expand_frac = 0.12,
                                        lbl_size_nb = 0.36,
                                        lbl_size_ego = 0.40,
                                        lbl_x_nb = 0.6,
                                        lbl_y_nb = -0.6,
                                        lbl_y_ego = 0.8) {
  i <- which(muni_proj$cod_dane == target_cod)[1]; stopifnot(!is.na(i))
  vec_i <- nb_full[[i]]
  
  ego_pt <- st_as_sf(data.frame(x = coords[i,1], y = coords[i,2],
                                label = muni_proj$MPIO_CNMBR[i]),
                     coords = c("x","y"), crs = st_crs(muni_proj))
  nb_pts <- st_as_sf(data.frame(x = coords[vec_i,1], y = coords[vec_i,2],
                                label = muni_proj$MPIO_CNMBR[vec_i]),
                     coords = c("x","y"), crs = st_crs(muni_proj))
  
  ego_lines <- lapply(vec_i, function(j)
    st_linestring(rbind(coords[i,], coords[j,]))) |> st_sfc(crs = st_crs(muni_proj)) |> st_as_sf()
  
  win <- st_buffer(ego_pt, dist = buffer_km * 1000)
  bb  <- st_bbox(win)
  dx  <- (bb$xmax - bb$xmin) * expand_frac
  dy  <- (bb$ymax - bb$ymin) * expand_frac
  bb_exp <- c(xmin = bb$xmin - dx, ymin = bb$ymin - dy,
              xmax = bb$xmax + dx, ymax = bb$ymax + dy)
  class(bb_exp) <- class(bb); attr(bb_exp, "crs") <- st_crs(win)
  muni_win <- st_crop(muni_proj, bb_exp)
  
  tm_shape(muni_win) + tm_borders(col = "grey80", lwd = 0.5) +
    tm_shape(ego_lines) + tm_lines(col = "darkorange", lwd = 2) +
    tm_shape(nb_pts) + tm_dots(col = "darkorange", size = 0.09) +
    tm_text("label", size = lbl_size_nb, xmod = lbl_x_nb, ymod = lbl_y_nb, shadow = TRUE) +
    tm_shape(ego_pt) + tm_dots(col = "#FF6666", size = 0.11) +
    tm_text("label", size = lbl_size_ego, ymod = lbl_y_ego, shadow = TRUE, fontface = "bold") +
    tm_layout(title = paste0("Ego-red K=4 — ", muni_proj$MPIO_CNMBR[i],
                             " (", muni_proj$DPTO_CNMBR[i], ")"),
              legend.show = FALSE, frame = FALSE)
}

tm_list <- lapply(cod_top10, function(cd)
  plot_ego_municipio_zoom_lbl(target_cod = cd, buffer_km = 42, expand_frac = 0.12))
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## 
## ── tmap v3 code detected ───────────────────────────────────────────────────────
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_text()`: migrate the layer options 'shadow' to 'options =
## opt_tm_text(<HERE>)'
## [v3->v4] `tm_layout()`: use `tm_title()` instead of `tm_layout(title = )`
mosaico_ego <- do.call(tmap::tmap_arrange, c(tm_list, list(ncol = 5)))
mosaico_ego
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.
## [plot mode] fit legend/component: Some legend items or map compoments do not
## fit well, and are therefore rescaled.
## ℹ Set the tmap option `component.autoscale = FALSE` to disable rescaling.

Los mozaicos de ego-redes muestran las relaciones espaciales locales de los diez municipios con mayor promedio del Índice de Incidencia del Conflicto Armado (IICA) durante el periodo 2017–2021. Cada recuadro representa un municipio “ego” (en rojo) y sus cuatro vecinos más cercanos (k = 4) conectados por líneas naranjas, permitiendo observar cómo se estructura territorialmente el conflicto en torno a los focos de mayor intensidad.

En la mayoría de los casos, los municipios con altos valores del IICA —como Tibú, El Tarra, Tarazá o Tumaco— se encuentran rodeados por vecinos que también presentan niveles elevados de violencia, lo que sugiere la existencia de conglomerados regionales de conflicto persistente. Esto evidencia que el fenómeno no se limita a un solo municipio, sino que se extiende en red hacia los territorios próximos, reafirmando la dependencia espacial positiva observada en los análisis globales.

En conjunto, los mozaicos permiten visualizar de forma sintética la conectividad territorial del conflicto armado, mostrando que los núcleos de alta violencia tienden a formar subregiones cohesionadas, donde los municipios se influyen mutuamente a través de factores geográficos, económicos y de control armado. Esta representación apoya la interpretación de que el conflicto opera bajo una lógica de difusión espacial, en la que las dinámicas locales se retroalimentan entre municipios contiguos.

Autocorrelación espacial global (Moran I)

moran_por_anio <- function(g_sf, anio, nb_global, nsim = 999) {
  colname <- paste0("IICA ", anio)
  vals <- g_sf[[colname]]
  sel <- !is.na(vals)
  n_sel <- sum(sel)
  if (n_sel < 10) {
    return(tibble(year = anio, moran_i = NA, p_perm = NA, n = n_sel))
  }
  nb_sub <- spdep::subset.nb(nb_global, sel)
  lw_sub <- nb2listw(nb_sub, style = "W", zero.policy = TRUE)
  set.seed(123)
  mc <- moran.mc(vals[sel], listw = lw_sub, nsim = nsim,
                 alternative = "two.sided", zero.policy = TRUE)
  tibble(year = anio,
         moran_i = as.numeric(mc$statistic),
         p_perm = as.numeric(mc$p.value),
         n = n_sel)
}

res_moran <- bind_rows(lapply(2017:2021, \(a) moran_por_anio(panel_wide, a, nb_full, 999)))

tabla_moran <- res_moran |>
  mutate(
    interpretacion = case_when(
      moran_i > 0 & p_perm < 0.05 ~ "Positiva y significativa",
      moran_i < 0 & p_perm < 0.05 ~ "Negativa y significativa",
      TRUE ~ "No significativa"
    ),
    `I de Moran` = round(moran_i, 3),
    `p (perm.)` = formatC(p_perm, format = "e", digits = 2)
  ) |>
  select(Año = year, `I de Moran`, `p (perm.)`, `N (mun.)` = n,
         Interpretación = interpretacion)
tabla_moran
## # A tibble: 5 × 5
##     Año `I de Moran` `p (perm.)` `N (mun.)` Interpretación          
##   <int>        <dbl> <chr>            <int> <chr>                   
## 1  2017        0.482 0.00e+00          1119 Positiva y significativa
## 2  2018        0.496 0.00e+00          1119 Positiva y significativa
## 3  2019        0.443 0.00e+00          1119 Positiva y significativa
## 4  2020        0.428 0.00e+00          1119 Positiva y significativa
## 5  2021        0.41  0.00e+00          1119 Positiva y significativa

El análisis de autocorrelación espacial mediante el Índice de Moran (I), este estadístico evalúa si los municipios con niveles similares de conflicto tienden a agruparse geográficamente (autocorrelación positiva) o si los municipios con valores contrastantes se encuentran próximos entre sí (autocorrelación negativa).

Los resultados obtenidos para el periodo 2017–2021 muestran de manera consistente una autocorrelación espacial positiva y estadísticamente significativa. El valor de Moran I oscila entre 0.48 y 0.41, con p < 0.001 en todos los años analizados. Estos valores se interpretan como evidencia de que el conflicto armado en Colombia no se distribuye de manera aleatoria, sino que presenta agrupamientos territoriales donde los municipios con alta incidencia del conflicto tienden a estar próximos entre sí, al igual que aquellos con baja incidencia. La magnitud del índice (moderadamente alta) indica una dependencia espacial, es decir, los niveles de conflicto en un municipio están correlacionados con los de sus vecinos.

Este patrón sugiere la existencia de “efectos de vecindad” o de “difusión territorial” del conflicto, donde la presencia de dinámicas violentas en un municipio puede influir en el comportamiento de los municipios circundantes. Esto es coherente con la literatura sobre geografía del conflicto (Anselin, 1995; Cederman, 2008), que plantea que los procesos de violencia armada tienden a reproducirse espacialmente a través de canales de interacción regional, tales como las rutas de movilidad, la continuidad territorial de grupos armados o la proximidad de economías ilegales.

En conjunto, los resultados del Índice de Moran global confirman que el conflicto armado colombiano presenta un patrón de concentración espacial persistente, estable a lo largo del tiempo y con una tendencia hacia la formación de núcleos regionales de alta intensidad.

Autocorrelación local (LISA)

lisa_por_anio <- function(g_sf, anio, nb_global, nsim = 999) {
  colname <- paste0("IICA ", anio)
  vals <- g_sf[[colname]]; sel <- !is.na(vals)
  if (sum(sel) < 10) return(tibble())
  nb_sub <- subset.nb(nb_global, sel)
  lw_sub <- nb2listw(nb_sub, style = "W", zero.policy = TRUE)
  zvals <- scale(vals[sel])[,1]
  set.seed(123)
  lmc <- localmoran_perm(zvals, lw_sub, nsim = nsim,
                         alternative = "two.sided", zero.policy = TRUE)
  wz <- lag.listw(lw_sub, zvals, zero.policy = TRUE)
  cluster <- case_when(
    zvals >= 0 & wz >= 0 ~ "Alto-Alto",
    zvals <= 0 & wz <= 0 ~ "Bajo-Bajo",
    zvals >= 0 & wz <= 0 ~ "Alto-Bajo",
    zvals <= 0 & wz >= 0 ~ "Bajo-Alto",
    TRUE ~ NA_character_)
  tibble(cod_dane = g_sf$cod_dane[sel],
         Ii = lmc[,"Ii"],
         p_value = lmc[,"Pr(z != E(Ii))"],
         cluster = cluster,
         sig95 = lmc[,"Pr(z != E(Ii))"] < 0.05,
         year = anio)
}
lisa_todos <- bind_rows(lapply(2017:2021, \(a) lisa_por_anio(panel_wide, a, nb_full, 999)))

Mientras el Índice de Moran global ofrece una medida resumen del grado general de dependencia espacial en el territorio, el análisis local mediante el Índice LISA (Local Indicators of Spatial Association) permite identificar dónde se concentran específicamente los patrones de agrupamiento o contraste. Esta herramienta, propuesta por Luc Anselin (1995), descompone la autocorrelación global en valores locales, asignando a cada municipio un estadístico de Moran individual que mide su grado de similitud con los municipios vecinos.

En términos prácticos, el LISA clasifica los municipios en cuatro tipos de clúster espacial: • Alto-Alto: municipios con valores altos de IICA rodeados de otros también altos, indicando zonas de conflicto concentrado. • Bajo-Bajo: municipios con valores bajos rodeados de otros similares, representando zonas estables o con baja incidencia de violencia. • Alto-Bajo: municipios con alto IICA rodeados de bajos, posibles puntos de transición o expansión del conflicto. • Bajo-Alto: municipios con bajo IICA rodeados de altos, que podrían reflejar bordes o amortiguadores del conflicto.

La aplicación del LISA al periodo 2017–2021, con 999 permutaciones y un nivel de significancia del 5%, reveló una estructura espacial heterogénea pero persistente en el tiempo. Predominan los clústeres Bajo-Bajo, que abarcan amplias zonas del centro y oriente del país, mientras que los Alto-Alto se concentran en regiones históricamente afectadas por la violencia, como el Catatumbo, el norte del Cauca, el Bajo Cauca antioqueño y el Pacífico nariñense. Estos focos de alta incidencia sugieren núcleos territoriales estables de conflicto armado, donde la intensidad y recurrencia de los eventos violentos mantienen una correlación significativa entre municipios vecinos.

La visualización de los resultados mediante mapas interactivos y gráficas comparativas por año facilita observar los cambios y persistencias en el tiempo. Las barras anuales muestran variaciones en la proporción de municipios dentro de cada clúster, evidenciando que, aunque hay fluctuaciones entre años, la estructura espacial del conflicto no se dispersa ni se reconfigura completamente, sino que mantiene una lógica de continuidad territorial. El mapa interactivo, por su parte, permite explorar dinámicamente los municipios según su clúster y condición PDET, mostrando cómo algunos de estos territorios priorizados aún permanecen dentro de las zonas de alta incidencia (Alto-Alto), lo que refuerza la hipótesis de persistencia territorial de la violencia armada.

En síntesis, los resultados del análisis LISA confirman que el conflicto armado en Colombia presenta patrones de agrupamiento local consistentes, donde los municipios con alta violencia tienden a formar bolsones regionales de conflicto.

Gráfica interactiva

niveles_cluster <- c("Alto-Alto","Bajo-Bajo","Alto-Bajo","Bajo-Alto")
pal_lisa <- c("Alto-Alto"="#ef3b2c","Bajo-Bajo"="#2b8cbe",
              "Alto-Bajo"="#fdae61","Bajo-Alto"="#9bd8e0")
dist_anual_abs <- lisa_todos |> mutate(cluster=factor(cluster,levels=niveles_cluster)) |> count(year,cluster,name="n")
p <- plot_ly(data=dist_anual_abs, x=~cluster, y=~n, frame=~year,
             type="bar", color=~cluster, colors=pal_lisa,
             text=~n, textposition="outside",
             textfont=list(size=12,color="black")) |>
  layout(title="Número de municipios por clúster LISA (selecciona el año)",
         xaxis=list(title="Clúster LISA",categoryorder="array",categoryarray=niveles_cluster),
         yaxis=list(title="Número de municipios"),
         showlegend=FALSE) |>
  animation_opts(frame=800,transition=300,easing="cubic-in-out") |>
  animation_slider(currentvalue=list(prefix="Año: ",font=list(size=14)))
p

Mapa interactivo

paleta_cluster <- colorFactor(palette=pal_lisa,levels=niveles_cluster)
lisa_pts_year <- function(anio,radio=7){
  iica_col<-paste0("IICA ",anio)
  lisa_y<-lisa_todos|>filter(year==anio)
  atrib_y<-panel_wide|>st_drop_geometry()|>select(cod_dane,all_of(iica_col),PDET)
  capa<-Shape|>st_transform(3116)|>
    left_join(lisa_y,by="cod_dane")|>
    left_join(atrib_y,by="cod_dane")|>
    mutate(cluster=factor(cluster,levels=niveles_cluster),
           PDET_flag=PDET%in%c("Sí","Si","SI","sí","1",1,TRUE,"TRUE"))
  pts<-st_point_on_surface(capa)
  pts84<-st_transform(pts,4326)
  coords<-st_coordinates(pts84)
  list(pts=pts84,coords=coords,iica_col=iica_col,radio=radio)
}
mapa_lisa_comparativo<-function(radio=7){
  anios<-2017:2021
  capas<-lapply(anios,lisa_pts_year,radio=radio);names(capas)<-paste0("Año ",anios)
  m<-leaflet(options=leafletOptions(zoomSnap=0.25,zoomDelta=0.25))%>%
    addProviderTiles(providers$CartoDB.Voyager)
  for(g in names(capas)){
    pts<-capas[[g]]$pts;coords<-capas[[g]]$coords;iica_col<-capas[[g]]$iica_col;r<-capas[[g]]$radio
    m<-m%>%
      addCircleMarkers(data=pts,lng=coords[,1],lat=coords[,2],group=g,
                       fillColor=~paleta_cluster(cluster),fillOpacity=0.95,radius=r,
                       stroke=TRUE,color="white",weight=2,
                       popup=~paste0("<b>Municipio:</b> ",MPIO_CNMBR," (",DPTO_CNMBR,")<br>",
                                     "<b>Clúster:</b> ",cluster,"<br>",
                                     "<b>",iica_col,":</b> ",round(get(iica_col),3),"<br>",
                                     "<b>LISA I:</b> ",round(Ii,3),"<br>",
                                     "<b>p-valor:</b> ",round(p_value,3),"<br>",
                                     "<b>PDET:</b> ",ifelse(PDET_flag,"Sí","No")))%>%
      addCircleMarkers(data=pts[pts$sig95,],lng=coords[pts$sig95,1],lat=coords[pts$sig95,2],
                       group=g,radius=r+1.5,fillOpacity=0,stroke=TRUE,color="black",weight=2)%>%
      addCircleMarkers(data=pts[pts$PDET_flag,],lng=coords[pts$PDET_flag,1],lat=coords[pts$PDET_flag,2],
                       group=g,radius=r+2.5,fillOpacity=0,stroke=TRUE,color="grey20",weight=2.5)
  }
  m%>%addLegend("bottomright",pal=paleta_cluster,values=niveles_cluster,
                title="Clúster LISA",opacity=1)%>%
    addLayersControl(overlayGroups=names(capas),
                     options=layersControlOptions(collapsed=FALSE))
}
mapa_lisa_comparativo(radio=7)
## Warning: st_point_on_surface assumes attributes are constant over geometries
## Warning: st_point_on_surface assumes attributes are constant over geometries
## Warning: st_point_on_surface assumes attributes are constant over geometries
## Warning: st_point_on_surface assumes attributes are constant over geometries
## Warning: st_point_on_surface assumes attributes are constant over geometries

Modelos espaciales: OLS, SAR, SEM

# 7.0 Variable indicadora PDET (lógica robusta)
panel_wide <- panel_wide |>
  dplyr::mutate(PDET_flag = PDET %in% c("Sí","Si","SI","sí","1",1,TRUE,"TRUE"))

# 7.1 OLS de referencia
ols <- lm(IICA_mean ~ PDET_flag, data = panel_wide)

# 7.2 Diagnóstico: autocorrelación en residuos OLS 
moran_res <- spdep::moran.test(residuals(ols), listw = W_full, zero.policy = TRUE)
lm_tests  <- spdep::lm.LMtests(ols, listw = W_full, test = "all", zero.policy = TRUE)
## Please update scripts to use lm.RStests in place of lm.LMtests
# 7.3 Modelos espaciales
# SAR (rezago en la variable dependiente)
sar <- spatialreg::lagsarlm(IICA_mean ~ PDET_flag, data = panel_wide,
                            listw = W_full, zero.policy = TRUE, method = "eigen")

# SEM (correlación en el término de error)
sem <- spatialreg::errorsarlm(IICA_mean ~ PDET_flag, data = panel_wide,
                              listw = W_full, zero.policy = TRUE, method = "eigen")

# 7.4 Resúmenes
sum_sar <- summary(sar)
sum_sem <- summary(sem)

print(moran_res)   # Moran en residuos OLS
## 
##  Moran I test under randomisation
## 
## data:  residuals(ols)  
## weights: W_full    
## 
## Moran I statistic standard deviate = 21.095, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.4116667741     -0.0008944544      0.0003824873
print(lm_tests)    # LM/Robustos
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = IICA_mean ~ PDET_flag, data = panel_wide)
## test weights: listw
## 
## RSerr = 430.98, df = 1, p-value < 2.2e-16
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = IICA_mean ~ PDET_flag, data = panel_wide)
## test weights: listw
## 
## RSlag = 476.97, df = 1, p-value < 2.2e-16
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = IICA_mean ~ PDET_flag, data = panel_wide)
## test weights: listw
## 
## adjRSerr = 2.9778, df = 1, p-value = 0.08442
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = IICA_mean ~ PDET_flag, data = panel_wide)
## test weights: listw
## 
## adjRSlag = 48.971, df = 1, p-value = 2.597e-12
## 
## 
##  Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
##  dependence
## 
## data:  
## model: lm(formula = IICA_mean ~ PDET_flag, data = panel_wide)
## test weights: listw
## 
## SARMA = 479.95, df = 2, p-value < 2.2e-16
summary(sar)       # coeficientes + rho, AIC
## 
## Call:spatialreg::lagsarlm(formula = IICA_mean ~ PDET_flag, data = panel_wide, 
##     listw = W_full, method = "eigen", zero.policy = TRUE)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.3165261 -0.0163639 -0.0083782  0.0052619  0.7404577 
## 
## Type: lag 
## Coefficients: (asymptotic standard errors) 
##                Estimate Std. Error z value  Pr(>|z|)
## (Intercept)   0.0100990  0.0024374  4.1433 3.423e-05
## PDET_flagTRUE 0.0753063  0.0064095 11.7493 < 2.2e-16
## 
## Rho: 0.69826, LR test value: 474.67, p-value: < 2.22e-16
## Asymptotic standard error: 0.022244
##     z-value: 31.391, p-value: < 2.22e-16
## Wald statistic: 985.37, p-value: < 2.22e-16
## 
## Log likelihood: 1326.999 for lag model
## ML residual variance (sigma squared): 0.0048119, (sigma: 0.069368)
## Number of observations: 1119 
## Number of parameters estimated: 4 
## AIC: -2646, (AIC for lm: -2173.3)
## LM test for residual autocorrelation
## test value: 0.85698, p-value: 0.35458
summary(sem)       # coeficientes + lambda, AIC
## 
## Call:
## spatialreg::errorsarlm(formula = IICA_mean ~ PDET_flag, data = panel_wide, 
##     listw = W_full, method = "eigen", zero.policy = TRUE)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.286754 -0.019984 -0.011724  0.003821  0.738385 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##                Estimate Std. Error z value  Pr(>|z|)
## (Intercept)   0.0528417  0.0080127  6.5948 4.259e-11
## PDET_flagTRUE 0.0989744  0.0096689 10.2363 < 2.2e-16
## 
## Lambda: 0.72938, LR test value: 429.73, p-value: < 2.22e-16
## Asymptotic standard error: 0.021873
##     z-value: 33.346, p-value: < 2.22e-16
## Wald statistic: 1112, p-value: < 2.22e-16
## 
## Log likelihood: 1304.532 for error model
## ML residual variance (sigma squared): 0.0049323, (sigma: 0.07023)
## Number of observations: 1119 
## Number of parameters estimated: 4 
## AIC: -2601.1, (AIC for lm: -2173.3)
AIC(ols, sar, sem) # comparación de ajuste
##     df       AIC
## ols  3 -2173.332
## sar  4 -2645.997
## sem  4 -2601.063
# 7.5 Tabla compacta de parámetros principales (sin paquetes extra)
tab_sar <- tibble::tibble(
  Modelo     = "SAR (lag)",
  Param      = c("(Intercepto)", "PDET_flag", "rho"),
  Estimador  = c(sum_sar$Coef["(Intercept)", 1],
                 sum_sar$Coef["PDET_flagTRUE", 1],
                 sum_sar$rho),
  `EE`       = c(sum_sar$Coef["(Intercept)", 2],
                 sum_sar$Coef["PDET_flagTRUE", 2],
                 sum_sar$rho.se),
  `p-valor`  = c(sum_sar$Coef["(Intercept)", 4],
                 sum_sar$Coef["PDET_flagTRUE", 4],
                 sum_sar$Wald1$p.value)
)

tab_sem <- tibble::tibble(
  Modelo     = "SEM (error)",
  Param      = c("(Intercepto)", "PDET_flag", "lambda"),
  Estimador  = c(sum_sem$Coef["(Intercept)", 1],
                 sum_sem$Coef["PDET_flagTRUE", 1],
                 sum_sem$lambda),
  `EE`       = c(sum_sem$Coef["(Intercept)", 2],
                 sum_sem$Coef["PDET_flagTRUE", 2],
                 sum_sem$lambda.se),
  `p-valor`  = c(sum_sem$Coef["(Intercept)", 4],
                 sum_sem$Coef["PDET_flagTRUE", 4],
                 sum_sem$Wald1$p.value)
)

tabla_modelos <- dplyr::bind_rows(tab_sar, tab_sem)
tabla_modelos
## # A tibble: 6 × 5
##   Modelo      Param        Estimador      EE `p-valor`
##   <chr>       <chr>            <dbl>   <dbl>     <dbl>
## 1 SAR (lag)   (Intercepto)    0.0101 0.00244  3.42e- 5
## 2 SAR (lag)   PDET_flag       0.0753 0.00641  0       
## 3 SAR (lag)   rho             0.698  0.0222   0       
## 4 SEM (error) (Intercepto)    0.0528 0.00801  4.26e-11
## 5 SEM (error) PDET_flag       0.0990 0.00967  0       
## 6 SEM (error) lambda          0.729  0.0219   0

El parámetro indica la presencia de una fuerte dependencia espacial positiva, lo que significa que los municipios tienden a reflejar el comportamiento de sus vecinos en términos de intensidad del conflicto. Por su parte, el coeficiente = 0.075 muestra que, manteniendo constante la influencia espacial, los municipios PDET presentan niveles promedio de conflicto aproximadamente un 7.5 % más altos que los municipios no PDET. Este resultado no implica una relación causal, sino que evidencia la persistencia de altos niveles de conflicto en los territorios priorizados, los cuales históricamente han concentrado la violencia y las condiciones estructurales que la perpetúan.

El modelo SEM (Spatial Error Model), que modela la autocorrelación en el término de error en lugar del rezago de la variable dependiente, presentó un ajuste inferior (AIC_{SEM} = -2601) y un coeficiente = 0.729 también significativo, pero de menor relevancia explicativa según las pruebas de diagnóstico. Esto confirma que el modelo SAR/SLM es el que mejor captura la estructura espacial del conflicto armado colombiano, al considerar que el nivel de conflicto en un municipio depende no solo de sus características internas, sino también de las condiciones presentes en los municipios vecinos.

Estos resultados reafirman la hipótesis de dependencia territorial del conflicto armado, donde las dinámicas violentas presentan un efecto de “contagio” o retroalimentación espacial. El valor elevado de y la significancia de demuestran que los esfuerzos institucionales en los municipios PDET se enfrentan a condiciones regionales que trascienden las fronteras administrativas, lo que sugiere la necesidad de políticas de intervención con enfoque territorial integrado.

Conclusiones

los resultados obtenidos a lo largo del estudio permite afirmar que el conflicto armado colombiano presenta una marcada dependencia espacial, tanto en su intensidad como en su persistencia territorial. Los análisis realizados —desde el Índice de Moran global y local (LISA) hasta los modelos econométricos espaciales (SAR y SEM)— evidencian que los niveles de violencia en un municipio no pueden explicarse de manera aislada, sino que están condicionados por las dinámicas de sus municipios vecinos. El conflicto, por tanto, actúa como un fenómeno relacional y regional, donde el espacio se convierte en un canal de propagación y retroalimentación de las dinámicas violentas.

El Índice de Moran global, con valores entre 0.41 y 0.48, mostró una autocorrelación positiva y significativa en todos los años analizados, confirmando la existencia de un patrón territorial estable. El análisis local (LISA) complementó esta evidencia al identificar zonas críticas de concentración —particularmente en el Catatumbo, el Pacífico nariñense, el norte del Cauca y el Bajo Cauca antioqueño— donde los municipios de alta violencia tienden a agruparse formando núcleos persistentes de conflicto armado. Estos resultados sustentan la hipótesis de continuidad territorial de la violencia, en la que los factores históricos, sociales y económicos se entrelazan con las interacciones espaciales para mantener la incidencia del conflicto en regiones específicas.

El modelo espacial autorregresivo (SAR) permitió formalizar esta dependencia al mostrar un coeficiente espacial de = 0.698 (p < 0.001), que indica una fuerte influencia del entorno geográfico en la magnitud del conflicto. Asimismo, el coeficiente de la variable PDET (= 0.075, p < 0.001) confirma que los municipios priorizados por las políticas de desarrollo territorial mantienen niveles significativamente más altos de conflicto, incluso después de controlar la influencia espacial. Esto no implica un efecto negativo de los programas, sino que demuestra que las intervenciones institucionales se desarrollan precisamente en los territorios de mayor vulnerabilidad estructural y persistencia histórica del conflicto.

Desde una perspectiva territorial, los resultados sugieren que el conflicto armado en Colombia opera como un sistema espacialmente interconectado, donde los cambios en un municipio pueden desencadenar efectos indirectos sobre los municipios circundantes. Esto tiene importantes implicaciones para las políticas de paz y desarrollo: los enfoques locales aislados resultan insuficientes si no se consideran las interdependencias regionales que sustentan las dinámicas de violencia.La evidencia empírica sugiere que las estrategias de intervención deben orientarse hacia regiones funcionales del conflicto, abordando de manera simultánea los municipios que comparten redes de interacción territorial y problemas estructurales comunes.

Asimismo, futuras líneas de investigación podrían incorporar variables complementarias —como la inversión pública, el acceso a servicios, la sustitución de cultivos o la presencia institucional— con el fin de explicar no solo la persistencia, sino también las causas de la propagación del conflicto. De esta manera, el trabajo avanza hacia un enfoque integral que combine evidencia empírica, análisis territorial y formulación de políticas basadas en el espacio, contribuyendo a la comprensión científica y aplicada del conflicto armado colombiano.