setBevertonHolt() does not respect erepro input

[patricksykes]

It is currently not possible to guarantee that setBevertonHolt() (and functions that rely on it, such as newTraitParams()) will respect the user's inputs for erepro and R_max. This is particularly an issue for erepro, which setBevertonHolt() is liable to set to unrealistic values greater than 1 despite the user specifying a value of 1 or less.

The documentation says

The values for erepro must be large enough to allow the required reproduction rate. If a smaller value is requested a warning is issued and the value is increased to the smallest possible value. The values for erepro should also be smaller than 1 to be physiologically sensible, but this is not enforced by the function.

However, I believe it should be possible to enforce this requirement.

[gustavdelius]

Hi Patrick. Thank you for raising this issue.

I think the issue is with newTraitParams() rather than with setBevertonHolt(). The setBevertonHolt() function is meant to be used when one already has a model with a desired steady state and should allow one to change the Beverton Holt curve without changing the steady state. Thus it does not enable the user to set erepro arbitrarily.

The real issue is that newTraitParams() by default imposes a reproduction level of 1/4, irrespective of whether such a reproduction level can actually be achieved, given the other parameters chosen by the user, without producing unrealistic values for erepro. Probably its behaviour should be changed so that it automatically uses a higher reproduction level if that is required to keep erepro realistic.

However, for now, you can work around this issue by specifying a higher reproduction level by hand via the reproduction_level argument to newTraitParams().

This may not address your true concern, so feel free to share more details.

[patricksykes]

Hi Gustav,

I see I must have been using the setBevertonHolt() function in the incorrect context - to set up the R_max parameter for a new ZooMizer run, rather than applying to an existing steady-state model. I still think though that the function should at least allow the user to specify that the value of erepro must be physically sensible, that is, less than or equal to 1. In my case, I specified erepro = 1 in the function call, and ended up with values ranging from about 5 to 200000.

For a little more context, I want to find reasonable erepro and R_max values for a range of chlorophyll inputs to ZooMizer. Having found that, as you said, newTraitParams() sets unrealistic values for erepro (even for quite high values of reproduction_level), I resolved to set erepro to 1 and look to tune the R_max using setBevertonHolt(). Is there another way to do so?

For now, I am setting the erepro and R_max values in ZooMizer by hand (to 1, and by an allometric relationship from the North Sea model: kappa * w_inf^(-1.5), respectively) but am finding that reproduction is higher than I would like for low chlorophyll and a too low at high chlorophyll relative to the biomass of the zooplankton resource. When I say reproduction is too high/low, I mean that fish eggs vastly outnumber the zooplankton spectrum for low chlorophyll (where the zooplankton spectrum is already relatively low) and the reproduction level is very low (on the order of 1e-10) and conversely at high chlorophyll, the reproduction level is sitting at or just below 1.

Do you have any suggestions on how to tune the R_max for different resource inputs, if setBevertonHolt() is not the right function to use here?

[gustavdelius]

Let us back up a bit and reflect on what the reproduction parameters are for. There is a two-dimensional space of parameter values. Itt is probably better to parametrise that no by erepro and R_max but by the rate of egg production and the reproduction level. The rate of egg production determines the abundance of your species. The reproduction level determines how sensitive the model is to perturbations. I calibrate a model in two steps. In the first step I get the steady state to agree with observations, in particular to have the abundances that are observed. For that I only need to tune the rate of egg production in that steady state. In the second step I get the model to have the desired sensitivity to changes in for example fishing pressure. For that I only need to change the reproduction level. If I do that using setBevertonHolt() then that second step does not mess up my first step because it does not affect the steady state.

Now I believe that you are working on the first step of getting an acceptable steady state. At that step I would set the reproduction level to 1 (which is equivalent to R_max = Inf). If in the process of getting your fish species to have the observed abundances you get an erepro > 1 then that means that there are problems elsewhere in your model. It means the species as modelled are just not viable. Perhaps they are not getting enough food or perhaps they experience too high mortality.

To answer your last question: don't tune R_max at all but leave it at Inf (set reproduction_level = 0). Use erepro to adjust the abundances of your species (with species_params(params)$erepro <- ). When an erepro gets above 1 for some species, given them more to eat or reduce their mortality.

[patricksykes]

Thank you for sharing your workflow for calibrating a model. I wonder how that differs in the case of ZooMizer where there is no data to fit to, rather we are interested in the kind of fish community that arises under different environmental conditions (and, in particular, differently-structured zooplankton communities). The current situation is that not all of the trait-based model species are surviving under all conditions. Often the planktivores lack food under low nutrient conditions, and in general a group or two will go extinct under most conditions I believe due to predation pressure and the competition for resources.

I will try your suggestion of first setting the reproduction level first to 1 (or rather 0.9999, since reproduction_level = 1 throws an error - is this expected behaviour?) and then adjusting the abundances with erepro.

[gustavdelius]

Hi Patrick. I am sorry, I wrote reproduction_level = 1 where I meant reproduction_level = 0. I have now edited the post to correct that.

Now I understand better the issue that you are facing. My previous comments would have been relevant only if you were intending to calibrate the model to an observed status quo and then investigate how changes in plankton would affect that status quo. However you are in the situation that you don't know the status quo but want some generic results and are frustrated by the fact that newTraitParams() does not set up good default values for that purpose. I am afraid I don't have good advice. But your scenario means that you will have no use for setBevertonHolt() (which should be used if you have a calibrates steady state that you want to preserve) but will want to set erepro and R_max directly.

[patricksykes]

Hi Gustav,

I have run the model again with reproduction_level = 0 and a higher value for alpha (the value used in Ken Andersen's book is 0.6). Combined, these allow the erepro to come out as less than 1. Ideally, I would like to have some density-dependence in reproduction though. Supposing I start with reproduction_level = 0 and obtain a steady state that's working for me, would that be the time to use setBevertonHolt(), specifying a reproduction_level greater than 0?

[gustavdelius]

Yes, but if you raise it too high (or lower R_max to far), then erepro will get above 1 again.

[gustavdelius]

I have transferred the issue of unrealistic erepro being produced by newTraitParams() to issue #249 and am closing this issue.