Hacer gráfica en R de precipitación a partir de archivos raster por departamento

En este tutorial explicare brevemente como realizar un gráfico tipo ojiva, en el cual se presente en el eje x los meses del año & en el eje y la precipitación acumulada para cada mes, a tres distintos departamentos de Colombia, estos son: Antioquía, Cundinamarca y Valle del Cauca. Los datos de entrada fueron tomados de la base de datos Worldclim (Hijmans et al., 2005) y extraídos para Colombia a partir del shape de límite departamental del SIGOT.

Los datos utilizados en este tutorial están disponibles en este enlace, y se describen a continuación:

1) 12 raster de precipitación para todo Colombia a resolución de 30 arcos de segundo (aproximadamente 1 km2)

2) Shapefile de los 32 departamentos de Colombia.

A continuación se describen los pasos, asumiendo que ya se tiene instalado R y RStudio en su computador.

1. Instalamos las siguientes librerías (importante tener conexión a Internet)

install.packages("raster")
install.packages("rgdal")
install.packages("ggplot2")
install.packages("dplyr")
install.packages("reshape2")

require("raster")
require("rgdal")
require("ggplot2")
require("dplyr")
require("reshape2")

2.  Llamamos nuestros datos en RStudio

myproj <- CRS("+proj=longlat +datum=WGS84") #proyección a utilizar
path   <- "E:/Blogs/_blogger/_graphClimate" #dirección donde está el directorio raíz de los datos
col_shp <- shapefile(paste0(path, "/_shp/LimiteDptal.shp")) # limite departamental
files  <- list.files(paste0(path, "/_climate/_colombia/_asc"), full.names = T, pattern = ".asc$") #archivos raster

 3. Escogemos en nuevos objetos los departamentos de interés

valleCauca_shp <- col_shp[col_shp$NOMBRE_DPT %in% "VALLE DEL CAUCA",]
cundinamarca_shp <- col_shp[col_shp$NOMBRE_DPT %in% "CUNDINAMARCA",]
antioquia_shp <- col_shp[col_shp$NOMBRE_DPT %in% "ANTIOQUIA",]

4. Archivos raster

files_prec <- grep("prec", files, value = T) %>%
              mixedsort() #seleccion de precipitación y ordenar datos
layers_prec <- lapply(files_prec, raster) #conversión de objetos a raster

5. Corte de archivos raster por departamento, se crean listas vacías donde se guardaran los raster para cada departamento.

prec_antioquia    <- list() #creación de lista vacía
prec_cundinamarca <- list() #creación de lista vacia
prec_valle        <- list() #creación de lista vacia

for(i in 1:length(layers_prec)){

  prec_antioquia[[i]] <- crop(layers_prec[[i]],    antioquia_shp)
  prec_antioquia[[i]] <- mask(prec_antioquia[[i]], antioquia_shp)

  prec_cundinamarca[[i]] <- crop(layers_prec[[i]], cundinamarca_shp)
  prec_cundinamarca[[i]] <- mask(prec_cundinamarca[[i]], cundinamarca_shp)

  prec_valle[[i]] <- crop(layers_prec[[i]], valleCauca_shp)
  prec_valle[[i]] <- mask(prec_valle[[i]],  valleCauca_shp)

}

6. Calculo de promedios de precipitación para cada mes en nuestros departamentos objeto de estudio. Se crean listas vacías donde se guardaran los datos de promedios para cada raster por departamento.

prec_mean_antioquia <- list(); prec_mean_cund <- list(); prec_mean_valle <- list()
prec_df_antioquia   <- NA; prec_df_cund <- NA; prec_df_valle <- NA

for(i in 1:length(prec_antioquia)){
  
  prec_mean_antioquia[[i]] <- mean(prec_antioquia[[i]][], na.rm = T)
  prec_df_antioquia <- t(as.data.frame((rbind(prec_mean_antioquia))))
  
  prec_mean_cund[[i]] <- mean(prec_cundinamarca[[i]][], na.rm = T)
  prec_df_cund      <- t(as.data.frame((rbind(prec_mean_cund))))
  
  prec_mean_valle[[i]]<- mean(prec_valle[[i]][], na.rm = T)
  prec_df_valle     <- t(as.data.frame((rbind(prec_mean_valle))))
  
}

7. Orden de datos para ggplot

df_prec <- data.frame(prec_df_antioquia, prec_df_cund, prec_df_valle)
df_prec <- df_prec %>%
      mutate(Month_numeric = 1:nrow(df_prec)) %>%
      mutate(prec_mean_antioquia = as.numeric(prec_mean_antioquia), 

      prec_mean_cund = as.numeric(prec_mean_cund),

      prec_mean_valle  = as.numeric(prec_mean_valle))

prec_df_rename <- rename(df_prec, Antioquia = prec_mean_antioquia, 
                 Cundinamarca = prec_mean_cund, Valle = prec_mean_valle)

prec_df_melt <- melt(prec_df_rename, id = 'Month_numeric')
prec_df_melt <- rename(prec_df_melt, Departamento = variable)

8. Creación de gráfico base

gg  <- ggplot(prec_df_melt, aes(Month_numeric, value,  fill = Departamento, colour = Departamento))+ geom_line() + xlab("Mes") + ylab("Precipitación (mm)") + ggtitle("Precipitación para cada mes") +
scale_x_continuous(expand = c(0,0), breaks = c(1:12), labels = c("Enero", "Febrero", "Marzo", "Abril", "Mayo", "Junio", "Julio", "Agosto", "Septiembre", "Octubre", "Noviembre", "Diciembre")) +        theme(axis.text.x = element_text(angle = 45, hjust = 1)) + 

theme(legend.position = "bottom")

9. Creación de gráfico base

gg <- gg + scale_colour_manual(labels = c("Antioquia", "Cundinamarca", "Valle del Cauca"), values = c("red", "green", "blue"))

10. Guardado de gráfico en buena resolución

ggsave(plot = gg, paste0(path, "/_png/precYear_dptos.png"), width = 9, height = 8, units = "in")
gg

El gráfico debe quedar así: