gadget-framework / gadget3

TMB-based gadget implemtation
GNU General Public License v2.0
8 stars 6 forks source link

g3_parameterized()-based defaults #118

Closed lentinj closed 11 months ago

lentinj commented 11 months ago

See https://github.com/gadget-framework/gadget3/pull/117 for more discussion - enough parametrized defaults so we can throw together a ~sensible model by just combining actions.

lentinj commented 11 months ago

@bthe says

Looks good, one thing that I can think of is that it might be worth reworking g3a_renewal_vonB so it would

a) use a more commonly used parametrisation of Linf, K and t_0 (instead of recl) b) adjust of the time of the initial conditions, as the mean length should be estimated at age - delta_t

@willbutler42 any suggestions on what this would look like? We might end up with 2 vonB functions though, which would be a bit annoying from a naming point of view.

bthe commented 11 months ago

I guess something like:

g3a_renewal_vonb <- function(
        Linf = g3_parameterized('Linf', by_stock = by_stock),
        K = g3_parameterized('K', by_stock = by_stock, scale = 0.001),
        t0 = g3_parameterized('t0', by_stock = by_stock),
        by_stock = TRUE) {
    f_substitute(
        quote( Linf * (1 - exp(-1 * K * (age - t0))) ),
        list(
            Linf = Linf,
            K = K,
            recl = t0))
}

would take care of my first bullet.

The second bullet could be hacked by adding deltat to t0, but it may be clearer to define an offset parameter to the function.

lentinj commented 11 months ago

Adding t0 = g3_parameterized('t0', by_stock = by_stock, offset = "deltat") should cover this (I think!).

But I think we'd need separate g3a_renewal_vonb_t0() & g3a_renewal_vonb_recl(). I guess the long-term ideal would be to switch the default g3a_renewal_vonb to being g3a_renewal_vonb_t0, but again that's a chunk of short-term pain.

bthe commented 11 months ago

Sound reasonable. I guess that you can go from g3a_renewal_vonb_t0 to g3a_renewal_vonb_recl by specifying t0 as (recage + log(1 - recl/Linf)/K)

lentinj commented 11 months ago

Okay, so I think the thing to do is:

So models that used g3a_renewal_vonb() explicitly won't break, which I think is more likely to be the case for older models. Models that don't specify mean_f OTOH will need updating.

We could eventually tidy up the alias if we felt it was worth it once the dust has settled.

lentinj commented 11 months ago

Also (sorry!) what do you both say to:

 g3a_renewal_normalparam <- function (
         stock,
         factor_f = g3_parameterized('rec',
             by_stock = by_stock,
             by_year = TRUE,
-            scale = g3_parameterized(name = 'rec.scalar', by_stock = by_stock),
+            scale = g3_parameterized(
+                name = 'rec.scalar',
+                by_stock = by_stock,
+                by_step = is.null(run_step)),
             ifmissing = NaN),
         mean_f = g3a_renewal_vonb(by_stock = by_stock),
         stddev_f = g3_parameterized('rec.sd', value = 10, by_stock = by_stock),
         alpha_f = g3_parameterized('walpha', by_stock = wgt_by_stock),
         beta_f = g3_parameterized('wbeta', by_stock = wgt_by_stock),
         by_stock = TRUE,
         wgt_by_stock = TRUE,
         run_age = quote(stock__minage),
         run_projection = FALSE,
         run_step = 1,
         run_f = NULL,
         run_at = 8) {

...i.e. if you set run_step = NULL, continuous recruitment, then you automatically get per-step rec.scalar parameters.

This would be needed for the anchovy continuous-renewal-apart-from-in-winter, but seems fairly generally applicable.

bthe commented 11 months ago

Looks fine to me. If one wants have different recruitment by year and step I assume you just do:

factor_f = g3_parameterized('rec',
             by_stock = by_stock,
             by_year = TRUE,
             scale = g3_parameterized(name = 'rec.scalar', by_stock = by_stock),
             by_step = TRUE)
lentinj commented 11 months ago

Exactly. Or we could make your suggestion the default if you specify continuous recruitment (although it would result in a fair few additional parameters).

Either would be useful I'm sure, I guess the question is which would "do the right thing" in more cases.

bthe commented 11 months ago

Exactly. Or we could make your suggestion the default if you specify continuous recruitment (although it would result in a fair few additional parameters).

Either would be useful I'm sure, I guess the question is which would "do the right thing" in more cases.

I think the timestep specific scalar is more appropriate, unless you have too much data lying around (I wish). But I think that it would good to have this as an example in the help page for renewal.