Figure 1: Example of some of the different environmental factors,species traits and processes that could be included in a simulation.
metaRange is a framework to build avariety of different process based species distribution models that caninclude a (basically) arbitrary number of environmental factors,processes, species and species interactions. The common denominator forall models built with metaRange is that they are raster and populationbased.
Each metaRangesimulation object contains oneenvironment that holds and manages all the environmentalfactors that may influence the simulation (as raster data) and one ormorespecies that are simulated in the environment.
Each species object has two main characteristics:traitsandprocesses. Speciestraits can be (somewhatarbitrary) pieces of data that describe and store information about thespecies whileprocesses are functions that describe how thespecies interacts with itself, time, climate and other species. Duringeach time step of the simulation, the processes of the species areexecuted in a user defined order and can access and modify the speciestraits.
While the models built with metaRange can be quite variable in theirstructure, they are all based on population dynamics. This means that inmost cases atrait will not be a single number but a matrixthat has the same size as the raster data in the environment, where eachvalue represents the trait value for one population in the correspondinggrid cell of the environment. Based on this, mostprocesseswill describe population (and meta population) dynamics and notindividual based mechanisms. Figure 1 shows an example of some of thedifferent environmental factors, species traits and processes that couldbe included in a simulation.
Figure 1: Example of some of the different environmental factors,species traits and processes that could be included in a simulation.
process priority queue to determine in whichorder the processes are executed within one time step.Figure 2: Overview of the different components of a simulation and howthey interact with each other. Note that the number of species as wellas the number of traits and processes per species is not limited, butonly a selection is shown for simplicity.
Following is a simple example of how to set up a simulation withmetaRange, in which we only use a single species and oneenvironmental factor (habitat quality). At the end of this introductionwe will see how the abundance of the species changes in relation to thequality of the habitat each population occupies.
To start, we need to load the packages.
The first step when setting up a simulation is the loading of theenvironment in which the simulation will take place. This can either bereal world data or “theoretical” / generated data and may include forexample different climate variables, land cover or elevation.
The simulation expects this data as anSpatRasterDataset(SDS) which is a collection of different raster files thatall share the same extent and resolution. Each sub-dataset in this SDSrepresents one environmental variable and each layer represents one timestep of the simulation. In other words, metaRange does not simulate theenvironmental conditions itself, but expects the user to provide theenvironmental data for each time step.
To create such a dataset one can use the functionterra::sds(). One important note: Since each layerrepresents the condition in one time step, all the raster files that gointo theSDS need to have the same number of layers(i.e. the desired number of time steps the simulation should have).After theSDS is created, the individual sub-datasetsshould be named, since this is how the simulation will refer tothem.
To simplify this introduction, we use an example landscape consistingonly of habitat quality data, with 10 time steps (layers) that are allthe same (i.e. no environmental change). Luckily theterrapackage has a built-in demo that we can use for this purpose.
Now we can turn this raster with one layer into anSDSthat has multiple layer (one for each time step).
r<-rep(r,10)landscape<-sds(r)names(landscape)<-c("habitat_quality")landscape#> class : SpatRasterDataset#> subdatasets : 1#> dimensions : 90, 95 (nrow, ncol)#> nlyr : 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 : habitat_qualityBefore creating the simulation, it may be helpful to enable extensivereporting, which will print out a lot of information each time ametaRange function is called. This can be enabled or disabled at anytime (i.e. also while the simulation is running), but in order tohighlight what each function call in this tutorial does, we enable it atthe beginning of the setup.
After the landscape is loaded, the simulation can be created usingthecreate_simulation() function. The only requiredargument issource_environment which is the landscape /environmentSDS that was created in the first step. One canoptionally specify an ID for the simulation and a seed for the randomnumber generator.
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 : 1#> dimensions : 90, 95 (nrow, ncol)#> nlyr : 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 : habitat_quality#>#> created simulation: example_simulationIf you want to inspect the simulation object, you can either printit, to lists its fields and methods or use thesummary()function to get an overview of the simulation state.
sim#> metaRangeSimulation object#> Fields:#> $ID#> $globals#> $environment#> $number_time_steps#> $time_step_layer#> $current_time_step#> $queue#> $processes#> $seed#> Species: none#> Methods:#> $species_names()#> $add_globals()#> $add_species()#> $add_traits()#> $add_process()#> $begin()#> $exit()#> $set_current_time_step()#> $set_time_layer_mapping()#> $print()#> $summary()summary(sim)#> ID: example_simulation#> Environment:#> Fields:#> $current ==== the environment at the current time step#> classes : all -> matrix#> number : 1#> names : habitat_quality#> $sourceSDS == the source raster data of the environment#> class : SpatRasterDataset#> subdatasets : 1#> dimensions : 90, 95 (nrow, ncol)#> nlyr : 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 : habitat_quality#> Time step layer mapping: 1 2 3 4 5 6 7 8 9 10#> Current time step: 1#> Seed: 1#> Species: 0#>#> Simulation level processes:#> NULL#> Gobal variables:#> NULL#> Queue:#> Remaining queue (this time step): 0#> NULL#> Future queue (next time step): 0#> NULLAfter the simulation is created, species can be added to it using theadd_species() function. At this point one has to switch tothe syntax of theR6 packagethat metaRange uses. This means thatadd_species() is amethod of thesimulation object and can be called using the$ operator (i.e. by indexing the simulation object andcalling a function that is stored inside of it). The only requiredargument isname which is the name of the species that willbe added.
This species can now be accessed by using the$ operatoragain.
After species are added to the simulation, traits can be assigned tothem using theadd_traits() method. The first argument isspecies which is a character vector of the species names towhich the trait should be assigned to. The second argument ispopulation_level, aTRUE/FALSE value, thatdecides if the trait should be stored with one value per population(i.e. as a matrix of the same size as the landscape) or not (i.e. onlyone value per species). All following arguments can be supplied in theform oftrait_name = trait_value.
For now we only add three traits:abundance (number ofindividuals in each population),reproduction_rate (howfast the populations can reproduce) andcarrying_capacity(maximum number of individuals per grid cell). Note: Traits alwaysrepresent the “current” state of a species. This means that theabundance we use as input here represents the initial state of thesimulation. Over the course of the simulation (i.e. in each time step)the traits can be updated and changed. In this example, the abundancewill change each time step while e.g. the reproduction rate stays thesame, but in other cases each trait might change with time.
We can check what traits a species has by printing them:
Or plotting them:
Note that the above plot is not very interesting, since the abundanceis the same for each population at the beginning of the simulation.
After the species and its traits are added, the processes thatdescribe how the species interacts with its environment can be added,using theadd_process() method. The arguments are:species which is a character vector of the species (names)that should receive the process,process_name which is ahuman readable name for the process andprocess_fun whichis the function that will be called when the process is executed.
One argument that might be confusing is theexecution_priority. This is a number that gives the processa priority “weight” and decides in which order the processes areexecuted within one time step. The smaller the number, the earlier theprocess will be executed (e.g. 1 gets executed before 2). In the casetwo (or more) processes have the same priority, it is assumed that theyare independent from each other and that their execution order does notmatter.
In this example we will only add a single process(reproduction) to the species, that is going to calculatethe abundance (for each population) in the next time step, depending onthe habitat quality. To do so, we can use a built-in functionricker_reproduction_model() that implements the “classic”Ricker reproduction model (Ricker, W.E. (1954)) [Ref. 1], whichdescribes the population dynamics of a species with non-overlappinggenerations in discrete time steps. This model features densitydependent growth and possibly overcompensatory dynamics (i.e. thepopulations can, depending on the reproduction rate, become larger thanthe carrying capacity, which inevitably leads to a decline in the nexttime step).
Note the use of theself keyword in the function. Inthis context,self refers to the species that the processis attached to. This means that the function can access the speciestraits and modify them and also access the environment (each speciesholds a reference to the simulation it was created in).
sim$add_process(species ="species_1",process_name ="reproduction",process_fun =function() {# use a ricker reproduction model# to calculate the new abundance# and let the carrying capacity# depend on the habitat qualityricker_reproduction_model( self$traits$abundance, self$traits$reproduction_rate, self$traits$carrying_capacity* self$sim$environment$current$habitat_quality )# print out the current mean abundanceprint(paste0("mean abundance: ",mean(self$traits$abundance)) ) },execution_priority =1)#> adding process: reproduction#> to species:#> [1] "species_1"#>After the species, traits and processes are added to the simulation,it can be executed via thebegin() method.
sim$begin()#> Starting simualtion.#> passed initial sanity checks.#> start of time step: 1#> |- species_1 : reproduction#> [1] "mean abundance: 84.1732579542268"#> |---- 0.00059 secs#> 10 % done | 0.053 secs remaining (estimate)#> start of time step: 2#> |- species_1 : reproduction#> [1] "mean abundance: 127.598064388471"#> |---- 0.00046 secs#> 20 % done | 0.086 secs remaining (estimate)#> start of time step: 3#> |- species_1 : reproduction#> [1] "mean abundance: 185.433176212344"#> |---- 0.00045 secs#> 30 % done | 0.036 secs remaining (estimate)#> start of time step: 4#> |- species_1 : reproduction#> [1] "mean abundance: 254.974925085775"#> |---- 0.00047 secs#> 40 % done | 0.051 secs remaining (estimate)#> start of time step: 5#> |- species_1 : reproduction#> [1] "mean abundance: 328.335965177599"#> |---- 0.00044 secs#> 50 % done | 0.025 secs remaining (estimate)#> start of time step: 6#> |- species_1 : reproduction#> [1] "mean abundance: 394.79908318955"#> |---- 0.00045 secs#> 60 % done | 0.019 secs remaining (estimate)#> start of time step: 7#> |- species_1 : reproduction#> [1] "mean abundance: 446.224976588256"#> |---- 0.00043 secs#> 70 % done | 0.015 secs remaining (estimate)#> start of time step: 8#> |- species_1 : reproduction#> [1] "mean abundance: 480.705159122178"#> |---- 0.00043 secs#> 80 % done | 0.0097 secs remaining (estimate)#> start of time step: 9#> |- species_1 : reproduction#> [1] "mean abundance: 501.356439544315"#> |---- 0.00042 secs#> 90 % done | 0.0047 secs remaining (estimate)#> start of time step: 10#> |- species_1 : reproduction#> [1] "mean abundance: 512.803082103058"#> |---- 5e-04 secs#> 100 % done | 0 secs remaining (estimate)#>#> Simulation: 'example_simulation' finished#> Exiting the Simulation#> Runtime: 0.065 secsTo investigate the results, you can use theplot()function.
To save results you can use thesave_species() function.This will save the (possibly specified) traits of a species, either as araster (.tif) or as a text (.csv) file, whatever is more appropriate forthe data. Note that this function doesnot save the speciesprocesses. One should keep a copy of the script that is used to run thesimulation to make it repeatable.