Title:	Dynamic Multi-Species Size Spectrum Modelling
Date:	2024-11-14
Type:	Package
Description:	A set of classes and methods to set up and run multi-species, trait based and community size spectrum ecological models, focused on the marine environment.
Maintainer:	Gustav Delius <gustav.delius@york.ac.uk>
Version:	2.5.3
License:	GPL-3
Imports:	assertthat, deSolve, dplyr, ggplot2 (≥ 3.4.0), ggrepel, grid, lubridate, methods, plotly, plyr, progress, Rcpp, reshape2, rlang, lifecycle
LinkingTo:	Rcpp
Depends:	R (≥ 3.1)
Suggests:	testthat (≥ 3.0.0), vdiffr, roxygen2, knitr, rmarkdown, pkgdown, covr, spelling
Collate:	'age_mat.R' 'helpers.R' 'MizerParams-class.R' 'MizerSim-class.R' 'reproduction.R' 'saveParams.R' 'species_params.R' 'setColours.R' 'setInteraction.R' 'setPredKernel.R' 'setSearchVolume.R' 'setMaxIntakeRate.R' 'setMetabolicRate.R' 'setMetadata.R' 'setExtMort.R' 'setExtEncounter.R' 'setReproduction.R' 'setResource.R' 'setFishing.R' 'setInitialValues.R' 'setBevertonHolt.R' 'upgrade.R' 'selectivity_funcs.R' 'pred_kernel_funcs.R' 'resource_dynamics.R' 'resource_semichemostat.R' 'resource_logistic.R' 'project.R' 'mizer-package.R' 'project_methods.R' 'rate_functions.R' 'summary_methods.R' 'plots.R' 'plotBiomassObservedVsModel.R' 'plotYieldObservedVsModel.R' 'animateSpectra.R' 'newMultispeciesParams.R' 'wrapper_functions.R' 'newSingleSpeciesParams.R' 'steady.R' 'extension.R' 'data.R' 'RcppExports.R' 'deprecated.R' 'get_initial_n.R' 'compareParams.R' 'customFunction.R' 'manipulate_species.R' 'calibrate.R' 'match.R' 'matchGrowth.R' 'steadySingleSpecies.R' 'defaults_edition.R' 'validSpeciesParams.R'
RoxygenNote:	7.3.2
Encoding:	UTF-8
LazyData:	true
URL:	https://sizespectrum.org/mizer/, https://github.com/sizespectrum/mizer
BugReports:	https://github.com/sizespectrum/mizer/issues
Language:	en-GB
RdMacros:	lifecycle
VignetteBuilder:	knitr
Config/testthat/edition:	3
NeedsCompilation:	yes
Packaged:	2024-10-16 22:10:15 UTC; gustav
Author:	Gustav Delius [cre, aut, cph], Finlay Scott [aut, cph], Julia Blanchard [aut, cph], Ken Andersen [aut, cph], Richard Southwell [ctb, cph]
Repository:	CRAN
Date/Publication:	2024-10-17 07:10:09 UTC

mizer: Multi-species size-based modelling in R

Description

The mizer package implements multi-species size-based modelling in R. It has been designed for modelling marine ecosystems.

Details

Using mizer is relatively simple. There are three main stages:

1. Setting the model parameters. This is done by creating an object of class MizerParams. This includes model parameters such as the life history parameters of each species, and the range of the size spectrum. There are several setup functions that help to create a MizerParams objects for particular types of models:

   • newSingleSpeciesParams()

   • newCommunityParams()

   • newTraitParams()

   • newMultispeciesParams()

2. Running a simulation. This is done by calling the project() function with the model parameters. This produces an object of MizerSim that contains the results of the simulation.

3. Exploring results. After a simulation has been run, the results can be explored using a range of plotting_functions, summary_functions and indicator_functions.

See the mizer website for full details of the principles behind mizer and how the package can be used to perform size-based modelling.

Author(s)

Maintainer: Gustav Delius gustav.delius@york.ac.uk (ORCID) [copyright holder]

Authors:

• Finlay Scott drfinlayscott@gmail.com [copyright holder]

• Julia Blanchard julia.blanchard@utas.edu.au (ORCID) [copyright holder]

• Ken Andersen kha@aqua.dtu.dk (ORCID) [copyright holder]

Other contributors:

• Richard Southwell richard.southwell@york.ac.uk [contributor, copyright holder]

See Also

Useful links:

• https://sizespectrum.org/mizer/

• https://github.com/sizespectrum/mizer

• Report bugs at https://github.com/sizespectrum/mizer/issues

Beverton Holt function to calculate density-dependent reproduction rate

Description

Takes the density-independent rates R_{di} of egg production (as calculated by getRDI()) and returns reduced, density-dependent reproduction rates R_{dd} given as

R_{dd} = R_{di} \frac{R_{max}}{R_{di} + R_{max}}

where R_{max} are the maximum possible reproduction rates that must be specified in a column in the species parameter dataframe. (All quantities in the above equation are species-specific but we dropped the species index for simplicity.)

Usage

BevertonHoltRDD(rdi, species_params, ...)

Arguments

Details

This is only one example of a density-dependence. You can write your own function based on this example, returning different density-dependent reproduction rates. Three other examples provided are RickerRDD(), SheperdRDD(), noRDD() and constantRDD(). For more explanation see setReproduction().

Value

Vector of density-dependent reproduction rates.

See Also

Other functions calculating density-dependent reproduction rate: RickerRDD(), SheperdRDD(), constantEggRDI(), constantRDD(), noRDD()

Alias for set_multispecies_model()

Description

An alias provided for backward compatibility with mizer version <= 1.0

Usage

MizerParams(
  species_params,
  interaction = matrix(1, nrow = nrow(species_params), ncol = nrow(species_params)),
  min_w_pp = 1e-10,
  min_w = 0.001,
  max_w = NULL,
  no_w = 100,
  n = 2/3,
  q = 0.8,
  f0 = 0.6,
  kappa = 1e+11,
  lambda = 2 + q - n,
  r_pp = 10,
  ...
)

Arguments

Value

A MizerParams object

A class to hold the parameters for a size based model.

Description

Although it is possible to build a MizerParams object by hand it is not recommended and several constructors are available. Dynamic simulations are performed using project() function on objects of this class. As a user you should never need to access the slots inside a MizerParams object directly.

Details

The MizerParams class is fairly complex with a large number of slots, many of which are multidimensional arrays. The dimensions of these arrays is strictly enforced so that MizerParams objects are consistent in terms of number of species and number of size classes.

The MizerParams class does not hold any dynamic information, e.g. abundances or harvest effort through time. These are held in MizerSim objects.

Slots

metadata

A list with metadata information. See setMetadata().

mizer_version

The package version of mizer (as returned by packageVersion("mizer")) that created or upgraded the model.

extensions

A named vector of strings where each name is the name of and extension package needed to run the model and each value is a string giving the information that the remotes package needs to install the correct version of the extension package, see https://remotes.r-lib.org/.

time_created

A POSIXct date-time object with the creation time.

time_modified

A POSIXct date-time object with the last modified time.

w

The size grid for the fish part of the spectrum. An increasing vector of weights (in grams) running from the smallest egg size to the largest maximum size.

dw

The widths (in grams) of the size bins

w_full

The size grid for the full size range including the resource spectrum. An increasing vector of weights (in grams) running from the smallest resource size to the largest maximum size of fish. The last entries of the vector have to be equal to the content of the w slot.

dw_full

The width of the size bins for the full spectrum. The last entries have to be equal to the content of the dw slot.

w_min_idx

A vector holding the index of the weight of the egg size of each species

maturity

An array (species x size) that holds the proportion of individuals of each species at size that are mature. This enters in the calculation of the spawning stock biomass with getSSB(). Set with setReproduction().

psi

An array (species x size) that holds the allocation to reproduction for each species at size, \psi_i(w). Changed with setReproduction().

intake_max

An array (species x size) that holds the maximum intake for each species at size. Changed with setMaxIntakeRate().

search_vol

An array (species x size) that holds the search volume for each species at size. Changed with setSearchVolume().

metab

An array (species x size) that holds the metabolism for each species at size. Changed with setMetabolicRate().

mu_b

An array (species x size) that holds the external mortality rate \mu_{ext.i}(w). Changed with setExtMort().

ext_encounter

An array (species x size) that holds the external encounter rate E_{ext.i}(w). Changed with setExtEncounter().

pred_kernel

An array (species x predator size x prey size) that holds the predation coefficient of each predator at size on each prey size. If this is NA then the following two slots will be used. Changed with setPredKernel().

ft_pred_kernel_e

An array (species x log of predator/prey size ratio) that holds the Fourier transform of the feeding kernel in a form appropriate for evaluating the encounter rate integral. If this is NA then the pred_kernel will be used to calculate the available energy integral. Changed with setPredKernel().

ft_pred_kernel_p

An array (species x log of predator/prey size ratio) that holds the Fourier transform of the feeding kernel in a form appropriate for evaluating the predation mortality integral. If this is NA then the pred_kernel will be used to calculate the integral. Changed with setPredKernel().

rr_pp

A vector the same length as the w_full slot. The size specific growth rate of the resource spectrum.

cc_pp

A vector the same length as the w_full slot. The size specific carrying capacity of the resource spectrum.

resource_dynamics

Name of the function for projecting the resource abundance density by one timestep.

other_dynamics

A named list of functions for projecting the values of other dynamical components of the ecosystem that may be modelled by a mizer extensions you have installed. The names of the list entries are the names of those components.

other_encounter

A named list of functions for calculating the contribution to the encounter rate from each other dynamical component.

other_mort

A named list of functions for calculating the contribution to the mortality rate from each other dynamical components.

other_params

A list containing the parameters needed by any mizer extensions you may have installed to model other dynamical components of the ecosystem.

rates_funcs

A named list with the names of the functions that should be used to calculate the rates needed by project(). By default this will be set to the names of the built-in rate functions.

sc

The community abundance of the scaling community

species_params

A data.frame to hold the species specific parameters. See species_params() for details.

given_species_params

A data.frame to hold the species parameters that were given explicitly rather than obtained by default calculations.

gear_params

Data frame with parameters for gear selectivity. See setFishing() for details.

interaction

The species specific interaction matrix, \theta_{ij}. Changed with setInteraction().

selectivity

An array (gear x species x w) that holds the selectivity of each gear for species and size, S_{g,i,w}. Changed with setFishing().

catchability

An array (gear x species) that holds the catchability of each species by each gear, Q_{g,i}. Changed with setFishing().

initial_effort

A vector containing the initial fishing effort for each gear. Changed with setFishing().

initial_n

An array (species x size) that holds the initial abundance of each species at each weight.

initial_n_pp

A vector the same length as the w_full slot that describes the initial resource abundance at each weight.

initial_n_other

A list with the initial abundances of all other ecosystem components. Has length zero if there are no other components.

resource_params

List with parameters for resource.

A

Abundance multipliers.

linecolour

A named vector of colour values, named by species. Used to give consistent colours in plots.

linetype

A named vector of linetypes, named by species. Used to give consistent line types in plots.

ft_mask

An array (species x w_full) with zeros for weights larger than the maximum weight of each species. Used to efficiently minimize wrap-around errors in Fourier transform calculations.

See Also

project() MizerSim() emptyParams() newMultispeciesParams() newCommunityParams() newTraitParams()

Constructor for the MizerSim class

Description

A constructor for the MizerSim class. This is used by project() to create MizerSim objects of the right dimensions. It is not necessary for users to use this constructor.

Usage

MizerSim(params, t_dimnames = NA, t_max = 100, t_save = 1)

Arguments

Value

An object of type MizerSim

A class to hold the results of a simulation

Description

A class that holds the results of projecting a MizerParams object through time using project().

Details

A new MizerSim object can be created with the MizerSim() constructor, but you will never have to do that because the object is created automatically by project() when needed.

As a user you should never have to access the slots of a MizerSim object directly. Instead there are a range of functions to extract the information. N() and NResource() return arrays with the saved abundances of the species and the resource population at size respectively. getEffort() returns the fishing effort of each gear through time. getTimes() returns the vector of times at which simulation results were stored and idxFinalT() returns the index with which to access specifically the value at the final time in the arrays returned by the other functions. getParams() returns the MizerParams object that was passed to project(). There are also several summary_functions and plotting_functions available to explore the contents of a MizerSim object.

The arrays all have named dimensions. The names of the time dimension denote the time in years. The names of the w dimension are weights in grams rounded to three significant figures. The names of the sp dimension are the same as the species name in the order specified in the species_params data frame. The names of the gear dimension are the names of the gears, in the same order as specified when setting up the MizerParams object.

Extensions of mizer can use the n_other slot to store the abundances of other ecosystem components and these extensions should provide their own functions for accessing that information.

The MizerSim class has changed since previous versions of mizer. To use a MizerSim object created by a previous version, you need to upgrade it with upgradeSim().

Slots

params

An object of type MizerParams.

n

Three-dimensional array (time x species x size) that stores the projected community number densities.

n_pp

An array (time x size) that stores the projected resource number densities.

n_other

A list array (time x component) that stores the projected values for other ecosystem components.

effort

An array (time x gear) that stores the fishing effort by time and gear.

Time series of size spectra

Description

Fetch the simulation results for the size spectra over time.

Usage

N(sim)

NResource(sim)

Arguments

Value

For N(): A three-dimensional array (time x species x size) with the number density of consumers

For NResource(): An array (time x size) with the number density of resource

Examples

str(N(NS_sim))
str(NResource(NS_sim))

Time series of other components

Description

Fetch the simulation results for other components over time.

Usage

NOther(sim)

Arguments

Value

A list array (time x component) that stores the projected values for other ecosystem components.

Example interaction matrix for the North Sea example

Description

The interaction coefficient between predator and prey species in the North Sea.

Usage

NS_interaction

Format

A 12 x 12 matrix.

Source

Blanchard et al.

Examples

params <- MizerParams(NS_species_params_gears,
                      interaction = NS_interaction)

Example MizerParams object for the North Sea example

Description

A MizerParams object created from the NS_species_params_gears species parameters and the inter interaction matrix together with an initial condition corresponding to the steady state obtained from fishing with an effort effort = c(Industrial = 0, Pelagic = 1, Beam = 0.5, Otter = 0.5).

Usage

NS_params

Format

A MizerParams object

Source

Blanchard et al.

See Also

Other example parameter objects: NS_sim

Examples

sim = project(NS_params, effort = c(Industrial = 0, Pelagic = 1, 
                                    Beam = 0.5, Otter = 0.5))
plot(sim)

Example MizerSim object for the North Sea example

Description

A MizerSim object containing a simulation with historical fishing mortalities from the North Sea, as created in the tutorial "A Multi-Species Model of the North Sea".

Usage

NS_sim

Format

A MizerSim object

Source

https://sizespectrum.org/mizer/articles/a_multispecies_model_of_the_north_sea.html

See Also

Other example parameter objects: NS_params

Examples

plotBiomass(NS_sim)

Example species parameter set based on the North Sea

Description

This data set is based on species in the North Sea (Blanchard et al.). It is a data.frame that contains all the necessary information to be used by the MizerParams() constructor. As there is no gear column, each species is assumed to be fished by a separate gear.

Usage

NS_species_params

Format

A data frame with 12 rows and 7 columns. Each row is a species.

species

Name of the species

w_max

Maximum size.

w_mat

Size at maturity

beta

Size preference ratio

sigma

Width of the size-preference

R_max

Maximum reproduction rate

k_vb

The von Bertalanffy k parameter

w_inf

The von Bertalanffy asymptotic size

Source

Blanchard et al.

Examples

params <- MizerParams(NS_species_params)

Example species parameter set based on the North Sea with different gears

Description

This data set is based on species in the North Sea (Blanchard et al.). It is similar to the data set NS_species_params except that this one has an additional column specifying the fishing gear that operates on each species.

Usage

NS_species_params_gears

Format

A data frame with 12 rows and 8 columns. Each row is a species.

species

Name of the species

w_max

Maximum size.

w_mat

Size at maturity

beta

Size preference ratio

sigma

Width of the size-preference

R_max

Maximum reproduction rate

k_vb

The von Bertalanffy k parameter

w_inf

The von Bertalanffy asymptotic size

gear

Name of the fishing gear

Source

Blanchard et al.

Examples

params <- MizerParams(NS_species_params_gears)

Ricker function to calculate density-dependent reproduction rate

Description

Takes the density-independent rates R_{di} of egg production and returns reduced, density-dependent rates R_{dd} given as

R_{dd} = R_{di} \exp(- b R_{di})

Usage

RickerRDD(rdi, species_params, ...)

Arguments

Value

Vector of density-dependent reproduction rates.

See Also

Other functions calculating density-dependent reproduction rate: BevertonHoltRDD(), SheperdRDD(), constantEggRDI(), constantRDD(), noRDD()

Sheperd function to calculate density-dependent reproduction rate

Description

Takes the density-independent rates R_{di} of egg production and returns reduced, density-dependent rates R_{dd} given as

R_{dd} = \frac{R_{di}}{1+(b\ R_{di})^c}

Usage

SheperdRDD(rdi, species_params, ...)

Arguments

Details

With b = 1/R_{max} and c = 1 this reduces to the Beverton-Holt reproduction rate, see BevertonHoltRDD().

Value

Vector of density-dependent reproduction rates.

See Also

Other functions calculating density-dependent reproduction rate: BevertonHoltRDD(), RickerRDD(), constantEggRDI(), constantRDD(), noRDD()

Add new species

Description

Takes a MizerParams object and adds additional species with given parameters to the ecosystem. It sets the initial values for these new species to their steady-state solution in the given initial state of the existing ecosystem. This will be close to the true steady state if the abundances of the new species are sufficiently low. Hence the abundances of the new species are set so that they are at most 1/100th of the resource power law. Their reproductive efficiencies are set so as to keep them at that low level.

Usage

addSpecies(
  params,
  species_params,
  gear_params = data.frame(),
  initial_effort,
  interaction
)

Arguments

Details

The resulting MizerParams object will use the same size grid where possible, but if one of the new species needs a larger range of w (either because a new species has an egg size smaller than those of existing species or a maximum size larger than those of existing species) then the grid will be expanded and all arrays will be enlarged accordingly.

If any of the rate arrays of the existing species had been set by the user to values other than those calculated as default from the species parameters, then these will be preserved. Only the rates for the new species will be calculated from their species parameters.

After adding the new species, the background species are not retuned and the system is not run to steady state. This could be done with steady(). The new species will have a reproduction level of 1/4, this can then be changed with setBevertonHolt()

Value

An object of type MizerParams

See Also

removeSpecies()

Examples

params <- newTraitParams()
species_params <- data.frame(
    species = "Mullet",
    w_max = 173,
    w_mat = 15,
    beta = 283,
    sigma = 1.8,
    h = 30,
    a = 0.0085,
    b = 3.11
)
params <- addSpecies(params, species_params)
plotSpectra(params)

Calculate age at maturity

Description

Uses the growth rate and the size at maturity to calculate the age at maturity

Usage

age_mat(params)

Arguments

Details

Using that by definition of the growth rate g(w) = dw/dt we have that

\mathrm{age_{mat}} = \int_0^{w_{mat}.}\frac{dw}{g(w)}

Value

A named vector. The names are the species names and the values are the ages at maturity.

Examples

age_mat(NS_params)

Calculate age at maturity from von Bertalanffy growth parameters

Description

This is not a good way to determine the age at maturity because the von Bertalanffy growth curve is not reliable for larvae and juveniles. However this was used in previous versions of mizer and is supplied for backwards compatibility.

Usage

age_mat_vB(object)

Arguments

Details

Uses the age at maturity that is implied by the von Bertalanffy growth curve specified by the w_inf, k_vb, t0, a and b parameters in the species_params data frame.

If any of k_vb is missing for a species, the function returns NA for that species. Default values of b = 3 and t0 = 0 are used if these are missing. If w_inf is missing, w_max is used instead.

Value

A named vector. The names are the species names and the values are the ages at maturity.

Animation of the abundance spectra

Description

Usage

animateSpectra(
  sim,
  species = NULL,
  time_range,
  wlim = c(NA, NA),
  ylim = c(NA, NA),
  power = 1,
  total = FALSE,
  resource = TRUE
)

Arguments

Value

A plotly object

See Also

Other plotting functions: plot,MizerParams,missing-method, plot,MizerSim,missing-method, plotBiomass(), plotDiet(), plotFMort(), plotFeedingLevel(), plotGrowthCurves(), plotPredMort(), plotSpectra(), plotYield(), plotYieldGear(), plotting_functions

Examples

animateSpectra(NS_sim, power = 2, wlim = c(0.1, NA), time_range = 1997:2007)

Box predation kernel

Description

A predation kernel where the predator/prey mass ratio is uniformly distributed on an interval.

Usage

box_pred_kernel(ppmr, ppmr_min, ppmr_max)

Arguments

Details

Writing the predator mass as w and the prey mass as w_p, the feeding kernel is 1 if w/w_p is between ppmr_min and ppmr_max and zero otherwise. The parameters need to be given in the species parameter dataframe in the columns ppmr_min and ppmr_max.

Value

A vector giving the value of the predation kernel at each of the predator/prey mass ratios in the ppmr argument.

See Also

setPredKernel()

Other predation kernel: lognormal_pred_kernel(), power_law_pred_kernel(), truncated_lognormal_pred_kernel()

Examples

params <- NS_params
# Set all required paramters before changing kernel type
species_params(params)$ppmr_max <- 4000
species_params(params)$ppmr_min <- 200
species_params(params)$pred_kernel_type <- "box"
plot(w_full(params), getPredKernel(params)["Cod", 10, ], type="l", log="x")

Calculate selectivity from gear parameters

Description

This function calculates the selectivity for each gear, species and size from the gear parameters. It is called by setFishing() when the selectivity is not set by the user.

Usage

calc_selectivity(params)

Arguments

Value

An array (gear x species x size) with the selectivity values

Examples

params <- NS_params
str(calc_selectivity(params))
calc_selectivity(params)["Pelagic", "Herring", ]

Calibrate the model scale to match total observed biomass

Description

Given a MizerParams object params for which biomass observations are available for at least some species via the biomass_observed column in the species_params data frame, this function returns an updated MizerParams object which is rescaled with scaleModel() so that the total biomass in the model agrees with the total observed biomass.

Usage

calibrateBiomass(params)

Arguments

Details

Biomass observations usually only include individuals above a certain size. This size should be specified in a biomass_cutoff column of the species parameter data frame. If this is missing, it is assumed that all sizes are included in the observed biomass, i.e., it includes larval biomass.

After using this function the total biomass in the model will match the total biomass, summed over all species. However the biomasses of the individual species will not match observations yet, with some species having biomasses that are too high and others too low. So after this function you may want to use matchBiomasses(). This is described in the blog post at https://bit.ly/2YqXESV.

If you have observations of the yearly yield instead of biomasses, you can use calibrateYield() instead of this function.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$biomass_observed <- 
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
species_params(params)$biomass_cutoff <- 10
params2 <- calibrateBiomass(params)
plotBiomassObservedVsModel(params2)

Calibrate the model scale to match total observed number

Description

Given a MizerParams object params for which number observations are available for at least some species via the number_observed column in the species_params data frame, this function returns an updated MizerParams object which is rescaled with scaleModel() so that the total number in the model agrees with the total observed number.

Usage

calibrateNumber(params)

Arguments

Details

Number observations usually only include individuals above a certain size. This size should be specified in a number_cutoff column of the species parameter data frame. If this is missing, it is assumed that all sizes are included in the observed number, i.e., it includes larval number.

After using this function the total number in the model will match the total number, summed over all species. However the numbers of the individual species will not match observations yet, with some species having numbers that are too high and others too low. So after this function you may want to use matchNumbers(). This is described in the blog post at https://bit.ly/2YqXESV.

If you have observations of the yearly yield instead of numbers, you can use calibrateYield() instead of this function.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$number_observed <-
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
species_params(params)$number_cutoff <- 10
params2 <- calibrateNumber(params)

Calibrate the model scale to match total observed yield

Description

Usage

calibrateYield(params)

Arguments

Details

This function has been deprecated and will be removed in the future unless you have a use case for it. If you do have a use case for it, please let the developers know by creating an issue at https://github.com/sizespectrum/mizer/issues.

Given a MizerParams object params for which yield observations are available for at least some species via the yield_observed column in the species_params data frame, this function returns an updated MizerParams object which is rescaled with scaleModel() so that the total yield in the model agrees with the total observed yield.

After using this function the total yield in the model will match the total observed yield, summed over all species. However the yields of the individual species will not match observations yet, with some species having yields that are too high and others too low. So after this function you may want to use matchYields().

If you have observations of species biomasses instead of yields, you can use calibrateBiomass() instead of this function.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$yield_observed <-
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
gear_params(params)$catchability <-
    c(1.3, 0.065, 0.31, 0.18, 0.98, 0.24, 0.37, 0.46, 0.18, 0.30, 0.27, 0.39)
params2 <- calibrateYield(params)
plotYieldObservedVsModel(params2)

Compare two MizerParams objects and print out differences

Description

Usage

compareParams(params1, params2)

Arguments

Value

String describing the differences

Examples

params1 <- NS_params
params2 <- params1
species_params(params2)$w_mat[1] <- 10
compareParams(params1, params2)

Alias for validSpeciesParams()

Description

An alias provided for backward compatibility with mizer version <= 2.5.2

Usage

completeSpeciesParams(species_params)

Arguments

Details

validGivenSpeciesParams() checks the validity of the given species parameter It throws an error if

• the species column does not exist or contains duplicates

• the maximum size is not specified for all species

If a weight-based parameter is missing but the corresponding length-based parameter is given, as well as the a and b parameters for length-weight conversion, then the weight-based parameters are added. If both length and weight are given, then weight is used and a warning is issued if the two are inconsistent.

If a w_inf column is given but no w_max then the value from w_inf is used. This is for backwards compatibility. But note that the von Bertalanffy parameter w_inf is not the maximum size of the largest individual, but the asymptotic size of an average individual.

Some inconsistencies in the size parameters are resolved as follows:

• Any w_mat that is not smaller than w_max is set to w_max / 4.

• Any w_mat25 that is not smaller than w_mat is set to NA.

• Any w_min that is not smaller than w_mat is set to 0.001 or w_mat /10, whichever is smaller.

• Any w_repro_max that is not larger than w_mat is set to 4 * w_mat.

The row names of the returned data frame will be the species names. If species_params was provided as a tibble it is converted back to an ordinary data frame.

The function tests for some typical misspellings of parameter names, like wrong capitalisation or missing underscores and issues a warning if it detects such a name.

validSpeciesParams() first calls validateGivenSpeciesParams() but then goes further by adding default values for species parameters that were not provided. The function sets default values if any of the following species parameters are missing or NA:

• w_repro_max is set to w_max

• w_mat is set to w_max/4

• w_min is set to 0.001

• alpha is set to 0.6

• interaction_resource is set to 1

• n is set to 3/4

Note that the species parameters returned by these functions are not guaranteed to produce a viable model. More checks of the parameters are performed by the individual rate-setting functions (see setParams() for the list of these functions).

Value

For validSpeciesParams(): A valid species parameter data frame with additional parameters with default values.

For validGivenSpeciesParams(): A valid species parameter data frame without additional parameters.

See Also

species_params(), validGearParams(), validParams(), validSim()

Choose egg production to keep egg density constant

Description

The new egg production is set to compensate for the loss of individuals from the smallest size class through growth and mortality. The result should not be modified by density dependence, so this should be used together with the noRDD() function, see example.

Usage

constantEggRDI(params, n, e_growth, mort, ...)

Arguments

Value

Vector with the value for each species

See Also

Other functions calculating density-dependent reproduction rate: BevertonHoltRDD(), RickerRDD(), SheperdRDD(), constantRDD(), noRDD()

Examples

# choose an example params object
params <- NS_params
# We set the reproduction rate functions
params <- setRateFunction(params, "RDI", "constantEggRDI")
params <- setRateFunction(params, "RDD", "noRDD")
# Now the egg density should stay fixed no matter how we fish
sim <- project(params, effort = 10, progress_bar = FALSE)
# To check that indeed the egg densities have not changed, we first construct
# the indices for addressing the egg densities
no_sp <- nrow(params@species_params)
idx <- (params@w_min_idx - 1) * no_sp + (1:no_sp)
# Now we can check equality between egg densities at the start and the end
all.equal(finalN(sim)[idx], initialN(params)[idx])

Give constant reproduction rate

Description

Simply returns the value from species_params$constant_reproduction.

Usage

constantRDD(rdi, species_params, ...)

Arguments

Value

Vector species_params$constant_reproduction

See Also

Other functions calculating density-dependent reproduction rate: BevertonHoltRDD(), RickerRDD(), SheperdRDD(), constantEggRDI(), noRDD()

Helper function to keep other components constant

Description

Helper function to keep other components constant

Usage

constant_other(params, n_other, component, ...)

Arguments

Value

The current value of the component

Replace a mizer function with a custom version

Description

This function allows you to make arbitrary changes to how mizer works by allowing you to replace any mizer function with your own version. You should do this only as a last resort, when you find that you can not use the standard mizer extension mechanism to achieve your goal.

Usage

customFunction(name, fun)

Arguments

Details

If the function you need to overwrite is one of the mizer rate functions, then you should use setRateFunction() instead of this function. Similarly you should use ⁠resource_dynamics()<-⁠ to change the resource dynamics and setReproduction() to change the density-dependence in reproduction. You should also investigate whether you can achieve your goal by introducing additional ecosystem components with setComponent().

If you find that your goal really does require you to overwrite a mizer function, please also create an issue on the mizer issue tracker at https://github.com/sizespectrum/mizer/issues to describe your goal, because it will be interesting to the mizer community and may motivate future improvements to the mizer functionality.

Note that customFunction() only overwrites the function used by the mizer code. It does not overwrite the function that is exported by mizer. This will become clear when you run the code in the Examples section.

This function does not in any way check that your replacement function is compatible with mizer. Calling this function can totally break mizer. However you can always undo the effect by reloading mizer with

detach(package:mizer, unload = TRUE)
library(mizer)

Value

No return value, called for side effects

Examples

## Not run: 
fake_project <- function(...) "Fake"
customFunction("project", fake_project)
mizer::project(NS_params) # This will print "Fake"
project(NS_params) # This will still use the old project() function
# To undo the effect:
customFunction("project", project)
mizer::project(NS_params) # This will again use the old project()

## End(Not run)

Set defaults for predation kernel parameters

Description

If the predation kernel type has not been specified for a species, then it is set to "lognormal" and the default values are set for the parameters beta and sigma.

Usage

default_pred_kernel_params(object)

Arguments

Value

The object with updated columns in the species params data frame.

Default editions

Description

Function to set and get which edition of default choices is being used.

Usage

defaults_edition(edition = NULL)

Arguments

Details

The mizer functions for creating new models make a lot of choices for default values for parameters that are not provided by the user. Sometimes we find better ways to choose the defaults and update mizer accordingly. When we do this, we will increase the edition number.

If you call defaults_edition() without an argument it returns the currently active edition. Otherwise it sets the active edition to the given value.

Users who want their existing code for creating models not to change behaviour when run with future versions of mizer should explicitly set the desired defaults edition at the top of their code.

The most recent edition is edition 2. It will become the default in the next release. The current default is edition 1. The following defaults are changed in edition 2:

• catchability = 0.3 instead of 1

• initial effort = 1 instead of 0

Value

The current edition number.

Check whether two objects are different

Description

Check whether two objects are numerically different, ignoring all attributes.

Usage

different(a, b)

Arguments

Details

We use this helper function in particular to see if a new value for a slot in MizerParams is different from the existing value in order to give the appropriate messages.

Value

TRUE or FALSE

Measure distance between current and previous state in terms of RDI

Description

This function can be used in projectToSteady() to decide when sufficient convergence to steady state has been achieved.

Usage

distanceMaxRelRDI(params, current, previous)

Arguments

Value

The largest absolute relative change in rdi: max(abs((current_rdi - previous_rdi) / previous_rdi))

See Also

Other distance functions: distanceSSLogN()

Measure distance between current and previous state in terms of fish abundances

Description

Calculates the sum squared difference between log(N) in current and previous state. This function can be used in projectToSteady() to decide when sufficient convergence to steady state has been achieved.

Usage

distanceSSLogN(params, current, previous)

Arguments

Value

The sum of squares of the difference in the logs of the (nonzero) fish abundances n: sum((log(current$n) - log(previous$n))^2)

See Also

Other distance functions: distanceMaxRelRDI()

Length based double-sigmoid selectivity function

Description

A hump-shaped selectivity function with a sigmoidal rise and an independent sigmoidal drop-off. This drop-off is what distinguishes this from the function sigmoid_length() and it is intended to model the escape of large individuals from the fishing gear.

Usage

double_sigmoid_length(w, l25, l50, l50_right, l25_right, species_params, ...)

Arguments

Details

The selectivity is obtained as the product of two sigmoidal curves, one rising and one dropping. The sigmoidal rise is based on the two parameters l25 and l50 which determine the length at which 25% and 50% of the stock is selected respectively. The sigmoidal drop-off is based on the two parameters l50_right and l25_right which determine the length at which the selectivity curve has dropped back to 50% and 25% respectively. The selectivity is given by the function

S(l) = \frac{1}{1 + \exp\left(\log(3)\frac{l50 -l}{l50 - l25}\right)}\frac{1}{1 + \exp\left(\log(3)\frac{l50_{right} -l}{l50_{right} - l25_{right}}\right)}

As the size-based model is weight based, and this selectivity function is length based, it uses the length-weight parameters a and b to convert between length and weight.

l = \left(\frac{w}{a}\right)^{1/b}

Value

Vector of selectivities at the given sizes.

See Also

gear_params() for setting the selectivity parameters.

Other selectivity functions: knife_edge(), sigmoid_length(), sigmoid_weight()

Create empty MizerParams object of the right size

Description

An internal function. Sets up a valid MizerParams object with all the slots initialised and given dimension names, but with some slots left empty. This function is to be used by other functions to set up full parameter objects.

Usage

emptyParams(
  species_params,
  gear_params = data.frame(),
  no_w = 100,
  min_w = 0.001,
  max_w = NA,
  min_w_pp = 1e-12
)

Arguments

Value

An empty but valid MizerParams object

Size grid

A size grid is created so that the log-sizes are equally spaced. The spacing is chosen so that there will be no_w fish size bins, with the smallest starting at min_w and the largest starting at max_w. For the resource spectrum there is a larger set of bins containing additional bins below min_w, with the same log size. The number of extra bins is such that min_w_pp comes to lie within the smallest bin.

Changes to species params

The species_params slot of the returned MizerParams object may differ from the data frame supplied as argument to this function because default values are set for missing parameters.

See Also

See newMultispeciesParams() for a function that fills the slots left empty by this function.

Size spectra at end of simulation

Description

Size spectra at end of simulation

Usage

finalN(sim)

finalNResource(sim)

idxFinalT(sim)

Arguments

Value

For finalN(): An array (species x size) holding the consumer number densities at the end of the simulation

For finalNResource(): A vector holding the resource number densities at the end of the simulation for all size classes

For idxFinalT(): An integer giving the index for extracting the results for the final time step

Examples

str(finalN(NS_sim))

# This could also be obtained using `N()` and `idxFinalT()`
identical(N(NS_sim)[idxFinalT(NS_sim), , ], finalN(NS_sim))
str(finalNResource(NS_sim))
idx <- idxFinalT(NS_sim)
idx
# This coincides with
length(getTimes(NS_sim))
# and corresponds to the final time
getTimes(NS_sim)[idx]
# We can use this index to extract the result at the final time
identical(N(NS_sim)[idx, , ], finalN(NS_sim))
identical(NResource(NS_sim)[idx, ], finalNResource(NS_sim))

Values of other ecosystem components at end of simulation

Description

Values of other ecosystem components at end of simulation

Usage

finalNOther(sim)

Arguments

Value

A named list holding the values of other ecosystem components at the end of the simulation

Gear parameters

Description

These functions allow you to get or set the gear parameters stored in a MizerParams object. These are used by setFishing() to set up the selectivity and catchability and thus together with the fishing effort determine the fishing mortality.

Usage

gear_params(params)

gear_params(params) <- value

Arguments

Details

The gear_params data has one row for each gear-species pair and one column for each parameter that determines how that gear interacts with that species. The columns are:

• species The name of the species

• gear The name of the gear

• catchability A number specifying how strongly this gear selects this species.

• sel_func The name of the function that calculates the selectivity curve.

• One column for each selectivity parameter needed by the selectivity functions.

For the details see setFishing().

There can optionally also be a column yield_observed that allows you to specify for each gear and species the total annual fisheries yield.

The fishing effort, which is also needed to determine the fishing mortality exerted by a gear is not set via the gear_params data frame but is set with initial_effort() or is specified when calling project().

If you change a gear parameter, this will be used to recalculate the selectivity and catchability arrays by calling setFishing(), unless you have previously set these by hand.

⁠gear_params<-⁠ automatically sets the row names to contain the species name and the gear name, separated by a comma and a space. The last example below illustrates how this facilitates changing an individual gear parameter.

Value

Data frame with gear parameters

See Also

validGearParams()

Other functions for setting parameters: setExtEncounter(), setExtMort(), setFishing(), setInitialValues(), setInteraction(), setMaxIntakeRate(), setMetabolicRate(), setParams(), setPredKernel(), setReproduction(), setSearchVolume(), species_params()

Examples

params <- NS_params

# gears set up in example
gear_params(params)

# setting totally different gears
gear_params(params) <- data.frame(
    gear = c("gear1", "gear2", "gear1"),
    species = c("Cod", "Cod", "Haddock"),
    catchability = c(0.5, 2, 1),
    sel_fun = c("sigmoid_weight", "knife_edge", "sigmoid_weight"),
    sigmoidal_weight = c(1000, NA, 800),
    sigmoidal_sigma = c(100, NA, 100),
    knife_edge_size = c(NA, 1000, NA)
    )
gear_params(params)

# changing an individual entry
gear_params(params)["Cod, gear1", "catchability"] <- 0.8

Calculate the total biomass of each species within a size range at each time step.

Description

Calculates the total biomass through time within user defined size limits. The default option is to use the whole size range. You can specify minimum and maximum weight or length range for the species. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used).

Usage

getBiomass(object, ...)

Arguments

Value

If called with a MizerParams object, a vector with the biomass in grams for each species in the model. If called with a MizerSim object, an array (time x species) containing the biomass in grams at each time step for all species.

See Also

Other summary functions: getDiet(), getGrowthCurves(), getN(), getSSB(), getYield(), getYieldGear()

Examples

biomass <- getBiomass(NS_sim)
biomass["1972", "Herring"]
biomass <- getBiomass(NS_sim, min_w = 10, max_w = 1000)
biomass["1972", "Herring"]

Calculate the slope of the community abundance

Description

Calculates the slope of the community abundance through time by performing a linear regression on the logged total numerical abundance at weight and logged weights (natural logs, not log to base 10, are used). You can specify minimum and maximum weight or length range for the species. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used). You can also specify the species to be used in the calculation.

Usage

getCommunitySlope(sim, species = NULL, biomass = TRUE, ...)

Arguments

Value

A data.frame with four columns: time step, slope, intercept and the coefficient of determination R^2.

See Also

Other functions for calculating indicators: getMeanMaxWeight(), getMeanWeight(), getProportionOfLargeFish()

Examples

# Slope based on biomass, using all species and sizes
slope_biomass <- getCommunitySlope(NS_sim)
slope_biomass[1, ] # in 1976
slope_biomass[idxFinalT(NS_sim), ] # in 2010

# Slope based on numbers, using all species and sizes
slope_numbers <- getCommunitySlope(NS_sim, biomass = FALSE)
slope_numbers[1, ] # in 1976

# Slope based on biomass, using all species and sizes between 10g and 1000g
slope_biomass <- getCommunitySlope(NS_sim, min_w = 10, max_w = 1000)
slope_biomass[1, ] # in 1976

# Slope based on biomass, using only demersal species and 
# sizes between 10g and 1000g
dem_species <- c("Dab","Whiting", "Sole", "Gurnard", "Plaice",
                 "Haddock", "Cod", "Saithe")
slope_biomass <- getCommunitySlope(NS_sim, species = dem_species, 
                                   min_w = 10, max_w = 1000)
slope_biomass[1, ] # in 1976

Get information about other ecosystem components

Description

Get information about other ecosystem components

Usage

getComponent(params, component)

Arguments

Value

A list with the entries initial_value, dynamics_fun, encounter_fun, mort_fun, component_params for the requested component. If the requested component does not exist, NULL is returned. If no component argument is given, then a list of lists for all components is returned.

Get critical feeding level

Description

The critical feeding level is the feeding level at which the food intake is just high enough to cover the metabolic costs, with nothing left over for growth or reproduction.

Usage

getCriticalFeedingLevel(params)

Arguments

Value

A matrix (species x size) with the critical feeding level

Examples

str(getFeedingLevel(NS_params))

Get diet of predator at size, resolved by prey species

Description

Calculates the rate at which a predator of a particular species and size consumes biomass of each prey species, resource, and other components of the ecosystem. Returns either the rates in grams/year or the proportion of the total consumption rate.

Usage

getDiet(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  proportion = TRUE
)

Arguments

Details

The rates D_{ij}(w) at which a predator of species i and size w consumes biomass from prey species j are calculated from the predation kernel \phi_i(w, w_p), the search volume \gamma_i(w), the feeding level f_i(w), the species interaction matrix \theta_{ij} and the prey abundance density N_j(w_p):

D_{ij}(w, w_p) = (1-f_i(w)) \gamma_i(w) \theta_{ij} \int N_j(w_p) \phi_i(w, w_p) w_p dw_p.

The prey index j runs over all species and the resource.

Extra columns are added for the external encounter rate and for any extra ecosystem components in your model for which you have defined an encounter rate function. These encounter rates are multiplied by 1-f_i(w) to give the rate of consumption of biomass from these extra components.

This function performs the same integration as getEncounter() but does not aggregate over prey species, and multiplies by 1-f_i(w) to get the consumed biomass rather than the available biomass. Outside the range of sizes for a predator species the returned rate is zero.

Value

An array (predator species x predator size x (prey species + resource + other components). Dimnames are "prey", "w", and "predator".

See Also

plotDiet()

Other summary functions: getBiomass(), getGrowthCurves(), getN(), getSSB(), getYield(), getYieldGear()

Examples

diet <- getDiet(NS_params)
str(diet)

Get energy rate available for growth

Description

Calculates the energy rate g_i(w) (grams/year) available by species and size for growth after metabolism, movement and reproduction have been accounted for.

Usage

getEGrowth(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A two dimensional array (prey species x prey size)

Your own growth rate function

By default getEGrowth() calls mizerEGrowth(). However you can replace this with your own alternative growth rate function. If your function is called "myEGrowth" then you register it in a MizerParams object params with

params <- setRateFunction(params, "EGrowth", "myEGrowth")

Your function will then be called instead of mizerEGrowth(), with the same arguments.

See Also

getERepro(), getEReproAndGrowth()

Other rate functions: getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the energy at a particular time step
growth <- getEGrowth(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], t = 15)
# Growth rate at this time for Sprat of size 2g
growth["Sprat", "2"]

Get energy rate available for reproduction

Description

Calculates the energy rate (grams/year) available for reproduction after growth and metabolism have been accounted for.

Usage

getERepro(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A two dimensional array (prey species x prey size) holding

\psi_i(w)E_{r.i}(w)

where E_{r.i}(w) is the rate at which energy becomes available for growth and reproduction, calculated with getEReproAndGrowth(), and \psi_i(w) is the proportion of this energy that is used for reproduction. This proportion is taken from the params object and is set with setReproduction().

Your own reproduction rate function

By default getERepro() calls mizerERepro(). However you can replace this with your own alternative reproduction rate function. If your function is called "myERepro" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ERepro", "myERepro")

Your function will then be called instead of mizerERepro(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the rate at a particular time step
erepro <- getERepro(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], t = 15)
# Rate at this time for Sprat of size 2g
erepro["Sprat", "2"]

Get energy rate available for reproduction and growth

Description

Calculates the energy rate E_{r.i}(w) (grams/year) available for reproduction and growth after metabolism and movement have been accounted for.

Usage

getEReproAndGrowth(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A two dimensional array (species x size) holding

E_{r.i}(w) = \max(0, \alpha_i\, (1 - {\tt feeding\_level}_i(w))\, {\tt encounter}_i(w) - {\tt metab}_i(w)).

Due to the form of the feeding level, calculated by getFeedingLevel(), this can also be expressed as

E_{r.i}(w) = \max(0, \alpha_i\, {\tt feeding\_level}_i(w)\, h_i(w) - {\tt metab}_i(w))

where h_i is the maximum intake rate, set with setMaxIntakeRate(). The assimilation rate \alpha_i is taken from the species parameter data frame in params. The metabolic rate metab is taken from params and set with setMetabolicRate().

The return value can be negative, which means that the energy intake does not cover the cost of metabolism and movement.

Your own energy rate function

By default getEReproAndGrowth() calls mizerEReproAndGrowth(). However you can replace this with your own alternative energy rate function. If your function is called "myEReproAndGrowth" then you register it in a MizerParams object params with

params <- setRateFunction(params, "EReproAndGrowth", "myEReproAndGrowth")

Your function will then be called instead of mizerEReproAndGrowth(), with the same arguments.

See Also

The part of this energy rate that is invested into growth is calculated with getEGrowth() and the part that is invested into reproduction is calculated with getERepro().

Other rate functions: getEGrowth(), getERepro(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the energy at a particular time step
e <- getEReproAndGrowth(params, n = N(sim)[15, , ], 
                        n_pp = NResource(sim)[15, ], t = 15)
# Rate at this time for Sprat of size 2g
e["Sprat", "2"]

Alias for getERepro()

Description

An alias provided for backward compatibility with mizer version <= 1.0

Usage

getESpawning(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A two dimensional array (prey species x prey size) holding

\psi_i(w)E_{r.i}(w)

where E_{r.i}(w) is the rate at which energy becomes available for growth and reproduction, calculated with getEReproAndGrowth(), and \psi_i(w) is the proportion of this energy that is used for reproduction. This proportion is taken from the params object and is set with setReproduction().

Your own reproduction rate function

By default getERepro() calls mizerERepro(). However you can replace this with your own alternative reproduction rate function. If your function is called "myERepro" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ERepro", "myERepro")

Your function will then be called instead of mizerERepro(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the rate at a particular time step
erepro <- getERepro(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], t = 15)
# Rate at this time for Sprat of size 2g
erepro["Sprat", "2"]

Fishing effort used in simulation

Description

Note that the array returned may not be exactly the same as the effort argument that was passed in to project(). This is because only the saved effort is stored (the frequency of saving is determined by the argument t_save).

Usage

getEffort(sim)

Arguments

Value

An array (time x gear) that contains the fishing effort by time and gear.

Examples

str(getEffort(NS_sim))

Get encounter rate

Description

Returns the rate at which a predator of species i and weight w encounters food (grams/year).

Usage

getEncounter(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A named two dimensional array (predator species x predator size) with the encounter rates.

Predation encounter

The encounter rate E_i(w) at which a predator of species i and weight w encounters food has contributions from the encounter of fish prey and of resource. This is determined by summing over all prey species and the resource spectrum and then integrating over all prey sizes w_p, weighted by predation kernel \phi(w,w_p):

E_i(w) = \gamma_i(w) \int \left( \theta_{ip} N_R(w_p) + \sum_{j} \theta_{ij} N_j(w_p) \right) \phi_i(w,w_p) w_p \, dw_p.

Here N_j(w) is the abundance density of species j and N_R(w) is the abundance density of resource. The overall prefactor \gamma_i(w) determines the predation power of the predator. It could be interpreted as a search volume and is set with the setSearchVolume() function. The predation kernel \phi(w,w_p) is set with the setPredKernel() function. The species interaction matrix \theta_{ij} is set with setInteraction() and the resource interaction vector \theta_{ip} is taken from the interaction_resource column in params@species_params.

Details

The encounter rate is multiplied by 1-f_0 to obtain the consumption rate, where f_0 is the feeding level calculated with getFeedingLevel(). This is used by the project() function for performing simulations.

The function returns values also for sizes outside the size-range of the species. These values should not be used, as they are meaningless.

If your model contains additional components that you added with setComponent() and for which you specified an encounter_fun function then the encounters of these components will be included in the returned value.

Your own encounter function

By default getEncounter() calls mizerEncounter(). However you can replace this with your own alternative encounter function. If your function is called "myEncounter" then you register it in a MizerParams object params with

params <- setRateFunction(params, "Encounter", "myEncounter")

Your function will then be called instead of mizerEncounter(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

encounter <- getEncounter(NS_params)
str(encounter)

Get the total fishing mortality rate from all fishing gears by time, species and size.

Description

Calculates the total fishing mortality (in units 1/year) from all gears by species and size and possibly time. See setFishing() for details of how fishing gears are set up.

Usage

getFMort(object, effort, time_range, drop = TRUE)

Arguments

Details

The total fishing mortality is just the sum of the fishing mortalities imposed by each gear, \mu_{f.i}(w)=\sum_g F_{g,i,w}. The fishing mortality for each gear is obtained as catchability x selectivity x effort.

Value

An array. If the effort argument has a time dimension, or object is of class MizerSim, the output array has three dimensions (time x species x size). If the effort argument does not have a time dimension, the output array has two dimensions (species x size).

The effort argument is only used if a MizerParams object is passed in. The effort argument can be a two dimensional array (time x gear), a vector of length equal to the number of gears (each gear has a different effort that is constant in time), or a single numeric value (each gear has the same effort that is constant in time). The order of gears in the effort argument must be the same as in the MizerParams object.

If the object argument is of class MizerSim then the effort slot of the MizerSim object is used and the effort argument is not used.

Your own fishing mortality function

By default getFMort() calls mizerFMort(). However you can replace this with your own alternative fishing mortality function. If your function is called "myFMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "FMort", "myFMort")

Your function will then be called instead of mizerFMort(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Get the total fishing mortality in the initial state
F <- getFMort(params, effort = 1)
str(F)
# Get the initial total fishing mortality when effort is different
# between the four gears:
F <- getFMort(params, effort = c(0.5,1,1.5,0.75))
# Get the total fishing mortality when effort is different
# between the four gears and changes with time:
effort <- array(NA, dim = c(20,4))
effort[, 1] <- seq(from = 0, to = 1, length = 20)
effort[, 2] <- seq(from = 1, to = 0.5, length = 20)
effort[, 3] <- seq(from = 1, to = 2, length = 20)
effort[, 4] <- seq(from = 2, to = 1, length = 20)
F <- getFMort(params, effort = effort)
str(F)
# Get the total fishing mortality using the effort already held in a 
# MizerSim object.
sim <- project(params, t_max = 20, effort = 0.5)
F <- getFMort(sim)
F <- getFMort(sim, time_range = c(10, 20))

Get the fishing mortality by time, gear, species and size

Description

Calculates the fishing mortality rate F_{g,i,w} by gear, species and size and possibly time (in units 1/year).

Usage

getFMortGear(object, effort, time_range)

Arguments

Value

An array. If the effort argument has a time dimension, or a MizerSim is passed in, the output array has four dimensions (time x gear x species x size). If the effort argument does not have a time dimension (i.e. it is a vector or a single numeric), the output array has three dimensions (gear x species x size).

Note

Here: fishing mortality = catchability x selectivity x effort.

The effort argument is only used if a MizerParams object is passed in. The effort argument can be a two dimensional array (time x gear), a vector of length equal to the number of gears (each gear has a different effort that is constant in time), or a single numeric value (each gear has the same effort that is constant in time). The order of gears in the effort argument must be the same the same as in the MizerParams object. If the effort argument is not supplied, its value is taken from the ⁠@initial_effort⁠ slot in the params object.

If the object argument is of class MizerSim then the effort slot of the MizerSim object is used and the effort argument is not used.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <-NS_params
# Get the fishing mortality in initial state
F <- getFMortGear(params, effort = 1)
str(F)
# Get the initial fishing mortality when effort is different
# between the four gears:
F <- getFMortGear(params, effort = c(0.5, 1, 1.5, 0.75))
# Get the fishing mortality when effort is different
# between the four gears and changes with time:
effort <- array(NA, dim = c(20, 4))
effort[, 1] <- seq(from=0, to = 1, length = 20)
effort[, 2] <- seq(from=1, to = 0.5, length = 20)
effort[, 3] <- seq(from=1, to = 2, length = 20)
effort[, 4] <- seq(from=2, to = 1, length = 20)
F <- getFMortGear(params, effort = effort)
str(F)
# Get the fishing mortality using the effort already held in a MizerSim object.
sim <- project(params, t_max = 20, effort = 0.5)
F <- getFMortGear(sim)
F <- getFMortGear(sim, time_range = c(10, 20))

Get feeding level

Description

Returns the feeding level. By default this function uses mizerFeedingLevel() to calculate the feeding level, but this can be overruled via setRateFunction().

Usage

getFeedingLevel(object, n, n_pp, n_other, time_range, drop = FALSE, ...)

Arguments

Value

If a MizerParams object is passed in, the function returns a two dimensional array (predator species x predator size) based on the abundances also passed in. If a MizerSim object is passed in, the function returns a three dimensional array (time step x predator species x predator size) with the feeding level calculated at every time step in the simulation. If drop = TRUE then the dimension of length 1 will be removed from the returned array.

Feeding level

The feeding level f_i(w) is the proportion of its maximum intake rate at which the predator is actually taking in fish. It is calculated from the encounter rate E_i and the maximum intake rate h_i(w) as

f_i(w) = \frac{E_i(w)}{E_i(w)+h_i(w)}.

The encounter rate E_i is passed as an argument or calculated with getEncounter(). The maximum intake rate h_i(w) is taken from the params object, and is set with setMaxIntakeRate(). As a consequence of the above expression for the feeding level, 1-f_i(w) is the proportion of the food available to it that the predator actually consumes.

Your own feeding level function

By default getFeedingLevel() calls mizerFeedingLevel(). However you can replace this with your own alternative feeding level function. If your function is called "myFeedingLevel" then you register it in a MizerParams object params with

params <- setRateFunction(params, "FeedingLevel", "myFeedingLevel")

Your function will then be called instead of mizerFeedingLevel(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Get initial feeding level
fl <- getFeedingLevel(params)
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the feeding level at all saved time steps
fl <- getFeedingLevel(sim)
# Get the feeding level for years 15 - 20
fl <- getFeedingLevel(sim, time_range = c(15, 20))

Get growth curves giving weight as a function of age

Description

Get growth curves giving weight as a function of age

Usage

getGrowthCurves(object, species = NULL, max_age = 20, percentage = FALSE)

Arguments

Value

An array (species x age) containing the weight in grams.

See Also

Other summary functions: getBiomass(), getDiet(), getN(), getSSB(), getYield(), getYieldGear()

Examples

growth_curves <- getGrowthCurves(NS_params, species = c("Cod", "Haddock"))
str(growth_curves)

library(ggplot2)
ggplot(melt(growth_curves)) +
  geom_line(aes(Age, value)) +
  facet_wrap(~ Species, scales = "free") +
  ylab("Size[g]") + xlab("Age[years]")

Deprecated function to get interaction matrix

Description

You should now use interaction_matrix() instead.

Usage

getInteraction(params)

Arguments

Alias for getPredMort()

Description

An alias provided for backward compatibility with mizer version <= 1.0

Usage

getM2(object, n, n_pp, n_other, time_range, drop = TRUE, ...)

Arguments

Value

If a MizerParams object is passed in, the function returns a two dimensional array (prey species x prey size) based on the abundances also passed in. If a MizerSim object is passed in, the function returns a three dimensional array (time step x prey species x prey size) with the predation mortality calculated at every time step in the simulation. Dimensions may be dropped if they have length 1 unless drop = FALSE.

Your own predation mortality function

By default getPredMort() calls mizerPredMort(). However you can replace this with your own alternative predation mortality function. If your function is called "myPredMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "PredMort", "myPredMort")

Your function will then be called instead of mizerPredMort(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Predation mortality in initial state
M2 <- getPredMort(params)
str(M2)
# With constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get predation mortality at one time step
M2 <- getPredMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ])
# Get predation mortality at all saved time steps
M2 <- getPredMort(sim)
str(M2)
# Get predation mortality over the years 15 - 20
M2 <- getPredMort(sim, time_range = c(15, 20))

Alias for getResourceMort()

Description

An alias provided for backward compatibility with mizer version <= 1.0

Usage

getM2Background(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A vector of mortality rate by resource size.

Your own resource mortality function

By default getResourceMort() calls mizerResourceMort(). However you can replace this with your own alternative resource mortality function. If your function is called "myResourceMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ResourceMort", "myResourceMort")

Your function will then be called instead of mizerResourceMort(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates()

Examples

params <- NS_params
# With constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get resource mortality at one time step
getResourceMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ])

Calculate the mean maximum weight of the community

Description

Calculates the mean maximum weight of the community through time. This can be calculated by numbers or biomass. The calculation is the sum of the w_max * abundance of each species, divided by the total abundance community, where abundance is either in biomass or numbers. You can specify minimum and maximum weight or length range for the species. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used). You can also specify the species to be used in the calculation.

Usage

getMeanMaxWeight(sim, species = NULL, measure = "both", ...)

Arguments

Value

Depends on the measure argument. If measure = “both” then you get a matrix with two columns, one with values by numbers, the other with values by biomass at each saved time step. If measure = “numbers” or “biomass” you get a vector of the respective values at each saved time step.

See Also

Other functions for calculating indicators: getCommunitySlope(), getMeanWeight(), getProportionOfLargeFish()

Examples

mmw <- getMeanMaxWeight(NS_sim)
years <- c("1967", "2010")
mmw[years, ]
getMeanMaxWeight(NS_sim, species=c("Herring","Sprat","N.pout"))[years, ]
getMeanMaxWeight(NS_sim, min_w = 10, max_w = 5000)[years, ]

Calculate the mean weight of the community

Description

Calculates the mean weight of the community through time. This is simply the total biomass of the community divided by the abundance in numbers. You can specify minimum and maximum weight or length range for the species. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used). You can also specify the species to be used in the calculation.

Usage

getMeanWeight(sim, species = NULL, ...)

Arguments

Value

A vector containing the mean weight of the community through time

See Also

Other functions for calculating indicators: getCommunitySlope(), getMeanMaxWeight(), getProportionOfLargeFish()

Examples

mean_weight <- getMeanWeight(NS_sim)
years <- c("1967", "2010")
mean_weight[years]
getMeanWeight(NS_sim, species = c("Herring", "Sprat", "N.pout"))[years]
getMeanWeight(NS_sim, min_w = 10, max_w = 5000)[years]

Get total mortality rate

Description

Calculates the total mortality rate \mu_i(w) (in units 1/year) on each species by size from predation mortality, background mortality and fishing mortality for a single time step.

Usage

getMort(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  effort = getInitialEffort(params),
  t = 0,
  ...
)

Arguments

Details

If your model contains additional components that you added with setComponent() and for which you specified a mort_fun function then the mortality inflicted by these components will be included in the returned value.

Value

A two dimensional array (prey species x prey size).

Your own mortality function

By default getMort() calls mizerMort(). However you can replace this with your own alternative mortality function. If your function is called "myMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "Mort", "myMort")

Your function will then be called instead of mizerMort(), with the same arguments.

See Also

getPredMort(), getFMort()

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the total mortality at a particular time step
mort <- getMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], 
                t = 15, effort = 0.5)
# Mortality rate at this time for Sprat of size 2g
mort["Sprat", "2"]

Calculate the number of individuals within a size range

Description

Calculates the number of individuals within user-defined size limits. The default option is to use the whole size range. You can specify minimum and maximum weight or lengths for the species. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used)

Usage

getN(object, ...)

Arguments

Value

If called with a MizerParams object, a vector with the numbers for each species in the model. If called with a MizerSim object, an array (time x species) containing the numbers at each time step for all species.

See Also

Other summary functions: getBiomass(), getDiet(), getGrowthCurves(), getSSB(), getYield(), getYieldGear()

Examples

numbers <- getN(NS_sim)
numbers["1972", "Herring"]
# The above gave a huge number, because that included all the larvae.
# The number of Herrings between 10g and 1kg is much smaller.
numbers <- getN(NS_sim, min_w = 10, max_w = 1000)
numbers["1972", "Herring"]

Extract the parameter object underlying a simulation

Description

Extract the parameter object underlying a simulation

Usage

getParams(sim)

Arguments

Value

The MizerParams object that was used to run the simulation

Examples

# This will be identical to the params object that was used to create the
# simulation
sim <- project(NS_params, t_max = 1)
identical(getParams(sim), NS_params)

Get available energy

Description

This is deprecated and is no longer used by the mizer project() method. Calculates the amount E_{a,i}(w) of food exposed to each predator as a function of predator size.

Usage

getPhiPrey(object, n, n_pp, ...)

Arguments

Value

A two dimensional array (predator species x predator size)

See Also

project()

Examples

params <-  NS_params
sim <- project(params, t_max = 20, effort = 0.5)
n <- sim@n[21,,]
n_pp <- sim@n_pp[21,]
getPhiPrey(params,n,n_pp)
# ->
getEncounter(params) / getSearchVolume(params)

Get total predation mortality rate

Description

Calculates the total predation mortality rate \mu_{p,i}(w_p) (in units of 1/year) on each prey species by prey size:

\mu_{p.i}(w_p) = \sum_j {\tt pred\_rate}_j(w_p)\, \theta_{ji}.

The predation rate pred_rate is returned by getPredRate().

Usage

getPredMort(object, n, n_pp, n_other, time_range, drop = TRUE, ...)

Arguments

Value

If a MizerParams object is passed in, the function returns a two dimensional array (prey species x prey size) based on the abundances also passed in. If a MizerSim object is passed in, the function returns a three dimensional array (time step x prey species x prey size) with the predation mortality calculated at every time step in the simulation. Dimensions may be dropped if they have length 1 unless drop = FALSE.

Your own predation mortality function

By default getPredMort() calls mizerPredMort(). However you can replace this with your own alternative predation mortality function. If your function is called "myPredMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "PredMort", "myPredMort")

Your function will then be called instead of mizerPredMort(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Predation mortality in initial state
M2 <- getPredMort(params)
str(M2)
# With constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get predation mortality at one time step
M2 <- getPredMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ])
# Get predation mortality at all saved time steps
M2 <- getPredMort(sim)
str(M2)
# Get predation mortality over the years 15 - 20
M2 <- getPredMort(sim, time_range = c(15, 20))

Get predation rate

Description

Calculates the potential rate (in units 1/year) at which a prey individual of a given size w is killed by predators from species j. In formulas

{\tt pred\_rate}_j(w_p) = \int \phi_j(w,w_p) (1-f_j(w)) \gamma_j(w) N_j(w) \, dw.

This potential rate is used in getPredMort() to calculate the realised predation mortality rate on the prey individual.

Usage

getPredRate(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A two dimensional array (predator species x prey size), where the prey size runs over fish community plus resource spectrum.

Your own predation rate function

By default getPredRate() calls mizerPredRate(). However you can replace this with your own alternative predation rate function. If your function is called "myPredRate" then you register it in a MizerParams object params with

params <- setRateFunction(params, "PredRate", "myPredRate")

Your function will then be called instead of mizerPredRate(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Predation rate in initial state
pred_rate <- getPredRate(params)
str(pred_rate)
# With constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the feeding level at one time step
pred_rate <- getPredRate(params, n = N(sim)[15, , ], 
                         n_pp = NResource(sim)[15, ], t = 15)

Calculate the proportion of large fish

Description

Calculates the proportion of large fish through time in the MizerSim class within user defined size limits. The default option is to use the whole size range. You can specify minimum and maximum size ranges for the species and also the threshold size for large fish. Sizes can be expressed as weight or size. Lengths take precedence over weights (i.e. if both min_l and min_w are supplied, only min_l will be used). You can also specify the species to be used in the calculation. This function can be used to calculate the Large Fish Index. The proportion is based on either abundance or biomass.

Usage

getProportionOfLargeFish(
  sim,
  species = NULL,
  threshold_w = 100,
  threshold_l = NULL,
  biomass_proportion = TRUE,
  ...
)

Arguments

Value

A vector containing the proportion of large fish through time

See Also

Other functions for calculating indicators: getCommunitySlope(), getMeanMaxWeight(), getMeanWeight()

Examples

lfi <- getProportionOfLargeFish(NS_sim, min_w = 10, max_w = 5000, 
                                threshold_w = 500)
years <- c("1972", "2010")
lfi[years]
getProportionOfLargeFish(NS_sim)[years]
getProportionOfLargeFish(NS_sim, species=c("Herring","Sprat","N.pout"))[years]
getProportionOfLargeFish(NS_sim, min_w = 10, max_w = 5000)[years]
getProportionOfLargeFish(NS_sim, min_w = 10, max_w = 5000,
    threshold_w = 500, biomass_proportion = FALSE)[years]

Get density dependent reproduction rate

Description

Calculates the density dependent rate of egg production R_i (units 1/year) for each species. This is the flux entering the smallest size class of each species. The density dependent rate is the density independent rate obtained with getRDI() after it has been put through the density dependence function. This is the Beverton-Holt function BevertonHoltRDD() by default, but this can be changed. See setReproduction() for more details.

Usage

getRDD(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  rdi = getRDI(params, n = n, n_pp = n_pp, n_other = n_other, t = t),
  ...
)

Arguments

Value

A numeric vector the length of the number of species.

See Also

getRDI()

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the rate at a particular time step
getRDD(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], t = 15)

Get density independent rate of egg production

Description

Calculates the density-independent rate of total egg production R_{di} (units 1/year) before density dependence, by species.

Usage

getRDI(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Details

This rate is obtained by taking the per capita rate E_r(w)\psi(w) at which energy is invested in reproduction, as calculated by getERepro(), multiplying it by the number of individualsN(w) and integrating over all sizes w and then multiplying by the reproductive efficiency \epsilon and dividing by the egg size w_min, and by a factor of two to account for the two sexes:

R_{di} = \frac{\epsilon}{2 w_{min}} \int N(w) E_r(w) \psi(w) \, dw

Used by getRDD() to calculate the actual, density dependent rate. See setReproduction() for more details.

Value

A numeric vector the length of the number of species.

Your own reproduction function

By default getRDI() calls mizerRDI(). However you can replace this with your own alternative reproduction function. If your function is called "myRDI" then you register it in a MizerParams object params with

params <- setRateFunction(params, "RDI", "myRDI")

Your function will then be called instead of mizerRDI(), with the same arguments. For an example of an alternative reproduction function see constantEggRDI().

See Also

getRDD()

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the density-independent reproduction rate at a particular time step
getRDI(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], t = 15)

Get all rates

Description

Calls other rate functions in sequence and collects the results in a list.

Usage

getRates(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  effort,
  t = 0,
  ...
)

Arguments

Details

By default this function returns a list with the following components:

• encounter from mizerEncounter()

• feeding_level from mizerFeedingLevel()

• e from mizerEReproAndGrowth()

• e_repro from mizerERepro()

• e_growth from mizerEGrowth()

• pred_rate from mizerPredRate()

• pred_mort from mizerPredMort()

• f_mort from mizerFMort()

• mort from mizerMort()

• rdi from mizerRDI()

• rdd from BevertonHoltRDD()

• resource_mort from mizerResourceMort()

However you can replace any of these rate functions by your own rate function if you wish, see setRateFunction() for details.

Value

List of rates.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getResourceMort()

Examples

rates <- getRates(NS_params)
names(rates)
identical(rates$encounter, getEncounter(NS_params))

Get reproduction level

Description

The reproduction level is the ratio between the density-dependent reproduction rate and the maximal reproduction rate.

Usage

getReproductionLevel(params)

Arguments

Value

A named vector with the reproduction level for each species.

Examples

getReproductionLevel(NS_params)

# The reproduction level can be changed without changing the steady state:
params <- setBevertonHolt(NS_params, reproduction_level = 0.9)
getReproductionLevel(params)

# The result is the ratio of RDD and R_max
identical(getRDD(params) / species_params(params)$R_max,
          getReproductionLevel(params))

Deprecated functions for getting resource parameters

Description

Use resource_dynamics(), resource_level(), resource_rate() and resource_capacity() instead.

Usage

getResourceDynamics(params)

getResourceLevel(params)

getResourceRate(params)

getResourceCapacity(params)

Arguments

Value

Name of the function that determines the resource dynamics

Get predation mortality rate for resource

Description

Calculates the predation mortality rate \mu_p(w) on the resource spectrum by resource size (in units 1/year).

Usage

getResourceMort(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  t = 0,
  ...
)

Arguments

Value

A vector of mortality rate by resource size.

Your own resource mortality function

By default getResourceMort() calls mizerResourceMort(). However you can replace this with your own alternative resource mortality function. If your function is called "myResourceMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ResourceMort", "myResourceMort")

Your function will then be called instead of mizerResourceMort(), with the same arguments.

See Also

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getMort(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates()

Examples

params <- NS_params
# With constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get resource mortality at one time step
getResourceMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ])

Calculate the SSB of species

Description

Calculates the spawning stock biomass (SSB) through time of the species in the MizerSim class. SSB is calculated as the total mass of all mature individuals.

Usage

getSSB(object)

Arguments

Value

If called with a MizerParams object, a vector with the SSB in grams for each species in the model. If called with a MizerSim object, an array (time x species) containing the SSB in grams at each time step for all species.

See Also

Other summary functions: getBiomass(), getDiet(), getGrowthCurves(), getN(), getYield(), getYieldGear()

Examples

ssb <- getSSB(NS_sim)
ssb[c("1972", "2010"), c("Herring", "Cod")]

Times for which simulation results are available

Description

Times for which simulation results are available

Usage

getTimes(sim)

Arguments

Value

A numeric vectors of the times (in years) at which simulation results have been stored in the MizerSim object.

Examples

getTimes(NS_sim)

Calculate the rate at which biomass of each species is fished

Description

This yield rate is given in grams per year. It is calculated at each time step saved in the MizerSim object.

Usage

getYield(object)

Arguments

Details

The yield rate y_i(t) for species i at time t is defined as

y_i(t)=\int\mu_{f.i}(w, t)N_i(w, t)w dw

where \mu_{f.i}(w, t) is the fishing mortality of an individual of species i and weight w at time t and N_i(w, t) is the abundance density of such individuals. The factor of w converts the abundance density into a biomass density and the integral aggregates the contribution from all sizes.

The total catch in a time period from t_1 to t_2 is the integral of the yield rate over that period:

C = \int_{t_1}^{t2}y_i(t)dt

In practice, as the yield rate is only available at the saved times, one can only approximate this integral by averaging over the available yield rates during the time period and multiplying by the time period. The less the yield changes between the saved values, the more accurate this approximation is. So the approximation can be improved by saving simulation results at smaller intervals, using the t_save argument to project(). But this is only a concern if abundances change quickly during the time period of interest.

Value

If called with a MizerParams object, a vector with the yield rate in grams per year for each species in the model. If called with a MizerSim object, an array (time x species) containing the yield rate at each time step for all species.

See Also

getYieldGear()

Other summary functions: getBiomass(), getDiet(), getGrowthCurves(), getN(), getSSB(), getYieldGear()

Examples

yield <- getYield(NS_sim)
yield[c("1972", "2010"), c("Herring", "Cod")]

# Running simulation for another year, saving intermediate time steps
params <- setInitialValues(getParams(NS_sim), NS_sim)
sim <- project(params, t_save = 0.1, t_max = 1, 
               t_start = 2010, progress_bar = FALSE)
# The yield rate for Herring decreases during the year
getYield(sim)[, "Herring"]
# We get the total catch in the year by averaging over the year
sum(getYield(sim)[1:10, "Herring"] / 10)

Calculate the rate at which biomass of each species is fished by each gear

Description

This yield rate is given in grams per year. It is calculated at each time step saved in the MizerSim object.

Usage

getYieldGear(object)

Arguments

Details

For details of how the yield rate is defined see the help page of getYield().

Value

If called with a MizerParams object, an array (gear x species) with the yield rate in grams per year from each gear for each species in the model. If called with a MizerSim object, an array (time x gear x species) containing the yield rate at each time step.

See Also

getYield()

Other summary functions: getBiomass(), getDiet(), getGrowthCurves(), getN(), getSSB(), getYield()

Examples

yield <- getYieldGear(NS_sim)
yield["1972", "Herring", "Herring"]
# (In this example MizerSim object each species was set up with its own gear)

Alias for getMort()

Description

An alias provided for backward compatibility with mizer version <= 1.0

Usage

getZ(
  params,
  n = initialN(params),
  n_pp = initialNResource(params),
  n_other = initialNOther(params),
  effort = getInitialEffort(params),
  t = 0,
  ...
)

Arguments

Details

If your model contains additional components that you added with setComponent() and for which you specified a mort_fun function then the mortality inflicted by these components will be included in the returned value.

Value

A two dimensional array (prey species x prey size).

Your own mortality function

By default getMort() calls mizerMort(). However you can replace this with your own alternative mortality function. If your function is called "myMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "Mort", "myMort")

Your function will then be called instead of mizerMort(), with the same arguments.

See Also

getPredMort(), getFMort()

Other rate functions: getEGrowth(), getERepro(), getEReproAndGrowth(), getEncounter(), getFMort(), getFMortGear(), getFeedingLevel(), getPredMort(), getPredRate(), getRDD(), getRDI(), getRates(), getResourceMort()

Examples

params <- NS_params
# Project with constant fishing effort for all gears for 20 time steps
sim <- project(params, t_max = 20, effort = 0.5)
# Get the total mortality at a particular time step
mort <- getMort(params, n = N(sim)[15, , ], n_pp = NResource(sim)[15, ], 
                t = 15, effort = 0.5)
# Mortality rate at this time for Sprat of size 2g
mort["Sprat", "2"]

Get default value for f0

Description

Fills in any missing values for f0 so that if the prey abundance was described by the power law \kappa w^{-\lambda} then the encounter rate coming from the given gamma parameter would lead to the feeding level f_0. This is thus doing the inverse of get_gamma_default(). Only for internal use.

Usage

get_f0_default(params)

Arguments

Details

For species for which no value for gamma is specified in the species parameter data frame, the f0 values is kept as provided in the species parameter data frame or it is set to 0.6 if it is not provided.

Value

A vector with the values of f0 for all species

See Also

Other functions calculating defaults: get_gamma_default(), get_h_default(), get_ks_default()

Get default value for gamma

Description

Fills in any missing values for gamma so that fish feeding on a resource spectrum described by the power law \kappa w^{-\lambda} achieve a feeding level f_0. Only for internal use.

Usage

get_gamma_default(params)

Arguments

Value

A vector with the values of gamma for all species

See Also

Other functions calculating defaults: get_f0_default(), get_h_default(), get_ks_default()

Get default value for h

Description

Sets h so that the species reaches maturity size w_mat at the maturity age age_mat if it feeds at feeding level f0.

Usage

get_h_default(params)

Arguments

Details

If age_mat is missing in the species parameter data frame, then it is calculated from the von Bertalanffy growth curve parameters k_vb and (optionally t0) taken from the species parameter data frame. This is not reliable and a warning is issued.

If no growth information is given at all for a species, the default is set to h = 30.

Value

A vector with the values of h for all species

See Also

Other functions calculating defaults: get_f0_default(), get_gamma_default(), get_ks_default()

Calculate initial population abundances

Description

This function uses the model parameters and other parameters to calculate initial values for the species number densities. These initial abundances are currently quite arbitrary and not close to the steady state. We intend to improve this in the future.

Usage

get_initial_n(params, n0_mult = NULL, a = 0.35)

Arguments

Value

A matrix (species x size) of population abundances.

Examples

init_n <- get_initial_n(NS_params)

Get default value for ks

Description

Fills in any missing values for ks so that the critical feeding level needed to sustain the species is as specified in the fc column in the species parameter data frame. If that column is not provided the default critical feeding level f_c = 0.2 is used.

Usage

get_ks_default(params)

Arguments

Value

A vector with the values of ks for all species

See Also

Other functions calculating defaults: get_f0_default(), get_gamma_default(), get_h_default()

Get values from feeding kernel function

Description

This involves finding the feeding kernel function for each species, using the pred_kernel_type parameter in the species_params data frame, checking that it is valid and all its arguments are contained in the species_params data frame, and then calling this function with the ppmr vector.

Usage

get_phi(species_params, ppmr)

Arguments

Value

An array (species x ppmr) with the values of the predation kernel function

Determine reproduction rate needed for initial egg abundance

Description

Determine reproduction rate needed for initial egg abundance

Usage

get_required_reproduction(params)

Arguments

Value

A vector of reproduction rates for all species

Get size range array

Description

Helper function that returns an array (species x size) of boolean values indicating whether that size bin is within the size limits specified by the arguments. Either the size limits can be the same for all species or they can be specified as vectors with one value for each species in the model.

Usage

get_size_range_array(
  params,
  min_w = min(params@w),
  max_w = max(params@w),
  min_l = NULL,
  max_l = NULL,
  ...
)

Arguments

Value

Boolean array (species x size)

Length to weight conversion

If min_l is specified there is no need to specify min_w and so on. However, if a length is specified (minimum or maximum) then it is necessary for the species parameter data.frame to include the parameters a and b that determine the relation between length l and weight w by

w = a l^b.

It is possible to mix length and weight constraints, e.g. by supplying a minimum weight and a maximum length, but this must be done the same for all species. The default values are the minimum and maximum weights of the spectrum, i.e., the full range of the size spectrum is used.

Get_time_elements

Description

Internal function to get the array element references of the time dimension for the time based slots of a MizerSim object.

Usage

get_time_elements(sim, time_range, slot_name = "n")

Arguments

Value

Named boolean vector indicating for each time whether it is included in the range or not.

Description of indicator functions

Description

Mizer provides a range of functions to calculate indicators from a MizerSim object.

Details

A list of available indicator functions for MizerSim objects is given in the table below

See Also

summary_functions, plotting_functions

Initial values for fish spectra

Description

Values used as starting values for simulations with project().

Usage

initialN(params) <- value

initialN(object)

Arguments

Value

A matrix with dimensions species x size holding the initial number densities for the fish spectra.

See Also

initialNResource(), initialNOther()

Examples

# Doubling abundance of Cod in the initial state of the North Sea model
params <- NS_params
initialN(params)["Cod", ] <- 2 * initialN(params)["Cod", ]
# Calculating the corresponding initial biomass
biomass <- initialN(params)["Cod", ] * dw(NS_params) * w(NS_params)
# Of course this initial state will no longer be a steady state
params <- steady(params)

Initial values for other ecosystem components

Description

Values used as starting values for simulations with project().

Usage

initialNOther(params) <- value

initialNOther(object)

Arguments

Value

A named list with the initial values of other ecosystem components

See Also

initialNResource(), initialN()

Initial value for resource spectrum

Description

Value used as starting value for simulations with project().

Usage

initialNResource(params) <- value

initialNResource(object)

Arguments

Value

A vector with the initial number densities for the resource spectrum

See Also

initialN(), initialNOther()

Examples

# Doubling resource abundance in the initial state of the North Sea model
params <- NS_params
initialNResource(params) <- 2 * initialNResource(params)
# Of course this initial state will no longer be a steady state
params <- steady(params)

Initial fishing effort

Description

The fishing effort is a named vector, specifying for each fishing gear the effort invested into fishing with that gear. The effort value for each gear is multiplied by the catchability and the selectivity to determine the fishing mortality imposed by that gear, see setFishing() for more details. The initial effort you have set can be overruled when running a simulation by providing an effort argument to project() which allows you to specify a time-varying effort.

Usage

initial_effort(params)

initial_effort(params) <- value

validEffortVector(effort, params)

Arguments

Details

A valid effort vector is a named vector with one effort value for each gear. However you can also supply the effort value in different ways:

• a scalar, which is then replicated for each gear

• an unnamed vector, which is then assumed to be in the same order as the gears in the params object

• a named vector in which the gear names have a different order than in the params object. This is then sorted correctly.

• a named vector which only supplies values for some of the gears. The effort for the other gears is then set to zero.

These conversions are done by the function validEffortVector().

An effort argument will lead to an error if it is either

• unnamed and of the wrong length

• named but where some names do not match any of the gears

• not numeric

Value

Effort vector

Alias for NS_interaction

Description

An alias provided for backward compatibility with mizer version <= 2.3

Usage

inter

Format

A 12 x 12 matrix.

Source

Blanchard et al.

Weight based knife-edge selectivity function

Description

A knife-edge selectivity function where weights greater or equal to knife_edge_size are fully selected and no fish smaller than this size are selected.

Usage

knife_edge(w, knife_edge_size, ...)

Arguments

Value

Vector of selectivities at the given sizes.

See Also

gear_params() for setting the knife_edge_size parameter.

Other selectivity functions: double_sigmoid_length(), sigmoid_length(), sigmoid_weight()

Length-weight conversion

Description

For each species, convert between length and weight using the relationship

w_i = a_i l_i^{b_i}

or

l_i = (w_i / a_i)^{1/b_i}

where a and b are taken from the species parameter data frame and i is the species index.

Usage

l2w(l, params)

w2l(w, params)

Arguments

Details

This is useful for converting a length-based species parameter to a weight-based species parameter.

If any a or b parameters are missing the default values a = 0.01 and b = 3 are used for the missing values.

Value

A vector with one entry for each species. l2w() returns a vector of weights in grams and w2l() returns a vector of lengths in cm.

Helper function to produce nice breaks on logarithmic axes

Description

This is needed when the logarithmic y-axis spans less than one order of magnitude, in which case the ggplot2 default produces no ticks.

Usage

log_breaks(n = 6)

Arguments

Details

Thanks to Heather Turner at https://stackoverflow.com/questions/14255533/pretty-ticks-for-log-normal-scale-using-ggplot2-dynamic-not-manual

Value

A function that can be used as the break argument in calls to scale_y_continuous() or scale_x_continuous()

Lognormal predation kernel

Description

This is the most commonly-used predation kernel. The log of the predator/prey mass ratio is normally distributed.

Usage

lognormal_pred_kernel(ppmr, beta, sigma)

Arguments

Details

Writing the predator mass as w and the prey mass as w_p, the feeding kernel is given as

\phi_i(w, w_p) = \exp \left[ \frac{-(\ln(w / w_p / \beta_i))^2}{2\sigma_i^2} \right]

if w/w_p is larger than 1 and zero otherwise. Here \beta_i is the preferred predator-prey mass ratio and \sigma_i determines the width of the kernel. These two parameters need to be given in the species parameter dataframe in the columns beta and sigma.

This function is called from setPredKernel() to set up the predation kernel slots in a MizerParams object.

Value

A vector giving the value of the predation kernel at each of the predator/prey mass ratios in the ppmr argument.

See Also

setPredKernel()

Other predation kernel: box_pred_kernel(), power_law_pred_kernel(), truncated_lognormal_pred_kernel()

Examples

params <- NS_params
plot(w_full(params), getPredKernel(params)["Cod", 10, ], type="l", log="x")
# The restriction that the kernel is zero for w/w_p < 1 is more
# noticeable for larger sigma
species_params(params)$sigma <- 4
plot(w_full(params), getPredKernel(params)["Cod", 10, ], type="l", log="x")

Match biomasses to observations

Description

The function adjusts the abundances of the species in the model so that their biomasses match with observations.

Usage

matchBiomasses(params, species = NULL)

Arguments

Details

The function works by multiplying for each species the abundance density at all sizes by the same factor. This will of course not give a steady state solution, even if the initial abundance densities were at steady state. So after using this function you may want to use steady() to run the model to steady state, after which of course the biomasses will no longer match exactly. You could then iterate this process. This is described in the blog post at https://bit.ly/2YqXESV.

Before you can use this function you will need to have added a biomass_observed column to your model which gives the observed biomass in grams. For species for which you have no observed biomass, you should set the value in the biomass_observed column to 0 or NA.

Biomass observations usually only include individuals above a certain size. This size should be specified in a biomass_cutoff column of the species parameter data frame. If this is missing, it is assumed that all sizes are included in the observed biomass, i.e., it includes larval biomass.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$biomass_observed <- 
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
species_params(params)$biomass_cutoff <- 10
params <- calibrateBiomass(params)
params <- matchBiomasses(params)
plotBiomassObservedVsModel(params)

Adjust model to produce observed growth

Description

Scales the search volume, the maximum consumption rate, the metabolic rate and the external encounter rate all by the same factor in order to achieve a growth rate that allows individuals to reach their maturity size by their maturity age while keeping the feeding level and the critical feeding level unchanged. Then recalculates the size spectra using steadySingleSpecies().

Usage

matchGrowth(params, species = NULL, keep = c("egg", "biomass", "number"))

Arguments

Details

Maturity size and age are taken from the w_mat and age_mat columns in the species_params data frame. If age_mat is missing, mizer calculates it from the von Bertalanffy growth curve parameters using age_mat_vB(). If those are not available either for a species, the growth rate for that species will not be changed.

Value

A modified MizerParams object with rescaled search volume, maximum consumption rate and metabolic rate and rescaled species parameters gamma,h, ks and k.

Match numbers to observations

Description

The function adjusts the numbers of the species in the model so that their numbers match with observations.

Usage

matchNumbers(params, species = NULL)

Arguments

Details

The function works by multiplying for each species the number density at all sizes by the same factor. This will of course not give a steady state solution, even if the initial number densities were at steady state. So after using this function you may want to use steady() to run the model to steady state, after which of course the numbers will no longer match exactly. You could then iterate this process. This is described in the blog post at https://bit.ly/2YqXESV.

Before you can use this function you will need to have added a number_observed column to your model which gives the observed number in grams. For species for which you have no observed number, you should set the value in the number_observed column to 0 or NA.

Number observations usually only include individuals above a certain size. This size should be specified in a number_cutoff column of the species parameter data frame. If this is missing, it is assumed that all sizes are included in the observed number, i.e., it includes larval number.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$number_observed <-
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
species_params(params)$number_cutoff <- 10
params <- calibrateNumber(params)
params <- matchNumbers(params)

Match yields to observations

Description

This function has been deprecated and will be removed in the future unless you have a use case for it. If you do have a use case for it, please let the developers know by creating an issue at https://github.com/sizespectrum/mizer/issues.

Usage

matchYields(params, species = NULL)

Arguments

Details

If you want to match the yields to observations, you should use the matchYield() function from the mizerExperimental package instead, which adjusts the catchability to match the yield rather than by adjusting the biomass.

The function adjusts the abundances of the species in the model so that their yearly yields under the given fishing mortalities match with observations.

The function works by multiplying for each species the abundance density at all sizes by the same factor. This will of course not give a steady state solution, even if the initial abundance densities were at steady state. So after using this function you may want to use steady() to run the model to steady state, after which of course the yields will no longer match exactly. You could then iterate this process. This is described in the blog post at https://bit.ly/2YqXESV.

Before you can use this function you will need to have added a yield_observed column to your model which gives the observed yields in grams per year. For species for which you have no observed biomass, you should set the value in the yield_observed column to 0 or NA.

Value

A MizerParams object

Examples

params <- NS_params
species_params(params)$yield_observed <- 
    c(0.8, 61, 12, 35, 1.6, 20, 10, 7.6, 135, 60, 30, 78)
gear_params(params)$catchability <-
    c(1.3, 0.065, 0.31, 0.18, 0.98, 0.24, 0.37, 0.46, 0.18, 0.30, 0.27, 0.39)
params <- calibrateYield(params)
params <- matchYields(params)
plotYieldObservedVsModel(params)

Get energy rate available for growth needed to project standard mizer model

Description

Calculates the energy rate g_i(w) (grams/year) available by species and size for growth after metabolism, movement and reproduction have been accounted for. Used by project() for performing simulations. You would not usually call this function directly but instead use getEGrowth(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerEGrowth(params, n, n_pp, n_other, t, e_repro, e, ...)

Arguments

Value

A two dimensional array (species x size) with the growth rates.

Your own growth rate function

By default getEGrowth() calls mizerEGrowth(). However you can replace this with your own alternative growth rate function. If your function is called "myEGrowth" then you register it in a MizerParams object params with

params <- setRateFunction(params, "EGrowth", "myEGrowth")

Your function will then be called instead of mizerEGrowth(), with the same arguments.

See Also

Other mizer rate functions: mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get energy rate available for reproduction needed to project standard mizer model

Description

Calculates the energy rate (grams/year) available for reproduction after growth and metabolism have been accounted for. You would not usually call this function directly but instead use getERepro(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerERepro(params, n, n_pp, n_other, t, e, ...)

Arguments

Value

A two dimensional array (species x size) holding

\psi_i(w)E_{r.i}(w)

where E_{r.i}(w) is the rate at which energy becomes available for growth and reproduction, calculated with mizerEReproAndGrowth(), and \psi_i(w) is the proportion of this energy that is used for reproduction. This proportion is taken from the params object and is set with setReproduction().

Your own reproduction rate function

By default getERepro() calls mizerERepro(). However you can replace this with your own alternative reproduction rate function. If your function is called "myERepro" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ERepro", "myERepro")

Your function will then be called instead of mizerERepro(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get energy rate available for reproduction and growth needed to project standard mizer model

Description

Calculates the energy rate E_{r.i}(w) (grams/year) available to an individual of species i and size w for reproduction and growth after metabolism and movement have been accounted for. You would not usually call this function directly but instead use getEReproAndGrowth(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerEReproAndGrowth(
  params,
  n,
  n_pp,
  n_other,
  t,
  encounter,
  feeding_level,
  ...
)

Arguments

Value

A two dimensional array (species x size) holding

E_{r.i}(w) = \max(0, \alpha_i\, (1 - {\tt feeding\_level}_i(w))\, {\tt encounter}_i(w) - {\tt metab}_i(w)).

Due to the form of the feeding level, calculated by getFeedingLevel(), this can also be expressed as

E_{r.i}(w) = \max(0, \alpha_i\, {\tt feeding\_level}_i(w)\, h_i(w) - {\tt metab}_i(w))

where h_i is the maximum intake rate, set with setMaxIntakeRate(). The assimilation rate \alpha_i is taken from the species parameter data frame in params. The metabolic rate metab is taken from params and set with setMetabolicRate().

The return value can be negative, which means that the energy intake does not cover the cost of metabolism and movement.

Your own energy rate function

By default getEReproAndGrowth() calls mizerEReproAndGrowth(). However you can replace this with your own alternative energy rate function. If your function is called "myEReproAndGrowth" then you register it in a MizerParams object params with

params <- setRateFunction(params, "EReproAndGrowth", "myEReproAndGrowth")

Your function will then be called instead of mizerEReproAndGrowth(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get encounter rate needed to project standard mizer model

Description

Calculates the rate E_i(w) at which a predator of species i and weight w encounters food (grams/year). You would not usually call this function directly but instead use getEncounter(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerEncounter(params, n, n_pp, n_other, t, ...)

Arguments

Value

A named two dimensional array (predator species x predator size) with the encounter rates.

Predation encounter

The encounter rate E_i(w) at which a predator of species i and weight w encounters food has contributions from the encounter of fish prey and of resource. This is determined by summing over all prey species and the resource spectrum and then integrating over all prey sizes w_p, weighted by predation kernel \phi(w,w_p):

E_i(w) = \gamma_i(w) \int \left( \theta_{ip} N_R(w_p) + \sum_{j} \theta_{ij} N_j(w_p) \right) \phi_i(w,w_p) w_p \, dw_p.

Here N_j(w) is the abundance density of species j and N_R(w) is the abundance density of resource. The overall prefactor \gamma_i(w) determines the predation power of the predator. It could be interpreted as a search volume and is set with the setSearchVolume() function. The predation kernel \phi(w,w_p) is set with the setPredKernel() function. The species interaction matrix \theta_{ij} is set with setInteraction() and the resource interaction vector \theta_{ip} is taken from the interaction_resource column in params@species_params.

Details

The encounter rate is multiplied by 1-f_0 to obtain the consumption rate, where f_0 is the feeding level calculated with getFeedingLevel(). This is used by the project() function for performing simulations.

The function returns values also for sizes outside the size-range of the species. These values should not be used, as they are meaningless.

If your model contains additional components that you added with setComponent() and for which you specified an encounter_fun function then the encounters of these components will be included in the returned value.

Your own encounter function

By default getEncounter() calls mizerEncounter(). However you can replace this with your own alternative encounter function. If your function is called "myEncounter" then you register it in a MizerParams object params with

params <- setRateFunction(params, "Encounter", "myEncounter")

Your function will then be called instead of mizerEncounter(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get the total fishing mortality rate from all fishing gears

Description

Calculates the total fishing mortality (in units 1/year) from all gears by species and size. The total fishing mortality is just the sum of the fishing mortalities imposed by each gear, \mu_{f.i}(w)=\sum_g F_{g,i,w}. You would not usually call this function directly but instead use getFMort(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerFMort(params, n, n_pp, n_other, t, effort, e_growth, pred_mort, ...)

Arguments

Value

An array (species x size) with the fishing mortality.

Your own fishing mortality function

By default getFMort() calls mizerFMort(). However you can replace this with your own alternative fishing mortality function. If your function is called "myFMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "FMort", "myFMort")

Your function will then be called instead of mizerFMort(), with the same arguments.

Note

Here: fishing mortality = catchability x selectivity x effort.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get the fishing mortality needed to project standard mizer model

Description

Calculates the fishing mortality rate F_{g,i,w} by gear, species and size. This is a helper function for mizerFMort().

Usage

mizerFMortGear(params, effort)

Arguments

Value

An three dimensional array (gear x species x size) with the fishing mortality

Note

Here: fishing mortality = catchability x selectivity x effort.

See Also

setFishing()

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get feeding level needed to project standard mizer model

Description

You would not usually call this function directly but instead use getFeedingLevel(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerFeedingLevel(params, n, n_pp, n_other, t, encounter, ...)

Arguments

Value

A two dimensional array (predator species x predator size) with the feeding level.

Feeding level

The feeding level f_i(w) is the proportion of its maximum intake rate at which the predator is actually taking in fish. It is calculated from the encounter rate E_i and the maximum intake rate h_i(w) as

f_i(w) = \frac{E_i(w)}{E_i(w)+h_i(w)}.

The encounter rate E_i is passed as an argument or calculated with getEncounter(). The maximum intake rate h_i(w) is taken from the params object, and is set with setMaxIntakeRate(). As a consequence of the above expression for the feeding level, 1-f_i(w) is the proportion of the food available to it that the predator actually consumes.

Your own feeding level function

By default getFeedingLevel() calls mizerFeedingLevel(). However you can replace this with your own alternative feeding level function. If your function is called "myFeedingLevel" then you register it in a MizerParams object params with

params <- setRateFunction(params, "FeedingLevel", "myFeedingLevel")

Your function will then be called instead of mizerFeedingLevel(), with the same arguments.

See Also

The feeding level is used in mizerEReproAndGrowth() and in mizerPredRate().

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get total mortality rate needed to project standard mizer model

Description

Calculates the total mortality rate \mu_i(w) (in units 1/year) on each species by size from predation mortality, background mortality and fishing mortality. You would not usually call this function directly but instead use getMort(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerMort(params, n, n_pp, n_other, t, f_mort, pred_mort, ...)

Arguments

Details

If your model contains additional components that you added with setComponent() and for which you specified a mort_fun function then the mortality inflicted by these components will be included in the returned value.

Value

A named two dimensional array (species x size) with the total mortality rates.

Your own mortality function

By default getMort() calls mizerMort(). However you can replace this with your own alternative mortality function. If your function is called "myMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "Mort", "myMort")

Your function will then be called instead of mizerMort(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get total predation mortality rate needed to project standard mizer model

Description

Calculates the total predation mortality rate \mu_{p,i}(w_p) (in units of 1/year) on each prey species by prey size:

\mu_{p.i}(w_p) = \sum_j {\tt pred\_rate}_j(w_p)\, \theta_{ji}.

You would not usually call this function directly but instead use getPredMort(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerPredMort(params, n, n_pp, n_other, t, pred_rate, ...)

Arguments

Value

A two dimensional array (prey species x prey size) with the predation mortality

Your own predation mortality function

By default getPredMort() calls mizerPredMort(). However you can replace this with your own alternative predation mortality function. If your function is called "myPredMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "PredMort", "myPredMort")

Your function will then be called instead of mizerPredMort(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredRate(), mizerRDI(), mizerRates(), mizerResourceMort()

Get predation rate needed to project standard mizer model

Description

Calculates the potential rate (in units 1/year) at which a prey individual of a given size w is killed by predators from species j. In formulas

{\tt pred\_rate}_j(w_p) = \int \phi_j(w,w_p) (1-f_j(w)) \gamma_j(w) N_j(w) \, dw.

This potential rate is used in the function mizerPredMort() to calculate the realised predation mortality rate on the prey individual. You would not usually call this function directly but instead use getPredRate(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerPredRate(params, n, n_pp, n_other, t, feeding_level, ...)

Arguments

Value

A named two dimensional array (predator species x prey size) with the predation rate, where the prey size runs over fish community plus resource spectrum.

Your own predation rate function

By default getPredRate() calls mizerPredRate(). However you can replace this with your own alternative predation rate function. If your function is called "myPredRate" then you register it in a MizerParams object params with

params <- setRateFunction(params, "PredRate", "myPredRate")

Your function will then be called instead of mizerPredRate(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerRDI(), mizerRates(), mizerResourceMort()

Get density-independent rate of reproduction needed to project standard mizer model

Description

Calculates the density-independent rate of total egg production R_{di} (units 1/year) before density dependence, by species. You would not usually call this function directly but instead use getRDI(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerRDI(params, n, n_pp, n_other, t, e_growth, mort, e_repro, ...)

Arguments

Details

This rate is obtained by taking the per capita rate E_r(w)\psi(w) at which energy is invested in reproduction, as calculated by getERepro(), multiplying it by the number of individualsN(w) and integrating over all sizes w and then multiplying by the reproductive efficiency \epsilon and dividing by the egg size w_min, and by a factor of two to account for the two sexes:

R_{di} = \frac{\epsilon}{2 w_{min}} \int N(w) E_r(w) \psi(w) \, dw

Used by getRDD() to calculate the actual, density dependent rate. See setReproduction() for more details.

Value

A numeric vector with the rate of egg production for each species.

Your own reproduction function

By default getRDI() calls mizerRDI(). However you can replace this with your own alternative reproduction function. If your function is called "myRDI" then you register it in a MizerParams object params with

params <- setRateFunction(params, "RDI", "myRDI")

Your function will then be called instead of mizerRDI(), with the same arguments. For an example of an alternative reproduction function see constantEggRDI().

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRates(), mizerResourceMort()

Get all rates needed to project standard mizer model

Description

Calls other rate functions in sequence and collects the results in a list.

Usage

mizerRates(params, n, n_pp, n_other, t = 0, effort, rates_fns, ...)

Arguments

Details

By default this function returns a list with the following components:

• encounter from mizerEncounter()

• feeding_level from mizerFeedingLevel()

• e from mizerEReproAndGrowth()

• e_repro from mizerERepro()

• e_growth from mizerEGrowth()

• pred_rate from mizerPredRate()

• pred_mort from mizerPredMort()

• f_mort from mizerFMort()

• mort from mizerMort()

• rdi from mizerRDI()

• rdd from BevertonHoltRDD()

• resource_mort from mizerResourceMort()

However you can replace any of these rate functions by your own rate function if you wish, see setRateFunction() for details.

Value

List of rates.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerResourceMort()

Get predation mortality rate for resource needed to project standard mizer model

Description

Calculates the predation mortality rate \mu_p(w) on the resource spectrum by resource size (in units 1/year). You would not usually call this function directly but instead use getResourceMort(), which then calls this function unless an alternative function has been registered, see below.

Usage

mizerResourceMort(params, n, n_pp, n_other, t, pred_rate, ...)

Arguments

Value

A vector of mortality rate by resource size.

Your own resource mortality function

By default getResourceMort() calls mizerResourceMort(). However you can replace this with your own alternative resource mortality function. If your function is called "myResourceMort" then you register it in a MizerParams object params with

params <- setRateFunction(params, "ResourceMort", "myResourceMort")

Your function will then be called instead of mizerResourceMort(), with the same arguments.

See Also

Other mizer rate functions: mizerEGrowth(), mizerERepro(), mizerEReproAndGrowth(), mizerEncounter(), mizerFMort(), mizerFMortGear(), mizerFeedingLevel(), mizerMort(), mizerPredMort(), mizerPredRate(), mizerRDI(), mizerRates()

Determine whether a MizerParams or MizerSim object needs to be upgraded

Description

Looks at the mizer version that was used to last update the object and returns TRUE if changes since that version require an upgrade of the object. You would not usually have to call this function. Upgrades are initiated automatically by validParams and validSim when necessary.

Usage

needs_upgrading(object)

Arguments

Value

TRUE or FALSE

Set up parameters for a community-type model

Description

This functions creates a MizerParams object describing a community-type model. The function has many arguments, all of which have default values.

Usage

newCommunityParams(
  max_w = 1e+06,
  min_w = 0.001,
  no_w = 100,
  min_w_pp = 1e-10,
  z0 = 0.1,
  alpha = 0.2,
  f0 = 0.7,
  h = 10,
  gamma = NA,
  beta = 100,
  sigma = 2,
  n = 2/3,
  kappa = 1000,
  lambda = 2.05,
  r_pp = 10,
  knife_edge_size = 1000,
  reproduction
)

Arguments

Details

A community model has several features that distinguish it from a multi-species model:

• Species identities of individuals are ignored. All are aggregated into a single community.

• The resource spectrum only extends to the start of the community spectrum.

• Reproductive rate is constant, independent of the energy invested in reproduction, which is set to 0.

• Standard metabolism is turned off (the parameter ks is set to 0). Consequently, the growth rate is now determined solely by the assimilated food

Fishing selectivity is modelled as a knife-edge function with one parameter, knife_edge_size, which determines the size at which species are selected.

The resulting MizerParams object can be projected forward using project() like any other MizerParams object. When projecting the community model it may be necessary to keep a small time step size dt of around 0.1 to avoid any instabilities with the solver. You can check for these numerical instabilities by plotting the biomass or abundance through time after the projection.

Value

An object of type MizerParams

References

K. H. Andersen,J. E. Beyer and P. Lundberg, 2009, Trophic and individual efficiencies of size-structured communities, Proceedings of the Royal Society, 276, 109-114

See Also

Other functions for setting up models: newMultispeciesParams(), newSingleSpeciesParams(), newTraitParams()

Examples

params <- newCommunityParams()
sim <- project(params, t_max = 10)
plotBiomass(sim)
plotSpectra(sim, power = 2)

# More satiation. More mortality
params <- newCommunityParams(f0 = 0.8, z0 = 0.4)
sim <- project(params, t_max = 10)
plotBiomass(sim)
plotSpectra(sim, power = 2)

Set up parameters for a general multispecies model

Description

Sets up a multi-species size spectrum model by filling all slots in the MizerParams object based on user-provided or default parameters. There is a long list of arguments, but almost all of them have sensible default values. The only required argument is the species_params data frame. All arguments are described in more details in the sections below the list.

Usage

newMultispeciesParams(
  species_params,
  interaction = NULL,
  no_w = 100,
  min_w = 0.001,
  max_w = NA,
  min_w_pp = NA,
  pred_kernel = NULL,
  search_vol = NULL,
  intake_max = NULL,
  metab = NULL,
  p = 0.7,
  ext_mort = NULL,
  z0pre = 0.6,
  z0exp = n - 1,
  ext_encounter = NULL,
  maturity = NULL,
  repro_prop = NULL,
  RDD = "BevertonHoltRDD",
  kappa = 1e+11,
  n = 2/3,
  resource_rate = 10,
  resource_capacity = kappa,
  lambda = 2.05,
  w_pp_cutoff = 10,
  resource_dynamics = "resource_semichemostat",
  gear_params = NULL,
  selectivity = NULL,
  catchability = NULL,
  initial_effort = NULL,
  info_level = 3,
  z0 = deprecated(),
  r_pp = deprecated()
)

Arguments

Value

An object of type MizerParams

Species parameters

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_max column with the maximum sizes of the species in grams. (You could alternatively specify the maximum length in cm in an l_max column.)

The ⁠species_params dataframe⁠ also needs to contain the parameters needed by any predation kernel function (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. So the species_params data frame saved in the returned MizerParams object will differ from the one you supply because it will have the missing species parameters filled in with default values.

If you are not happy with any of the species parameter values used you can always change them later with species_params<-().

All the parameters will be mentioned in the following sections.

Setting initial values

The initial values for the species number densities are set using the function get_initial_n(). These are quite arbitrary and not very close to the steady state abundances. We intend to improve this in the future.

The initial resource number density N_R(w) is set to a power law with coefficient kappa (\kappa) and exponent -lambda (-\lambda):

N_R(w) = \kappa\, w^{-\lambda}

for all w less than w_pp_cutoff and zero for larger sizes.

Size grid

A size grid is created so that the log-sizes are equally spaced. The spacing is chosen so that there will be no_w fish size bins, with the smallest starting at min_w and the largest starting at max_w. For the resource spectrum there is a larger set of bins containing additional bins below min_w, with the same log size. The number of extra bins is such that min_w_pp comes to lie within the smallest bin.

Units in mizer

Mizer uses grams to measure weight, centimetres to measure lengths, and years to measure time.

Mizer is agnostic about whether abundances are given as

1. numbers per area,

2. numbers per volume or

3. 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.

Your 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 in plots, it will choose the units appropriate for case 3. So for example in plotBiomass() it gives the unit as grams.

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 function scaleModel() is useful to change the units you are using.

Setting interaction matrix

You do not need to specify an interaction matrix. If you do not, then the predator-prey interactions are purely determined by the size of predator and prey and totally independent of the species of predator and prey.

The interaction matrix \theta_{ij} modifies the interaction of each pair of species in the model. This can be used for example to allow for different spatial overlap among the species. The values in the interaction matrix are used to scale the encountered food and predation mortality (see on the website the section on predator-prey encounter rate and on predation mortality). The first index refers to the predator species and the second to the prey species.

The interaction matrix is used when calculating the food encounter rate in getEncounter() and the predation mortality rate in getPredMort(). Its entries are dimensionless numbers. If all the values in the interaction matrix are equal then predator-prey interactions are determined entirely by size-preference.

This function checks that the supplied interaction matrix is valid and then stores it in the interaction slot of the params object.

The order of the columns and rows of the interaction argument should be the same as the order in the species params data frame in the params object. If you supply a named array then the function will check the order and warn if it is different. One way of creating your own interaction matrix is to enter the data using a spreadsheet program and saving it as a .csv file. The data can then be read into R using the command read.csv().

The interaction of the species with the resource are set via a column interaction_resource in the species_params data frame. By default this column is set to all 1s.

Setting predation kernel

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 for 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 function box_pred_kernel(), and "power_law", implemented by the function power_law_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. Except for beta and sigma, 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 for pred_kernel_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).

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.

Setting search volume

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 chosen 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).

Setting maximum intake rate

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.

If h_i is set to Inf, fish of species i will consume all encountered food.

Setting metabolic rate

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) = k_s w^p + k w,

where k_s 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 k_s, 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

k_s = 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.

Setting external mortality rate

The external mortality is all the mortality that is not due to fishing or predation by predators included in the model. The external mortality could be due to predation by predators that are not explicitly included in the model (e.g. mammals or seabirds) or due to other causes like illness. It is a rate with units 1/year.

The ext_mort argument allows you to specify an external 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 setExtMort().

If the ext_mort argument is not supplied, then the external mortality is assumed to depend only on the species, not on the size of the individual: \mu_{ext.i}(w) = z_{0.i}. The value of the constant z_0 for each species is taken from the z0 column of the species parameter data frame, if that column exists. Otherwise it is calculated as

z_{0.i} = {\tt z0pre}_i\, w_{inf}^{\tt z0exp}.

Setting external encounter rate

The external encounter rate is the rate at which a predator encounters food that is not explicitly modelled. It is a rate with units mass/year.

The ext_encounter argument allows you to specify an external encounter rate that depends on species and body size. You can see an example of this in the Examples section of the help page for setExtEncounter().

Setting reproduction

For each species and at each size, the proportion \psi 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. There is a size w_repro_max at which all the energy is invested into reproduction and therefore all growth stops. There can be no fish larger than w_repro_max. If you have not specified the w_repro_max column in the species parameter data frame, then the maximum size w_max is used instead.

Maturity ogive

If the the proportion of individuals that are mature is not supplied via the maturity argument, 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, although each species has its own maturity ogive.) 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 maximum weight in the w_max 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_{mat25}).

The sigmoidal function given above would strictly reach 1 only asymptotically. Mizer instead sets the function equal to 1 already at a size taken from the w_repro_max column in the species parameter data frame, if it exists, or otherwise from the w_max column. Also, for computational simplicity, any proportion smaller than 1e-8 is set to 0.

Investment into reproduction

If the the energy available to a mature individual that is invested into reproduction is not supplied via the repro_prop argument, it is set to the allometric form

{\tt repro\_prop}(w) = \left(\frac{w}{w_{mat_max}}\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 maximum sizes are taken from the w_repro_max column in the species parameter data frame, if it exists, or otherwise from the w_max column.

The total proportion of energy invested into reproduction of an individual of size w is then

\psi(w) = {\tt maturity}(w){\tt repro\_prop}(w)

Reproductive efficiency

The reproductive efficiency \epsilon, i.e., the proportion of energy allocated to reproduction that results in egg biomass, is set through 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():

R_{di} = \frac{\epsilon}{2 w_{min}} \int N(w) E_r(w) \psi(w) \, dw

Density dependence

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, to encode all the density dependence in the stock-recruitment relationship that is not already included in the other two sources of density dependence, mizer puts the the density-independent rate of egg production through a density-dependence function. The result is returned by getRDD(). The name of the density-dependence function is specified by the RDD argument. The default is the Beverton-Holt function BevertonHoltRDD(), 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 RickerRDD() and SheperdRDD() and you can easily use these as models for writing your own functions.

Setting fishing

Gears

In mizer, fishing mortality is imposed on species by fishing gears. The total per-capita fishing mortality (1/year) 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.

Selectivity

The selectivity at size of each gear for each species is saved as a three dimensional array (gear x species x size). Each entry has a range between 0 (that gear is not selecting that species at that size) to 1 (that gear is selecting all individuals of that species of that size). This three dimensional array can be specified explicitly via the selectivity argument, but usually mizer calculates it from the gear_params slot of the MizerParams object.

To allow the calculation of the selectivity array, the gear_params slot must be a data frame with one row for each gear-species combination. So if for example a gear can select three species, then that gear contributes three rows to the gear_params data frame, one for each species it can select. The data frame must have columns gear, holding the name of the gear, species, holding the name of the species, and sel_func, holding the name of the function that calculates the selectivity curve. Some selectivity functions are included in the package: knife_edge(), sigmoid_length(), double_sigmoid_length(), and sigmoid_weight(). Users are able to write their own size-based selectivity function. The first argument to the function must be w and the function must return a vector of the selectivity (between 0 and 1) at size.

Each selectivity function may have parameters. Values for these parameters must be included as columns in the gear 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() that a 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 gear 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.

The most commonly-used selectivity function is sigmoid_length(). It has a smooth transition from 0 to 1 at a certain size. The sigmoid_length() function has the two parameters l50 and l25 that are the lengths in cm at which 50% or 25% of the fish are selected by the gear. If you choose this selectivity function then the l50 and l25 columns must be included in the gear parameters data.frame.

In case each species is only selected by one gear, the columns of the gear_params data frame can alternatively be provided as columns of the species_params data frame, if this is more convenient for the user to set up. Mizer will then copy these columns over to create the gear_params data frame when it creates the MizerParams object. However changing these columns in the species parameter data frame later will not update the gear_params data frame.

Catchability

Catchability is used as an additional factor 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.

Catchability is stored as a two dimensional array (gear x species). This can either be provided explicitly via the catchability argument, or the information can be provided via a catchability column in the gear_params data frame.

In the case where each species is selected by only a single gear, the catchability column can also be provided in the species_params data frame. Mizer will then copy this over to the gear_params data frame when the MizerParams object is created.

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.

Setting resource dynamics

You would usually set the resource dynamics only after having finished the calibration of the steady state. Then setting the resource dynamics with this function will preserve that steady state, unless you explicitly choose to set balance = FALSE. Your choice of the resource dynamics only affects the dynamics around the steady state. The higher the resource rate or the lower the resource capacity the less sensitive the model will be to changes in the competition for resource.

The resource_dynamics argument allows you to choose the resource dynamics function. By default, mizer uses a semichemostat model to describe the resource dynamics in each size class independently. This semichemostat dynamics is implemented by the function resource_semichemostat(). You can change that to use a logistic model implemented by resource_logistic() or you can use resource_constant() which keeps the resource constant or you can write your own function.

Both the resource_semichemostat() and the resource_logistic() dynamics are parametrised in terms of a size-dependent rate r_R(w) and a size-dependent capacity c_R. The help pages of these functions give the details.

The resource_rate argument can be a vector (with the same length as w_full(params)) specifying the intrinsic resource growth rate for each size class. Alternatively it can be a single number, which is then used as the coefficient in a power law: then the intrinsic growth rate r_R(w) at size w is set to

r_R(w) = r_R w^{n-1}.

The power-law exponent n is taken from the n argument.

The resource_capacity argument can be a vector specifying the intrinsic resource carrying capacity for each size class. Alternatively it can be a single number, which is then used as the coefficient in a truncated power law: then the intrinsic growth rate c_R(w) at size w is set to

c(w) = \kappa\, w^{-\lambda}

for all w less than w_pp_cutoff and zero for larger sizes. The power-law exponent \lambda is taken from the lambda argument.

The values for lambda, n and w_pp_cutoff are stored in a list in the resource_params slot of the MizerParams object so that they can be re-used automatically in the future. That list can be accessed with resource_params(). It also holds the coefficient kappa that describes the steady-state resource abundance.

See Also

Other functions for setting up models: newCommunityParams(), newSingleSpeciesParams(), newTraitParams()

Examples

params <- newMultispeciesParams(NS_species_params)

Set up parameters for a single species in a power-law background

Description

This functions creates a MizerParams object with a single species. This species is embedded in a fixed power-law community spectrum

N_c(w) = \kappa w^{-\lambda}

This community provides the food income for the species. Cannibalism is switched off. The predation mortality arises only from the predators in the power-law community and it is assumed that the predators in the community have the same feeding parameters as the foreground species. The function has many arguments, all of which have default values.

Usage

newSingleSpeciesParams(
  species_name = "Target species",
  w_max = 100,
  w_min = 0.001,
  eta = 10^(-0.6),
  w_mat = w_max * eta,
  no_w = log10(w_max/w_min) * 20 + 1,
  n = 3/4,
  p = n,
  lambda = 2.05,
  kappa = 0.005,
  alpha = 0.4,
  h = 30,
  beta = 100,
  sigma = 1.3,
  f0 = 0.6,
  fc = 0.25,
  ks = NA,
  gamma = NA,
  ext_mort_prop = 0,
  reproduction_level = 0,
  R_factor = deprecated(),
  w_inf = deprecated(),
  k_vb = deprecated()
)

Arguments