Relleno de datos de lluvia mensual de una estación con CHIRPS

Anteriormente, generamos una función para descargar los raster de lluvia estimada a partir de satélite (CHIRPS V 2.0), ahora lo que vamos hacer es extraer los datos para un punto en especifico o la ubicación de una estación climática, el periodo de tiempo de la serie de datos chirps debe ser la misma que la serie de datos de la estación (no importa que existan datos faltantes, solo especificar el datos faltante como -99.9) para que sean comparables.

Luego se creara un modelo de ajuste entre los datos estimados y los datos reales, esto con el fin de que los datos a rellenar sean acordes a la tendencia de los datos observados o pluviometricos.

La importancia de completar series de datos radica darle continuidad a la serie y que esto me permita realizar mejores análisis, crear modelos hidrológicos estables, balances hídricos, etc.

####################################################################

library(raster)

# Ruta del directorio de trabajo
setwd("D:/Documents/chirps/")

# Directorio de los datos de chirp V-2.0
ruta_chirps = "D:/Documents/chirps/raster/"

# Periodo de datos a analizar
fecha_inicial = c("1981/1/1")
fecha_final = c("2018/5/1")

# Coordenada geografica de las estacion (longitud, Latitud)
coor.x = c(-89.8106)
coor.y = c(15.6083)

# Convertir las coordenadas de las estaciones a puntos
ubicacion = data.frame(lon=coor.x, lat=coor.y)
coordinates(ubicacion) = c("lon", "lat")
proj4string(ubicacion) = CRS("+init=epsg:4326") # WGS 84

# Cargar los datos raster
lista_chirps = paste0(ruta_chirps,"chirps-v2.0.", format(seq(as.Date(fecha_inicial), as.Date(fecha_final), "month"), "%Y.%m"), ".tif")
ras_chirps = stack(lista_chirps)

# Extraer datos de lluvia de la estacion
chirps = as.vector((extract(ras_chirps, ubicacion)))

# Cargar datos pluviometricos (reales u observados) de la estacion
estacion = scan("estaciones/Cahabon_1_1981_5_2018.txt", sep = "\t", na.strings = "-99.9")
serie = estacion


# Modelo de ajuste de datos chirps a pluviometicos
model = lm(estacion~chirps)

# estadisticos del modelo de ajuste lineal
summary(model)

# coeficientes del modelo de ajuste lineal
alfa = coef(model)[[2]]
beta = coef(model)[[1]]

# Grafica
plot(estacion~chirps)
abline(b=1,a=1, col="red")

# posiciones de datos faltantes pluviometricos
faltantes = which(is.na(serie))

# Datos chirps ajustados para relleno
relleno = round(alfa*chirps[faltantes]+beta,1)

# Completacion de la serie real
serie[faltantes] = relleno

# Grafico de serie completa (negro) y serie incompleta (rojo)
plot(serie,type="l", xlab="Lluvia [mm]", main="Comparacion de serie completa e incompleta de lluvia")
lines(estacion, col="red")

# Crear tabla de fechas y serie completa
fecha = format(seq(as.Date(fecha_inicial), as.Date(fecha_final), "month"), "%Y.%m")
serie.completa = data.frame(fecha,serie)

write.csv(serie.completa, "Cahabon_rellenada.csv", row.names = FALSE)

#############################################################

Descargar
Archivo_R