# 03: Niche suitability

#### Fallert, S. and Cabral, J.S.

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.

# References

1. Yin, X., Kropff, M.J., McLaren, G., Visperas, R.M., (1995) A nonlinear model for crop development as a function of temperature, Agricultural and Forest Meteorology, Volume 77, Issues 1-2, Pages 1–16, doi:10.1016/0168-1923(95)02236-Q

2. Yan, W., Hunt, L.A. (1999) An Equation for Modelling the Temperature Response of Plants using only the Cardinal Temperatures, Annals of Botany, Volume 84, Issue 5, Pages 607–614, ISSN 0305-7364, doi:10.1006/anbo.1999.0955