library(dynamicSDM)In this tutorial, we will learn about functions within the dynamicSDM package for projecting intra- and inter-annual geographical distributions and abundances at high spatiotemporal resolution. This will involve extracting rasters for each explanatory variable across southern Africa at various intervals. Then we can project our dynamic SDMs onto these data to visualise species distribution patterns through space and time.
The dynamicSDM functions for extracting dynamic rasters require Google Earth Engine and Google Drive to be initialised. Fill in the code below with your Google account email, and run the code to check that rgee and googledrive have been correctly installed and authorised.
library(rgee)
rgee::ee_check()
library(googledrive)
googledrive::drive_user()
#user.email<-"your_google_email_here"Note: You will need internet connection for this tutorial. Raster extraction may take some time depending on your internet connection strength.
We will be extracting rasters, generating projection covariate data frames, exporting projection rasters and saving GIFs. Let’s first make separate folders in our project directory to save this various outputs too. This is especially important as this suite of dynamicSDM functions rely on unique naming of files in the local directory or drive folder.
variablenames <- c("eight_sum_prec", "year_sum_prec", "grass_crop_percentage")
project_directory <- file.path(tempdir(), "dynamicSDM_vignette")
# project_directory<-"your_path_here"
dir.create(project_directory, showWarnings = FALSE)
projection_directories <- file.path(project_directory, "projection")
dir.create(projection_directories, showWarnings = FALSE)
projectionrasters <- file.path(project_directory, "projectionrasters")
dir.create(projectionrasters, showWarnings = FALSE)
projection_covariates <- file.path(project_directory, "projectioncovariates")
dir.create(projection_covariates, showWarnings = FALSE)
projection_GIF <- file.path(project_directory, "projection_GIF")
dir.create(projection_GIF, showWarnings = FALSE)dynamic_proj_dates() generates all dates in a given period at the interval type and size specified. This can be used to create a vector of dates that you want to project a species distribution at.
# 5 day intervals
dynamic_proj_dates(startdate = "2018-01-01",
enddate = "2018-12-01",
interval = 5,
interval.level = "day")
#> [1] "2018-01-01" "2018-01-06" "2018-01-11" "2018-01-16" "2018-01-21"
#> [6] "2018-01-26" "2018-01-31" "2018-02-05" "2018-02-10" "2018-02-15"
#> [11] "2018-02-20" "2018-02-25" "2018-03-02" "2018-03-07" "2018-03-12"
#> [16] "2018-03-17" "2018-03-22" "2018-03-27" "2018-04-01" "2018-04-06"
#> [21] "2018-04-11" "2018-04-16" "2018-04-21" "2018-04-26" "2018-05-01"
#> [26] "2018-05-06" "2018-05-11" "2018-05-16" "2018-05-21" "2018-05-26"
#> [31] "2018-05-31" "2018-06-05" "2018-06-10" "2018-06-15" "2018-06-20"
#> [36] "2018-06-25" "2018-06-30" "2018-07-05" "2018-07-10" "2018-07-15"
#> [41] "2018-07-20" "2018-07-25" "2018-07-30" "2018-08-04" "2018-08-09"
#> [46] "2018-08-14" "2018-08-19" "2018-08-24" "2018-08-29" "2018-09-03"
#> [51] "2018-09-08" "2018-09-13" "2018-09-18" "2018-09-23" "2018-09-28"
#> [56] "2018-10-03" "2018-10-08" "2018-10-13" "2018-10-18" "2018-10-23"
#> [61] "2018-10-28" "2018-11-02" "2018-11-07" "2018-11-12" "2018-11-17"
#> [66] "2018-11-22" "2018-11-27"
# 2 week intervals
dynamic_proj_dates(startdate = "2018-01-01",
enddate = "2018-12-01",
interval = 2,
interval.level = "week")
#> [1] "2018-01-01" "2018-01-15" "2018-01-29" "2018-02-12" "2018-02-26"
#> [6] "2018-03-12" "2018-03-26" "2018-04-09" "2018-04-23" "2018-05-07"
#> [11] "2018-05-21" "2018-06-04" "2018-06-18" "2018-07-02" "2018-07-16"
#> [16] "2018-07-30" "2018-08-13" "2018-08-27" "2018-09-10" "2018-09-24"
#> [21] "2018-10-08" "2018-10-22" "2018-11-05" "2018-11-19"For this tutorial, we will use the following projection dates to explore the intra-annual distribution patterns of red-billed quelea in southern Africa.
projectiondates <- dynamic_proj_dates(startdate = "2018-01-01",
enddate = "2018-12-01",
interval = 3,
interval.level = "month")extract_dynamic_raster() and extract_buffered_raster() are sister functions to those used in stage 2 to extract explanatory variables. They require largely the same arguments input, but with arguments to specify the temporal (projection dates) and spatial extents for the rasters.
Run the following code to extract rasters for each explanatory variable and projection date.
data("sample_extent_data")
extract_dynamic_raster(dates = projectiondates,
datasetname = "UCSB-CHG/CHIRPS/DAILY",
bandname = "precipitation",
user.email = user.email,
spatial.res.metres = 5566,
GEE.math.fun = "sum",
temporal.direction = "prior",
temporal.res = 56,
spatial.ext = sample_extent_data,
varname = variablenames[1],
save.directory = projectionrasters)
extract_dynamic_raster(dates = projectiondates,
datasetname = "UCSB-CHG/CHIRPS/DAILY",
bandname = "precipitation",
user.email = user.email,
spatial.res.metres = 5566,
GEE.math.fun = "sum",
temporal.direction = "prior",
temporal.res = 364,
spatial.ext = sample_extent_data,
varname = variablenames[2],
save.directory = projectionrasters)
matrix<-dynamicSDM::get_moving_window(radial.distance = 10000,
spatial.res.degrees = 0.05,
spatial.ext = sample_extent_data)
extract_buffered_raster(dates=projectiondates,
datasetname = "MODIS/006/MCD12Q1",
bandname = "LC_Type5",
spatial.res.metres = 500,
GEE.math.fun = "sum",
moving.window.matrix = matrix,
user.email = user.email,
categories = c(6,7),
agg.factor = 12,
spatial.ext = sample_extent_data,
varname = variablenames[3],
save.directory = projectionrasters)Now we need to combine the rasters for each variable into single projection covariate data frames for each projection date.
dynamic_proj_covariates(dates = projectiondates,
varnames = variablenames,
local.directory = projectionrasters,
spatial.ext = sample_extent_data,
spatial.mask = sample_extent_data,
spatial.res.degrees = 0.05,
resample.method = c("bilinear","bilinear","ngb"),
cov.file.type = "csv",
prj="+proj=longlat +datum=WGS84",
save.directory = projection_covariates)To project distribution onto these covariates, our fitted boosted regression tree models from stage 3 need to be imported and input into dynamic_proj(). You can read these in from your local directory, or simply rerun the code here.
We choose to project distribution suitability (projection.method = proportional) but this function can also generate binary projections if given threshold (sdm.thresh) or abundance if given fitted species abundance models (sam.mod).
#sample_brt_models<- readRDS(paste0(project_directory, "/fitted_quelea_SDMs.rds"))
data("sample_explan_data")
sample_explan_data$weights <- (1 - sample_explan_data$REL_SAMP_EFFORT)
sample_brt_models <- brt_fit(sample_explan_data,
response.col = "presence.absence",
varnames = variablenames,
block.col = "BLOCK.CATS",
weights.col = "weights",
distribution = "bernoulli")
dynamic_proj(dates = projectiondates,
projection.method = c("proportional"),
local.directory = projection_covariates,
cov.file.type = "csv",
sdm.mod = sample_brt_models,
spatial.mask = sample_extent_data,
save.directory = projection_directories)Let’s read in and plot our projections!
terra::plot(terra::rast(list.files(projection_directories)[1]))
terra::plot(terra::rast(list.files(projection_directories)[2]))
terra::plot(terra::rast(list.files(projection_directories)[3]))
terra::plot(terra::rast(list.files(projection_directories)[4]))To better visualise our dynamic distribution projections through time, dynamic_proj_GIF() can generate animated GIFs from the projection rasters.
We will set our own colour scheme and add borders for the southern African countries.
cols1<- c("#F0F0F0","#40863A","#FBF357","#ED8E07","#cc0000","#660000")
border.countries<- c('South Africa', 'Botswana','Lesotho', 'Swaziland','Mozambique','Namibia'
,'Zimbabwe','Angola','Zambia','Malawi')
dynamic_proj_GIF(
dates = projectiondates,
projection.type = "proportional",
local.directory = projection_directories,
save.directory = projection_GIF,
width = 7,
height = 5,
colour.palette.custom = cols1,
borders = TRUE,
border.regions = border.countries,
border.colour = "grey50",
legend.max = 1,
legend.min = 0,
legend.name = "Distribution\n suitability",
file.name = "RBQ_proportional_GIF")This concludes the dynamicSDM tutorial. We hope that you have been able to utilise the package to model the dynamic distributions of our case study species, the red-billed quelea, and can now apply these functions to your own study species!
We have demonstrated a variety of options across our functions, but even more are available with the diverse and flexible arguments. Please see our package manual for full details of each dynamicSDM function manual.
Also, please visit the GitHub discussion page for dynamicSDM if you have any questions or are facing errors with any of our functions github. We also have a quick guide to common bugs that are thrown up when installing and running dynamicSDM.