A quick go through
Creating a model
The definition of ATN is based on a model object (formally a S4 class in R). The model object is initialised specifying a fixed set of parameters: the number of species, the number of basal species, species body masses, a matrix defining the trophic interactions and, for the version including the nutrient dynamics, the number of nutrients. The first thing to do is therefore to create the corresponding R variables. While one can use an empirical food web for its analysis, it is also possible to generate synthetic food webs using the niche model from Williams and Martinez (2000) or using allometric scaling as defined in Schneider et al. (2016).
Generating synthetic food webs (if needed)
ATNr has two functions to generate synthetic food webs,
create_niche_model() for the niche model (Williams and Martinez (2000)) and
create_Lmatrix() for the allometric scaling model (Schneider et al. (2016)). The niche model
requires information on the number of species and connectance of the
desired food web:
library(ATNr)
set.seed(123)
n_species <- 20 # number of species
conn <- 0.3 # connectance
fw <- create_niche_model(n_species, conn)
# The number of basal species can be calculated:
n_basal <- sum(colSums(fw) == 0)As the niche model does not rely on allometry, it is possible to
estimate species body masses from their trophic levels, which can be
calculated form the TroLev function of the package. For
instance:
TL = TroLev(fw) #trophic levels
masses <- 1e-2 * 10 ^ (TL - 1)The allometric scaling model generate links based on species body masses. Therefore, it requires as an input a vector containing the body mass of species, as well as a parameter informing on the number of basal species desired. It produces a so-called L matrix which formally quantifies the probability for a consumer to successfully attack and consumer an encountered resource:
n_species <- 20
n_basal <- 5
masses <- sort(10^runif(n_species, 2, 6)) #body mass of species
L <- create_Lmatrix(masses, n_basal)This L matrix can then be transformed into a binary food web:
fw <- L
fw[fw > 0] <- 1More details about the generative models and the the usage precaution around them can be found in the section “The food web generative functions”
Creating a specific ATN model
As soon as a food web is stored in a matrix, it is possible to create a model object that refers to the desired specific model
# initialisation of the model object. It is possible to create a ode corresponding to
# Schneider et al. 2016, Delmas et al. 2016 or Binzer et al. 2016:
# 1) Schneider et al. 2016
n_nutrients <- 3
model_unscaled_nuts <- create_model_Unscaled_nuts(n_species, n_basal, n_nutrients, masses, fw)
# 2) Delmas et al. 2016:
model_scaled <- create_model_Scaled(n_species, n_basal, masses, fw)
# 3) Binzer et al. 2016
model_unscaled <- create_model_Unscaled(n_species, n_basal, masses, fw)Once created, it is possible to access to the methods and attributes of the object to initialise or update them:
# updating the hill coefficient of consumers in the Unscaled_nuts model:
model_unscaled_nuts$q <- rep(1.4, model_unscaled_nuts$nb_s - model_unscaled_nuts$nb_s)
# Changing the assimilation efficiencies of all species to 0.5 in the Scaled model:
model_scaled$e = rep(0.5, model_scaled$nb_s)
# print the different fields that can be updated and their values:
# str(model_unscaled_nuts)It is important to keep in mind that some rules apply here:
The order of the species in the different fields must be consistent: the first species in the
$BMobject corresponds to the first species in the$fwobject and in the$eobject.The objects that are specific to a species type (i.e. basal species or consumers) are dimensioned accordingly: the handling time (
$h) sets the handling time of consumers on resources. Therefore, thehmatrix has a number of rows equal to the number of species and a number of columns equal to the number of consumers (as non consumer species do not have a handling time by definition). In that case, the first row correspond to the first species and the first column to the first consumer.The object describing the interactions between plants and nutrients (
$Kor$V) are matrices for which the number of rows equals to the number of nutrients and a number of columns matches the number of basal species (this point is specific to the Schneider model which is the only one including an explicit dynamics of the nutrient pool).
To run the population dynamics, all the parameters must be defined. It is possible to automatically load a by default parametrisation using the dedicated functions:
# for a model created by create_model_Unscaled_nuts():
model_unscaled_nuts <- initialise_default_Unscaled_nuts(model_unscaled_nuts, L)
# for a model created by create_model_Scaled():
model_scaled <- initialise_default_Scaled(model_scaled)
# for a model created by create_model_Unscaled():
model_unscaled <- initialise_default_Unscaled(model_unscaled)Importantly, for the unscaled with nutrients model of Schneider et al. (2016), the calculation of
consumption rates rely on the L matrix created above or, in case of
empirical networks, by a matrix that defines the probability of a
consumer to successfully attack and consume an encountered prey. The
default initialisation of the Unsacled_nuts and
Unscaled models can also include temperature effects (20° C
by default).
Running the population dynamics
Once all the parameters are properly defined, the ODEs can be
integrated by any solver. We present here a solution based on
lsoda from library deSolve (Soetaert, Petzoldt, and Setzer (2010)), but
other solutions exist (sundialr is also a possibility). The
package propose a direct wrapper to lsoda with the function
lsoda_wrapper:
biomasses <-runif(n_species, 2, 3) # starting biomasses
biomasses <- append(runif(3, 20, 30), biomasses) # nutrient concentration
# defining the desired integration time
times <- seq(0, 1500, 5)
sol <- lsoda_wrapper(times, biomasses, model_unscaled_nuts)To have more control of the integration, it is however possible to
not use the wrapper proposed in the package and directly work with the
lsoda function. Here is an example:
# running simulations for the Schneider model
model_unscaled_nuts$initialisations()
sol <- deSolve::lsoda(
biomasses,
times,
function(t, y, params) {
return(list(params$ODE(y, t)))
},
model_unscaled_nuts
)Not that the call of
model_unscaled_nuts$initialisation() is important here as
it pre-computes some variables to optimise code execution. This function
is normally internally called by lsoda_wrapper. In case of
an integration that does not rely on this wrapper function, the call to
$initialisation() is needed for ALL the model types.
The package also contains a simple function to plot the time series obtained: plot_odeweb. The colours only differentiate the species using their ranks in the food web matrix (from blue to red).
par(mar = c(4, 4, 1, 1))
plot_odeweb(sol, model_unscaled_nuts$nb_s)
