Hoy vamos a hacer algunos ejercicios para trabajar con datos bioclimáticos de algunas bases de datos. Se instalarán algunos paquetes necesarios de CRAN y GitHub.

# install.packages(c("raster","dismo","devtools","envirem"))
library(raster)
## Loading required package: sp
library(dismo)
library(envirem)
## Loading required package: palinsol
## Loading required package: gsl
## 
## Please see the vignette at https://envirem.github.io/ENVIREM_tutorial.html for a detailed walk-through of this package.

Descarga de datos directamente de Worldclim

# worldclim (10)
# 5 and 2.5 
clim <- getData("worldclim",var="bio",res=10)
clim
## class      : RasterStack 
## dimensions : 900, 2160, 1944000, 19  (nrow, ncol, ncell, nlayers)
## resolution : 0.1666667, 0.1666667  (x, y)
## extent     : -180, 180, -60, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## names      :  bio1,  bio2,  bio3,  bio4,  bio5,  bio6,  bio7,  bio8,  bio9, bio10, bio11, bio12, bio13, bio14, bio15, ... 
## min values :  -269,     9,     8,    72,   -59,  -547,    53,  -251,  -450,   -97,  -488,     0,     0,     0,     0, ... 
## max values :   314,   211,    95, 22673,   489,   258,   725,   375,   364,   380,   289,  9916,  2088,   652,   261, ...
names(clim)
##  [1] "bio1"  "bio2"  "bio3"  "bio4"  "bio5"  "bio6"  "bio7"  "bio8"  "bio9" 
## [10] "bio10" "bio11" "bio12" "bio13" "bio14" "bio15" "bio16" "bio17" "bio18"
## [19] "bio19"

Extraigo solo la temperatura promedio anual y la precipitación anual

tpa <- clim[[1]]
pa <- clim[[12]]

Graficar

par(mfrow=c(2,1))
plot(tpa, main="Temperatura promedio anual") 
plot(pa, main="Precipitación anual")

tpa
## class      : RasterLayer 
## dimensions : 900, 2160, 1944000  (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667  (x, y)
## extent     : -180, 180, -60, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## source     : /home/julian/Documents/Taller_R/Bioclimaticos/wc10/bio1.bil 
## names      : bio1 
## values     : -269, 314  (min, max)
pa
## class      : RasterLayer 
## dimensions : 900, 2160, 1944000  (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667  (x, y)
## extent     : -180, 180, -60, 90  (xmin, xmax, ymin, ymax)
## crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## source     : /home/julian/Documents/Taller_R/Bioclimaticos/wc10/bio12.bil 
## names      : bio12 
## values     : 0, 9916  (min, max)

Estadísticas a través de celdas (media, máximo, mínimo y desviación estándar)

mean.tpa <- cellStats(tpa, stat='mean')
mean.pa <- cellStats(pa, stat='mean')

max.tpa <- cellStats(tpa, stat='max')
max.pa <- cellStats(pa, stat='max')

min.tpa <- cellStats(tpa, stat='min')
min.pa <- cellStats(pa, stat='min')

sd.tpa <- cellStats(tpa, stat='sd')
sd.pa <- cellStats(pa, stat='sd')

Vamos a subir los datos de temperaturas y precipitación que bajaron con Oscar. Borrar primero el raster del 2019

# primero temperatura mínima mensual
tmin <- stack(list.files(path="/home/julian/Documents/Taller_R/WC/WC_Tmin", pattern='*.tif', full.names=TRUE))

# luego temperatura máxima mensual
tmax <- stack(list.files(path="/home/julian/Documents/Taller_R/WC/WC_Tmax", pattern='*.tif', full.names=TRUE))

# precipitación mensual
prec <- stack(list.files(path="/home/julian/Documents/Taller_R/WC/WC_Prec", pattern='*.tif', full.names=TRUE))

Calcular los promedios mensuales para todo el periodo (1961-2018). 696 rasters son 58 años de datos.

t <- 696/12

# Separar los años
t1 <- split(seq_len(696), rep(1:58, each=12))

# Agrupar cada mes de todos los años
t2 <- split(seq_len(696), rep(1:12, length.out=696))

t.tmin <- lapply(1:12, function(x) tmin[[t2[[x]]]])
t.tmax <- lapply(1:12, function(x) tmax[[t2[[x]]]])
t.prec <- lapply(1:12, function(x) prec[[t2[[x]]]])

# calcular el promedio de cada mes de todos los años
t.tmin.mean <- lapply(1:12, function(x) mean(t.tmin[[x]]))
t.tmax.mean <- lapply(1:12, function(x) mean(t.tmax[[x]]))
t.prec.mean <- lapply(1:12, function(x) mean(t.prec[[x]]))

# convierto la lista de rasters en un stack
t.tmin.mean2 <- stack(t.tmin.mean)
t.tmax.mean2 <- stack(t.tmin.mean)
t.prec.mean2 <- stack(t.prec.mean)

names(t.tmin.mean2) <- c("ene","feb","mar","abr","may","jun","jul","ago","sep","oct","nov","dic")
names(t.tmax.mean2) <- c("ene","feb","mar","abr","may","jun","jul","ago","sep","oct","nov","dic")
names(t.prec.mean2) <- c("ene","feb","mar","abr","may","jun","jul","ago","sep","oct","nov","dic")

# indexado manual
# tmin.ene <- tmin[[c(1,13,25,37,49,61,73,85,97,109,121,133,...,)]]

Ahora vamos a generar los bioclimáticos (19) para los 58 años

library(dismo)
bios1961_2018 <- biovars(t.prec.mean2, t.tmin.mean2, t.tmax.mean2)

# generar el mapa para la bio1
plot(bios1961_2018[[1]])

# Generar el mapa para la bio12
plot(bios1961_2018[[12]])

Calcular las diferencias entre los bioclimáticos generados y los descargados en la primera parte.

# remuestrear al pixel más grande
clim2 <- raster::resample(clim, bios1961_2018[[1]], method="bilinear")

# diferencia
dif.tpa <- bios1961_2018[[1]] - clim2[[1]]/10
plot(dif.tpa)

# histograma

dif.tpa2 <- rasterToPoints(dif.tpa)
hist(dif.tpa2[,3])

Eso es todo por ahora.