Autor: Martín Alalú

Accesibilidad a centros de salud

Con el objetivo de garantizar el acceso a la salud y descomprimir la creciente demanda hacia los grandes hospitales el Gobierno de la Ciudad Autónoma de Buenos Aires impulsó la creación de Centros Médicos Barriales y Centros de Salud (CeSAC) en la Ciudad. Con esta estrategia se buscó poder establecer un primer nivel de atención de la salud que, según datos del Ministerio de Salud porteño, potencialmente puede resolver el 80 % de los problemas de salud de la población.

Más allá de si esta estrategia logró o no descromprimir los grandes hospitales el propósito del presente trabajo será poder analizar la cobertura de los centros de salud a lo largo del territorio de la Ciudad y la accesibilidad a ellos por parte de la población. ¿La cobertura de los centros de salud está equilibrada entre los barrios? ¿Cuál es la distancia promedio a los centros de salud en los distintos barrios? ¿Es la misma entre los barrios del sur y del norte? ¿Cuál es la potencial demanda de cada uno de los centros de salud?

Para poder responder estas preguntas se usarán datos abiertos publicados en el portal de datos abiertos del Gobierno de la Ciudad, data.buenosaires.gob.ar.

Principalmente se usarán los datasets de: - Radios Censales - Barrios - Hospitales - Centros Médicos Barriales - CeSAC

library(sf)
library(osmdata)
library(tidyverse)
library(ggmap)
library(knitr)
library(leaflet)

#Levantamos los datos de BA Data

radios <- st_read("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/informacion-censal-por-radio/CABA_rc.geojson")

## Reading layer `CABA_rc' from data source `http://cdn.buenosaires.gob.ar/datosabiertos/datasets/informacion-censal-por-radio/CABA_rc.geojson' using driver `GeoJSON'
## Simple feature collection with 3554 features and 8 fields
## geometry type:  MULTIPOLYGON
## dimension:      XY
## bbox:           xmin: -58.53092 ymin: -34.70574 xmax: -58.33455 ymax: -34.528
## epsg (SRID):    4326
## proj4string:    +proj=longlat +datum=WGS84 +no_defs

hospitales<- read.csv("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/hospitales/hospitales.csv", fileEncoding = "UTF-8")
centros_barriales <- read.csv("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/centros-medicos-barriales/centros-medicos-barriales.csv",fileEncoding = "UTF-8")    
cesac <- read.csv("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/centros-de-salud-y-accion-comunitaria-cesac/centros-de-salud-y-accion-comunitaria.csv",fileEncoding = "UTF-8")
barrios <- st_read("http://cdn.buenosaires.gob.ar/datosabiertos/datasets/barrios/barrios.geojson")

## Reading layer `barrios_badata' from data source `http://cdn.buenosaires.gob.ar/datosabiertos/datasets/barrios/barrios.geojson' using driver `GeoJSON'
## Simple feature collection with 48 features and 4 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: -58.53152 ymin: -34.70529 xmax: -58.33515 ymax: -34.52649
## epsg (SRID):    4326
## proj4string:    +proj=longlat +datum=WGS84 +no_defs

#Filtramos los nulos, definimos los campos de coordenadas y la proyección, creamos una columna con la zona del barrio (norte-sur) y con el tipo de establecimiento. Nos quedamos con las columnas que nos interesan
 
#Para CESAC borramos la palabra "COMUNA" de la columna asi nos queda todo igual
cesac <- cesac %>% 
  mutate(tipo="cesac") %>%
  mutate(region = case_when(grepl("NUÑEZ",barrio) | 
                              grepl("SAAVEDRA",barrio) | 
                              grepl("COGHLAN",barrio) | 
                              grepl("BELGRANO",barrio) | 
                              grepl("VILLA URQUIZA",barrio) | 
                              grepl("VILLA PUEYRREDON",barrio) | 
                              grepl("PARQUE CHAS",barrio) | 
                              grepl("VILLA ORTUZAR",barrio) | 
                              grepl("COLEGIALES",barrio) | 
                              grepl("PALERMO",barrio) | 
                              grepl("RECOLETA",barrio) | 
                              grepl("RETIRO",barrio) | 
                              grepl("SAN NICOLAS",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("ALMAGRO",barrio) | 
                              grepl("PATERNAL",barrio) | 
                              grepl("VILLA DEL PARQUE",barrio) |
                              grepl("AGRONOMÍA",barrio) | 
                              grepl("VILLA DEVOTO",barrio) | 
                              grepl("VILLA REAL",barrio) | 
                              grepl("MONTE CASTRO",barrio) | 
                              grepl("VILLA SANTA RITA",barrio) | 
                              grepl("VILLA GRAL. MITRE",barrio) | 
                              grepl("CABALLITO",barrio) | 
                              grepl("FLORESTA",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("VELEZ SARSFIELD",barrio) | 
                              grepl("VILLA LURO",barrio) ~ "NORTE",
                            TRUE~"SUR")) %>% 
  mutate(barrio=str_to_upper(barrio)) %>% 
  mutate(comuna=gsub("COMUNA","",comuna)) %>% 
  filter(!is.na(long), !is.na(lat)) %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326) %>% 
  select(nombre,comuna,barrio,tipo, calle, altura,region)

centros_barriales <- centros_barriales %>% 
  mutate(tipo="centro_barrial") %>%
  mutate(barrio=str_to_upper(barrio)) %>%
  mutate(region = case_when(grepl("NUÑEZ",barrio) | 
                              grepl("SAAVEDRA",barrio) | 
                              grepl("COGHLAN",barrio) | 
                              grepl("BELGRANO",barrio) | 
                              grepl("VILLA URQUIZA",barrio) | 
                              grepl("VILLA PUEYRREDON",barrio) | 
                              grepl("PARQUE CHAS",barrio) | 
                              grepl("VILLA ORTUZAR",barrio) | 
                              grepl("COLEGIALES",barrio) | 
                              grepl("PALERMO",barrio) | 
                              grepl("RECOLETA",barrio) | 
                              grepl("RETIRO",barrio) | 
                              grepl("SAN NICOLAS",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("ALMAGRO",barrio) | 
                              grepl("PATERNAL",barrio) | 
                              grepl("VILLA DEL PARQUE",barrio) |
                              grepl("AGRONOMÍA",barrio) | 
                              grepl("VILLA DEVOTO",barrio) | 
                              grepl("VILLA REAL",barrio) | 
                              grepl("MONTE CASTRO",barrio) | 
                              grepl("VILLA SANTA RITA",barrio) | 
                              grepl("VILLA GRAL. MITRE",barrio) | 
                              grepl("CABALLITO",barrio) | 
                              grepl("FLORESTA",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("VELEZ SARSFIELD",barrio) | 
                              grepl("VILLA LURO",barrio) ~ "NORTE",
                            TRUE~"SUR")) %>% 
  filter(!is.na(long), !is.na(lat)) %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326) %>% 
  select(nombre,comuna,barrio,tipo,calle, altura, region)

hospitales <- hospitales %>%
  mutate(tipo="hospital") %>% 
  mutate(barrio=str_to_upper(barrio)) %>% 
  mutate(region = case_when(grepl("NUÑEZ",barrio) | 
                              grepl("SAAVEDRA",barrio) | 
                              grepl("COGHLAN",barrio) | 
                              grepl("BELGRANO",barrio) | 
                              grepl("VILLA URQUIZA",barrio) | 
                              grepl("VILLA PUEYRREDON",barrio) | 
                              grepl("PARQUE CHAS",barrio) | 
                              grepl("VILLA ORTUZAR",barrio) | 
                              grepl("COLEGIALES",barrio) | 
                              grepl("PALERMO",barrio) | 
                              grepl("RECOLETA",barrio) | 
                              grepl("RETIRO",barrio) | 
                              grepl("SAN NICOLAS",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("ALMAGRO",barrio) | 
                              grepl("PATERNAL",barrio) | 
                              grepl("VILLA DEL PARQUE",barrio) |
                              grepl("AGRONOMÍA",barrio) | 
                              grepl("VILLA DEVOTO",barrio) | 
                              grepl("VILLA REAL",barrio) | 
                              grepl("MONTE CASTRO",barrio) | 
                              grepl("VILLA SANTA RITA",barrio) | 
                              grepl("VILLA GRAL. MITRE",barrio) | 
                              grepl("CABALLITO",barrio) | 
                              grepl("FLORESTA",barrio) | 
                              grepl("BALVANERA",barrio) | 
                              grepl("VELEZ SARSFIELD",barrio) | 
                              grepl("VILLA LURO",barrio) ~ "NORTE",
                            TRUE~"SUR")) %>% 
  filter(!is.na(long), !is.na(lat)) %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326) %>% 
  select(nombre,comuna,barrio,tipo, calle, altura, region)

radios <- radios %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326)

#Unimos centros de salud y centros medicos barriales
establecimientos_salud <- rbind(cesac,centros_barriales,hospitales)

#Contamos cantidad
tipos_establecimientos <-  establecimientos_salud %>% 
  count(tipo)
  
cantidad_tipo <- ggplot(tipos_establecimientos) +
    geom_col(aes(x=reorder(tipo, n),
               y=n,
               fill=tipo))+
  coord_flip()+
  labs(title = "Tipos de establecimientos de salud",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       color = "Tipo",
       fill="Tipo",
       y="Cantidad")+
  theme_minimal()

print(cantidad_tipo)

# Distribución territorial

En la Ciudad hay 44 CeSACs y 34 Centros Médicos Barriales, mientras que hospitales hay 36. Es decir que con la creación de estas nuevas unidades de atención primaria se duplicó la cantidad de efectores de salud públicos en la Ciudad.

Barracas (14) es el barrio con mayor cantidad, seguido por Villa Lugano (11), Flores (8) y Parque Patricios (7). Sin embargo, y a pesar de esta expansión, aún hay barrios que no cuentan con ningún establecimiento de salud: CONSTITUCION, BOEDO, VILLA LURO, VILLA ORTUZAR, VILLA RIACHUELO, VERSALLES, PUERTO MADERO, MONSERRAT y SAN NICOLAS

#Cantidad establecimientos por barrio
establecimientos_barrios <- count(establecimientos_salud,barrio)
st_geometry(establecimientos_barrios) <- NULL

poblacion_barrios <- radios %>% 
    group_by(BARRIO) %>% 
    summarise(total = sum(POBLACION)) %>% 
    arrange(desc(total))

st_geometry(poblacion_barrios) <- NULL
salud_barrios <- left_join(barrios, establecimientos_barrios)
salud_barrios <- left_join(salud_barrios,poblacion_barrios,by = c("barrio" = "BARRIO"))
salud_barrios <- salud_barrios %>% 
  mutate(ratio=total/n)


#Mapeamos
plot(ggplot()+
  geom_sf(data=salud_barrios, aes(fill = n))+
  scale_fill_continuous(type = "viridis")+
    labs(title = "Cantidad de establecimientos de salud por barrio",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       fill="Establecimientos")+
    theme_void())

#Cantidad Cesac
cesac_centros_barriales <- filter(establecimientos_salud,tipo!="hospital")
cuenta_cesac_cmb <- count(cesac_centros_barriales,barrio)
st_geometry(cuenta_cesac_cmb) <- NULL
cesac_cmb_barrios <- left_join(barrios, cuenta_cesac_cmb)
cesac_cmb_barrios <- cesac_cmb_barrios %>% 
  mutate(tipo="cesac_cmb")

#Cantidad hospitales
cuenta_hospitales <- count(hospitales,barrio)
st_geometry(cuenta_hospitales) <- NULL
hospitales_barrios <- left_join(barrios, cuenta_hospitales)
hospitales_barrios <- hospitales_barrios %>% 
  mutate(tipo="hospitales")
  
#Unimos todo
hospitales_cesac_cmb_barrios <- rbind(cesac_cmb_barrios,hospitales_barrios)

#Mapeamos
plot(ggplot()+
  geom_sf(data=hospitales_cesac_cmb_barrios, aes(fill = n))+
  facet_wrap(~tipo)+
  scale_fill_continuous(type = "viridis")+
    labs(title = "Cantidad de establecimientos de salud por barrio según tipo",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       fill="Establecimientos")+
    theme_void())

barrios_sin_salud <- salud_barrios%>% 
  filter(is.na(n)==TRUE) %>% 
  select(barrio) 
st_geometry(barrios_sin_salud) <- NULL

barrios_sin_salud

##             barrio
## 1     CONSTITUCION
## 2            BOEDO
## 3       VILLA LURO
## 4        AGRONOMIA
## 5    VILLA ORTUZAR
## 6  VILLA RIACHUELO
## 7        VERSALLES
## 8    PUERTO MADERO
## 9        MONSERRAT
## 10     SAN NICOLAS

view(salud_barrios)

barrios_nuevos <- spread(hospitales_cesac_cmb_barrios,tipo,n)
barrios_nuevos <- barrios_nuevos %>% 
  filter(is.na(hospitales)==TRUE & is.na(cesac_cmb)==FALSE) %>% 
  select(barrio)
st_geometry(barrios_nuevos) <- NULL

barrios_nuevos

##               barrio
## 1          CHACARITA
## 2         COLEGIALES
## 3           FLORESTA
## 4          MATADEROS
## 5      NUEVA POMPEYA
## 6              NUÑEZ
## 7  PARQUE AVELLANEDA
## 8   PARQUE CHACABUCO
## 9             RETIRO
## 10          SAAVEDRA
## 11         SAN TELMO
## 12   VELEZ SARSFIELD
## 13      VILLA CRESPO
## 14  VILLA DEL PARQUE
## 15  VILLA PUEYRREDON
## 16        VILLA REAL
## 17  VILLA SANTA RITA
## 18     VILLA SOLDATI
## 19     VILLA URQUIZA

#Demanda potencial

Más allá de la cantidad neta de establecimientos de salud por barrio resulta interesante ver cuál es la potencial demanda que tienen, es decir la cantidad de efectores relativa a la población del barrio. Por ejemplo, Parque Patricios es el barrio con mejor accesibilidad en este sentido, con 1 efector por cada 5855 habitantes, seguido por Barracas con 6389, por el otro lado Belgrano, Villa Crespo y Almagro tienen una demanda potencial muy elevada teniendo en cuenta la cantidad de establecimientos de Salud.

options(scipen = 999)
#Mapeamos
plot(ggplot()+
  geom_sf(data=salud_barrios, aes(fill = ratio))+
  scale_fill_continuous(type = "viridis")+
    labs(title = "Cantidad de establecimientos de salud por población",
       subtitle = "Ciudad Autónoma de Buenos Aires")+
    theme_void())+
  labs(fill="Establecimientos / pob")

ggplot(salud_barrios) +
    geom_col(data=filter(salud_barrios,is.na(ratio)==FALSE),aes(x=reorder(barrio, ratio),
               y=ratio))+
    coord_flip()+
  labs(title="Cantidad de población por establecimiento de salud",
       subtitle="Ciudad Autónoma de Buenos Aires")

leaflet() %>% 
  setView(lng = -58.446216, lat = -34.615766, zoom = 14) %>% 
  addTiles() %>% 
  addPolygons(data=filter(barrios,barrio=="CABALLITO"), weight = 2, fillColor = "yellow") %>% 
  addMarkers(data = establecimientos_salud,popup = ~nombre) %>% 
  addProviderTiles(providers$CartoDB.Positron)

# Calculamos la densidad poblacional de cada radio
radios <-  radios %>% 
  mutate(densidad=POBLACION/AREA_KM2)

#Mapeamos radios censales y densidad poblacional
plot(ggplot() + 
  geom_sf(data = radios, aes(fill = densidad), color = NA) +
  scale_fill_viridis_c() +
  labs(title = "Densidad de población",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       fill = "hab/km2")+
  theme_void())

# Extraemos los centroides de los barrios y los radios censales. Con st_point_on_surface nos aseguramos que el centroide caiga en el polígono mismo.
radios_c <- st_point_on_surface(radios)
barrios_c <- st_point_on_surface(barrios)

plot(ggplot() + 
  geom_sf(data=radios,size=.0)+
  geom_sf(data=radios_c, aes(color = densidad), size=1, alpha=1)+
  scale_colour_viridis_c()+
  labs(title = "Densidad de población en centroides",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       color = "hab/km2"))

binpal <- colorBin("viridis", radios$densidad, 6, pretty = FALSE)

leaflet() %>% 
  setView(lng = -58.446216, lat = -34.615766, zoom = 12) %>% 
  addTiles() %>% 
  addPolygons(data=radios, stroke=FALSE, smoothFactor = 0.2, fillOpacity = .8,
    color = ~binpal(densidad),popup=paste(radios$RADIO_ID,"<br>",radios$BARRIO,"<br>",
                                          "Densidad:",radios$densidad)) %>% 
  addProviderTiles(providers$CartoDB.Positron)

#Traemos la información del establecimiento de salud más cercano a cada centroide del radio
cesac_cmb <- filter(establecimientos_salud,tipo!="hospital")

radios_c_cesac_cmb <- st_join(radios_c, cesac_cmb, st_nearest_feature)

#Agregamos una columna con la distancia y lo pasamos a numero
radios_c_cesac_cmb <- radios_c_cesac_cmb %>% 
  mutate(distancia=st_distance(radios_c, cesac_cmb[st_nearest_feature(radios_c, cesac_cmb),], by_element = TRUE)) %>% 
  mutate(distancia=as.numeric(distancia))

#Seteamos mismas coordenadas que df de establecimientos de salud
radios_c_cesac_cmb <- radios_c_cesac_cmb %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326)

#Joineamos la info de los centroides con el poligono de radio
radios_cesac_cmb <- st_join(radios, radios_c_cesac_cmb, st_contains)

#Mapeamos radios con distancia
plot(ggplot() + 
  geom_sf(data = radios_cesac_cmb, aes(fill = distancia), color = NA) +
  scale_fill_viridis_c() +
  labs(title = "Distancia a establecimiento de salud",
       subtitle = "Ciudad Autónoma de Buenos Aires",
       fill = "Distancia")+
  theme_void())

#Mapa interactivo
binpal <- colorBin("viridis", radios_cesac_cmb$distancia, 8, pretty = FALSE)

leaflet() %>% 
  setView(lng = -58.446216, lat = -34.615766, zoom = 12) %>% 
  addTiles() %>% 
  addPolygons(data=radios_cesac_cmb, stroke=FALSE, smoothFactor = 0.2, fillOpacity = .8,
    color = ~binpal(distancia),popup=paste("Distancia:",radios_cesac_cmb$distancia)) %>% 
  addProviderTiles(providers$CartoDB.Positron)

#Promedio de distancia
mean(radios_c_cesac_cmb$distancia)

## [1] 804.0117

ggplot(radios_c_cesac_cmb, aes(x=distancia)) + 
  geom_histogram()+
  labs(title = "Distancia a centros de salud en radios censales",
       subtitle = "Ciudad Autónoma de Buenos Aires")

radios_c_cesac_cmb %>% 
  group_by(region) %>% 
  summarise(promedio_totales = mean(distancia))

## Simple feature collection with 2 features and 2 fields
## geometry type:  MULTIPOINT
## dimension:      XY
## bbox:           xmin: -58.52916 ymin: -34.69602 xmax: -58.34681 ymax: -34.53745
## epsg (SRID):    4326
## proj4string:    +proj=longlat +datum=WGS84 +no_defs
## # A tibble: 2 x 3
##   region promedio_totales                                          geometry
## * <chr>             <dbl>                                  <MULTIPOINT [°]>
## 1 NORTE              897. (-58.52916 -34.61674, -58.52911 -34.62055, -58.5~
## 2 SUR                688. (-58.52854 -34.63144, -58.52844 -34.63825, -58.5~

distancia_barrios <- radios_c_cesac_cmb %>% 
  group_by(barrio) %>% 
  summarise(promedio_totales = mean(distancia)) %>% 
  arrange(promedio_totales)

#Distancia promedio de barrios
ggplot(distancia_barrios) +
    geom_col(aes(x=reorder(barrio, promedio_totales),
               y=promedio_totales))+
    coord_flip()+
  labs(title="Distancia promedio a centro de salud por barrio",
       subtitle="Ciudad Autónoma de Buenos Aires")

norte_sur <- radios_c_cesac_cmb %>% 
  group_by(region) %>% 
  summarise(promedio_totales = mean(distancia)) %>% 
  arrange(promedio_totales)


#Distancia promedio Norte vs. Sur
ggplot(norte_sur) +
    geom_col(aes(x=reorder(region, promedio_totales),
               y=promedio_totales,fill=region))+
    coord_flip()+
  labs(title="Distancia promedio a centro de salud según región",
       subtitle="Ciudad Autónoma de Buenos Aires")

st_geometry(radios_c_cesac_cmb) <- NULL
radios_c_cesac_cmb <- as.data.frame(radios_c_cesac_cmb)

#Ordeno por distancia y agrego columna con el porcentaje acumulado de población y la población acumulada
radios_c_cesac_cmb <- radios_c_cesac_cmb %>% 
  arrange(distancia) %>% 
  mutate(distancia=round(distancia, digits = 0)) %>% 
  mutate(cum_pob_porc = (cumsum(POBLACION)/sum(POBLACION)*100)) %>% 
  mutate(cum_pob=cumsum(POBLACION))

ggplot(data=radios_c_cesac_cmb, aes(x=distancia, y=cum_pob_porc)) +
    geom_line()+
    geom_point(data=filter(radios_c_cesac_cmb,RADIO_ID=="1_3_18" | RADIO_ID=="7_17_4"),aes(color=BARRIO),size=4)+
    labs(title="Porcentaje de población según distancia a centro de salud",
         subtitle="Ciudad Autónoma de Buenos Aires",
         color="Radio")

#Detectamos outliers
ggplot(radios_cesac_cmb, aes(x =comuna , y = distancia))+
  geom_boxplot()

#Creamos dataframe con los outliers
alpha <- 1.5
isOutlier <- function(values) {
  # Determine the UL and LL
  UL <- median(values, na.rm = TRUE) + alpha*IQR(values, na.rm = TRUE)
  LL <- median(values, na.rm = TRUE) - alpha*IQR(values, na.rm = TRUE)
  values > UL | values < LL  # Return a boolean vector that can be used as a filter later on. 
}

radios_cesac_cmb_outliers <- subset(radios_cesac_cmb, isOutlier(radios_cesac_cmb$distancia))
mean(radios_cesac_cmb_outliers$distancia)

## [1] 1929.483

#Mapeamos solo outliers
binpal2 <- colorBin("viridis", radios_cesac_cmb_outliers$distancia, 6, pretty = FALSE)

leaflet() %>% 
  setView(lng = -58.446216, lat = -34.615766, zoom = 12) %>% 
  addTiles() %>% 
  addPolygons(data=radios_cesac_cmb_outliers, stroke=FALSE, smoothFactor = 0.2, fillOpacity = .8,
    color = ~binpal(distancia),popup=paste("Distancia:",radios_cesac_cmb_outliers$distancia)) %>% 
  addProviderTiles(providers$CartoDB.Positron)