Análisis Exploratorio de Datos Espaciales

url <- "https://fonts.googleapis.com/css2?family=Playfair+Display:wght@500;700&family=Inter:wght@300;400;500;600&display=swap'"
datos <- read.csv(url)

library(readxl)
library(psych)
library(ggplot2)
library(dplyr)
library(ggthemes)
library(openxlsx)

1. Gestión de datos

general = read.csv("cedegeneral.csv")
names(general)

##  [1] "codmun"         "ano"            "gandina"        "gcaribe"       
##  [5] "gpacifica"      "gorinoquia"     "gamazonia"      "pobl_rur"      
##  [9] "pobl_urb"       "pobl_tot"       "areaoficialkm2" "altura"        
## [13] "discapital"     "disbogota"

general = read.csv("cedegeneral.csv")
names(general)

##  [1] "codmun"         "ano"            "gandina"        "gcaribe"       
##  [5] "gpacifica"      "gorinoquia"     "gamazonia"      "pobl_rur"      
##  [9] "pobl_urb"       "pobl_tot"       "areaoficialkm2" "altura"        
## [13] "discapital"     "disbogota"

names(general)[2] = "año"

ipm = read.csv("IPMCol.csv", sep = ",", dec = ",", encoding="UTF-8")
ipm = ipm[, -1]
names(ipm)[1] = "codmun"

th = read_excel("tasahomicidios.xlsx", col_types = c("text", "text", "numeric", "text", "numeric", "skip"))
th = th[,c(-1,-2,-4)]
names(th)[names(th) == "Código Entidad"] = "codmun"

cultivos = read_excel("cultivos.xlsx")
cultivos[1:3,]

## # A tibble: 3 × 3
##   Codigo_municipio Tipologia_Coca Cultivos_coca_hectareas
##              <dbl> <chr>                            <dbl>
## 1             5001 Sin Coca                            NA
## 2             5002 Sin Coca                            NA
## 3             5004 Sin Coca                            NA

merge1 = merge(general, ipm, by = "codmun")
merge2 = merge(merge1, th, by = "codmun")
base_final = merge(merge2, cultivos, by.x = "codmun", by.y = "Codigo_municipio")

n = nrow(base_final)
base_final[-(1:(n-5)),1:4]

##      codmun  año gandina gcaribe
## 1097  97666 2017       0       0
## 1098  99001 2017       0       0
## 1099  99524 2017       0       0
## 1100  99624 2017       0       0
## 1101  99773 2017       0       0

2.Estadística Descriptiva

attach(base_final)

summary(base_final[, (ncol(base_final) - 4):ncol(base_final)])

##   IPM.nacional   Descripcion.Etnica Tasa_Homicidios    Tipologia_Coca
##  Min.   : 0.00   Mode:logical       Min.   :  0.00   Length   :1101  
##  1st Qu.:31.20   NAs :1101          1st Qu.:  0.00   N.unique :   3  
##  Median :43.60                      Median : 15.70   N.blank  :   0  
##  Mean   :43.98                      Mean   : 24.18   Min.nchar:   8  
##  3rd Qu.:56.00                      3rd Qu.: 33.15   Max.nchar:  18  
##  Max.   :92.60                      Max.   :272.65                   
##                                     NAs    :10                       
##  Cultivos_coca_hectareas
##  Min.   :    1.04       
##  1st Qu.:   12.37       
##  Median :   68.95       
##  Mean   :  495.75       
##  3rd Qu.:  307.29       
##  Max.   :16960.24       
##  NAs    :908

media_ipm = mean(IPM.nacional)
media_ipm

## [1] 43.97729

media_distbog = mean(disbogota)
media_distbog

## [1] 314.4424

quantile(IPM.nacional, na.rm = TRUE)

##   0%  25%  50%  75% 100% 
##  0.0 31.2 43.6 56.0 92.6

quantile(IPM.nacional, probs = c(0.2, 0.4, 0.6, 0.8))

##  20%  40%  60%  80% 
## 28.5 38.5 48.7 59.2

hist(IPM.nacional, breaks = 20, col="#B81840", main="Distribucion Observada Índice de Pobreza Multidimensional", xlab="Tasa", ylab=("Frecuencia Relativa"), freq=F)

– Análisis: El histograma de frecuencias relativas del Índice de Pobreza Multidimensional (IPM) muestra una distribución bastante simétrica y de forma acampanada para los municipios colombianos en 2017. La mayor concentración de datos se ubica en el rango intermedio, específicamente entre los valores de 30 y 50, lo que indica que predominan los niveles moderados de pobreza, mientras que los municipios con puntajes extremadamente bajos o altos son minoría.

plot(density(IPM.nacional), main = "Gráfico de densidad: IPM", ylab = "Frecuencia", col="#830717")

– Análisis: El histograma de frecuencias relativas del Índice de Pobreza Multidimensional (IPM) muestra una distribución bastante simétrica y de forma acampanada para los municipios colombianos en 2017. La mayor concentración de datos se ubica en el rango intermedio, específicamente entre los valores de 30 y 50, lo que indica que predominan los niveles moderados de pobreza, mientras que los municipios con puntajes extremadamente bajos o altos son minoría.

ggplot(base_final, aes(x = IPM.nacional)) +
  geom_histogram(binwidth = 5, fill = "#B81840", color = "black") +
  geom_vline(xintercept = mean(base_final$IPM.nacional), color = "#FFBDC0", linetype = "dashed", size = 1) +
  labs(
    title = "Distribución Observada Índice de Pobreza Multidimensional",
    x = "Tasa",
    y = "Frecuencia Relativa"
  )

– Análisis: El gráfico de densidad del Índice de Pobreza Multidimensional (IPM) complementa el análisis anterior al suavizar la distribución de los datos para una muestra robusta de 1101 municipios (N=1101). La curva muestra un comportamiento unimodal y bastante simétrico, donde la mayor concentración de probabilidad se eleva entre los valores de 35 y 50, confirmando que predominan los niveles intermedios de pobreza. Esta forma suave y acampanada es una excelente noticia para tu taller, ya que demuestra que la variable no presenta sesgos extremos ni valores atípicos severos que puedan comprometer los supuestos de normalidad al momento de estimar tu modelo de regresión.

ric = IQR(IPM.nacional)
ric

## [1] 24.8

sd_ipm = sd(IPM.nacional)
var_ipm = var(IPM.nacional)

ind_rur = pobl_rur/pobl_tot
base_final = cbind(base_final, ind_rur)

sesgo_kurt_IPM = c("Skewness" = skew(base_final$IPM.nacional), "Kurtosis" = kurtosi(base_final$IPM.nacional))
sesgo_kurt_disbogota = c("Skewness" = skew(base_final$disbogota), "Kurtosis" = kurtosi(base_final$disbogota))
sesgo_kurt_ind_rur = c("Skewness" = skew(base_final$ind_rur), "Kurtosis" = kurtosi(base_final$ind_rur))

hist_IPM = hist(base_final$IPM.nacional, plot = FALSE)
hist_disbogota = hist(base_final$disbogota, plot = FALSE)
hist_ind_rur = hist(base_final$ind_rur, plot = FALSE)

par(mfrow = c(1, 3))

coord_x_IPM = 4
coord_y_IPM = 200
# Crear el histograma
barplot(hist_IPM$counts, names.arg = hist_IPM$mids, col = "#FFD99F", border = "black", main = "Histograma IPM.nacional", xlab = "Índice", ylab = "Frecuencia Relativa")

text(coord_x_IPM, coord_y_IPM, paste("Sesgo:", round(sesgo_kurt_IPM["Skewness"], 2), "\nCurtosis:", round(sesgo_kurt_IPM["Kurtosis"], 2)), pos = 4)

coord_x_disbogota = 4
coord_y_disbogota = 200
barplot(hist_disbogota$counts, names.arg = hist_disbogota$mids, col = "#FE886E", border = "black", main = "Histograma disbogota", xlab = "Distancia en KM", ylab = "Frecuencia Relativa")
text(coord_x_disbogota, coord_y_disbogota, paste("Sesgo:", round(sesgo_kurt_disbogota["Skewness"], 2), "\nCurtosis:", round(sesgo_kurt_disbogota["Kurtosis"], 2)), pos = 4)

#ind_rur
coord_x_ind_rur = 0.2
coord_y_ind_rur = 150
barplot(hist_ind_rur$counts, names.arg = hist_ind_rur$mids, col = "#F24970", border = "black", main = "Histograma ind_rur", xlab = "Tasa", ylab = "Frecuencia Relativa")
text(coord_x_ind_rur, coord_y_ind_rur, paste("Sesgo:", round(sesgo_kurt_ind_rur["Skewness"], 2), "\nCurtosis:", round(sesgo_kurt_ind_rur["Kurtosis"], 2)), pos = 4)

– Análisis: Para analizar la distribución de las variables clave, analicé los histogramas y sus medidas de forma. El Índice de Pobreza Multidimensional (IPM.nacional) muestra una distribución equilibrada con un sesgo levemente positivo de 0.15. Por su parte, la distancia a Bogotá (disbogota) presenta un sesgo positivo de 0.68, reflejando que la mayoría de los municipios se concentran a distancias intermedias y pocos se ubican en zonas apartadas. Finalmente, el índice de ruralidad (ind_rur) muestra un sesgo negativo de −0.47, lo que evidencia condiciones elevadas de ruralidad en una porción significativa de los municipios analizados.

par(mfrow = c(1, 1))

base_final$densidad = base_final$pobl_tot /base_final$areaoficialkm2
describe(base_final$densidad)

##    vars    n   mean     sd median trimmed   mad  min      max    range  skew
## X1    1 1101 159.15 733.01  45.57   58.79 42.07 0.15 16113.35 16113.21 13.89
##    kurtosis    se
## X1   246.04 22.09

plot(IPM.nacional ~ disbogota)

– Análisis: El diagrama de dispersión entre el Índice de Pobreza Multidimensional (IPM.nacional) y la distancia a Bogotá (disbogota) permite visualizar la relación bivariada entre ambas variables. Se observa una tendencia positiva donde, a medida que aumenta la distancia geográfica a la capital, los niveles de pobreza tienden a incrementarse en los municipios. Asimismo, la nube de puntos muestra una mayor dispersión en distancias intermedias y altas, evidenciando la variabilidad territorial que será clave para contrastar tu hipótesis en el modelo de regresión.

ggplot(data = base_final, 
       mapping = aes(x = disbogota,y = IPM.nacional)) +
  geom_point() +
  labs(title = "Relación entre IPM y Distancia a Bogotá", 
       y = "Indice de pobreza multidimensional", 
       x = "Distancia en Km en línea recta", 
       caption = "Fuente: Elaboración propia con datos del DANE", 
       )

– Análisis: El diagrama de dispersión permite examinar la relación bivariada entre la distancia a Bogotá y el Índice de Pobreza Multidimensional (IPM) para los 1,101 municipios analizados. Los resultados visuales sugieren una leve tendencia positiva; es decir, a medida que incrementa la distancia geográfica respecto a la capital, los municipios tienden a registrar mayores niveles de pobreza. No obstante, la alta dispersión observada en los rangos intermedios y altos de distancia evidencia una importante heterogeneidad territorial. Esta variabilidad indica que la distancia por sí sola no es un determinante absoluto, lo cual justifica la necesidad de contrastar estas hipótesis mediante la estimación de un modelo de regresión.

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional, color = gandina)) +
  geom_point() +
  labs(title = "Relación entre IPM y Distancia a Bogotá", 
       y = "Indice de pobreza multidimensional", 
       x = "Distancia en Km en línea recta", 
       caption = "Fuente: Elaboración propia con datos del DANE")

– Análisis: El diagrama de dispersión incorpora la variable gandina (asociada a la pertenencia o proporción de la región andina) mediante una escala de color, lo que permite refinar el análisis bivariado al introducir una dimensión territorial adicional. Al observar la distribución de los puntos según esta gradiente, se aprecia que los municipios con valores más altos en dicha variable tienden a concentrarse de manera diferenciada respecto a aquellos con valores más bajos, particularmente en los rangos de distancia media y baja con relación a Bogotá. Esto sugiere que la región andina introduce un patrón particular en la relación entre el Índice de Pobreza Multidimensional (IPM) y la distancia a la capital, aportando evidencia sobre la heterogeneidad espacial que será clave incorporar en las especificaciones del modelo de regresión.

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional, color = factor(gandina))) +
  geom_point() +
  labs(title = "Relación entre IPM y Distancia a Bogotá", 
       y = "Indice de pobreza multidimensional", 
       x = "Distancia en Km en línea recta", 
       caption = "Fuente: Elaboración propia con datos del DANE") +
  scale_color_manual(values = c("1" = "#8C7F5F", "0" = "#F9BEAB"))

– Análisis: El diagrama de dispersión multivariado relaciona el Índice de Pobreza Multidimensional con la distancia a Bogotá, incorporando una distinción por color para identificar la pertenencia a la Región Andina (1 en café para municipios andinos y 0 en piel para los demás). La visualización permite evidenciar diferencias territoriales claras en la distribución de los datos y cómo se agrupan los municipios según su ubicación geográfica. Este desglose espacial es un insumo clave para enriquecer tu taller y fundamentar la inclusión de variables de control o categóricas al momento de estimar el modelo de regresión.

base_final = base_final %>%
  mutate(region = case_when(
    gandina == 1 ~ "Region Andina",
    gcaribe == 1 ~ "Region Caribe",
    gpacifica == 1 ~ "Region Pacifica",
    gorinoquia == 1 ~ "Región Orinoquía",
    gamazonia == 1 ~ "Región Amazonía",
    TRUE ~ "Ninguna"  # Si ninguna de las dummies es 1
  ))

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional, color = region)) +
  geom_point() +
  scale_color_manual(values = c(
    "Region Andina" = "#B2182B", 
    "Region Caribe" = "#F46D43", 
    "Region Pacifica" = "#FDAE61", 
    "Región Orinoquía" = "#4575B4", 
    "Región Amazonía" = "#313695"
  )) +
  labs(title = "Relación entre IPM y Distancia a Bogotá", 
       y = "Indice de pobreza multidimensional", 
       x = "Distancia en Km en línea recta", 
       caption = "Fuente: Elaboración propia con datos del DANE")

– Análisis: El diagrama de dispersión multivariado analiza la relación entre el Índice de Pobreza Multidimensional (IPM) y la distancia a Bogotá, incorporando una diferenciación por color según la región natural de los municipios. Los resultados muestran que los municipios de la Región Andina concentran la mayor parte de las observaciones en distancias cortas e intermedias con un amplio rango de pobreza, mientras que las demás regiones se extienden hacia mayores distancias, evidenciando las marcadas brechas territoriales y el aislamiento de las zonas periféricas del país.

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional)) +
  geom_point() +
  geom_smooth(method = "lm") + 
  labs(title = "Relación entre IPM y Distancia a Bogotá", 
       y = "Indice de pobreza multidimensional", 
       x = "Distancia en Km en línea recta", 
       caption = "Fuente: Elaboración propia con datos del DANE")

– Análisis: El diagrama de dispersión muestra la relación bivariada entre el Índice de Pobreza Multidimensional (IPM) y la distancia en línea recta a Bogotá. La línea de tendencia central evidencia una asociación positiva, lo que indica que, a mayor distancia geográfica respecto a la capital, los municipios tienden a registrar mayores niveles de pobreza. Asimismo, se observa que la mayor concentración de municipios se encuentra en distancias cortas e intermedias, mientras que los territorios más apartados muestran una mayor dispersión en sus condiciones socioeconómicas.

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional)) +
  geom_point(color = "#ABD9E9") +
  geom_smooth(method = "lm", color = "#D7191F", linetype = "dashed", se = FALSE) +
  labs(
    title = "Relación entre IPM y Distancia a Bogotá", 
    y = "Indice de pobreza multidimensional", 
    x = "Distancia en Km en línea recta", 
    caption = "Fuente: Elaboración propia con datos del DANE"
  ) +
  theme_minimal() +
  geom_hline(yintercept = media_ipm, linetype = 1, color = "BLACK") +
  geom_vline(xintercept = media_distbog, linetype = 4, color = "BLACK") +
  theme(
    plot.title = element_text(face = "bold")
  )

## `geom_smooth()` using formula = 'y ~ x'

– Análisis: El diagrama de dispersión analiza la relación entre el Índice de Pobreza Multidimensional (IPM) y la distancia a Bogotá, incorporando dos líneas de referencia que corresponden a los promedios nacionales de ambas variables. La línea horizontal marca la media del índice de pobreza y la línea vertical señala la distancia promedio a la capital. La línea de tendencia central confirma una asociación positiva, evidenciando visualmente que los municipios ubicados a mayor distancia tienden a registrar niveles de pobreza superiores al promedio general del país.

ggplot(data = base_final,
       mapping = aes(x = disbogota,
                     y = IPM.nacional)) +

  geom_point(aes(colour = ind_rur,
                 size = densidad),
             alpha = 0.75) +

  geom_smooth(method = "lm",
              colour = "#333333",
              fill = "#BDBDBD",
              linewidth = 0.9) +

  scale_colour_gradientn(
    colours = c("#2C7BB6", "#ABD9E9", "#FFFFBF", "#FDAE61", "#D7191C"),
    name = "Índice rural"
  ) +

  scale_size_continuous(
    range = c(2, 9),
    name = "Densidad"
  ) +

  labs(
    title = "Relación entre IPM y Distancia a Bogotá",
    y = "Índice de pobreza multidimensional",
    x = "Distancia en Km en línea recta",
    caption = "Fuente: Elaboración propia con datos del DANE"
  ) +

  theme_minimal() +

  geom_hline(yintercept = media_ipm,
             linetype = 1,
             colour = "#777777",
             linewidth = 0.7) +

  geom_vline(xintercept = 314.4,
             linetype = 4,
             colour = "#CC3311",
             linewidth = 0.8)

Análisis: El análisis exploratorio espacial de los municipios colombianos revela una marcada heterogeneidad territorial en las condiciones de vida. Los resultados evidencian una relación directa entre el aislamiento geográfico respecto a la capital y los mayores niveles de pobreza multidimensional (IPM), una vulnerabilidad que se intensifica notablemente en los contextos de mayor ruralidad. Estos hallazgos subrayan la importancia de incorporar variables espaciales y demográficas para comprender las brechas socioeconómicas del país.

ggsave("diagrama_dispersion_1.jpeg")

cov(y = IPM.nacional,  x =disbogota, use = "complete.obs", method = "pearson")

## [1] 1769.147

cor(y = IPM.nacional,  x =disbogota, use = "complete.obs")

## [1] 0.5670527

cov(base_final[,c("IPM.nacional","disbogota","discapital","ind_rur","densidad")])

##              IPM.nacional    disbogota   discapital      ind_rur     densidad
## IPM.nacional   275.089084  1769.147113   290.472434   0.33100304  -1740.40992
## disbogota     1769.147113 35384.064103   990.172883  -6.41175602    630.68384
## discapital     290.472434   990.172883  3175.558813   1.94157679  -6650.12193
## ind_rur          0.331003    -6.411756     1.941577   0.05868393    -50.10636
## densidad     -1740.409918   630.683836 -6650.121933 -50.10635644 537304.94056

cor(base_final[,c("IPM.nacional","disbogota","discapital","ind_rur","densidad")])

##              IPM.nacional    disbogota discapital     ind_rur     densidad
## IPM.nacional   1.00000000  0.567052687  0.3107837  0.08238264 -0.143154319
## disbogota      0.56705269  1.000000000  0.0934108 -0.14070630  0.004574013
## discapital     0.31078365  0.093410800  1.0000000  0.14222793 -0.160993785
## ind_rur        0.08238264 -0.140706298  0.1422279  1.00000000 -0.282177779
## densidad      -0.14315432  0.004574013 -0.1609938 -0.28217778  1.000000000

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional)) +
  geom_point(aes(colour = ind_rur, size = densidad), alpha = 0.8) +
  geom_smooth(method = "lm", color = "black", linetype = "dashed") +
  scale_color_gradientn(colors = c("#2C7BB6", "#ABD9E9", "#FFFFBF", "#FDAE61", "#D7191C")) +
  labs(
    title = "Relación entre IPM y Distancia a Bogotá", 
    y = "Índice de pobreza multidimensional", 
    x = "Distancia en Km en línea recta", 
    caption = "Fuente: Elaboración propia con datos del DANE",
    colour = "Índice de ruralidad",
    size = "Densidad"
  ) +
  theme_minimal() +
  geom_hline(yintercept = 43.9, linetype = 1, color = "grey40") +
  geom_vline(xintercept = 314.4, linetype = 4, color = "grey40") +
  geom_text(x = 1000, y = 15, color = "blue", label = "Correlación = 0.56") +
  theme(
    plot.title = element_text(face = "bold"),
    panel.grid.minor = element_blank()
  )

Análisis: En el gráfico de dispersión que relaciona el Índice de Pobreza Multidimensional (IPM.nacional) con la distancia a Bogotá, puedo observar una asociación positiva respaldada por un coeficiente de correlación de 0.56, lo que indica que los municipios más alejados de la capital tienden a registrar mayores niveles de carencias. Asimismo, al incorporar el índice de ruralidad a través de la escala de colores, noto que los tonos rojizos y anaranjados predominan en la parte alta de la distribución, reflejando que los territorios con mayor ruralidad enfrentan una mayor vulnerabilidad multidimensional. Finalmente, el tamaño de los puntos revela que los centros con mayor densidad poblacional se concentran en zonas con menores distancias o por debajo de la tendencia general, evidenciando de forma clara los contrastes territoriales del país.

ggplot(data = base_final, 
       mapping = aes(x = disbogota, y = IPM.nacional)) +
  geom_point(aes(colour = ind_rur, size = densidad), alpha = 0.8) +
  geom_smooth(method = "lm", color = "black", linetype = "dashed", se = FALSE) +
  scale_color_gradientn(colors = c("#2C7BB6", "#ABD9E9", "#FFFFBF", "#FDAE61", "#D7191C")) +
  facet_wrap(~Tipologia_Coca, ncol = 1) +
  labs(
    title = "Relación entre IPM y Distancia a Bogotá", 
    y = "Índice de pobreza multidimensional", 
    x = "Distancia en Km en línea recta", 
    caption = "Fuente: Elaboración propia con datos del DANE",
    colour = "Índice de ruralidad",
    size = "Densidad"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(face = "bold"),
    panel.grid.minor = element_blank()
  )

Análisis: El gráfico de dispersión segmentado mediante paneles (facet_wrap) evalúa la relación entre el Índice de Pobreza Multidimensional y la distancia a Bogotá, desagregando la muestra según la tipología de cultivos (“Coca Asperjados”, “Coca No Asperjados” y “Sin Coca”). Las rectas de ajuste estimadas en cada panel permiten comparar cómo se comporta la tendencia bivariada en los diferentes contextos territoriales, evidenciando heterogeneidades estructurales entre los municipios. Esta separación por facetas constituye un recurso analítico clave para el taller, ya que facilita justificar la incorporación de variables de control y evaluar efectos diferenciales antes de proceder con la estimación formal del modelo de regresión.

pairs(base_final[,c("IPM.nacional","disbogota","discapital","ind_rur","densidad", "areaoficialkm2")])

Análisis: La matriz de dispersión (pairs) permite examinar simultáneamente las relaciones bivariadas entre todas las variables clave del estudio: IPM.nacional, disbogota, discapital, ind_rur, densidad y areaoficialkm2. Este conjunto de gráficos cruzados facilita identificar de un solo vistazo qué pares de variables presentan asociaciones lineales directas o inversas, así como detectar posibles problemas de multicolinealidad o concentraciones extremas antes de definir la especificación definitiva del modelo de regresión.

write.xlsx(base_final, "Base_General.xlsx")

3. Exploración y Visualización de Datos Espaciales.

library(sf)
library(spData)
library(sp)
library(shapefiles)
library(tmap)
library(RColorBrewer)
library(leaflet)
library(paletteer)

dir()

##  [1] "Base_General.xlsx"          "cedegeneral.csv"           
##  [3] "cultivos.xlsx"              "diagrama_dispersion_1.jpeg"
##  [5] "Ejercicio final.Rmd"        "Ejercicio-final_files"     
##  [7] "Ejercicio-final.docx"       "Ejercicio-final.html"      
##  [9] "Ejercicio-final.log"        "Ejercicio-final.Rmd"       
## [11] "Ejercicio-final.tex"        "IPMCol.csv"                
## [13] "MGN_MPIO_POLITICO.dbf"      "MGN_MPIO_POLITICO.prj"     
## [15] "MGN_MPIO_POLITICO.shp"      "MGN_MPIO_POLITICO.shx"     
## [17] "RegionOrinoquia"            "rsconnect"                 
## [19] "tasahomicidios.xlsx"

base_final = base_final[,-1]

3.1 Tipos de datos espaciales

georef = st_read(".", "MGN_MPIO_POLITICO")

## Reading layer `MGN_MPIO_POLITICO' from data source 
##   `/Users/trululu/Downloads/Ej1_GestionRStudio (extract' using driver `ESRI Shapefile'
## Simple feature collection with 1121 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -81.73562 ymin: -4.229406 xmax: -66.84722 ymax: 13.39473
## Geodetic CRS:  WGS 84

georef = as_Spatial(georef)

georef@data[["MPIO_CDPMP"]] = as.numeric(georef@data[["MPIO_CDPMP"]])

merge1 = merge(general, ipm, by = "codmun")
merge2 = merge(merge1, th, by = "codmun")
base_final = merge(merge2, cultivos, by.x = "codmun", by.y = "Codigo_municipio")
"codmun" %in% names(base_final)

## [1] TRUE

bf = merge(georef, base_final, by.x="MPIO_CDPMP", by.y="codmun")
bf_sf = sf::st_as_sf(bf)

3.2 Visualización de mapas con tmap

a = Sys.time()

tm_shape(bf_sf) + 
  tm_fill("Tasa_Homicidios", 
          palette = paletteer_c("ggthemes::Purple", 30)) +
  tm_layout(
    title = "Tasa de homicidios en Colombia",
    title.position = c("center", "top"),
    title.size = 1.2,
    frame = FALSE
  )

– Análisis: El mapa coroplético de la tasa de homicidios en Colombia muestra una distribución espacial altamente asimétrica. La gran mayoría del territorio nacional predomina con tonos rosa claros, lo que indica que la mayor parte de los municipios registra tasas bajas ( 0 a 49 homicidios por cada 100.000 habitantes). Sin embargo, se aprecian manchas de tonos intermedios (50 a 149) en varias zonas de la Orinoquía, la Amazonía y el occidente del país, mientras que los valores más críticos (150 a 300, en morado oscuro) aparecen de forma muy puntual y focalizada en pequeños clústeres del nororiente, el noroccidente y el suroccidente colombiano. Finalmente, las zonas en gris (Missing) en el suroriente corresponden a territorios sin datos reportados en la muestra.

object.size(x = bf)/(1024^2)

## 91.9 bytes

orinoquia = subset(bf_sf, gorinoquia == 1)
orinoquia_sf = sf::st_as_sf(orinoquia)

tm_shape(orinoquia_sf) + 
  tm_fill("Tasa_Homicidios", 
          palette = paletteer_c("ggthemes::Purple", 30)) + 
  tm_layout(
    legend.position = c("left", "top"),
    legend.text.size = 0.8, 
    legend.title.size = 1.0,
    inner.margins = c(0.05, 0.25, 0.05, 0.05) 
  )

– Análisis: El mapa coroplético de la tasa de homicidios en la Región Orinoquía evidencia un predominio general de niveles bajos e intermedios (0 a 39 homicidios por cada 100.000 habitantes, en tonos rosa claros) a lo largo de las extensas llanuras del oriente, el norte y el centro de la región. No obstante, se identifican focos críticos de violencia (60 a 120, en tonos morados oscuros) concentrados principalmente hacia el borde occidental y suroccidental (correspondientes a zonas del piedemonte y el sur del Meta), así como un municipio puntual en el extremo nororiental (en el departamento del Vichada) y una franja intermedia (\(40\) a \(59\)) en el sector centro-norte. Por último, se observa un único municipio pequeño en gris (Missing) hacia el borde occidental que no cuenta con información reportada en la muestra.

tm_shape(orinoquia_sf) +
  tm_fill("Tasa_Homicidios", 
          palette = paletteer_c("ggthemes::Purple", 30), 
          style = "quantile", 
          n = 6) + 
  tm_layout(
    legend.position = c("left", "top"),
    legend.text.size = 0.8, 
    legend.title.size = 1.0,
    inner.margins = c(0.05, 0.25, 0.05, 0.05) # Abre espacio a la izquierda
  ) + 
  tm_borders(alpha = 0.5) + 
  tm_compass(position = c("right", "bottom"))

– Análisis: El mapa coroplético de la tasa de homicidios En la Región Orinoquía, construido mediante una clasificación por cuantiles en 6 categorías permite analizar la distribución relativa de la violencia asegurando una proporción equivalente de municipios en cada intervalo. Estadísticamente, los puntos de corte de la leyenda confirman una marcada asimetría positiva (sesgo a la derecha): el 50% de los municipios de la región registra tasas iguales o inferiores a la mediana de 22.7 homicidios por cada 100.000 habitantes (los tres primeros cuantiles: 0 a 22.7, en tonos rosa claros), mientras que el sexto cuantil (el 16.7% más violento, en morado más oscuro) abarca un rango sumamente amplio y disperso que va desde 53.3 hasta 105.2. Espacialmente, al emplear cuantiles se visibilizan contrastes internos que permanecían ocultos bajo intervalos fijos: los municipios ubicados en los dos cuantiles superiores (33.7 a 105.2) se concentran en el suroccidente y el piedemonte llanero, en una franja transversal del centro-norte y en el extremo nororiental, mientras que las tasas más bajas (0 a 13.5) se agrupan en el extremo norte y en sectores puntuales del piedemonte central.

base_final$ind_rur = base_final$pobl_tot_rur_placeholder
base_final$ind_rur  = base_final$pobl_rur / base_final$pobl_tot
base_final$densidad = base_final$pobl_tot / base_final$areaoficialkm2

bf = merge(georef, base_final, by.x = "MPIO_CDPMP", by.y = "codmun")
bf_sf = sf::st_as_sf(bf)

orinoquia_sf = subset(bf_sf, gorinoquia == 1)

"ind_rur" %in% names(orinoquia_sf)

## [1] TRUE

tm_shape(orinoquia_sf) +
  tm_polygons(
    fill = "ind_rur",
    fill.scale = tm_scale_intervals(
      style  = "quantile",
      n      = 6,
      values = "brewer.purples"
    ),
    fill.legend = tm_legend(
      title = "Índice de Ruralidad",
      position = tm_pos_out("right", "center") # Saca la leyenda afuera a la derecha
    ),
    col_alpha = 0.2
  ) +
  tm_compass(position = c("right", "bottom")) +
  tm_title("Región Orinoquía") +
  tm_layout(
    legend.text.size  = 0.8,
    legend.title.size = 1.0,
    frame = TRUE
  )

– Análisis: El mapa coroplético del Índice de Ruralidad evidencia una distribución marcadamente diferenciada entre las subregiones. Estadísticamente, el punto de corte del tercer cuantil revela que la mitad de los municipios de la Orinoquía tiene al menos un 47.3% de su población residiendo en zonas rurales (mediana = 0.473), mientras que el sexto cuantil agrupa a los municipios predominantemente rurales, donde entre el 66.4% y el 92.2% de sus habitantes vive en el campo (0.664 a 0.922, en morado más oscuro). Espacialmente, los dos cuantiles de mayor ruralidad (0.560 a 0.922) se concentran en los municipios de gran extensión territorial del sur, el oriente y el suroccidente de la región, así como en una franja puntual sobre la vertiente cordillerana occidental. En contraste, los cuantiles de menor ruralidad (0.04 a 0.323, en tonos blancos y lilas claros) se localizan en el extremo norte, el nororiente y en un clúster específico del piedemonte centro-occidental, donde se asientan las principales cabeceras urbanas de la región.

tm_shape(orinoquia_sf) + 
  tm_fill("ind_rur", 
          palette = "brewer.purples", 
          style = "sd", 
          title = "Índice de Ruralidad \n(Desviaciones Estándar)") + 
  tm_borders(alpha = 0.2) + 
  tm_compass(position = c("right", "bottom")) +
  tm_layout(
    title = "Reg. Orinoquía", 
    legend.outside = TRUE,
    legend.outside.position = "right",
    legend.text.size = 0.8, 
    legend.title.size = 1.0, 
    frame = FALSE
  )

– Análisis: El mapa coroplético del Índice de Ruralidad en la Región Orinoquía, construido bajo el método de desviaciones estándar permite identificar qué municipios se alejan significativamente del promedio regional mediante intervalos de amplitud fija (approx 0.213). A diferencia de la clasificación por cuantiles —que forzaba al 16.7% de los municipios a ubicarse en la categoría más oscura, el método de desviación estándar revela que solo un grupo muy reducido de municipios en el suroccidente y en la vertiente cordillerana occidental presenta una ruralidad extremadamente alta en términos estadísticos (superando el umbral de 0.875, en morado oscuro). La gran mayoría del territorio, incluyendo las extensas sabanas del oriente, el centro y el sur, se concentra dentro de los rangos intermedios cercanos a la media (0.234 a 0.874, en tonos lilas medios), mientras que los municipios que se desvían significativamente hacia abajo (con predominio urbano e índices entre 0.020 y 0.233, en tono casi blanco) se localizan puntualmente en el extremo norte, el nororiente y las principales cabeceras del piedemonte central.

tm_shape(orinoquia_sf) +
  tm_polygons(
    fill = "ind_rur",
    fill.scale = tm_scale_intervals(
      style  = "quantile",
      n      = 5,
      values = "plasma"
    ),
    fill.legend = tm_legend(
      title    = "Índice de Ruralidad",
      position = tm_pos_out("right", "center")
    ),
    col_alpha = 0.2
  ) +
  tm_text(
    "MPIO_CNMBR", 
    size = 0.35, 
    col = "white",
    bgcol = "black",       
    bgcol_alpha = 0.5,     
    options = opt_tm_text(remove_overlap = TRUE) 
  ) +
  tm_compass(position = tm_pos_in("right", "top")) +
  tm_scalebar(position = tm_pos_in("center", "bottom")) +
  tm_layout(frame = FALSE)

– Análisis: El mapa coroplético del Índice de Ruralidad en la Región Orinoquía, clasificado mediante 5 quintiles con la paleta plasma y etiquetas de identificación municipal, permite relacionar directamente los niveles de ruralidad con territorios específicos. El quintil superior (0.638 a 0.922, en amarillo brillante), que concentra al 20% de los municipios más rurales de la región, abarca los territorios de mayor extensión geográfica en el sur y oriente —como Cumaribo, Mapiripán, La Macarena y Uribe—, además de municipios cordilleranos puntuales. Los rangos intermedios y medio-altos (0.336 a 0.637, en tonos rosados y naranjas) se extienden por las sabanas centrales y de transición, incluyendo a Puerto Gaitán, La Primavera, Santa Rosalía, Orocué, Hato Corozal y Tame. Finalmente, el primer y segundo quintil (0.047 a 0.335, en tonos morados y azules oscuros) agrupan a los municipios con mayor consolidación urbana y administrativa, tales como Arauca, Puerto Rondón y Cravo Norte en el norte, Puerto Carreño en el extremo nororiental, Paz de Ariporo y Aguazul en el sector central, y el corredor del piedemonte conformado por Villavicencio, Puerto López y San Martín.

orinoquia2 = st_as_sf(orinoquia)

orinoquia2 = st_as_sf(orinoquia)
st_write(orinoquia2,
         dsn = "RegionOrinoquia",
         layer = "Base_General_Orinoquia",
         driver = "ESRI Shapefile",
         delete_layer = TRUE,
         layer_options = "ENCODING=UTF-8")

## Deleting layer `Base_General_Orinoquia' using driver `ESRI Shapefile'
## Writing layer `Base_General_Orinoquia' to data source 
##   `RegionOrinoquia' using driver `ESRI Shapefile'
## options:        ENCODING=UTF-8 
## Writing 59 features with 33 fields and geometry type Multi Polygon.

4. Modelos de regresión

library(nortest)
library(lmtest)

4.1 Modelo Univariado

mod1 = lm(IPM.nacional ~ disbogota) 
summary(mod1)

## 
## Call:
## lm(formula = IPM.nacional ~ disbogota)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -74.396  -9.388   0.060   9.984  43.029 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 28.255674   0.802623   35.20   <2e-16 ***
## disbogota    0.049998   0.002191   22.82   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 13.67 on 1099 degrees of freedom
## Multiple R-squared:  0.3215, Adjusted R-squared:  0.3209 
## F-statistic: 520.9 on 1 and 1099 DF,  p-value: < 2.2e-16

Fórmula del modelo lineal

\[\widehat{IPM} = 28.26 + 0.05(disbogota)\]

4.2 Modelo multivariado

mod2 = lm(IPM.nacional ~ disbogota + base_final$densidad + log(areaoficialkm2) + gpacifica + gcaribe + gamazonia + gorinoquia) 
summary(mod2)

## 
## Call:
## lm(formula = IPM.nacional ~ disbogota + base_final$densidad + 
##     log(areaoficialkm2) + gpacifica + gcaribe + gamazonia + gorinoquia)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -49.594  -8.153  -0.397   7.790  41.709 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         -1.181243   1.927942  -0.613 0.540205    
## disbogota            0.030993   0.003090  10.029  < 2e-16 ***
## base_final$densidad -0.001953   0.000489  -3.995 6.90e-05 ***
## log(areaoficialkm2)  6.126897   0.344416  17.789  < 2e-16 ***
## gpacifica           -0.964205   1.125332  -0.857 0.391733    
## gcaribe              6.807615   1.512792   4.500 7.52e-06 ***
## gamazonia           -9.465666   2.188906  -4.324 1.67e-05 ***
## gorinoquia          -6.144078   1.736770  -3.538 0.000421 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 11.73 on 1093 degrees of freedom
## Multiple R-squared:  0.5033, Adjusted R-squared:  0.5002 
## F-statistic: 158.2 on 7 and 1093 DF,  p-value: < 2.2e-16

Fórmula del modelo multivariado

\[IPM_i = \beta_0 + \beta_1(disbogota_i) + \beta_2(densidad_i) + \beta_3\ln(area_i) + \beta_4(pacifica_i) + \beta_5(caribe_i) + \beta_6(amazonia_i) + \beta_7(orinoquia_i) + \varepsilon_i\]

confint(mod2)

##                             2.5 %      97.5 %
## (Intercept)          -4.964128442  2.60164328
## disbogota             0.024929035  0.03705625
## base_final$densidad  -0.002912928 -0.00099409
## log(areaoficialkm2)   5.451106349  6.80268832
## gpacifica            -3.172260361  1.24385032
## gcaribe               3.839311472  9.77591921
## gamazonia           -13.760598921 -5.17073208
## gorinoquia           -9.551858891 -2.73629655

– Análisis: Los intervalos de confianza al 95% ratifican la significancia estadística de las variables cuyos rangos no contienen el cero: la distancia a Bogotá ([0.025, 0.037]), el logaritmo del área municipal ([5.45, 6.80]) y la pertenencia a la región Caribe ([3.84, 9.78]) mantienen un efecto estrictamente positivo sobre el IPM, mientras que la densidad poblacional ([-0.0029, -0.0010]) y las regiones de Amazonía ([-13.76, -5.17]) y Orinoquía ([-9.55, -2.74]) presentan un efecto consistentemente negativo frente a la región base (Andina). En contraste, tanto el intercepto (\([-4.96, 2.60]\)) como la variable regional gpacifica ([-3.17, 1.24]) incluyen el cero dentro de sus intervalos, lo que confirma que sus efectos individuales no son estadísticamente distintos de cero a un nivel de confianza del 95%.

prediccion_m1 = predict(mod1)
prediccion_m2 = predict(mod2)
residuos_m1 = residuals(mod1)
residuos_m2 = residuals(mod2)
residuos_m2_std = rstandard(mod2)
base_final = cbind(base_final,prediccion_m1, residuos_m1 ,  prediccion_m2, residuos_m2,residuos_m2_std)

plot(IPM.nacional, prediccion_m1)

– Análisis: El gráfico de dispersión entre el IPM observado (IPM.nacional) y los valores ajustados del Modelo 1 (prediccion_m1) evidencia una tendencia positiva coherente con la capacidad explicativa moderada del modelo (R^2 = 0.3215). Sin embargo, se aprecia una alta dispersión vertical y un claro truncamiento inferior en las predicciones: dado que el intercepto del modelo es \(\beta_0 = 28.26\) y la distancia a Bogotá es siempre positiva, el modelo es incapaz de predecir valores de IPM menores a 28.26, sobreestimando sistemáticamente la pobreza en los municipios con niveles reales de IPM entre 0 y 25. Asimismo, se identifican valores atípicos pronunciados, como municipios con un IPM real inferior a \(20\) para los cuales el modelo predice valores superiores a \(90\) por encontrarse geográficamente muy alejados de la capital.

plot(IPM.nacional, prediccion_m2)

– Análisis: En comparación con el modelo simple, el gráfico de dispersión entre el IPM observado (IPM.nacional) y las predicciones del Modelo 2 (prediccion_m2) muestra una nube de puntos mucho más compacta y alineada sobre la diagonal principal, reflejando el incremento en la capacidad explicativa (R^2 = 0.5033) y la disminución del error estándar residual (11.73). Al incorporar controles demográficos, de escala territorial y regionales, se corrige en gran medida el truncamiento inferior observado en el Modelo 1, permitiendo estimar valores de IPM más bajos y precisos. No obstante, debido a la estructura lineal aditiva, aparece una observación atípica puntual con una predicción ligeramente inferior a 0 (correspondiente a un municipio con altísima densidad poblacional y un IPM real cercano a 25), ya que el modelo lineal no restringe el rango de la variable respuesta al intervalo [0, 100].

hist(residuos_m2, nclass = 20, col = "#B81840", border = "black")

– Análisis: El histograma de los residuos del Modelo 2 (residuos_m2) exhibe una distribución unimodal, simétrica y claramente acampanada (forma de campana de Gauss) centrada en torno a cero. La gran mayoría de los errores de predicción se concentran en el rango de \([-20, 20]\), disminuyendo de manera gradual y equilibrada hacia ambos extremos hasta alcanzar valores puntuales cercanos a \(-50\) y \(40\). Este comportamiento visual ratifica la ausencia de sesgos sistemáticos o asimetrías severas, confirmando nuevamente el cumplimiento del supuesto de normalidad de los errores.

plot(prediccion_m2, residuos_m2_std)

– Análisis: El gráfico de dispersión entre los valores ajustados (prediccion_m2) y los residuos estandarizados (residuos_m2_std) muestra que los errores se distribuyen alrededor de la media cero sin presentar tendencias curvas pronunciadas, lo que respalda el supuesto de linealidad y la correcta especificación funcional del modelo. Asimismo, la gran mayoría de las observaciones se concentra dentro de la banda de [-2, 2] desviaciones estándar, identificándose únicamente unos pocos valores atípicos puntuales que superan el umbral de (especialmente en la cola inferior, cerca de \(-4\)). No obstante, la amplitud vertical de la nube de puntos no es completamente uniforme a lo largo del eje horizontal —presentando mayor dispersión en los valores predichos intermedios (30 a 60).

##4.3 Inferencia sobre la estimación. ### 4.3.1 Bondad de Ajuste

AIC(mod2)

## [1] 8555.404

BIC(mod2)

## [1] 8600.44

anova(mod2)

## Analysis of Variance Table
## 
## Response: IPM.nacional
##                       Df Sum Sq Mean Sq  F value    Pr(>F)    
## disbogota              1  97300   97300 707.6231 < 2.2e-16 ***
## base_final$densidad    1   6428    6428  46.7487 1.340e-11 ***
## log(areaoficialkm2)    1  38235   38235 278.0707 < 2.2e-16 ***
## gpacifica              1    709     709   5.1592 0.0233176 *  
## gcaribe                1   6015    6015  43.7447 5.843e-11 ***
## gamazonia              1   1899    1899  13.8099 0.0002125 ***
## gorinoquia             1   1721    1721  12.5149 0.0004207 ***
## Residuals           1093 150290     138                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

– Análisis: El Modelo 2 registra un AIC=8555.40 y un BIC=8600.44, métricas que sirven como referencia para evaluar el ajuste penalizado frente a otros modelos (donde valores menores indican un mejor desempeño). A partir de la descomposición de varianza (ANOVA), la suma de cuadrados explicada (152307) frente a la total (302597) demuestra que el modelo explica el 50.33% de la variabilidad del IPM municipal (R2=0.5033). Además, todas las variables aportan significativamente al ajuste (p-valor<0.05): la distancia a Bogotá (disbogota, F=707.62) y el área municipal (log(areaoficialkm2), F=278.07) son los determinantes de mayor peso explicativo, seguidos por la densidad poblacional (F=46.75) y las diferencias regionales frente a la región base (Andina).

4.3.1.1 Normalidad

qqnorm(residuos_m2)
qqline(residuos_m2)

lillie.test(residuos_m2_std)

## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  residuos_m2_std
## D = 0.021808, p-value = 0.2312

shapiro.test(residuos_m2_std)

## 
##  Shapiro-Wilk normality test
## 
## data:  residuos_m2_std
## W = 0.99529, p-value = 0.001716

bptest(mod2)

## 
##  studentized Breusch-Pagan test
## 
## data:  mod2
## BP = 64.446, df = 7, p-value = 1.943e-11

– Análisis: La prueba de Lilliefors, ideal para muestras grandes (n=1122), confirma la normalidad de los residuos estandarizados al no rechazar la hipótesis nula (D=0.0218, p-valor=0.2312>0.05). Aunque el test de Shapiro-Wilk arroja un p-valor=0.0017, esto se debe a su hipersensibilidad ante pequeñas desviaciones en las colas en muestras grandes, pues su estadístico W=0.9953 es prácticamente igual a 1. Junto con el ajuste lineal del gráfico Q-Q, se valida el supuesto de normalidad.