The Ecological Model and Components

EcologicalNetworksDynamics represents an ecological network as a julia value of type Model.

m = default_model(Foodweb([:a => :b, :b => :c]))
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 17 components:
  - Species: 3 (:a, :b, :c)
  - Foodweb: 2 links
  - Body masses: [100.0, 10.0, 1.0]
  - Metabolic classes: [:invertebrate, :invertebrate, :producer]
  - GrowthRate: [·, ·, 1.0]
  - Carrying capacity: [·, ·, 1.0]
  - ProducersCompetition: 1.0
  - LogisticGrowth
  - Efficiency: 0.45 to 0.85.
  - MaximumConsumption: [8.0, 8.0, ·]
  - Hill exponent: 2.0
  - Consumers preferences: 1.0
  - Intra-specific interference: [·, ·, ·]
  - Half-saturation density: [0.5, 0.5, ·]
  - BioenergeticResponse
  - Metabolism: [0.09929551852928711, 0.1765751761097696, 0.0]
  - Mortality: [0.0, 0.0, 0.0]

Values of this type essentially describe a graph, with various nodes compartments representing e.g. species or nutrients and various edges compartments representing e.g. trophic links, or facilitation links. In addition to the network topology, the model also holds data describing the model further, and brought by the various models components. There are three possible "levels" for this data:

  • Graph-level data describe properties of the whole system. e.g. temperature, hill-exponent etc. These are typically scalar values.
  • Node-level data describe properties of particular nodes in the graph: e.g. species body mass, nutrients turnover etc. These are typically vector values.
  • Edge-level data describe properties of particular links: e.g. trophic links efficiency, half-saturation of a producer-to-nutrient links etc. These are typically matrix values.

Model Properties

The data held by the model can be accessed via the various model properties, with functions named get_<X>:

get_hill_exponent(m) # Graph-level data (a number).
2.0
get_body_masses(m) # Node-level data (a vector with one value per species).
3-element EcologicalNetworksDynamics.BodyMasses:
 100.0
  10.0
   1.0
get_efficiency(m) # Edge-level data (a matrix with one value per species interaction).
3×3 EcologicalNetworksDynamics.EfficiencyRates:
 0.0  0.85  0.0
 0.0  0.0   0.45
 0.0  0.0   0.0

Alternately, the data can also be accessed via julia's m.<x> property accessor:

m.hill_exponent # Same as get_hill_exponent(m).
m.body_masses   # Same as get_body_masses(m).
m.efficiency    # Same as get_efficiency(m).

Some data can be modified this way, either with set_<x>!(m, value). But not all:

# Okay: this is terminal data.
set_hill_exponent!(m, 2.1)
m.hill_exponent = 2.1 # (alternate syntax for the same operation)

# Not okay: this could make the rest of the model data inconsistent.
m.species_richness = 4
ERROR: In property 'species_richness' of '<inner parms>': This property is read-only.

If you need a model with different values for read-only data, you need to build a new model with the values you desire.

m = default_model(Foodweb([:a => :b, :b => [:c, :d]])) # Re-construct with a 4th species.
m.species_richness # Now the value is what you want.
4

The full list of available model properties can be queried with:

properties(m)

Model Components

The Model value is very flexible and can represent a variety of different networks. It is made from the combination of various components.

Empty Model and the add! Method.

When you start from a default_model, you typically obtain a full-fledged value, with all the components required to simulate the dynamics. Alternately, you can start from an empty model:

m = Model()
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 0 component.

In this situation, you need to add the components one by one. But this gives you full control over the model content.

An empty model cannot be simulated, because the data required for simulation is missing from it.

simulate(m, 0.5)
ERROR: In method 'simulate' for '<inner parms>': Requires a component 'Mortality'.

Also, an empty model cannot be queried for data, because there is no data inside:

m.richness
ERROR: In property 'richness' of '<inner parms>': Component 'Species' is required to read this property.

The most basic way to add a Species component to your model is to use the add! function:

add!(m, Species(3))
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 1 component:
  - Species: 3 (:s1, :s2, :s3)

Now that the Species component has been added, the related properties can be queried from the model:

m.richness
3
m.species_names
3-element EcologicalNetworksDynamics.SpeciesNames:
 :s1
 :s2
 :s3

But the other properties cannot be queried, because the associated components are still missing:

m.trophic_links
ERROR: In property 'trophic_links' of '<inner parms>': Component 'Foodweb' is required to read this property.

Before we add the missing Foodweb component, let us explain that the component addition we did above actually happened in two stages.

Blueprints Expand into Components.

To add a component to a model, we first need to create a blueprint for the component. A blueprint is a julia value containing all the data needed to construct a component.

sp = Species(3) # This is a blueprint, useful to later expand into a model component.
blueprint for Species:
  names: 3-element Vector{Symbol}:
 :s1
 :s2
 :s3

When you call the add! function, you feed it with a model and a blueprint. The blueprint is read and expanded into a component within the given model:

m = Model() # Empty model.
add!(m, sp) # Expand blueprint `sp` into a `Species` component within `m`.
m           # The result is a model with 1 component inside.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 1 component:
  - Species: 3 (:s1, :s2, :s3)

As we have seen before: once it has been expanded into the model, you cannot always edit the component data directly. For instance, the following does not work:

m.species_names[2] = "rhino"
ERROR: View error (EcologicalNetworksDynamics.SpeciesNames): This view into graph nodes data is read-only.

However, you can always edit the blueprint, then re-expand it later into other models.

sp.names[2] = :rhino    # Edit one species name within the blueprint.
push!(sp.names, :ficus) # Append a new species to the blueprint.
m2 = Model(sp)          # Create a new model from the modified blueprint.
m2                      # This new model contains the alternate data.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 1 component:
  - Species: 4 (:s1, :rhino, :s3, :ficus)

Blueprints can get sophisticated. For instance, here are various ways to create blueprints for a Foodweb component.

fw = Foodweb(:niche, S = 5, C = 0.2) # From a random model.
fw = Foodweb([0 1 0; 1 0 1; 0 0 0])  # From an adjacency matrix.
fw = Foodweb([1 => 2, 2 => 3])       # From an adjacency list.

If you want to test the corresponding Foodweb component, but you don't want to loose the original model, you can keep a safe copy of it before you actually expand the blueprint:

base = copy(m) # Keep a safe, basic, incomplete version of the model.
add!(m, fw)    # Expand the foodweb into a new component within `m`: `base` remains unchanged.

A shorter way to do so is to directly use julia's + operator, which always leaves the original model unchanged and creates an augmented copy of it:

m = base + fw # Create a new model `m` with a Foodweb inside, leaving model `base` unchanged.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 2 components:
  - Species: 3 (:s1, :s2, :s3)
  - Foodweb: 2 links

Separating blueprints creation from final components expansion gives you flexibility when creating your models. Blueprints can either be thrown after use, or kept around to be modified and reused without limits.

Model Constraints.

Of course, you cannot expand blueprints into components that would yield inconsistent models:

base = Model(Species(3)) # A model a with 3-species compartment.
m = base + Foodweb([0 1; 0 0]) # An adjacency matrix with only 2×2 values.
ERROR:
Could not expand blueprint for 'Foodweb' to system for '<inner parms>':
  Invalid size for parameter 'A': expected (3, 3), got (2, 2).

Components cannot be removed from a model, because it could lead to inconsistent model values. Components cannot either be duplicated or replaced within a model:

m = Model(Foodweb(:niche, S = 5, C = 0.2))
m += Foodweb([:a => :b]) # Nope: already added.
ERROR:
Could not expand blueprint for 'Foodweb' to system for '<inner parms>':
  cannot add components twice.

If you ever feel like you need to "change a component" or "remove a component" from a model, the correct way to do so is to construct a new model from the blueprints and/or the other base models you have kept around.

Components also require each other: you cannot specify trophic links efficiency in your model without having first specified what trophic links are:

m = Model(Species(3))
m += Efficiency(4)
ERROR:
Could not expand blueprint for 'Efficiency' to system for '<inner parms>':
  missing required component 'Foodweb'.

Blueprint Nesting (advanced).

To help you not hit the above problem too often, some blueprints take advantage of the fact that they contain the information needed to also expand into some of the components they require. Conceptually, they embed smaller blueprints within them.

For instance, the following blueprint for a foodweb contains enough information to expand into both a Foodweb component, and the associated Species component if needed:

fw = Foodweb([1 => 2, 2 => 3]) # Species nodes can be inferred from this blueprint..
m = Model(fw) # .. a blank model given only this blueprint becomes equiped with the 2 components.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 2 components:
  - Species: 3 (:s1, :s2, :s3)
  - Foodweb: 2 links

So it is not an error to expand the Foodweb component into a model not already having a Species compartment. We say that the Foodweb blueprint implies a Species blueprint.

If you need more species in your model than appear in your foodweb blueprint, you can still explicitly expand the Species blueprint before you add the foodweb:

m = Model(Species(5), Foodweb([1 => 2, 2 => 3])) # A model with 2 isolated species.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 2 components:
  - Species: 5 (:s1, :s2, :s3, :s4, :s5)
  - Foodweb: 2 links

Some blueprints, on the other hand, explicitly bring other blueprints. For instance, the LinearResponse brings both ConsumptionRate and ConsumersPreference sub-blueprints:

lin = LinearResponse()
blueprint for LinearResponse:
  alpha: blueprint for ConsumptionRate:
    alpha: 1.0
  w: blueprint for ConsumersPreferences:
    w: :homogeneous

So a model given this single blueprint can expand with 3 additional components.

m += lin
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 5 components:
  - Species: 5 (:s1, :s2, :s3, :s4, :s5)
  - Foodweb: 2 links
  - Consumption rate: [1.0, 1.0, 0.0, 0.0, 0.0]
  - Consumers preferences: 1.0
  - LinearResponse

The difference with "implication" though, is that the sub-blueprints "brought" do conflict with existing components:

m = Model(fw, ConsumptionRate(2)) # This model already has a consumption rate.
m += lin # So it is an error to bring another consumption rate with this blueprint.
ERROR:
Could not expand blueprint for 'LinearResponse' to system for '<inner parms>':
  blueprint also brings 'ConsumptionRate', which is already in the system.

This protects you from obtaining a model value with ambiguous consumption rates.

To prevent the sub-blueprint ConsumptionRate from being brought, you need to explicitly remove it from the blueprint containing it:

lin.alpha = nothing # Remove the brought sub-blueprint.
lin = LinearResponse(alpha = nothing) # Or create directly without the brought sub-blueprint.
m += lin # Consistent model obtained.
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 5 components:
  - Species: 3 (:s1, :s2, :s3)
  - Foodweb: 2 links
  - Consumption rate: [2.0, 2.0, 0.0]
  - Consumers preferences: 1.0
  - LinearResponse

Using the Default Model.

Building a model from scratch can be tedious, because numerous components are required for the eventual simulation to take place.

Here is how you could do it with only temporary blueprints immediately dismissed:

m = Model(
  Foodweb([:a => :b, :b => :c]),
  BodyMass(1),
  MetabolicClass(:all_invertebrates),
  BioenergeticResponse(),
  LogisticGrowth(),
  Metabolism(:Miele2019),
  Mortality(0),
)

Here is how you could do it with blueprints that you would keep around to later reassemble into other models:

# Basic blueprints saved into variables for later edition.
fw = Foodweb([:a => :b, :b => :c])
bm = BodyMass(1)
mc = MetabolicClass(:all_invertebrates)
be = BioenergeticResponse()
lg = LogisticGrowth()
mb = Metabolism(:Miele2019)
mt = Mortality(0)

# One model with all the associated components.
m = Model() + fw + bm + mc + be + lg + mb + mt

If this is too tedious, you can use the default_model function instead to automatically create a model with all (or most) components required for simulation. The only mandatory argument to default_model is a Foodweb blueprint:

fw = Foodweb([:a => :b, :b => :c])
m = default_model(fw)
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 17 components:
  - Species: 3 (:a, :b, :c)
  - Foodweb: 2 links
  - Body masses: [100.0, 10.0, 1.0]
  - Metabolic classes: [:invertebrate, :invertebrate, :producer]
  - GrowthRate: [·, ·, 1.0]
  - Carrying capacity: [·, ·, 1.0]
  - ProducersCompetition: 1.0
  - LogisticGrowth
  - Efficiency: 0.45 to 0.85.
  - MaximumConsumption: [8.0, 8.0, ·]
  - Hill exponent: 2.0
  - Consumers preferences: 1.0
  - Intra-specific interference: [·, ·, ·]
  - Half-saturation density: [0.5, 0.5, ·]
  - BioenergeticResponse
  - Metabolism: [0.09929551852928711, 0.1765751761097696, 0.0]
  - Mortality: [0.0, 0.0, 0.0]

But you can feed other blueprints into it to fine-tweak only the parameters you want to modify.

m = default_model(fw, BodyMass(Z = 1.5), Efficiency(2))
(m.body_masses, m.efficiency)
([2.25, 1.5, 1.0], [0.0 2.0 0.0; 0.0 0.0 2.0; 0.0 0.0 0.0])

The function default_model tries hard to figure the default model you expect based on the few blueprints you input. For instance, it assumes that you need a different type of functional response if you input a Temperature component, and temperature-dependent allometry rates:

m = default_model(fw, Temperature(220))
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 18 components:
  - Species: 3 (:a, :b, :c)
  - Foodweb: 2 links
  - Temperature: 220.0K
  - Body masses: [100.0, 10.0, 1.0]
  - Metabolic classes: [:invertebrate, :invertebrate, :producer]
  - GrowthRate: [·, ·, 2.4457111255767084e-12]
  - Carrying capacity: [·, ·, 34342.667662374224]
  - ProducersCompetition: 1.0
  - LogisticGrowth
  - Efficiency: 0.45 to 0.85.
  - Hill exponent: 2.0
  - Consumers preferences: 1.0
  - Intra-specific interference: [·, ·, ·]
  - Handling time: 170440.9238100917 to 178473.55707771532.
  - Attack rates: 6.894034926202549e-9 to 2.4460958976605684e-8.
  - ClassicResponse
  - Metabolism: [1.7880778943692315e-12, 3.650786484958455e-12, 0.0]
  - Mortality: [0.0, 0.0, 0.0]

Or if you wish to explicitly represent Nutrients as a separate nodes compartment in your ecological network:

m = default_model(fw, Nutrients.Nodes(2))
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 20 components:
  - Species: 3 (:a, :b, :c)
  - Foodweb: 2 links
  - Nutrients: 2 (:n1, :n2)
  - Body masses: [100.0, 10.0, 1.0]
  - Metabolic classes: [:invertebrate, :invertebrate, :producer]
  - GrowthRate: [·, ·, 1.0]
  - Nutrients turnover: [0.25, 0.25]
  - Nutrients supply: [4.0, 4.0]
  - Nutrients concentration: 0.5
  - Nutrients half-saturation: 0.15
  - NutrientIntake
  - Efficiency: 0.45 to 0.85.
  - MaximumConsumption: [8.0, 8.0, ·]
  - Hill exponent: 2.0
  - Consumers preferences: 1.0
  - Intra-specific interference: [·, ·, ·]
  - Half-saturation density: [0.5, 0.5, ·]
  - BioenergeticResponse
  - Metabolism: [0.09929551852928711, 0.1765751761097696, 0.0]
  - Mortality: [0.0, 0.0, 0.0]

But the function will not choose between two similar blueprints if you bring both, even implicitly:

m = default_model(
  fw,
  BodyMass(Z = 1.5),      # <- Customize body mass.
  ClassicResponse(e = 2), # <- This blueprint also brings a BodyMassy
)
Model (alias for EcologicalNetworksDynamics.Framework.System{<inner parms>}) with 17 components:
  - Species: 3 (:a, :b, :c)
  - Foodweb: 2 links
  - Body masses: [2.25, 1.5, 1.0]
  - Metabolic classes: [:invertebrate, :invertebrate, :producer]
  - GrowthRate: [·, ·, 1.0]
  - Carrying capacity: [·, ·, 1.0]
  - ProducersCompetition: 1.0
  - LogisticGrowth
  - Efficiency: 2.0
  - Hill exponent: 2.0
  - Consumers preferences: 1.0
  - Intra-specific interference: [·, ·, ·]
  - Handling time: 0.155544054398016 to 0.24694341535807782.
  - Attack rates: 60.00826508049385 to 76.53601152370896.
  - ClassicResponse
  - Metabolism: [0.25637992641130597, 0.2837310291334913, 0.0]
  - Mortality: [0.0, 0.0, 0.0]

In this situation, either stop implicitly bringing BodyMass with ClassicResponse(e=2, M=nothing), or directly move you custom body mass input into the larger blueprint:

m = default_model(
  fw,
  ClassicResponse(e = 2, M = (; Z = 1.5)),
)