The question of how the distribution of species is influenced by environmental conditions is an important topic of study in ecology. Because of that, we will show in this tutorial how to implement a simple niche suitability model in metaRange. For this, we will create a landscape with two environmental variables (temperature and precipitation) and then add two similar species to it that only differ in their environmental preferences. At the end we can run the simulation and compare how this difference affects the distribution of the species.

Setup

First we load the packages and create the example landscape.

Code

library(metaRange)
library(terra)
set_verbosity(2)

# find the example raster file
raster_file <- system.file("ex/elev.tif", package = "terra")

# load it
r <- rast(raster_file)

# adjust the values
temperature <- scale(r, center = FALSE, scale = TRUE) * 10 + 273.15
precipitation <- r * 2

Now we can (again) turn the raster into an SDS with one layer per time step.

Code

temperature <- rep(temperature, 10)
precipitation <- rep(precipitation, 10)
landscape <- sds(temperature, precipitation)
names(landscape) <- c("temperature", "precipitation")

Code

terra::plot(
    landscape$temperature[[1]],
    col = hcl.colors(100, "RdYlBu", rev = TRUE),
    main = "Temperature [K]"
)

Figure 1: The temperature of the example landscape. Only the first layer of 10 identical ones is shown.

Code

terra::plot(
    landscape$precipitation[[1]],
    col = hcl.colors(100, "Earth"),
    main = "Precipitation [mm]"
)

Figure 2: The precipitation of the example landscape. Only the first layer of 10 identical ones is shown.

Creating the simulation and species

Code

sim <- create_simulation(
    source_environment = landscape,
    ID = "example_simulation",
    seed = 1
)
#> number of time steps: 10
#> time step layer mapping: 1, 2, 3, 4, 5, 6, 7, 8, 9, 10
#> added environment
#> class       : SpatRasterDataset 
#> subdatasets : 2 
#> dimensions  : 90, 95 (nrow, ncol)
#> nlyr        : 10, 10 
#> resolution  : 0.008333333, 0.008333333  (x, y)
#> extent      : 5.741667, 6.533333, 49.44167, 50.19167  (xmin, xmax, ymin, ymax)
#> coord. ref. : lon/lat WGS 84 (EPSG:4326) 
#> source(s)   : memory 
#> names       : temperature, precipitation
#> 
#> created simulation: example_simulation
sim$add_species(name = "species_1")
#> adding species
#> name: species_1
sim$add_species(name = "species_2")
#> adding species
#> name: species_2

Adding traits to species

To start, we add the same traits as in the previous tutorials and also add a trait called climate_suitability, where we will store the information about how suitable the environment in each cell is for the population that lives there.

Code

sim$add_traits(
    species = c("species_1", "species_2"),
    population_level = TRUE,
    abundance = 500,
    climate_suitability = 1,
    reproduction_rate = 0.3,
    carrying_capacity = 1000
)
#> adding traits:
#> [1] "abundance"           "climate_suitability" "reproduction_rate"  
#> [4] "carrying_capacity"
#> 
#> to species:
#> [1] "species_1" "species_2"
#>

Contrary to the above, some traits may not require to be stored at the population level. In this example, this could be the case for the environmental preferences of a species. If we assume that they are a property of the species as a whole (i.e. the same for all populations), we can set the parameter population_level to FALSE and the traits will be added as they are, without further processing.

As mentioned in the introduction paragraph, we will give both species different environmental preferences for the two environmental variables in the simulation environment (temperature & precipitation).

Note that the names of the traits are arbitrary and can be chosen by the user and that there is no predetermined connection between e.g. “min_temperature” and the temperature variable in the environment. To establish these connections, the user needs to add processes to the species that access the correct traits and use them in a sensible way (This is why meaningful trait names are important).

Code

sim$add_traits(
    species = "species_1",
    population_level = FALSE,
    max_temperature = 300,     # Kelvin
    optimal_temperature = 288, # Kelvin
    min_temperature = 280,     # Kelvin
    max_precipitation = 1000,    # mm
    optimal_precipitation = 700, # mm
    min_precipitation = 200      # mm
)
#> adding traits:
#> [1] "max_temperature"       "optimal_temperature"   "min_temperature"      
#> [4] "max_precipitation"     "optimal_precipitation" "min_precipitation"
#> 
#> to species:
#> [1] "species_1"
#> 
sim$add_traits(
    species = "species_2",
    population_level = FALSE,
    max_temperature = 290,
    optimal_temperature = 285,
    min_temperature = 270,
    max_precipitation = 1000,
    optimal_precipitation = 500,
    min_precipitation = 0
)
#> adding traits:
#> [1] "max_temperature"       "optimal_temperature"   "min_temperature"      
#> [4] "max_precipitation"     "optimal_precipitation" "min_precipitation"
#> to species:
#> [1] "species_2"
#>

Adding processes

Calculate the suitability

To calculate the suitability, we use the metaRange function calculate_suitability() that was adapted from a equation published by Yin et al. in 1995 (Ref. 1) and simplified by Yan and Hunt in 1999 (eq:4 in Ref. 2). The equation takes the three cardinal values of an environmental niche (minimum tolerable value, optimal vale and maximum tolerable value) and constructs a suitability curve based on a beta distribution.

Code

min_value <- 270
opt_value <- 285
max_value <- 290
x <- seq(min_value, max_value, length.out = 100)
y <- calculate_suitability(max_value, opt_value, min_value, x)
plot(x, y, type = "l", xlab = "Temperature [K]", ylab = "Suitability")

Figure 3: Example suitability curve for the temperature niche of species 2.

In the following code we add a process to both species that calculates the environmental suitability for precipitation and temperature and then multiplies the values to create a joint suitability over the two niches. Note that one could also define a custom function to calculate the suitability, if this built-in function does not adequately describe the ecology of the target species.

Suitability

Code

sim$add_process(
    species = c("species_1", "species_2"),
    process_name = "calculate_suitability",
    process_fun = function() {
        self$traits$climate_suitability <-
            calculate_suitability(
                self$traits$max_temperature,
                self$traits$optimal_temperature,
                self$traits$min_temperature,
                self$sim$environment$current$temperature
            ) *
            calculate_suitability(
                self$traits$max_precipitation,
                self$traits$optimal_precipitation,
                self$traits$min_precipitation,
                self$sim$environment$current$precipitation
            )
    },
    execution_priority = 1
)
#> adding process: calculate_suitability
#> to species:
#> [1] "species_1" "species_2"
#>

Reproduction

As in the previous tutorials, we use a Ricker reproduction model to calculate the new abundance of the species, but this time we let both the carrying capacity and the reproduction rate depend on the niche suitability of the environment.

Code

sim$add_process(
    species = c("species_1", "species_2"),
    process_name = "reproduction",
    process_fun = function() {
        self$traits$abundance <-
            ricker_reproduction_model(
                self$traits$abundance,
                self$traits$reproduction_rate * self$traits$climate_suitability,
                self$traits$carrying_capacity * self$traits$climate_suitability
            )
    },
    execution_priority = 2
)
#> adding process: reproduction
#> to species:
#> [1] "species_1" "species_2"
#>

Results

Now, we can execute the simulation and compare the results.

Code

set_verbosity(1)
sim$begin()
#> Starting simualtion.
#> start of time step: 1
#>  10 % done | 0.017 secs remaining (estimate)
#> start of time step: 2
#>  20 % done | 0.13 secs remaining (estimate)
#> start of time step: 3
#>  30 % done | 0.079 secs remaining (estimate)
#> start of time step: 4
#>  40 % done | 0.069 secs remaining (estimate)
#> start of time step: 5
#>  50 % done | 0.058 secs remaining (estimate)
#> start of time step: 6
#>  60 % done | 0.045 secs remaining (estimate)
#> start of time step: 7
#>  70 % done | 0.033 secs remaining (estimate)
#> start of time step: 8
#>  80 % done | 0.024 secs remaining (estimate)
#> start of time step: 9
#>  90 % done | 0.017 secs remaining (estimate)
#> start of time step: 10
#> 100 % done | 0 secs remaining (estimate)
#> 
#> Simulation: 'example_simulation' finished
#> Exiting the Simulation
#> Runtime: 0.12 secs

Code

# define a nice color palette
plot_cols <- hcl.colors(100, "Purple-Yellow", rev = TRUE)
plot(
    sim,
    obj = "species_1",
    name = "abundance",
    main = "Species 1: abundance",
    col = plot_cols
)

Figure 6: The resulting abundance distribution of species 1 after 10 simulation time steps.

Code

plot(
    sim$species_2,
    trait = "abundance",
    main = "Species 2: abundance",
    col = plot_cols
)

Figure 7: The resulting abundance distribution of species 2 after 10 simulation time steps.

03: Niche suitability

Fallert, S. and Cabral, J.S.

Setup

Creating the simulation and species

Adding traits to species

Adding processes

Calculate the suitability

Suitability

Reproduction

Results

References