In drfinlayscott/mizer: Dynamic Multi-Species Size Spectrum Modelling

library(mizer)

The previous sections have used wrapper functions to set up MizerParams objects that are appropriate for single-species, community- and trait-based models. We now turn our attention to multispecies, or species-specific, models. These are potentially more complicated than the community and trait-based models and use the full power of the mizer package.

In multispecies type models multiple species are resolved. However, unlike in the trait-based model which also resolves multiple species, the species parameters will be those of real-world species. There are several advantages to this approach. As well as investigating the community as a whole (as was done for the community and trait-based models), we are able to investigate the dynamics of individual species. This means that species specific management rules can be tested and species specific metrics, such as yield, can be compared to reference levels.

A multispecies model can take more effort to set up. For example, each species will have different life-history parameters; there may be multiple gear types with different selectivities targeting different groups of species; the fishing effort of each gear may change with time instead of just being constant (which has been the case in the simulations we have looked at so far); the interactions between the species needs to be considered.

In later sections we build up a multispecies model for the North Sea. To effectively use mizer for a multispecies model we are going to have to take a closer look at the MizerParams class and the project() function. This will all be done in the context of examples so hopefully everything will be clear.

We also take a closer look at some of the summary plots and analyses that can be performed, for example, calculating a range of size-based indicators.

Setting up a multispecies model {#sec:MizerParams}

Overview

The MizerParams class is used for storing model parameters. We have already met the MizerParams class when we looked at community and trait-based models. However, to set up a multispecies model we will need to specify many more parameters.This is probably the most complicated part of using the mizer package, so we will take it slowly.

A MizerParams object stores the:

• life-history parameters of the species in the community, such as maximum size $w_{max}$;
• size-based biological parameters for the species, such as the search volume;
• density-dependent reproduction functions and parameters of each species;
• an interaction matrix to describe the spatial overlap of pairs of species;
• parameters relating to the growth and dynamics of the resource spectrum;
• fishing gear parameters: selectivity and catchability.

Note that the MizerParams class does not store any parameters that can vary through time, such as fishing effort or population abundance. These are stored in the MizerSim class which we will come to later in the section on running a simulation.

Although the MizerParams class contains a lot of information, it is relatively straightforward to set up and use. Objects of class MizerParams are created using the constructor method newMultispeciesParams() (this constructor method was called MizerParams() in previous version of mizer). This constructor method can take many arguments. However, creation is simplified because many of the arguments have default values.

In the rest of this section we look at the main arguments to the newMultispeciesParams() function. To help understand how the constructor is used and how the MizerParams class relates to the equations given in the model description section, there is an example section where we create example parameter objects using data that comes with the mizer package.

The species parameters {#sec:species_parameters_dataframe}

Although many of the arguments used when creating a MizerParams object are optional, there is one argument that must be supplied by the user: the species specific parameters. These are stored in a single data.frame object. The data.frame is arranged species by parameter, so each column is a parameter and each row has the parameters for one of the species in the model. Although it is possible to create the data.frame by hand in R, it is probably easier to create the data externally as a .csv file (perhaps using a suitable open source spreadsheet such as LibreOffice) and then read the data into R.

For each species in the model community there are certain parameters that are essential and that do not have default values. The user must provide values for these parameters. There are also some essential parameters that have default values, such as the selectivity function parameters, and some that are calculated internally using default relationships if not explicitly provided. These defaults are used if the parameters are not found in the data.frame.

The essential columns of the species parameters data.frame that have no default values are: species, the names of the species in the community and w_max, the maximum mass of the species.

The gear parameters {#sec:fishing_gear}

In mizer, fishing mortality is imposed on species by fishing gears. The total fishing mortality is obtained by summing over the mortality from all gears, \begin{equation} % {#eq:muf} \mu_{f.i}(w) = \sum_g F_{g,i}(w), \end{equation} where the fishing mortality $F_{g,i}(w)$ imposed by gear $g$ on species $i$ at size $w$ is calculated as: \begin{equation} % {#eq:sel} F_{g,i}(w) = S_{g,i}(w) Q_{g,i} E_{g} \end{equation} where $S$ is the selectivity by species, gear and size, $Q$ is the catchability by species and gear and $E$ is the fishing effort by gear. The selectivity at size has a range between 0 (not selected at that size) to 1 (fully selected at that size). Catchability is used as an additional scalar to make the link between gear selectivity, fishing effort and fishing mortality. For example, it can be set so that an effort of 1 gives a desired fishing mortality. In this way effort can then be specified relative to a 'base effort', e.g. the effort in a particular year.

Selectivity and catchability are stored as arrays in the MizerParams object. However the user does not have to create these arrays by hand if they provide a data frame with the necessary information. In particular the selectivity can be calculate by specifying functions for the selectivity curves. Mizer provides a range of such selectivity functions and the user just needs to specify their parameters for each gear and each species in the gear_params data frame. All the details can be found on the help page for setFishing().

Fishing effort is not stored in the MizerParams object. Instead, effort is set when the simulation is run and can vary through time (see the section on running a simulation).

Example of making MizerParams objects {#sec:params_example}

As mentioned in the preceding sections, an object of MizerParams is created by using the newMultispeciesParams() constructor method.

The first step is to prepare the species specific parameter data.frame. As mentioned above, one way of doing this is to use a spreadsheet and save it as a .csv file. We will use this approach here. An example .csv file has been included in the package. This contains the species parameters for a multispecies North Sea model. The location of the file can be found by running

params_location <- system.file("extdata", "NS_species_params.csv",
                               package = "mizer")

This file can be opened with most spreadsheets or a text editor for you to inspect. This can be loaded into R with

species_params <- read.csv(params_location)

This reads the .csv file into R in the form of a data.frame. You can check this with the class:

class(species_params)

Let's have a look at the data frame:

species_params

You can see that there are $r nrow(species_params)$ species and $r ncol(species_params)$ columns of parameters: species, w_max,w_mat,beta,sigma,R_max and k_vb.

Of these parameters, species and w_max are essential and have no default values (as described in the section on species parameters). w_max is the maximum size of the species, w_mat is its maturity size, and beta and sigma are parameters of the predation kernel (by default mizer uses a log-normal predation kernel). R_max is a parameter introducing additional density dependence into reproduction parameter using a Beverton-Holt type function (see setReproduction() for details). The final column, k_vb, will be used to calculate values for h and then gamma. This column is only essential here because the h and gamma are not included in the data.frame. It would also have been possible to include h and gamma columns in the data.frame and not include the k_vb column.

The values of the non-essential species specific parameters, like for example alpha, k, ks, z0, w_min and erepro, were not included in the data.frame. This means that the default values will be automatically used when we create the MizerParams object.

For this example we will not set up gear parameters. There are no columns describing the fishing selectivity. There is no sel_func column to determine the selectivity function. This means that the default selectivity function, knife_edge, will be used. As mentioned in the section on fishing gears, this function also needs another argument, knife_edge_size. This is not present in the data.frame and so it will be set to the default value of w_mat. Also, there is no catchability column so a default value for catchability of 1 will be used for all gears and species.

To create the MizerParams object we pass the species parameter data.frame into the newMultispeciesParams() constructor method:

params <- newMultispeciesParams(species_params)

We have just created a MizerParams object:

class(params)

The MizerParams object also stores a copy of the species parameter data frame that we provided. We can look at it with species_params():

species_params(params)

We can see that this returns the original species data.frame (with w_max and so on), plus any default values that may not have been included in the original data.frame. For example, we can see that there are now columns for alpha and h and gamma etc.

Also note how the default fishing gears have been set up. Even though we did not provide a gear parameter data frame, the MizerParams object has one that we can access with

gear_params(params)

All species are caught by a gear called "knife_edge_gear". The selectivity function for each fishing gear has been set in the sel_func column to the default function, knife_edge(). A catchability column has been added with a default value of 1 for each of the species that the gear catches. An example of setting the catchability by hand can be seen in the section on the North Sea.

There is a summary() method for MizerParams objects which prints a useful summary of the model parameters:

summary(params)

As well as giving a summary of the species in the model and what gear is fishing what species, it gives a summary of the size structure of the community. For example there are $r length(w(params))$ size classes in the community, ranging from $r signif(min(w(params)), 3)$ to $r signif(max(w(params)), 3)$ . These values are controlled by the arguments no_w, min_w and max_w respectively. For example, if we wanted 200 size classes in the model we would use:

params200 <- newMultispeciesParams(species_params, no_w = 200)
summary(params200)

Setting the interaction matrix

So far we have created a MizerParams object by passing in only the species parameter data.frame argument. We did not specify an interaction matrix. The interaction matrix \eqn{\theta_{ij}} describes the interaction of each pair of species in the model. This can be viewed as a proxy for spatial interaction e.g. to model predator-prey interaction that is not size based. The values in the interaction matrix are used to scale the encountered food in [getEncounter()] and the predation mortality rate in [getPredMort()] (see the section on predator-prey encounter rate and on predation mortality).

The entries of the interaction matrix are dimensionless numbers taking values are between 0 (species do not overlap and therefore do not interact with each other) to 1 (species overlap perfectly). By default mizer sets all values to 1, implying that all species fully interact with each other, i.e. the species are spread homogeneously across the model area.

getInteraction(params)

For the North Sea this is not the case and so the model would be improved by also including an interaction matrix which describes the spatial overlap between species.

An example interaction matrix for the North Sea has been included in mizer as a .csv file. The location of the file can be found by running:

inter_location <- system.file("extdata", "NS_interaction.csv",
                              package = "mizer")

Take a look at it in a spreadsheet if you want. As mentioned above, to read this file into R we can make use of the read.csv() function. However, this time we want the first column of the .csv file to be the row names. We therefore use an additional argument to the read.csv() function: row.names.

inter <- read.csv(inter_location, row.names = 1)
inter

We can set the interaction matrix in our existing MizerParams object params with the setInteraction() function:

params <- setInteraction(params, interaction = inter)

Alternatively, instead of changing the interaction matrix in the existing MizerParams object, we could have created a new object from scratch with our interaction matrix by passing it to newMultispeciesParams():

params_new <- newMultispeciesParams(species_params, interaction = inter)

Note that the first argument must be the species parameters data.frame. The remaining arguments can be in any order but should be named. We are using the default values for all other parameters.

We now have all we need to start running projections. Before we get to that though, we'll take a quick look at how different fishing gears can be set up.

Setting different gears {#sec:setting_different_gears}

In the above example, each species is caught by the same gear (named "knife_edge_gear"). This is the default when no gear information is provided.

gear_params(params)

Here, we look at an example where we set up four different gears: Industrial, Pelagic, Beam and Otter trawl, that catch different combinations of species. We can achieve that by only changing the gear column in the gear_params data frame.

gear_params(params)$gear <- c("Industrial", "Industrial", "Industrial",
                              "Pelagic", "Beam", "Otter",
                              "Beam", "Otter", "Beam",
                              "Otter", "Otter", "Otter")

You can see the result by calling summary() on the params object.

summary(params)

In this example the same gear now catches multiple stocks. For example, the Industrial gear catches Sprat, Sandeel and Norway Pout. Why would we want to set up the gears like this? In the next section on running a multispecies model we will see that to project the model through time you can specify the fishing effort for each gear through time. By setting the gears up in this way you can run different management scenarios of changing the efforts of the fishing gears rather than on individual species. It also means that after a simulation has been run you can examine the catches by gear.

Setting to steady state

Once the MizerParams object has been properly set up, it may be the case that one wishes put the system in steady state. Sometimes this can be done simply by running the model using project() until it reaches steady state. However, this method is not guaranteed to work, and there is a function called steady() that is more reliable. The function steady() must be supplied with a MizerParams object. It takes that MizerParams object, looks at the initial system state, computes the levels of reproduction of the different species, hold them fixed, and evolves the system until a steady state is reached (or more precisely, until the amount that the population abundances change during a time-step is below some small tolerance level). After this, the reproductive efficiency of each species is altered so that when the reproduction dynamics are turned back on (i.e., when we stop holding recruitment levels fixed), the values of the reproduction levels which we held the system fixed at will be realized. The steady function is not sure to converge, and the way it re-tunes the reproductive efficiency values may not be realistic, but the idea is to alter the other parameters in the system until steady() does arrive at a steady state with sensible reproductive efficiency values.

Now that we know how to create a multispecies model we shall discuss how to run a multispecies model.

drfinlayscott/mizer documentation built on April 13, 2024, 9:16 a.m.