Trabajo Aplicado: Análisis espacial de datos geoespaciales
Luis Alfonso Salsavilca Cayetano
2025-07-05
Librerias
Objetivos
Analizar la distribución de la cantidad de centros de salud por distrito en el Perú mediante un enfoque areal, aplicando un modelo de regresión espacial (SAR) que considere la dependencia entre distritos vecinos y variables socioespaciales seleccionadas, con énfasis comparativo en Lima Metropolitana.
Analizar la distribución espacial de la tasa poblacional de centros de salud por distrito en el Perú mediante un modelo de regresión espacial (SAR), identificando los principales factores socioespaciales asociados y evaluando su desempeño predictivo con énfasis en Lima Metropolitana.
Analizar la distribución espacial de la densidad territorial de centros de salud por distrito en el Perú mediante un modelo de regresión espacial (SAR), identificando los principales factores socioespaciales asociados y evaluando su desempeño predictivo con énfasis en Lima Metropolitana.
Antecedentes
El acceso equitativo a los servicios de salud es fundamental para garantizar el derecho a la salud y promover el desarrollo sostenible. En países como el Perú, la distribución territorial de centros de salud no solo obedece a criterios poblacionales, sino también a desigualdades estructurales de tipo socioeconómico y geográfico, especialmente evidentes en entornos urbanos como Lima Metropolitana.
A pesar de los esfuerzos por descentralizar la atención sanitaria, persisten brechas en la cobertura, muchas veces invisibilizadas por análisis tradicionales que no consideran la dimensión espacial del territorio. En este sentido, los modelos de regresión espacial, como el modelo SAR (Spatial Autoregressive Model), ofrecen herramientas más adecuadas para identificar patrones y relaciones entre variables socioespaciales, incorporando la influencia de los distritos vecinos en la distribución de servicios.
Este estudio aplica un enfoque areal para analizar tres dimensiones de la infraestructura sanitaria distrital: la cantidad de centros de salud, su tasa en relación con la población, y su densidad territorial por superficie. El uso de modelos SAR permite evaluar estas dimensiones tanto a nivel nacional como en Lima Metropolitana, integrando factores como pobreza multidimensional, servicios básicos y características demográficas, con el fin de entender mejor las dinámicas territoriales que condicionan el acceso a la salud.
Resultados (salidas de R e interpretaciones)
Leer shapefile de centros de salud
shp_peru <- st_read("centros_salud.shp")## Reading layer `centros_salud' from data source
## `F:\U-TECNOLOGIA EMERGENTE\Mineria de Datos espaciales\centros de salud\centros_salud.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 11930 features and 47 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -81.31063 ymin: -18.33542 xmax: -68.65998 ymax: -0.09737184
## Geodetic CRS: WGS 84Explorar estructura de datos
names(shp_peru)## [1] "OBJECTID" "Instituci" "Código_Ú" "Nombre_del" "Clasificac"
## [6] "Tipo" "Departamen" "Provincia" "Distrito" "UBIGEO"
## [11] "Dirección" "Código_DI" "Código_Re" "Código_Mi" "DISA"
## [16] "Red" "Microrred" "Código_UE" "Unidad_Eje" "Categoria"
## [21] "Teléfono" "Tipo_Doc_C" "Nro_Doc_Ca" "Horario" "Inicio_de_"
## [26] "Director_M" "Estado" "Situación" "Condición" "Inspecció"
## [31] "NORTE" "ESTE" "COTA" "CAMAS" "RUC"
## [36] "x" "y" "NOMDEP" "CCDD" "NOMPROV"
## [41] "CCPP" "NOMDIST" "CCDI" "IDPROV" "TIPO_SERVI"
## [46] "contacto" "descargar" "geometry"shp_peru %>% distinct(NOMDEP) %>% nrow() # Departamentos## [1] 26shp_peru %>% distinct(NOMDEP, NOMPROV) %>% nrow() # Provincias## [1] 197shp_peru %>% distinct(UBIGEO) %>% nrow() # Distritos## [1] 1858Mapa de puntos de centros de salud
ggplot() +
geom_sf(data = shp_peru, color = "darkblue", size = 0.5, alpha = 0.7) +
theme_minimal() +
labs(
title = "Centros de Salud del Perú",
subtitle = "Fuente: INEI - Geogpsperu"
)
Mapa de Peru con sus Limites
peru_mapa <- gadm(country = "PER", level = 1, path = tempdir()) %>%
st_as_sf()
ggplot() +
geom_sf(data = peru_mapa, fill = "white", color = "gray60") +
geom_sf(data = shp_peru, color = "red", size = 0.5, alpha = 0.4) +
theme_minimal() +
labs(
title = "Centros de salud en el Perú",
)
Agregar límites distritales y contar centros por distritos
distritos_peru <- gadm("PER", level = 3, path = tempdir()) %>%
st_as_sf()
shp_peru <- st_transform(shp_peru, crs = st_crs(distritos_peru))
centros_distrito <- st_join(distritos_peru, shp_peru) %>%
group_by(UBIGEO, geometry) %>%
summarise(n_centros = n(), .groups = "drop")Shapefile de centros de salud con la cantidad de centros medicos por distrito n_centros
# Leer centros de salud
shp_peru <- st_read("centros_salud.shp")## Reading layer `centros_salud' from data source
## `F:\U-TECNOLOGIA EMERGENTE\Mineria de Datos espaciales\centros de salud\centros_salud.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 11930 features and 47 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -81.31063 ymin: -18.33542 xmax: -68.65998 ymax: -0.09737184
## Geodetic CRS: WGS 84# Descargar distritos desde GADM (nivel 3)
distritos_peru <- gadm("PER", level = 3, path = tempdir()) %>%
st_as_sf()
# Igualar sistemas de coordenadas
shp_peru <- st_transform(shp_peru, crs = st_crs(distritos_peru))
# Contar centros por distrito
centros_distrito <- st_join(distritos_peru, shp_peru) %>%
group_by(UBIGEO) %>%
summarise(n_centros = n(), .groups = "drop")Mapa de Peru con datos areales
library(tmap)
tm_shape(centros_distrito) +
tm_polygons("n_centros", palette = "Blues", title = "N° de centros de salud") +
tm_layout(main.title = "Distribución de centros de salud por distrito")
## Integracion de Variable predictoras
Integracion 1: Altitud promedio (cota_media) por distrito
# Asegúrate de tener 'COTA' en formato numérico
shp_peru$COTA <- as.numeric(shp_peru$COTA)
# Promedio de altitud por distrito
elevacion_distrito <- shp_peru %>%
st_drop_geometry() %>%
group_by(UBIGEO) %>%
summarise(cota_media = mean(COTA, na.rm = TRUE))centros_datos <- centros_distrito %>%
left_join(elevacion_distrito, by = "UBIGEO")
head(centros_datos)## Simple feature collection with 6 features and 3 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -78.03814 ymin: -6.987491 xmax: -77.60135 ymax: -5.818956
## Geodetic CRS: WGS 84
## # A tibble: 6 × 4
## UBIGEO n_centros geometry cota_media
## <chr> <int> <POLYGON [°]> <dbl>
## 1 010101 52 ((-77.91036 -6.179178, -77.91026 -6.181372, -77.9… 2344.
## 2 010102 1 ((-77.73849 -6.027038, -77.73689 -6.028538, -77.7… 2822
## 3 010103 3 ((-77.87025 -6.873619, -77.87189 -6.87369, -77.87… 1772.
## 4 010104 1 ((-77.71581 -6.251809, -77.7157 -6.252254, -77.71… 2155
## 5 010105 4 ((-77.78757 -6.12462, -77.78617 -6.126541, -77.78… 2434.
## 6 010106 5 ((-77.99744 -6.911234, -77.99992 -6.913089, -78.0… 2680.Integracion 2: Casos totales de anemia por distrito a diciembre del 2024
knitr::include_graphics("imagenes/imagen_anemia.jpg")
Fuente: https://www.minsa.gob.pe/reunis/data/Indicadores_Multisectoriales_Anemia.asp
Clasificacion de los indicadores:
- Niños de 4 meses que inician suplementación de hierro.
- Niños entre 6 y 8 meses con tamizaje de anemia.
- Niños entre 6а 11 meses sin anemia, suplementados.
- Niños de 6 a 11 meses con anemia, reciben tratamiento
oportuno.
- Niños de 4 a 5 meses con visita domiciliaria.
- Niños de 6 a 11 meses con anemia y visita domiciliaria.
- Madre de niños de 6 a 8 meses, asisten a sesion demostrativa de alimentos.
Se consideraran los casos totales de anemia sumando las 7 clasificaciones de los indicadores registrados por el MINSA y agrupamos en una variable num_total .
# Leer hoja "Multis" del archivo Excel
anemia_minsa <- read_excel("Indicadores_Multisectorial_Diciembre2024.xlsx", sheet = "Multis")# Asegurar que todos estén en formato uniforme (por si hay espacios o mayúsculas)
anemia_minsa$mes <- str_trim(tolower(anemia_minsa$mes)) # ejemplo: "Diciembre" → "diciembre"
# Filtrar solo las filas correspondientes al mes de diciembre
data_dic <- anemia_minsa %>% filter(mes == "diciembre")
# Verificar que se filtró correctamente
print("Casos totales a diciembre del 2024: ")## [1] "Casos totales a diciembre del 2024: "print(nrow(data_dic))## [1] 5656# Agrupar por ubigeo y sumar todos los numX y denX
resultados_dic <- data_dic %>%
group_by(ubigeo) %>%
summarise(
num_total = sum(num1, num2, num3, num4, num5, num6, num7, na.rm = TRUE),
) %>%
ungroup()# 1. Crear columna UBIGEO con ceros a la izquierda en resultados_dic
resultados_dic <- resultados_dic %>%
mutate(UBIGEO = str_pad(as.character(ubigeo), width = 6, pad = "0"))
# 2. Asegurar que UBIGEO en centros_datos sea también texto
centros_datos <- centros_datos %>%
mutate(UBIGEO = as.character(UBIGEO))
# 3. Unir resultados_dic (con solo num_total) al dataset principal
centros_datos <- centros_datos %>%
left_join(resultados_dic[, c("UBIGEO", "num_total")], by = "UBIGEO")head(centros_datos)## Simple feature collection with 6 features and 4 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -78.03814 ymin: -6.987491 xmax: -77.60135 ymax: -5.818956
## Geodetic CRS: WGS 84
## # A tibble: 6 × 5
## UBIGEO n_centros geometry cota_media num_total
## <chr> <int> <POLYGON [°]> <dbl> <dbl>
## 1 010101 52 ((-77.91036 -6.179178, -77.91026 -6.181… 2344. 493
## 2 010102 1 ((-77.73849 -6.027038, -77.73689 -6.028… 2822 4
## 3 010103 3 ((-77.87025 -6.873619, -77.87189 -6.873… 1772. 9
## 4 010104 1 ((-77.71581 -6.251809, -77.7157 -6.2522… 2155 10
## 5 010105 4 ((-77.78757 -6.12462, -77.78617 -6.1265… 2434. 14
## 6 010106 5 ((-77.99744 -6.911234, -77.99992 -6.913… 2680. 30Integracion 3: Datos de los distrito de la plataforma nacional georeferenciados
knitr::include_graphics("imagenes/imagen_geo_distritos.jpg")
Fuente: https://visor.geoperu.gob.pe/
Vamos a integrar 32 variables de los distritos del Peru que estan relacionados con las condiciones demograficas, economicas, laboral, condiciones de vivienda, etc. que se detallaran en el tablero final de variables.
distritos <- readxl::read_excel("GeoPeru-peru_distritos.xlsx",
sheet = "Perú Distritos",
.name_repair = "universal")distritos$UBIGEO <- stringr::str_pad(as.character(distritos$cod_dist),
width = 6,
pad = "0")vars <- c("UBIGEO","total_pers","almenosu_1","nbi1_abs","nbi2_abs",
"nbi3_abs","nbi4_abs","nbi5_abs","nbi1_porc","nbi2_porc",
"nbi3_porc","nbi4_porc","nbi5_porc","num_hog","sup_tot",
"pvs_agua_r","pvs_sh","pvs_aelec","phs_inter","phs_tcelu",
"p_lees_no","p_dni","p_no_docum","pea","peao","pead","pei",
"pgr_quin1","pgr_quin2","pgr_quin3","pgr_quin4","pgr_quin5")
distritos_final <- distritos[, vars[vars %in% names(distritos)]]Tablero final con todas las variables de estudio
head(distritos_final)## # A tibble: 6 × 32
## UBIGEO total_pers almenosu_1 nbi1_abs nbi2_abs nbi3_abs nbi4_abs nbi5_abs
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 030220 4041 22.3 29 429 131 47 364
## 2 030415 2260 31.1 1 239 447 0 73
## 3 030206 708 20.2 27 72 19 5 36
## 4 030404 1803 60.7 45 166 880 92 132
## 5 030405 683 9.4 0 35 10 5 14
## 6 030304 1660 49.3 172 320 568 16 62
## # ℹ 24 more variables: nbi1_porc <dbl>, nbi2_porc <dbl>, nbi3_porc <dbl>,
## # nbi4_porc <dbl>, nbi5_porc <dbl>, num_hog <dbl>, sup_tot <dbl>,
## # pvs_agua_r <dbl>, pvs_sh <dbl>, pvs_aelec <dbl>, phs_inter <dbl>,
## # phs_tcelu <dbl>, p_lees_no <dbl>, p_dni <dbl>, p_no_docum <dbl>, pea <dbl>,
## # peao <dbl>, pead <dbl>, pei <dbl>, pgr_quin1 <dbl>, pgr_quin2 <dbl>,
## # pgr_quin3 <dbl>, pgr_quin4 <dbl>, pgr_quin5 <dbl># Verifica que ambas columnas UBIGEO sean de tipo carácter
centros_datos$UBIGEO <- as.character(centros_datos$UBIGEO)
distritos_final$UBIGEO <- as.character(distritos_final$UBIGEO)
# Unir usando left_join por UBIGEO
centros_datos_full <- centros_datos %>%
left_join(distritos_final, by = "UBIGEO")head(centros_datos_full)## Simple feature collection with 6 features and 35 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -78.03814 ymin: -6.987491 xmax: -77.60135 ymax: -5.818956
## Geodetic CRS: WGS 84
## # A tibble: 6 × 36
## UBIGEO n_centros geometry cota_media num_total total_pers
## <chr> <int> <POLYGON [°]> <dbl> <dbl> <dbl>
## 1 010101 52 ((-77.91036 -6.179178, -77.9… 2344. 493 30812
## 2 010102 1 ((-77.73849 -6.027038, -77.7… 2822 4 259
## 3 010103 3 ((-77.87025 -6.873619, -77.8… 1772. 9 1077
## 4 010104 1 ((-77.71581 -6.251809, -77.7… 2155 10 642
## 5 010105 4 ((-77.78757 -6.12462, -77.78… 2434. 14 585
## 6 010106 5 ((-77.99744 -6.911234, -77.9… 2680. 30 1730
## # ℹ 30 more variables: almenosu_1 <dbl>, nbi1_abs <dbl>, nbi2_abs <dbl>,
## # nbi3_abs <dbl>, nbi4_abs <dbl>, nbi5_abs <dbl>, nbi1_porc <dbl>,
## # nbi2_porc <dbl>, nbi3_porc <dbl>, nbi4_porc <dbl>, nbi5_porc <dbl>,
## # num_hog <dbl>, sup_tot <dbl>, pvs_agua_r <dbl>, pvs_sh <dbl>,
## # pvs_aelec <dbl>, phs_inter <dbl>, phs_tcelu <dbl>, p_lees_no <dbl>,
## # p_dni <dbl>, p_no_docum <dbl>, pea <dbl>, peao <dbl>, pead <dbl>,
## # pei <dbl>, pgr_quin1 <dbl>, pgr_quin2 <dbl>, pgr_quin3 <dbl>, …Integracion 4: Nombres de los DEPARTAMENTOS, PROVINCIA y DISTRITOS del Peru al tablero final
# Leer el shapefile si aún no está cargado
shp_peru <- st_read("F:/U-TECNOLOGIA EMERGENTE/Mineria de Datos espaciales/centros de salud/centros_salud.shp")## Reading layer `centros_salud' from data source
## `F:\U-TECNOLOGIA EMERGENTE\Mineria de Datos espaciales\centros de salud\centros_salud.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 11930 features and 47 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -81.31063 ymin: -18.33542 xmax: -68.65998 ymax: -0.09737184
## Geodetic CRS: WGS 84# Extraer UBIGEO y nombres únicos de ubicación
nombres_distritos <- shp_peru %>%
st_drop_geometry() %>% # quitamos la geometría para trabajar como data.frame
dplyr::select(UBIGEO, NOMDEP, NOMPROV, NOMDIST) %>%
distinct()centros_datos_full <- centros_datos_full %>%
left_join(nombres_distritos, by = "UBIGEO")Descripción de las variables integradas
(centros_datos_full)
| Variable | Descripción | Tipo |
|---|---|---|
UBIGEO |
Código único del distrito (6 dígitos) | chr |
n_centros |
Número de centros de salud en el distrito | int |
geometry |
Geometría espacial del distrito (objeto sf tipo
POLYGON) |
sfc |
cota_media |
Altura media del distrito (en metros sobre el nivel del mar) | dbl |
num_total |
Niños intervenidos por anemia con DNI (Dic. 2024) | dbl |
total_pers |
Total de personas en el distrito | dbl |
almenosu_1 |
% de personas con al menos una necesidad básica insatisfecha (NBI) | dbl |
nbi1_abs |
Personas con NBI por vivienda inadecuada | dbl |
nbi2_abs |
Personas con NBI por falta de servicios básicos | dbl |
nbi3_abs |
Personas con NBI por hacinamiento | dbl |
nbi4_abs |
Personas con NBI por dependencia económica | dbl |
nbi5_abs |
Personas con NBI por asistencia escolar | dbl |
nbi1_porc |
% de población con NBI por vivienda inadecuada | dbl |
nbi2_porc |
% de población con NBI por falta de servicios | dbl |
nbi3_porc |
% de población con NBI por hacinamiento | dbl |
nbi4_porc |
% de población con NBI por dependencia económica | dbl |
nbi5_porc |
% de población con NBI por asistencia escolar | dbl |
num_hog |
Número de hogares en el distrito | dbl |
sup_tot |
Superficie total del distrito (en km²) | dbl |
pvs_agua_r |
% de viviendas con acceso restringido a agua potable | dbl |
pvs_sh |
% de viviendas sin servicios higiénicos adecuados | dbl |
pvs_aelec |
% de viviendas sin acceso a electricidad | dbl |
phs_inter |
% de hogares con acceso a internet | dbl |
phs_tcelu |
% de hogares con acceso a telefonía celular | dbl |
p_lees_no |
% de personas que no saben leer | dbl |
p_dni |
% de personas con DNI | dbl |
p_no_docum |
% de personas sin ningún documento de identidad | dbl |
pea |
Población económicamente activa | dbl |
peao |
PEA ocupada | dbl |
pead |
PEA desocupada | dbl |
pei |
Población en edad de ir al colegio | dbl |
pgr_quin1 |
% de población en el quintil 1 de ingreso (más bajo) | dbl |
pgr_quin2 |
% de población en el quintil 2 de ingreso | dbl |
pgr_quin3 |
% de población en el quintil 3 de ingreso | dbl |
pgr_quin4 |
% de población en el quintil 4 de ingreso | dbl |
pgr_quin5 |
% de población en el quintil 5 de ingreso (más alto) | dbl |
NOMDEP |
Nombre del departamento | chr |
NOMPROV |
Nombre de la provincia | chr |
NOMDIST |
Nombre del distrito |
Preprocesamiento de valores N.A
Variable cota_media:
Reemplazar NA de las alturas con el promedio provincial.
cota_provincia <- centros_datos_full %>%
group_by(NOMPROV) %>%
summarise(cota_media_prom = mean(cota_media, na.rm = TRUE))cota_provincia <- cota_provincia %>% st_drop_geometry()# Unir al dataset original
centros_datos_full <- centros_datos_full %>%
left_join(cota_provincia, by = "NOMPROV") %>%
mutate(
# Reemplazar NA de cota_media por promedio provincial
cota_media = ifelse(is.na(cota_media), cota_media_prom, cota_media),
# Eliminar columna auxiliar
cota_media_prom = NULL
)sum(is.na(centros_datos_full$cota_media)) # debe ser 0## [1] 4Existen 4 valores faltantes.
Completamos la cota_media con 106 m.s.n.m en los 4 distritos de la provincia de Putumayo que son los unicos faltantes.
centros_datos_full <- centros_datos_full %>%
mutate(
cota_media = case_when(
UBIGEO %in% c("160801", "160802", "160803", "160804") ~ 106,
TRUE ~ cota_media
)
)Variable num_total de casos de anemia por distrito:
centros_datos_full <- centros_datos_full %>%
mutate(num_total = ifelse(is.na(num_total), 0, num_total)) # reemplazamos NA por 0sum(is.na(centros_datos_full$num_total)) # debe ser 0## [1] 0Completamos los valores faltantes con “0” debido a que no se registraron casos en varios de los distritos
Data completa sin NAs
head(centros_datos_full)## Simple feature collection with 6 features and 38 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -78.03814 ymin: -6.987491 xmax: -77.60135 ymax: -5.818956
## Geodetic CRS: WGS 84
## # A tibble: 6 × 39
## UBIGEO n_centros geometry cota_media num_total total_pers
## <chr> <int> <POLYGON [°]> <dbl> <dbl> <dbl>
## 1 010101 52 ((-77.91036 -6.179178, -77.9… 2344. 493 30812
## 2 010102 1 ((-77.73849 -6.027038, -77.7… 2822 4 259
## 3 010103 3 ((-77.87025 -6.873619, -77.8… 1772. 9 1077
## 4 010104 1 ((-77.71581 -6.251809, -77.7… 2155 10 642
## 5 010105 4 ((-77.78757 -6.12462, -77.78… 2434. 14 585
## 6 010106 5 ((-77.99744 -6.911234, -77.9… 2680. 30 1730
## # ℹ 33 more variables: almenosu_1 <dbl>, nbi1_abs <dbl>, nbi2_abs <dbl>,
## # nbi3_abs <dbl>, nbi4_abs <dbl>, nbi5_abs <dbl>, nbi1_porc <dbl>,
## # nbi2_porc <dbl>, nbi3_porc <dbl>, nbi4_porc <dbl>, nbi5_porc <dbl>,
## # num_hog <dbl>, sup_tot <dbl>, pvs_agua_r <dbl>, pvs_sh <dbl>,
## # pvs_aelec <dbl>, phs_inter <dbl>, phs_tcelu <dbl>, p_lees_no <dbl>,
## # p_dni <dbl>, p_no_docum <dbl>, pea <dbl>, peao <dbl>, pead <dbl>,
## # pei <dbl>, pgr_quin1 <dbl>, pgr_quin2 <dbl>, pgr_quin3 <dbl>, …Índice de Moran
# Crear vecinos espaciales a partir de los polígonos
nb <- poly2nb(centros_datos_full, snap = 1e-6)
# Detectar distritos sin vecinos
sin_vecinos <- which(card(nb) == 0)
# Crear lista de pesos espaciales (matriz W)
lw <- nb2listw(nb, style = "W", zero.policy = TRUE)
# Test de Moran con todos los distritos
cat("Índice de Moran (con todos los distritos):\n")## Índice de Moran (con todos los distritos):print(moran.test(centros_datos_full$n_centros, lw, zero.policy = TRUE))##
## Moran I test under randomisation
##
## data: centros_datos_full$n_centros
## weights: lw
##
## Moran I statistic standard deviate = 33.906, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3890083919 -0.0005382131 0.0001319963# Test sin distritos sin vecinos (solo si existen)
if (length(sin_vecinos) > 0) {
centros_filtrado <- centros_datos_full[-sin_vecinos, ]
nb_filtrado <- poly2nb(centros_filtrado)
lw_filtrado <- nb2listw(nb_filtrado, style = "W")
cat("\nÍndice de Moran (sin distritos sin vecinos):\n")
print(moran.test(centros_filtrado$n_centros, lw_filtrado))
} Valor de Moran’s I: 0.389 El valor 0.389 sugiere una moderada autocorrelación espacial positiva.
p-valor < 2.2e-16 Rechazamos la hipótesis nula de aleatoriedad.
Existe una autocorrelación espacial significativa y positiva del número de centros de salud entre distritos. Es decir, los distritos que tienen muchos centros de salud tienden a estar cerca de otros distritos que también tienen muchos, y viceversa.
variograma
# Asegúrate que centros_datos_full es un sf con geometría POLYGON
# Extraer los centroides (puntos representativos por distrito)
centros_points <- st_centroid(centros_datos_full)
# Convertir a objeto Spatial para usar gstat
centros_points_sp <- as(centros_points, "Spatial")# Definir objeto gstat con la variable a interpolar
gstat_model <- gstat(formula = n_centros ~ 1, data = centros_points_sp)# Calcular el variograma empírico
variograma_emp <- variogram(gstat_model)
# Visualizar
plot(variograma_emp, main = "Variograma empírico de n_centros")
# Ajustar modelo teórico al variograma empírico
variograma_teorico <- fit.variogram(variograma_emp, model = vgm(psill = 150, model = "Sph", nugget = 30, range = 400))
# Visualizar ajuste
plot(variograma_emp, variograma_teorico, main = "Ajuste del modelo esférico al variograma empírico")
El variograma empírico obtenido describe cómo varía el número de centros de salud entre distritos en función de la distancia. Si bien esta herramienta se usa comúnmente en datos continuos o puntuales (como temperatura o contaminación), en este caso estamos aplicándola sobre datos areales, es decir, valores agregados por unidades geográficas (distritos).
En este tipo de datos, la interpretación del variograma debe realizarse con cautela, ya que no se trata de una variable medida en puntos específicos, sino de un valor resumen por área. Aun así, el variograma puede ofrecer una visión
Cabe destacar que el uso de kriging, una técnica de interpolación espacial basada en el variograma, asume variables continuas en el espacio. Por tanto, su aplicación directa sobre datos areales no busca interpolar dentro de cada distrito, sino entender si existe una continuidad espacial que permita detectar regiones con valores similares. Esta interpretación es válida siempre que se comprenda que el análisis se realiza entre áreas, no dentro de ellas
Analisis areal para la variable predictora Cantidad de centros de salud por distrito
Modelo completo para el Perú
library(spdep)
modelo_sar <- lagsarlm(
n_centros ~ num_total + total_pers + almenosu_1 + sup_tot + cota_media,
data = centros_datos_full,
listw = lw # lista de vecinos con style = "W"
)
summary(modelo_sar)##
## Call:lagsarlm(formula = n_centros ~ num_total + total_pers + almenosu_1 +
## sup_tot + cota_media, data = centros_datos_full, listw = lw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -68.23952 -2.01011 -0.79559 0.78205 103.20320
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.2974e+00 6.4572e-01 3.5578 0.0003739
## num_total 2.6383e-03 7.6444e-04 3.4513 0.0005578
## total_pers 1.4630e-04 7.3703e-06 19.8501 < 2.2e-16
## almenosu_1 -2.7640e-02 1.0476e-02 -2.6384 0.0083287
## sup_tot 6.5594e-04 1.0227e-04 6.4135 1.422e-10
## cota_media 1.2966e-04 1.4309e-04 0.9061 0.3648724
##
## Rho: 0.20774, LR test value: 77.123, p-value: < 2.22e-16
## Asymptotic standard error: 0.027491
## z-value: 7.5569, p-value: 4.13e-14
## Wald statistic: 57.107, p-value: 4.13e-14
##
## Log likelihood: -6388.731 for lag model
## ML residual variance (sigma squared): 56.475, (sigma: 7.515)
## Number of observations: 1858
## Number of parameters estimated: 8
## AIC: 12793, (AIC for lm: 12869)
## LM test for residual autocorrelation
## test value: 0.02474, p-value: 0.87502El modelo SAR logra explicar bien la variabilidad espacial en el número de centros de salud por distrito, identificando efectos significativos de factores como población, extensión territorial y pobreza.
# Filtrar filas con UBIGEO válido (eliminar la última fila basura)
centros_datos_full <- centros_datos_full %>%
filter(!is.na(UBIGEO))Mapa predicho de todo el Perú
# Redondear los valores predichos
centros_datos_full$pred_sar <- round(fitted(modelo_sar))
# Categorizar manualmente para clases más interpretables
centros_datos_full$pred_cat <- cut(
centros_datos_full$pred_sar,
breaks = c(0, 2, 5, 10, 20, 50, 100, Inf),
labels = c("1-2", "3-5", "6-10", "11-20", "21-50", "51-100", "101+"),
right = FALSE
)
# Crear mapa con clases definidas
library(tmap)
tm_shape(centros_datos_full %>% filter(!is.na(pred_cat))) +
tm_polygons("pred_cat",
palette = "YlOrRd",
title = "Centros de salud\n(predicción SAR)",
colorNA = "transparent",
textNA = "") +
tm_layout(
main.title = "Predicción de centros de salud por distrito (modelo SAR)",
legend.outside = TRUE
)
Interpretación:
- Muestra la cantidad estimada por el modelo en base a las variables y
la estructura espacial.
- Reproduce correctamente patrones regionales: alta densidad en Lima,
sierra central y costa norte.
- El uso de clases como 1–2, 3–5, 6–10, etc., permite visualizar mejor agrupamientos regionales.
Mapa real de todo el Perú
library(tmap)
tm_shape(centros_distrito) +
tm_polygons("n_centros",
palette = "Blues",
title = "N° de centros de salud",
colorNA = "transparent", # Oculta polígonos con NA
textNA = "") + # Oculta la etiqueta "Missing"
tm_layout(
main.title = "Distribucion de centros de salud por distrito",
legend.outside = TRUE
)
Interpretación:
- Muestra la cantidad observada de centros por distrito.
- Zonas costeras y Lima destacan con mayores concentraciones (hasta
140).
- Muchas zonas rurales tienen entre 1 y 20 centros.
Distritos de Lima Metropolitana
Cantidad de centros en LM
# Filtrar solo Lima Metropolitana
lima_metropolitana <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA", !is.na(pred_cat))
# Mapa solo para Lima Metropolitana
tm_shape(lima_metropolitana) +
tm_polygons("pred_cat",
palette = "YlOrRd",
title = "Centros de salud\n(predicción SAR)",
colorNA = "transparent",
textNA = "") +
tm_layout(
main.title = "Predicción SAR - Lima Metropolitana",
legend.outside = TRUE
)
# Filtrar solo Lima Metropolitana
centros_distrito_lima <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA")
# Mapa real para Lima Metropolitana
library(tmap)
tm_shape(centros_distrito_lima) +
tm_polygons("n_centros",
palette = "Blues",
title = "N° de centros de salud",
colorNA = "transparent",
textNA = "") +
tm_layout(
main.title = "Distribución de centros de salud por distrito (Lima Metropolitana)",
legend.outside = TRUE
)
Interpretación:
Mapa observado: muestra la distribución real de centros de salud, concentrados fuertemente en distritos como Surco, San Juan de Lurigancho, Los Olivos y Miraflores.
Mapa predicho: reproduce parcialmente este patrón, pero subestima fuertemente en distritos con alta concentración.
library(dplyr)
tabla_lima <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA") %>%
st_drop_geometry() %>%
dplyr::select(NOMDIST, n_centros, pred_sar) %>%
mutate(
diferencia = pred_sar - n_centros,
error_absoluto = abs(diferencia)
) %>%
arrange(desc(error_absoluto))print(tabla_lima) # o usar print(tabla_lima)## # A tibble: 43 × 5
## NOMDIST n_centros pred_sar diferencia error_absoluto
## <chr> <int> <dbl> <dbl> <dbl>
## 1 SANTIAGO DE SURCO 136 63 -73 73
## 2 SAN JUAN DE LURIGANCHO 113 181 68 68
## 3 MIRAFLORES 86 27 -59 59
## 4 SAN MARTIN DE PORRES 56 114 58 58
## 5 LOS OLIVOS 113 61 -52 52
## 6 SAN ISIDRO 73 21 -52 52
## 7 SAN BORJA 73 29 -44 44
## 8 SANTA ANITA 85 44 -41 41
## 9 CARABAYLLO 26 66 40 40
## 10 JESUS MARIA 63 23 -40 40
## # ℹ 33 more rowsInterpretación:
De los 43 distritos de Lima Metropolitana:
- Distritos con mayor subestimación:
– Santiago de Surco (−73 centros)
–San Juan de Lurigancho (−68)
– Miraflores (−59)
– San Martín de Porres (−58)
- Predicciones más precisas (error = 0):
– Chorrillos, Pueblo Libre, Chaclacayo
- Error absoluto medio (estimado visualmente): alrededor de 20–30 centros
Análisis areal para la variable predictora tasa poblacional
centros_datos_full <- centros_datos_full %>%
mutate(
tasa_poblacional = (n_centros / total_pers) * 10000 # Centros por 10,000 personas
)library(spdep)
# Asegúrate de que la geometría esté bien definida
centros_vecinos <- poly2nb(centros_datos_full) # vecinos
lw <- nb2listw(centros_vecinos, style = "W") # lista de pesoslibrary(dplyr)
# Convertir el sf en data frame base (sin geometría)
df <- as.data.frame(st_drop_geometry(centros_datos_full))
# Seleccionar columnas predictoras quitando las no deseadas
predictoras <- df %>%
dplyr::select(-any_of(c(
"UBIGEO", "n_centros", "tasa_poblacional", "densidad_territorial",
"NOMDEP", "NOMPROV", "NOMDIST", "pred_sar", "pred_cat"
))) %>%
colnames()# Crear fórmula dinámica para el modelo SAR
formula_sar <- as.formula(
paste("tasa_poblacional ~", paste(predictoras, collapse = " + ")))Modelo completo
modelo_sar_completo <- lagsarlm(
formula = formula_sar,
data = centros_datos_full,
listw = lw)
summary(modelo_sar_completo)##
## Call:lagsarlm(formula = formula_sar, data = centros_datos_full, listw = lw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -23.37983 -4.27031 -0.81195 2.74387 53.02265
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## (1 not defined because of singularities)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -9.8516e+00 2.7384e+01 -0.3598 0.7190324
## cota_media 6.7242e-04 2.1070e-04 3.1913 0.0014162
## num_total 5.9765e-03 1.5503e-03 3.8551 0.0001157
## total_pers 2.6950e-04 1.4167e-04 1.9024 0.0571246
## almenosu_1 6.3262e-03 4.8763e-02 0.1297 0.8967758
## nbi1_abs -1.0062e-04 8.2922e-05 -1.2134 0.2249856
## nbi2_abs -4.8742e-04 2.2010e-04 -2.2146 0.0267879
## nbi3_abs -4.7054e-04 1.6781e-04 -2.8040 0.0050466
## nbi4_abs -1.1607e-03 1.1527e-03 -1.0069 0.3139612
## nbi5_abs -6.1045e-04 5.1715e-04 -1.1804 0.2378365
## nbi1_porc 4.6645e-02 3.2224e-02 1.4475 0.1477535
## nbi2_porc 6.2211e-02 3.8612e-02 1.6112 0.1071407
## nbi3_porc 3.7240e-02 3.4682e-02 1.0738 0.2829248
## nbi4_porc 2.7575e-01 1.2481e-01 2.2094 0.0271441
## nbi5_porc 8.6847e-02 6.8631e-02 1.2654 0.2057220
## num_hog -8.8435e-04 2.4343e-04 -3.6329 0.0002803
## sup_tot 2.3022e-04 1.3160e-04 1.7493 0.0802323
## pvs_agua_r -1.3396e-02 1.1300e-02 -1.1855 0.2358189
## pvs_sh 2.1770e-02 1.1401e-02 1.9094 0.0562059
## pvs_aelec -1.2947e-02 1.7088e-02 -0.7576 0.4486639
## phs_inter 9.0347e-02 2.7879e-02 3.2406 0.0011927
## phs_tcelu 6.4711e-02 1.8142e-02 3.5669 0.0003612
## p_lees_no -1.3027e-01 5.5049e-02 -2.3665 0.0179562
## p_dni 2.9477e-01 2.7809e-01 1.0600 0.2891627
## p_no_docum 1.3105e+00 5.6526e-01 2.3184 0.0204293
## pea -1.4993e-03 9.2615e-04 -1.6188 0.1054882
## peao 1.6810e-03 9.3581e-04 1.7963 0.0724543
## pead NA NA NA NA
## pei -1.1753e-04 2.1297e-04 -0.5518 0.5810585
## pgr_quin1 1.2564e-01 1.8450e-01 0.6810 0.4958796
## pgr_quin2 -5.9110e-01 1.9523e-01 -3.0277 0.0024639
## pgr_quin3 -9.0415e-01 1.7367e-01 -5.2060 1.929e-07
## pgr_quin4 -4.9891e-01 1.7403e-01 -2.8668 0.0041467
## pgr_quin5 -1.3996e+00 1.6140e-01 -8.6720 < 2.2e-16
##
## Rho: 0.49694, LR test value: 294.79, p-value: < 2.22e-16
## Asymptotic standard error: 0.027274
## z-value: 18.22, p-value: < 2.22e-16
## Wald statistic: 331.98, p-value: < 2.22e-16
##
## Log likelihood: -6431.596 for lag model
## ML residual variance (sigma squared): 57.451, (sigma: 7.5796)
## Number of observations: 1858
## Number of parameters estimated: 35
## AIC: 12933, (AIC for lm: 13226)
## LM test for residual autocorrelation
## test value: 0.0012815, p-value: 0.97144Interpretación:
El modelo SAR con todas las variables explica adecuadamente la distribución espacial de la tasa poblacional de centros de salud por distrito. Se identificaron factores clave como altitud, intervenciones en anemia, pobreza estructural y conectividad digital. La autocorrelación espacial fue significativa y correctamente modelada.
Modelo simplificado
formula_tasa_opt <- tasa_poblacional ~ num_total + total_pers + nbi2_abs + nbi3_abs +
nbi4_abs + num_hog + sup_tot + phs_inter + phs_tcelu +
p_lees_no + p_no_docum + pgr_quin2 + pgr_quin3 +
pgr_quin4 + pgr_quin5
modelo_tasa_opt <- lagsarlm(
formula = formula_tasa_opt,
data = centros_datos_full,
listw = lw
)
summary(modelo_tasa_opt)##
## Call:lagsarlm(formula = formula_tasa_opt, data = centros_datos_full,
## listw = lw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -22.16815 -4.42751 -0.93485 2.71128 52.93357
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5375e+01 2.4401e+00 6.3007 2.962e-10
## num_total 1.9438e-03 1.1259e-03 1.7265 0.0842496
## total_pers 2.3507e-04 6.8666e-05 3.4234 0.0006185
## nbi2_abs -2.6501e-04 1.8808e-04 -1.4091 0.1588201
## nbi3_abs -4.5497e-04 1.3961e-04 -3.2589 0.0011185
## nbi4_abs -1.1889e-03 9.8647e-04 -1.2052 0.2281104
## num_hog -7.2768e-04 2.2667e-04 -3.2103 0.0013260
## sup_tot 2.1044e-04 1.2732e-04 1.6528 0.0983735
## phs_inter 1.4240e-01 2.3770e-02 5.9910 2.085e-09
## phs_tcelu 1.0090e-01 1.5830e-02 6.3742 1.840e-10
## p_lees_no -7.5660e-02 4.5790e-02 -1.6523 0.0984643
## p_no_docum 9.9128e-01 3.0967e-01 3.2011 0.0013691
## pgr_quin2 -4.0705e-01 1.4543e-01 -2.7990 0.0051265
## pgr_quin3 -8.4365e-01 1.7009e-01 -4.9599 7.053e-07
## pgr_quin4 -5.3731e-01 1.7117e-01 -3.1391 0.0016945
## pgr_quin5 -1.2117e+00 1.4807e-01 -8.1835 2.220e-16
##
## Rho: 0.52895, LR test value: 358.04, p-value: < 2.22e-16
## Asymptotic standard error: 0.026143
## z-value: 20.233, p-value: < 2.22e-16
## Wald statistic: 409.38, p-value: < 2.22e-16
##
## Log likelihood: -6461.142 for lag model
## ML residual variance (sigma squared): 58.998, (sigma: 7.681)
## Number of observations: 1858
## Number of parameters estimated: 18
## AIC: 12958, (AIC for lm: 13314)
## LM test for residual autocorrelation
## test value: 0.15297, p-value: 0.69571Interpretación:
- Rho = 0.528, altamente significativo → evidencia fuerte de
autocorrelación espacial correctamente modelada.
- LM test de residuos p = 0.69 → sin autocorrelación residual, buen ajuste.
Gráfico de “tasa poblacional” para Lima Metropolitana
centros_datos_full$pred_tasa <- predict(modelo_tasa_opt)
centros_datos_full$pred_tasa <- pmax(centros_datos_full$pred_tasa, 0)
# Mapa predicho
lima_metropolitana_pred_tasa <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA", !is.na(pred_tasa))
tm_shape(lima_metropolitana_pred_tasa) +
tm_polygons("pred_tasa",
palette = "YlOrRd",
title = "Centros de salud\npor 10,000 hab. (Predicción)",
style = "quantile",
colorNA = "transparent", textNA = "") +
tm_layout(
main.title = "Predicción SAR de tasa poblacional\nLima Metropolitana",
legend.outside = TRUE
)
# Mapa real de tasa poblacional
lima_metropolitana_tasa_real <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA", !is.na(tasa_poblacional))
tm_shape(lima_metropolitana_tasa_real) +
tm_polygons("tasa_poblacional",
palette = "Blues",
title = "Centros de salud\npor 10,000 hab. (Real)",
style = "quantile",
colorNA = "transparent", textNA = "") +
tm_layout(
main.title = "Densidad poblacional real de centros de salud\nLima Metropolitana",
legend.outside = TRUE
)
Interpretación:
Predicción SAR (Mapa rojo-amarillo)
– El modelo predice mayor tasa de centros de salud por habitante en distritos periféricos y vulnerables como Villa El Salvador, Pachacámac o Punta Hermosa.
– Los distritos del centro urbano aparecen con tasas más bajas, posiblemente por alta densidad poblacional.Observado (Mapa azul)
– Muestra mayor densidad poblacional de centros en distritos como Jesús María, Breña, San Luis o El Agustino, es decir, zonas urbanas centrales.
– Las zonas periféricas tienen menor tasa observada, lo que indica que el modelo sobrestima la cobertura en periferia y subestima en zonas densamente urbanas.
Análisis areal para la variable predictora densidad territorial
# Redefinir densidad por cada 10,000 km²
centros_datos_full$densidad_territorial_10k <- (centros_datos_full$n_centros / centros_datos_full$sup_tot) * 10000# Crear fórmula igual que antes
formula_sar_densidad <- as.formula(
paste("densidad_territorial_10k ~", paste(predictoras, collapse = " + "))
)Modelo completo
# Ajustar modelo SAR
modelo_sar_densidad <- lagsarlm(
formula = formula_sar_densidad,
data = centros_datos_full,
listw = lw
)
summary(modelo_sar_densidad)##
## Call:lagsarlm(formula = formula_sar_densidad, data = centros_datos_full,
## listw = lw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -49007.95 -1466.35 260.51 1709.53 187378.16
##
## Type: lag
## Coefficients: (numerical Hessian approximate standard errors)
## (1 not defined because of singularities)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 7.6431e+04 1.3832e+04 5.5259 3.279e-08
## cota_media 7.1746e-01 2.1534e-01 3.3317 0.0008632
## num_total -3.8537e+00 1.5735e+00 -2.4492 0.0143188
## total_pers -9.1397e-01 1.4517e-01 -6.2959 3.057e-10
## almenosu_1 -2.2787e+01 NaN NaN NaN
## nbi1_abs -6.3323e-03 NaN NaN NaN
## nbi2_abs 1.9843e-01 2.0855e-01 0.9515 0.3413580
## nbi3_abs -2.3440e-03 NaN NaN NaN
## nbi4_abs 2.7145e+00 1.1232e+00 2.4167 0.0156611
## nbi5_abs 9.0936e-01 4.8424e-01 1.8779 0.0603958
## nbi1_porc 4.7948e+01 1.3078e+01 3.6665 0.0002459
## nbi2_porc 4.7856e+01 2.1453e+01 2.2308 0.0256952
## nbi3_porc 6.5882e+00 NaN NaN NaN
## nbi4_porc -1.5542e+02 1.2656e+02 -1.2280 0.2194419
## nbi5_porc 1.2478e+02 5.1321e+01 2.4314 0.0150389
## num_hog 9.2313e-01 2.5600e-01 3.6059 0.0003110
## sup_tot -6.9258e-02 1.3925e-02 -4.9737 6.570e-07
## pvs_agua_r -1.2921e+01 1.1935e+01 -1.0825 0.2790145
## pvs_sh 4.2194e+01 1.1888e+01 3.5494 0.0003861
## pvs_aelec -1.4397e+01 1.9743e+01 -0.7292 0.4658860
## phs_inter -5.7818e+02 3.1643e+01 -18.2720 < 2.2e-16
## phs_tcelu 6.7226e+01 1.9339e+01 3.4761 0.0005087
## p_lees_no -2.8519e-01 NaN NaN NaN
## p_dni -2.1942e+02 1.3778e+02 -1.5925 0.1112693
## p_no_docum -1.1741e+02 NaN NaN NaN
## pea -2.4717e+00 9.6834e-01 -2.5525 0.0106950
## peao 3.4577e+00 9.7960e-01 3.5297 0.0004160
## pead NA NA NA NA
## pei 6.4965e-01 2.2599e-01 2.8746 0.0040448
## pgr_quin1 -9.9512e+01 1.4285e+02 -0.6966 0.4860425
## pgr_quin2 -1.6376e+02 1.4869e+02 -1.1014 0.2707394
## pgr_quin3 2.9199e+01 NaN NaN NaN
## pgr_quin4 -7.3276e+01 1.1679e+02 -0.6274 0.5303869
## pgr_quin5 -4.4466e+02 1.2828e+02 -3.4662 0.0005278
##
## Rho: 0.26543, LR test value: 61.702, p-value: 3.9968e-15
## Approximate (numerical Hessian) standard error: 0.032728
## z-value: 8.1101, p-value: 4.4409e-16
## Wald statistic: 65.773, p-value: 5.5511e-16
##
## Log likelihood: -19408.82 for lag model
## ML residual variance (sigma squared): 68711000, (sigma: 8289.2)
## Number of observations: 1858
## Number of parameters estimated: 35
## AIC: 38888, (AIC for lm: 38947)Interpretación:
Rho = 0.265, p < 0.001 → existe autocorrelación espacial moderada, correctamente modelada.
Wald y LR tests significativos → el modelo es globalmente significativo.
AIC = 38888, ligeramente mejor que el modelo lineal equivalente.
Algunos coeficientes no se estiman (NaN) por colinealidad o redundancia, lo cual sugiere revisar multicolinealidad en siguientes pasos.
Modelo simplificado
modelo_densidad_final <- lagsarlm(
densidad_territorial_10k ~ cota_media + num_total + total_pers + nbi4_abs +
nbi1_porc + nbi2_porc + nbi5_porc + num_hog + sup_tot +
pvs_sh + phs_inter + phs_tcelu + pea + peao + pei + pgr_quin5,
data = centros_datos_full,
listw = lw
)
summary(modelo_densidad_final)##
## Call:lagsarlm(formula = densidad_territorial_10k ~ cota_media + num_total +
## total_pers + nbi4_abs + nbi1_porc + nbi2_porc + nbi5_porc +
## num_hog + sup_tot + pvs_sh + phs_inter + phs_tcelu + pea +
## peao + pei + pgr_quin5, data = centros_datos_full, listw = lw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -47762.93 -1346.28 224.25 1654.13 187803.30
##
## Type: lag
## Coefficients: (numerical Hessian approximate standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 52282.75495 2921.13873 17.8981 < 2.2e-16
## cota_media 0.81370 0.20126 4.0431 5.276e-05
## num_total -2.05987 1.32374 -1.5561 0.1196841
## total_pers -0.92645 0.14151 -6.5468 5.880e-11
## nbi4_abs 3.22272 0.88933 3.6238 0.0002904
## nbi1_porc 27.98084 14.06405 1.9895 0.0466429
## nbi2_porc 26.53856 25.94755 1.0228 0.3064133
## nbi5_porc 94.62662 41.60572 2.2744 0.0229440
## num_hog 0.96881 0.25273 3.8333 0.0001264
## sup_tot -0.13460 0.12583 -1.0697 0.2847513
## pvs_sh 28.62894 10.22152 2.8008 0.0050968
## phs_inter -575.90846 31.36498 -18.3615 < 2.2e-16
## phs_tcelu 58.12719 17.59143 3.3043 0.0009522
## pea -2.82420 0.99644 -2.8343 0.0045926
## peao 3.75988 1.01685 3.6976 0.0002177
## pei 0.78569 0.20078 3.9132 9.110e-05
## pgr_quin5 -537.68448 129.40258 -4.1551 3.251e-05
##
## Rho: 0.27216, LR test value: 65.425, p-value: 5.5511e-16
## Approximate (numerical Hessian) standard error: 0.03266
## z-value: 8.3333, p-value: < 2.22e-16
## Wald statistic: 69.443, p-value: < 2.22e-16
##
## Log likelihood: -19415.77 for lag model
## ML residual variance (sigma squared): 69194000, (sigma: 8318.3)
## Number of observations: 1858
## Number of parameters estimated: 19
## AIC: 38870, (AIC for lm: 38933)El parámetro de autocorrelación espacial (rho = 0.27) fue estadísticamente significativo (p < 0.01), lo que confirma que existen dependencias espaciales que justifican el uso del modelo SAR. Además, el AIC disminuyó respecto al modelo completo, lo que indica una mejora en la eficiencia del modelo con un menor número de variables.
Este modelo permite comprender mejor cómo se distribuyen los centros de salud considerando las características del territorio, y proporciona una base sólida para orientar políticas de equidad en el acceso a la salud en zonas desatendidas.
centros_datos_full$pred_densidad <- fitted(modelo_densidad_final)Gráfico de “densidad_territorial” para Lima Metropolitana
# Reemplazar valores negativos con cero en pred_densidad
centros_datos_full$pred_densidad <- pmax(centros_datos_full$pred_densidad, 0)
# Filtrar Lima Metropolitana
lima_metropolitana_dens <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA", !is.na(pred_densidad))
# Mapa de predicción SAR corregido
tm_shape(lima_metropolitana_dens) +
tm_polygons("pred_densidad",
palette = "YlGnBu",
title = "Centros de salud\npor 10,000 km² (Predicción)",
style = "quantile",
colorNA = "transparent", textNA = "") +
tm_layout(
main.title = "Predicción SAR de densidad territorial\nLima Metropolitana",
legend.outside = TRUE
)
# Filtrar Lima Metropolitana con datos reales
lima_metropolitana_real <- centros_datos_full %>%
filter(NOMDEP == "LIMA", NOMPROV == "LIMA", !is.na(densidad_territorial_10k))
# Mapa real
tm_shape(lima_metropolitana_real) +
tm_polygons("densidad_territorial_10k",
palette = "Blues",
title = "Centros de salud\npor 10,000 km² (Real)",
style = "quantile",
colorNA = "transparent", textNA = "") +
tm_layout(
main.title = "Densidad territorial real de centros de salud\nLima Metropolitana",
legend.outside = TRUE
)
Interpretación:
El modelo predice una mayor concentración territorial de centros de salud en los distritos costeros y céntricos de Lima Metropolitana, destacando áreas como Lima Cercado, San Juan de Lurigancho, y Villa El Salvador, donde se observa una alta densidad predicha (> 34,000 por 10,000 km²).
La distribución real muestra un patrón similar, con los distritos urbanos más densamente poblados también mostrando una alta densidad territorial de centros de salud. Sin embargo, algunos distritos como Chaclacayo o Cieneguilla presentan valores observados más bajos que lo predicho.
Conclusiones
El modelo SAR aplicado a datos areales permitió explicar parcialmente la variación territorial en la cantidad de centros de salud por distrito, considerando factores como población total, intervenciones por anemia, superficie del distrito y pobreza. A nivel nacional, el modelo capturó patrones espaciales generales; sin embargo, en Lima Metropolitana presentó limitaciones al subestimar valores en distritos con alta concentración de servicios, lo que sugiere que variables adicionales específicas del contexto urbano podrían mejorar su desempeño.
El modelo SAR optimizado permitió identificar que variables como el acceso a tecnologías (internet y celular), la pobreza multidimensional, el número de hogares y los quintiles de ingreso están significativamente asociadas a la tasa de centros de salud por habitante. Aunque el modelo mostró un buen ajuste espacial general, en Lima Metropolitana tiende a sobreestimar la cobertura en zonas periféricas y subestimar en distritos céntricos y densamente urbanizados, lo que sugiere la necesidad de incluir variables propias del entorno urbano para una mejor predicción local.
El modelo SAR permitió capturar adecuadamente los patrones espaciales de la densidad territorial de centros de salud en Lima Metropolitana, revelando que variables como la altitud media, la población total, el acceso a servicios básicos y los niveles de pobreza explican significativamente la distribución. La comparación entre los valores predichos y observados muestra un buen ajuste en zonas urbanas centrales, aunque persisten ciertas desviaciones en distritos periféricos, lo que sugiere oportunidades de mejora en la asignación territorial de infraestructura sanitaria.
Bibliografía
Ministerio de Salud del Perú. (s.f.). Indicadores multisectoriales – Anemia. REUNIS. https://www.minsa.gob.pe/reunis/data/Indicadores_Multisectoriales_Anemia.asp
GeoPerú. (s.f.). Visor geográfico nacional. https://visor.geoperu.gob.pe/
Gobierno del Perú. (s.f.). Plataforma Nacional de Datos Abiertos – Anemia. https://www.datosabiertos.gob.pe/dataset/anemia
Gobierno del Perú. (s.f.). Sistema de Información REUNIS. https://www.minsa.gob.pe/reunis/
Gobierno del Perú. (s.f.). Acciones del gobierno frente a la anemia. https://www.gob.pe/740