Change model parameters without using @

Mizer provides dedicated functions for changing model parameters. Using them protects you from pitfalls arising from manipulations with the @ symbol.

Gustav Delius true
08-08-2021

Introduction

Often, after creating a mizer model, you will want to make changes to the model parameters, either to improve the model or to investigate the effect of changes. For example you might want to study the community consequence of changes in the physiology due to warmer waters, or whatever your research question is.

Mizer makes this very easy. However, it is important to use the right syntax. In particular, one should avoid the use of the @ symbol, as I will explain in this blog post. So this post serves as a reminder of how one can make changes to a model after it has been created.

change

As you will know, all the parameters describing a mizer model are contained in an object of class MizerParams. You may either have created this object yourself, for example with the newMultispeciesParams() function, or received it from a colleague. Such an object has many slots holding the various pieces of information about the model.

There are often two ways of accessing these slots. One involves the use of the @ symbol, the other uses a dedicated mizer function. If you look at the mizer code, you will see a lot of @ symbols. However I will discuss the potential pitfalls when using the @ symbol and advocate for the use of dedicated functions that mizer provides for the purpose of changing model parameters.

This blog post was written with mizer version 2.2.1.9001

library(mizer)
packageVersion("mizer")

[1] '2.2.1.9001'

Throughout this post we will use the NS_params MizerParams object that comes as an example with the mizer package.

params <- NS_params

Species parameters

Let’s start with the species_params data frame. This holds species-specific parameters that are used by mizer to calculate physiological rates according to specific assumptions about the size-dependence of these rates. (The details do not concern us here, but you can find all the information on the help page of species_params(). These size-dependent rates are then stored as large arrays in the slots of the MizerParams object, to be efficiently used during model projections. But also the species parameters are stored, and can therefore be used to recalculate the rates when some species parameters change.

You can get the species parameters out of a MizerParams object in two ways, both giving identical results:

# Using function
sp <- species_params(params)
# Using @
spa <- params@species_params
# These are identical
identical(sp, spa)

[1] TRUE

The fact that these are identical is not surprising if you look at the code for the species_params() function:

species_params

function(params) {
    params@species_params
}
<bytecode: 0x55ea4c495ba8>
<environment: namespace:mizer>

The function actually only contains one statement, which accesses the species_params slot using the @ notation.

The difference between the two notations becomes apparent only when you want to make a change to the species parameters. Let’s assume we want to change the h parameter for the first species in our example model. It currently has the value

species_params(params)$h[1]

[1] 18.20276

As you know, this parameter h is used to calculate the maximum intake rate \(h(w)\) as a power-law function of size: \(h(w) = h w^n\). So for example the maximum intake rate at the smallest size is

getMaxIntakeRate(params)[1, 1]

[1] 0.1820276

Now let’s increase the value of the h parameter using the @ notation:

params@species_params$h[1] <- 20

If we look at the maximum intake rate, we see that it has not changed:

getMaxIntakeRate(params)[1, 1]

[1] 0.1820276

All we have done is change the value in the species parameter data frame, but this did not trigger a recalculation of the maximum intake rate. We should instead have used

species_params(params)$h[1] <- 20

Take a look at the syntax, which is really a bit weird when compared to other programming languages. If you want to dig deeper into this, a good place to look is https://adv-r.hadley.nz/functions.html#replacement-functions.

This way of changing \(h\) does indeed change the maximum intake rate:

getMaxIntakeRate(params)[1, 1]

[1] 0.2

The reason is clear if we look at the code:

`species_params<-`

function(params, value) {
    value <- validSpeciesParams(value)
    params@species_params <- value
    suppressMessages(setParams(params))
}
<bytecode: 0x55ea4fdb2030>
<environment: namespace:mizer>

So three things actually happen when you change a species parameter via the setter function:

  1. Your new value is checked for validity
  2. Your new value is saved in the species_params slot
  3. The other slots in the MizerParams object are updated by calling setParams().

Gear parameters

Similar comments apply to the gear parameters. The gear parameters are used by mizer to set up the catchability and selectivity arrays. You can find more details on the help page for setting fishing. You can get the current gear parameters in two equivalent ways:

# Using function
gp <- gear_params(params)
# Using @
gpa <- params@gear_params
# These are identical
identical(gp, gpa)

[1] TRUE

But if you want to actually change the selectivity or catchability by changing the gear parameters you need to use the functional form. Here is the current gear_params data frame in the example model:

gear_params(params)

         gear species   sel_func knife_edge_size catchability
1  Industrial   Sprat knife_edge              13            1
2  Industrial Sandeel knife_edge               4            1
3  Industrial  N.pout knife_edge              23            1
4     Pelagic Herring knife_edge              99            1
5        Beam     Dab knife_edge              21            1
6       Otter Whiting knife_edge              75            1
7        Beam    Sole knife_edge              78            1
8       Otter Gurnard knife_edge              39            1
9        Beam  Plaice knife_edge             105            1
10      Otter Haddock knife_edge             165            1
11      Otter     Cod knife_edge            1606            1
12      Otter  Saithe knife_edge            1076            1

If, for example, we want to reduce the catchability of Sprat with the Industrial gear to 0.8 we would do

gear_params(params)$catchability[1] <- 0.8

Note that changing gear parameters in the species_params data frame will not have the desired effect. You need to change them in the gear_params data frame.

Resource parameters

Not surprisingly, the same applies to the resource parameters. These are used to set up the size-dependent carrying capacity and replenishment rate for the resource. You should access them with

resource_params(params)

$kappa
[1] 1e+11

$lambda
[1] 2.133333

$r_pp
[1] 10

$n
[1] 0.6666667

$w_pp_cutoff
[1] 9.820907

and change them with, for example,

resource_params(params)$r_pp <- 4

Rate arrays

When you call newMultispeciesParams(), then mizer uses the information in species_params, gear_params and resource_params to set up various arrays that will later make it much faster to run simulations of the model. If you are not happy with how mizer fills these arrays, you can also change them directly. And again you can do that either with @ notation or without.

Let’s take the example of the maximum intake rate \(h(w)\) that we already discussed earlier. This is stored in the intake_max slot of the MizerParams object, as a two-dimensional array, with one row for each species and one column for each size bin. You can get at this way in the two equivalent ways:

identical(getMaxIntakeRate(params), params@intake_max)

[1] TRUE

You may wonder, why the function is called getMaxIntakeRate() rather than intake_max(), and I am wondering too. Naming things is difficult, and I now think that I made a bad choice when choosing those names. It gets worse when we now look at the syntax for changing intake_max.

We can either do it with the @ notation, for example

params@intake_max <- 2 * params@intake_max

or we can do it without @ notation

params <- setMaxIntakeRate(params, intake_max = 2 * getMaxIntakeRate(params))

This is ugly, and in future versions of mizer I think we will also allow

intake_max(params) <- 2 * intake_max(params)

But for now we are stuck with the functions that are all listed on the [help page] for setParams(). Note that these functions do not modify the params object in place, but create a new MizerParams object, which we then have to assign to a variable.

Again there are benefits in avoiding directly accessing the slot with @. These are:

  1. The new value you assign is checked for validity. If you make a mistake in an assignment using @ you will not get any warning and instead will run into mysterious and cryptic error messages later.

  2. The new value gets protected from automatically being overwritten when you make changes to other parameters.

We illustrate the second point with a simple example. Let’s set one entry in the metab slot, which holds the size-dependent metabolic rates, to a particular value.

params@metab[1, 1] <- 2
params@metab[1, 1]

[1] 2

Let us then make a totally unrelated change, say by changing the reproductive efficiency of the 5th species and look again at our entry in metab.

species_params(params)$erepro[5] <- 0.1
params@metab[1, 1]

[1] 0.02891793

The reason for this is that the change in the species parameter has triggered a recalculation of the rate arrays from the species parameters, overwriting our manual change.

Now let’s try the same with the proper way of changing the metabolic rate.

metab <- getMetabolicRate(params)
metab[1, 1] <- 2
params <- setMetabolicRate(params, metab = metab)

Now this will not get overwritten when some other parameter changes.

species_params(params)$erepro[5] <- 0.1
params@metab[1, 1]

[1] 2

The way that was done internally is by attaching a comment to the metab slot.

comment(params@metab)

[1] "set manually"

You could have chosen a more informative comment, for example

comment <- "Just changed the [1, 1] entry for test purposes."
params <- setMetabolicRate(params, metab = metab,
                           comment_metab = comment)

All that matters is that there is a comment. This also tells us how we can un-protect a slot so that it can be auto-computed from the species parameters again:

comment(params@metab) <- NULL

species_params(params)$erepro[5] <- 0.1
params@metab[1, 1]

[1] 0.02891793

Summary

We have seen how to change species parameters, gear parameters, resource parameters or other slots in a MizerParams object by using the appropriate functions. You can find a complete list of these functions in the mizer reference pages. We discussed how this avoids the pitfalls that arise when accessing slots directly with the @ notation. I hope this will be useful to you when you explore your own mizer model.

This blog post was motivated by a question by Leslie Garay-Narváez. Please keep the questions coming.

Corrections

If you see mistakes or want to suggest changes, please create an issue on the source repository.

Reuse

Text and figures are licensed under Creative Commons Attribution CC BY 4.0. Source code is available at https://github.com/sizespectrum/MizerBlog/, unless otherwise noted. The figures that have been reused from other sources don't fall under this license and can be recognized by a note in their caption: "Figure from ...".

Citation

For attribution, please cite this work as

Delius (2021, Aug. 8). mizer blog: Change model parameters without using @. Retrieved from https://blog.mizer.sizespectrum.org/posts/2021-08-08-change-model-parameters-without-using/

BibTeX citation

@misc{delius2021change,
  author = {Delius, Gustav},
  title = {mizer blog: Change model parameters without using @},
  url = {https://blog.mizer.sizespectrum.org/posts/2021-08-08-change-model-parameters-without-using/},
  year = {2021}
}