Contexto geográfico

Incorporación de la zona municipal de Torreón en el polígono urbano Uso de consultas para la lectura filtrada del Shapefile del contexto estatal filtrado municipal

Fuentes de datos

https://inegi.org.mx/contenidos/productos/prod_serv/contenidos/espanol/bvinegi/productos/geografia/marcogeo/794551132173_s.zip

Calidad de aire https://www.kaggle.com/datasets/elianaj/mexico-air-quality-dataset

ruta<-"C:/Users/cguer/Documents/Claudia/Midropbox/Investigacion y escritos/Marcogeoestadisticonacional/mg_2023_integrado/conjunto_de_datos/"
infile <- "00ent.shp"
nomarchi<-paste0(ruta,infile)
datageom<- read_sf(nomarchi,options = "ENCODING=latin1")
datageom= st_transform(datageom, crs = 4326)

datageom <- datageom %>% 
   filter(CVE_ENT %in% c("21"))


ruta<-"C:/Users/cguer/Documents/Claudia/Midropbox/Investigacion y escritos/Marcogeoestadisticonacional/mg_2023_integrado/conjunto_de_datos/"
infile <- "00mun.shp"
nomarchi<-paste0(ruta,infile)
datageom<- read_sf(nomarchi,options = "ENCODING=latin1")
datageom= st_transform(datageom, crs = 4326)

datageom <- datageom %>% 
   filter(CVE_ENT %in% c("21"))%>% 
   filter(CVE_MUN %in% c("114"))

#datageo_geom <- st_geometry(datageom)
#polygon = st_cast(datageo_geom, "POLYGON") 

#plot(polygon[c(2,3,4)])
#datageom<-polygon[c(2,3,4)] 
#datageom<-polygon 

ruta<-"C:/Users/cguer/Documents/Claudia/Midropbox/Investigacion y escritos/karamanis/data/"
infile <- "stations_rsinaica.CSV"
nomarchi<-paste0(ruta,infile)

#*stations_hourly.CSV
#*stations_daily.CSV
estaciones<-read.csv(nomarchi)
estaciones <- estaciones[which(!is.na(estaciones$lat) & !is.na(estaciones$lon)),]
estaciones_sf <- st_as_sf(estaciones, crs = 4326,
  coords = c("lon", "lat"))

estaciones_sf<- estaciones_sf%>% 
   filter(network_code %in% c("PUE"))|> dplyr::filter(station_id<=164) 



ggplot() + geom_sf(data = datageom) + 
    geom_sf(data = estaciones_sf, color="red")

Información por hora

Se selecciona NO2, se filtra a los datos posteriores al 01 de enero de 2021, se eliminan las estaciones sin datos o con valores de cero Se grafica el promedio por estación

#(µg/m3)
library(hrbrthemes)
library(gganimate)

infile <- "calairepuebla.CSV"
nomarchi<-paste0(ruta,infile)
medidasestacioneshor<-read.csv(nomarchi)
fecha <- as.POSIXct(medidasestacioneshor$fecha, format = "%d/%m/%Y %H:%M")
PM10 <- c(medidasestacioneshor$pm10aguasanta)
PM2.5<- c(medidasestacioneshor$pm2.5aguasanta)
pest1<- data.frame(fecha = fecha, PM10 = as.numeric(PM10),PM2.5=as.numeric(PM2.5),station_id = 161)

fecha <- as.POSIXct(medidasestacioneshor$fecha, format = "%d/%m/%Y %H:%M")
PM10 <- c(medidasestacioneshor$pm10bine)
PM2.5 <- c(medidasestacioneshor$pm2.5bine)
pest2<- data.frame(fecha = fecha, PM10 = as.numeric(PM10),PM2.5=as.numeric(PM2.5),station_id = 163)

fecha <- as.POSIXct(medidasestacioneshor$fecha, format = "%d/%m/%Y %H:%M")
PM10 <- c(medidasestacioneshor$pm10ninfas)
PM2.5 <- c(medidasestacioneshor$pm2.5ninfas)

pest3<- data.frame(fecha = fecha, PM10 = as.numeric(PM10),PM2.5=as.numeric(PM2.5), station_id = 162)
fecha <- as.POSIXct(medidasestacioneshor$fecha, format = "%d/%m/%Y %H:%M")
PM10 <- c(medidasestacioneshor$pm10utp)
PM2.5 <- c(medidasestacioneshor$pm2.5utp)
pest4<- data.frame(fecha = fecha, PM10 = as.numeric(PM10),PM2.5=as.numeric(PM2.5), station_id = 164)


medidashor<-rbind(pest1,pest2,pest3,pest4)

medidashor2 <- medidashor[, c("fecha", "station_id", "PM2.5")]


library(dplyr)

dfs<-medidashor2  %>%
  mutate(date = lubridate::ymd_hms(fecha),
         date_dia = format(date, "%Y-%m")) %>%
    group_by(date_dia,station_id)%>%
  summarize(promedio = mean(PM2.5,na.rm=TRUE))
dfs<-dfs %>% filter(promedio>0)

#dfs<-dfs %>% filter(date_dia == '2020-01-01')
#dfs<-dfs %>% filter(date_dia >= '2024-01-01')


mediasPM2.5D<-inner_join(dfs,estaciones_sf)


mediasPM2.5D_sf <- st_as_sf(mediasPM2.5D, crs = st_crs(4326))
#mediasPM2.5D_sf$date2 = as.POSIXct(mediasPM2.5D_sf$date_dia,format = "%Y-%m-%d")

ggplot() + geom_sf(data = datageom) + 
    geom_sf(data = mediasPM2.5D_sf, aes(color = promedio))

map2.5<-mediasPM2.5D_sf |>
  ggplot() +
  geom_sf(data = datageom) + 
  geom_sf(
    aes(color = promedio),
    linewidth = 0.1) +
  scale_color_viridis_c(
    #labels = mean_PM2.5,
    # = c(0, 100),
    option = "C",
    guide = guide_colorsteps(show.limits = FALSE)
  )+
  theme_ipsum_inter(grid = TRUE) +
  theme(
    legend.position = "top",
    plot.title = element_text(
      face = "bold",
      hjust = 0.5
    ),
    axis.text.x = element_blank(),
    axis.text.y = element_blank(),
    legend.key.width = unit(2, "cm"),
    legend.text = element_text(
      size = 12
    ))

map2.5 + facet_wrap(vars(date_dia))

Animación

map2.5<-mediasPM2.5D_sf |>
  ggplot() +
  geom_sf(data = datageom) + 
  geom_sf(
    aes(color = promedio),
    linewidth = 0.1) +
  scale_color_viridis_c(
    #labels = mean_PM2.5,
    # = c(0, 100),
    option = "C",
    guide = guide_colorsteps(show.limits = FALSE)
  )+
  theme_ipsum_inter(grid = TRUE) +
  theme(
    legend.position = "top",
    plot.title = element_text(
      face = "bold",
      hjust = 0.5
    ),
    axis.text.x = element_blank(),
    axis.text.y = element_blank(),
    legend.key.width = unit(2, "cm"),
    legend.text = element_text(
      size = 12
    ))+
  transition_manual(date_dia) +
  labs(title = "Puebla PM2.5 (µg/m3) promedio mensual  {current_frame}")

#animate(map2.5,renderer = gifski_renderer(),nframes = 500,fps=2, width = 800, height = 800)

#anim_save('pm25puebla.gif')

if (knitr:::is_latex_output()) {
  knitr::asis_output('\\url{....}')
} else {
  knitr::include_graphics("pm25puebla.gif")
}