R/MizerParams-class.R
newMultispeciesParams.RdSets up a multi-species size spectrum model by filling all slots in the
MizerParams object based on user-provided or default
parameters. It does this by creating an empty MizerParams object with
emptyParams and then filling the slots by passing its arguments
to setParams. There is a long list of arguments, but almost
all of them have sensible default values. All arguments are described in more
details in the sections below the list.
newMultispeciesParams( species_params, interaction = matrix(1, nrow = nrow(species_params), ncol = nrow(species_params)), no_w = 100, min_w = 0.001, max_w = NA, min_w_pp = NA, no_w_pp = NA, n = 2/3, p = 0.7, f0 = 0.6, pred_kernel = NULL, search_vol = NULL, intake_max = NULL, metab = NULL, mu_b = NULL, z0pre = 0.6, z0exp = n - 1, maturity = NULL, repro_prop = NULL, srr = "srrBevertonHolt", rr_pp = NULL, cc_pp = NULL, r_pp = 10, kappa = 1e+11, lambda = 2.05, w_pp_cutoff = 10, plankton_dynamics = "plankton_semichemostat", initial_effort = NULL )
| species_params | A data frame of species-specific parameter values. |
|---|---|
| interaction | Interaction matrix of the species (predator by prey). Entries should be numbers between 0 and 1. See "Setting interactions" section below. |
| no_w | The number of size bins in the consumer spectrum. |
| min_w | Sets the size of the eggs of all species for which this is not
given in the |
| max_w | The largest size of the consumer spectrum. By default this is set to the largest w_inf specified in the species_params data frame. |
| min_w_pp | The smallest size of the plankton spectrum. By default this is set to the smallest value at which any of the consumers can feed. |
| no_w_pp | Obsolete argument that is no longer used because the number of plankton size bins is determined because all size bins have to be logarithmically equally spaced. |
| n | Scaling exponent of the intrinsic growth rate |
| pred_kernel | Optional. An array (species x predator size x prey size) that holds the predation coefficient of each predator at size on each prey size. If not supplied, a default is set as described in section "Setting predation kernel". |
| search_vol | Optional. An array (species x size) holding the search volume for each species at size. If not supplied, a default is set as described in the section "Setting search volume". |
| intake_max | Optional. An array (species x size) holding the maximum intake rate for each species at size. If not supplied, a default is set as described in the section "Setting maximum intake rate". |
| metab | Optional. An array (species x size) holding the metabolic rate for each species at size. If not supplied, a default is set as described in the section "Setting metabolic rate". |
| mu_b | Optional. An array (species x size) holding the background mortality rate. |
| z0pre | If |
| z0exp | If |
| maturity | Optional. An array (species x size) that holds the proportion of individuals of each species at size that are mature. If not supplied, a default is set as described in the section "Setting reproduction". |
| repro_prop | Optional. An array (species x size) that holds the proportion of consumed energy that a mature individual allocates to reproduction for each species at size. If not supplied, a default is set as described in the section "Setting reproduction". |
| srr | The name of the stock recruitment function. Defaults to
|
| rr_pp | Optional. Vector of plankton intrinsic growth rates |
| cc_pp | Optional. Vector of plankton intrinsic carrying capacity |
| r_pp | Coefficient of the intrinsic growth rate |
| kappa | Coefficient of the intrinsic carrying capacity |
| lambda | Scaling exponent of the intrinsic carrying capacity |
| w_pp_cutoff | The upper cut off size of the plankton spectrum. Default is 10 g. |
| plankton_dynamics | Function that determines plankton dynamics by calculating the plankton spectrum at the next time step from the current state. |
| initial_effort | Optional. A number or a named numeric vector specifying the fishing effort. If a number, the same effort is used for all gears. If a vector, must be named by gear. |
An object of type MizerParams
The only essential argument is a data frame that contains the species parameters. The data frame is arranged species by parameter, so each column of the parameter data frame is a parameter and each row has the values of the parameters for one of the species in the model.
There are two essential columns that must be included in the species
parameter data.frame and that do not have default values: the
species column that should hold strings with the names of the
species and the w_inf column with the asymptotic sizes of the species.
The species_params dataframe also needs to contain the parameters needed by any predation kernel function or size selectivity function. This will be mentioned in the appropriate sections below.
For all other species parameters, mizer will calculate default values if they
are not included in the species parameter data frame. They will be
automatically added when the MizerParams object is created. For these
parameters you can also specify values for only some species and leave the
other entries as NA and the missing values will be set to the defaults.
All the parameters will be mentioned in the following sections.
Mizer uses grams to measure weight, centimetres to measure lengths, and years to measure time.
Mizer is agnostic about whether abundances are given as
numbers per area,
numbers per volume or
total numbers for the entire study area.
You should make the choice most convenient for your application and then stick with it. If you make choice 1 or 2 you will also have to choose a unit for area or volume. Your choice will then determine the units for some of the parameters. This will be mentioned when the parameters are discussed in the sections below.
You choice will also affect the units of the quantities you may want to
calculate with the model. For example, the yield will be in grams/year/m^2 in
case 1 if you choose m^2 as your measure of area, in grams/year/m^3 in case 2
if you choose m^3 as your unit of volume, or simply grams/year in case 3. The
same comment applies for other measures, like total biomass, which will be
grams/area in case 1, grams/volume in case 2 or simply grams in case 3. When
mizer puts units on axes, for example in plotBiomass, it will simply
put grams, as appropriate for case 3.
You can convert between these choices. For example, if you use case 1, you need to multiply with the area of the ecosystem to get the total quantity. If you work with case 2, you need to multiply by both area and the thickness of the productive layer. In that respect, case 2 is a bit cumbersome.
The species interaction matrix \(\theta_{ij}\), is used when calculating the
food encounter rate in getEncounter and the predation mortality rate in
getPredMort. Its entries are dimensionless numbers between
0 and 1 that characterise the strength at which predator species \(i\)
predates on prey species \(j\).
This function checks that the supplied interaction
matrix is valid and then stores it in the interaction slot of the
params object before returning that object.
The order of the columns and rows of the interaction argument should be the
same as the order in the species params dataframe in the params object.
If you supply a named array then the function will check the order and warn
if it is different.
The interaction of the species with the plankton are set via a column
interaction_p in the species_params data frame. Again the entries
have to be numbers between 0 and 1. By default this column is set to all
1s.
Kernel dependent on predator to prey size ratio
If the pred_kernel argument is not supplied, then this function sets a
predation kernel that depends only on the ratio of predator mass to prey
mass, not on the two masses independently. The shape of that kernel is then
determined by the pred_kernel_type column in species_params.
The default pred_kernel_type is "lognormal". This will call the function
lognormal_pred_kernel to calculate the predation kernel.
An alternative pred_kernel type is "box", implemented by the functions
box_pred_kernel. These functions require certain species
parameters in the species_params data frame. For the lognormal kernel these
are beta and sigma, for the box kernel they are
ppmr_min and ppmr_max. They are explained in the help pages
for the kernel functions. No defaults are set for these parameters. If they
are missing from the species_params data frame then mizer will issue an
error message.
You can use any other string as the type. If for example you choose "my" then
you need to define a function my_pred_kernel that you can model on the
existing functions like lognormal_pred_kernel.
When using a kernel that depends on the predator/prey size ratio only, mizer
does not need to store the entire three dimensional array in the MizerParams
object. Such an array can be very big when there is a large number of size
bins. Instead, mizer only needs to store two two-dimensional arrays that hold
Fourier transforms of the feeding kernel function that allow the encounter
rate and the predation rate to be calculated very efficiently. However, if
you need the full three-dimensional array you can calculate it with the
getPredKernel function.
Kernel dependent on both predator and prey size
If you want to work with a feeding kernel that depends on predator mass and prey mass independently, you can specify the full feeding kernel as a three-dimensional array (predator species x predator size x prey size). The dimensions are thus (no_sp, no_w, no_w_full).
You should use this option only if a kernel dependent only on the predator/prey mass ratio is not appropriate. Using a kernel dependent on predator/prey mass ratio only allows mizer to use fast Fourier transform methods to significantly reduce the running time of simulations.
The order of the predator species in pred_kernel should be the same
as the order in the species params dataframe in the `params` object. If you
supply a named array then the function will check the order and warn if it is
different.
The search volume \(\gamma_i(w)\) of an individual of species \(i\)
and weight \(w\) multiplies the predation kernel when
calculating the encounter rate in getEncounter and the
predation rate in getPredRate.
The name "search volume" is a bit misleading, because \(\gamma_i(w)\) does not have units of volume. It is simply a parameter that determines the rate of predation. Its units depend on your choice, see section "Units in mizer". If you have chose to work with total abundances, then it is a rate with units 1/year. If you have chosen to work with abundances per m^2 then it has units of m^2/year. If you have chosen to work with abundances per m^3 then it has units of m^3/year.
If the search_vol argument is not supplied, then the search volume is
set to
$$\gamma_i(w) = \gamma_i w^q_i.$$
The values of \(\gamma_i\) (the search volume at 1g) and \(q_i\) (the
allometric exponent of the search volume) are taken from the gamma and
q columns in the species parameter dataframe. If the gamma
column is not supplied in the species parameter dataframe, a default is
calculated by the get_gamma_default function. Note that only
for predators of size \(w = 1\) gram is the value of the species parameter
\(\gamma_i\) the same as the value of the search volume \(\gamma_i(w)\).
The maximum intake rate \(h_i(w)\) of an individual of species \(i\) and
weight \(w\) determines the feeding level, calculated with
getFeedingLevel. It is measured in grams/year.
If the intake_max argument is not supplied, then the maximum intake
rate is set to $$h_i(w) = h_i w^n_i.$$
The values of \(h_i\) (the maximum intake rate of an individual of size 1
gram) and \(n_i\) (the allometric exponent for the intake rate) are taken
from the h and n columns in the species parameter dataframe. If
the h column is not supplied in the species parameter dataframe, it is
calculated by the get_h_default function, using f0 and
the k_vb column, if they are supplied.
If \(h_i\) is set to Inf, fish will consume all encountered food.
The metabolic rate is subtracted from the energy income rate to calculate
the rate at which energy is available for growth and reproduction, see
getEReproAndGrowth. It is measured in grams/year.
If the metab argument is not supplied, then for each species the
metabolic rate \(k(w)\) for an individual of size \(w\) is set to
$$k(w) = ks w^p + k w,$$
where \(ks w^p\) represents the rate of standard metabolism and \(k w\)
is the rate at which energy is expended on activity and movement. The values
of \(ks\), \(p\) and \(k\) are taken from the ks, p and
k columns in the species parameter dataframe. If any of these
parameters are not supplied, the defaults are \(k = 0\), \(p = n\) and
$$ks = f_c h \alpha w_{mat}^{n-p},$$
where \(f_c\) is the critical feeding level taken from the fc column
in the species parameter data frame. If the critical feeding level is not
specified, a default of \(f_c = 0.2\) is used.
The background mortality is all the mortality that is not due to either predation or fishing. It is a rate with units 1/year.
The mu_b argument allows you to specify a background mortality rate
that depends on species and body size. You can see an example of this in
the Examples section of the help page for setBMort.
If the mu_b argument is not supplied, then the background mortality
is assumed to depend only on the asymptotic size of the species, not on the
size of the individual: \(\mu_{b.i}(w) = z_{0.i}\). The value of the
constant \(z_0\) for each species is taken from the z0 column of the
species_params data frame, if that column exists. Otherwise it is calculated
as
$$z_{0.i} = {\tt z0pre}_i\, w_{inf}^{\tt z0exp}.$$
For each species and at each size, the proportion of the available energy
that is invested into reproduction is the product of two factors: the
proportion maturity of individuals that are mature and the proportion
repro_prop of the energy available to a mature individual that is
invested into reproduction.
If the maturity argument is not supplied, then it is set to a sigmoidal
maturity ogive that changes from 0 to 1 at around the maturity size:
$${\tt maturity}(w) = \left[1+\left(\frac{w}{w_{mat}}\right)^{-U}\right]^{-1}.$$
(To avoid clutter, we are not showing the species index in the equations.)
The maturity weights are taken from the w_mat column of the
species_params data frame. Any missing maturity weights are set to 1/4 of the
asymptotic weight in the w_inf column.
The exponent \(U\) determines the steepness of the maturity ogive. By
default it is chosen as \(U = 10\), however this can be overridden by
including a column w_mat25 in the species parameter dataframe that
specifies the weight at which 25% of individuals are mature, which sets
\(U = \log(3) / \log(w_{mat} / w_{25}).\)
The sigmoidal function given above would strictly reach 1 only asymptotically.
Mizer instead sets the function equal to 1 already at the species'
maximum size, taken from the compulsory w_inf column in the
species_params data frame.
If the repro_prop argument is not supplied, it is set to the
allometric form
$${\tt repro\_prop}(w) = \left(\frac{w}{w_{inf}}\right)^{m-n}.$$
Here \(n\) is the scaling exponent of the energy income rate. Hence
the exponent \(m\) determines the scaling of the investment into
reproduction for mature individuals. By default it is chosen to be
\(m = 1\) so that the rate at which energy is invested into reproduction
scales linearly with the size. This default can be overridden by including a
column m in the species parameter dataframe. The asymptotic sizes
are taken from the compulsory w_inf column in the species_params
data frame.
The reproductive efficiency, i.e., the proportion of energy allocated to
reproduction that results in egg biomass, is set from the erepro
column in the species_params data frame. If that is not provided, the default
is set to 1 (which you will want to override). The offspring biomass divided
by the egg biomass gives the rate of egg production, returned by
getRDI.
The stock-recruitment relationship is an emergent phenomenon in mizer, with several sources of density dependence. Firstly, the amount of energy invested into reproduction depends on the energy income of the spawners, which is density-dependent due to competition for prey. Secondly, the proportion of larvae that grow up to recruitment size depends on the larval mortality, which depends on the density of predators, and on larval growth rate, which depends on density of prey.
Finally, the proportion of eggs that are viable and hatch to larvae can be density dependent. Somewhat misleadingly, mizer refers to this relationship between the number of eggs and the number of hatched larvae as the stock-recruitment relationship, even though it is only one part of the full stock-recruitment relationship. However it is the only part that can be set independently, while the other parts are already determined by the predation parameters and other model parameters. Thus in practice this part of the density dependence is used to encode all the density dependence that is not already included in the other two sources of density dependence.
To calculate the density-dependent rate of larvae production, mizer puts the
the density-independent rate of egg production through a "stock-recruitment"
function. The result is returned by getRDD. The name of the
stock-recruitment function is specified by the srr argument. The
default is the Beverton-Holt function srrBevertonHolt, which
requires an R_max column in the species_params data frame giving the
maximum egg production rate. If this column does not exist, it is initialised
to Inf, leading to no density-dependence. Other functions provided by
mizer are srrRicker and srrSheperd and you can
easily use these as models for writing your own functions.
Gears
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,
$$\mu_{f.i}(w) = \sum_g F_{g,i}(w),$$
where the fishing mortality \(F_{g,i}(w)\) imposed by gear \(g\) on
species \(i\) at size \(w\) is calculated as:
$$F_{g,i}(w) = S_{g,i}(w) Q_{g,i} E_{g},$$
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.
At the moment a species can only be selected by one fishing gear, although
each gear can select more than one species (this is a limitation with the
current package that will be developed in future releases). The gear
selecting each species can be specified in the gear column in the
species_params data frame. If no gear is specified, the default gear is
"knife_edge_gear".
Selectivity
The selectivity at size of each gear has a range between 0 (not selected at
that size) to 1 (fully selected at that size). It is given by a selectivity
function. The name of the selectivity function is given by the sel_func
column in the species parameters data frame. Some selectivity functions are
included in the package: knife_edge(), sigmoid_length(),
double_sigmoid_length(), and sigmoid_weight(). New functions can be defined
by the user. Each
gear has the same selectivity function for all the species it selects, but
the parameter values for each species may be different, e.g. the lengths of
species that a gear selects may be different.
Each selectivity function has a range of parameters. Values for these
parameters must be included as columns in the species parameters data.frame.
The names of the columns must exactly match the names of the corresponding
arguments of the selectivity function. For example, the default selectivity
function is knife_edge() which has sudden change of selectivity from 0 to 1
at a certain size. In its help page you can see that the knife_edge()
function has arguments w and knife_edge_size The first argument, w, is
size (the function calculates selectivity at size). All selectivity functions
must have w as the first argument. The values for the other arguments must
be found in the species parameters data.frame. So for the knife_edge()
function there should be a knife_edge_size column. Because knife_edge()
is the default selectivity function, the knife_edge_size argument has a
default value = w_mat.
Catchability
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.
Because currently mizer only allows one gear to select each species, the
catchability matrix \(Q_{g,i}\) reduces to a catchability vector
\(Q_{i}\), and this is given as a column catchability in the species
parameter data frame. If it is not specified, it defaults to 1.
Effort
The initial fishing effort is stored in the MizerParams object. If it is
not supplied, it is set to zero. The initial effort can be overruled when
the simulation is run with project(), where it is also possible to specify
an effort that varies through time.
By default, mizer uses a semichemostat model to describe the plankton
dynamics in each size class independently. This semichemostat dynamics is implemented
by the function plankton_semichemostat. You can change the
plankton dynamics by writing your own function, modelled on
plankton_semichemostat, and then passing the name of your
function in the plankton_dynamics argument.
The rr_pp argument is a vector specifying the intrinsic plankton
growth rate for each size class. If it is not supplied, then the intrinsic growth
rate \(r_p(w)\) at size \(w\)
is set to $$r_p(w) = r_p\, w^{n-1}.$$
The values of \(r_p\) and \(n\) are taken from the r_pp
and n arguments.
The cc_pp argument is a vector specifying the intrinsic plankton
carrying capacity for each size class. If it is not supplied, then the intrinsic carrying
capacity \(c_p(w)\) at size \(w\)
is set to $$c_p(w) = \kappa\, w^{-\lambda}$$
for all \(w\) less than w_pp_cutoff and zero for larger sizes.
The values of \(\kappa\) and \(\lambda\) are taken from the kappa
and lambda arguments.
project, MizerSim,
newCommunityParams, newTraitParams
Other functions for setting up models:
newCommunityParams(),
newSheldonParams(),
newTraitParams()