Add spawning to a g3 model

g3a_spawn_recruitment_fecundity(
        p0 = g3_parameterized('spawn.p0', value = 1, by_stock = by_stock),
        p1 = g3_parameterized('spawn.p1', value = 1, by_stock = by_stock),
        p2 = g3_parameterized('spawn.p2', value = 1, by_stock = by_stock),
        p3 = g3_parameterized('spawn.p3', value = 1, by_stock = by_stock),
        p4 = g3_parameterized('spawn.p4', value = 1, by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn_recruitment_simplessb(
        mu = g3_parameterized('spawn.mu', by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn_recruitment_ricker(
        mu = g3_parameterized('spawn.mu', by_stock = by_stock),
        lambda = g3_parameterized('spawn.lambda', by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn_recruitment_bevertonholt(
        mu = g3_parameterized('spawn.mu', by_stock = by_stock),
        lambda = g3_parameterized('spawn.lambda', by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn_recruitment_bevertonholt_ss3(
        # Steepness parameter
        h = g3_parameterized('spawn.h', lower = 0.1, upper = 1, value = 0.5,
            by_stock = by_stock ),
        # Recruitment deviates
        R = g3_parameterized('spawn.R', by_year = TRUE, exponentiate = TRUE,
            # Unfished equilibrium recruitment
            scale = "spawn.R0",
            by_stock = by_stock),
        # Unfished equilibrium spawning biomass (corresponding to R0)
        B0 = g3_parameterized('spawn.B0', by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn_recruitment_hockeystick(
        r0 = g3_parameterized('spawn.r0', by_stock = by_stock),
        blim = g3_parameterized('spawn.blim', value = 1, by_stock = by_stock),
        by_stock = TRUE )

g3a_spawn(
        stock,
        recruitment_f,
        proportion_f = 1,
        mortality_f = 0,
        weightloss_f = 0,
        weightloss_args = list(),
        output_stocks = list(),
        output_ratios = rep(1 / length(output_stocks), times = length(output_stocks)),
        mean_f = g3a_renewal_vonb_t0(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_step = NULL,
        run_f = ~TRUE,
        run_at = g3_action_order$spawn,
        recruit_at = g3_action_order$renewal)

Arguments

p0,p1,p2,p3,p4

Substituted into g3a_spawn_recruitment_fecundity formula, see below.

mu,lambda,h,R,B0,r0,blim

Substituted into g3a_spawn_recruitment_* formula, see below.

stock

The mature g3_stock that will spawn in this action.

recruitment_f

A list of formula generated by one of the g3a_spawn_recruitment_* functions, containing

s

Formula run for each subset of stock

r

Final formula for calculating number of recruits for spawning action

proportion_f

formula generated by one of the g3_suitability_* functions, describing the proportion of stock that will spawn at this timestep.

mortality_f

formula generated by one of the g3_suitability_* functions, describing the proportion of spawning stock that will die during spawning.

weightloss_f

formula generated by one of the g3_suitability_* functions, describing the overall weight loss during spawning.

DEPRECATED: Use weightloss_args for new models.

weightloss_args

list of options to pass to g3a_weightloss, e.g. rel_loss & abs_loss. If not empty, a weightloss action will be included as part of spawning.

output_stocks

List of g3_stocks that will be spawned into.

output_ratios

Vector of proportions for how to distribute into output_stocks, summing to 1, default evenly spread.

mean_f,stddev_f,alpha_f,beta_f

formula substituted into stock structure calculations, see g3a_renewal_normalparam for details.

run_step

Which step to perform renewal in, or NULL for continuous spawning. Adds cur_step == (run_step) into default run_f.

run_f

formula specifying a condition for running this action, default always runs.

run_at

Integer order that spawning actions will be run within model, see g3_action_order.

recruit_at

Integer order that recruitment from spawning will be run within model, see g3_action_order.

by_stock,wgt_by_stock

Controls how parameters are grouped, see g3_parameterized

Details

To restrict spawning to a particular step in a year, or a particular area, use run_f. For example:

cur_step == 1

Spawning will happen on first step of every year

cur_step == 1 && cur_year >= 1990

Spawning will happen on first step of every year after 1990

cur_step == 2 && area = 1

Spawning will happen on second step of every year, in the first area

The action will define the following stock instance variables for each given stock and output_stock:

stock__spawnprop

Proportion of (stock) that are spawning in this spawning event

stock__spawningnum

Numbers of (stock) that are spawning in this spawning event

output_stock__spawnednum

Numbers of (output_stock) that will be produced in this spawning event

Value

g3a_spawn_recruitment_fecundity

A pair of formula objects: $$ S = l ^{p_{1}} a^{p_{2}} (p N_{al})^{p_{3}} W_{al}^{p_{4}} $$ $$ R = p_{0} S $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

\(p_{0..4}\)

Arguments provided to function

g3a_spawn_recruitment_simplessb

A pair of formula objects: $$ S = N_{al} p W_{al} $$ $$ R = \mu S $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

Argument provided to function

g3a_spawn_recruitment_ricker

A pair of formula objects: $$ S = N_{al} p W_{al} $$ $$ R = \mu S e^{-\lambda S} $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

Argument provided to function

Argument provided to function

g3a_spawn_recruitment_bevertonholt

A pair of formula objects: $$ S = N_{al} p W_{al} $$ $$ R = \frac{\mu S}{\lambda + S} $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

Argument provided to function

Argument provided to function

g3a_spawn_recruitment_bevertonholt_ss3

A modified beverton-holt implementation from SS3 returning a pair of formula objects: $$ S = N_{al} p W_{al} $$ $$ R = \frac{ 4 h R_0 s R }{ B_0 (1 - h) + S (5 h - 1) } $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

\(h\)

Steepness parameter, by default provided by the srr_h parameter

\(R_0\)

Unfished equilibrium recruitment, by default provided by the R0 parameter

\(R\)

Recruitment deviates, by default provided by the R parameter table

\(B_0\)

Unfished equilibrium spawning biomass (corresponding to R0), by default provided by the B0 parameter

Argument provided to function

g3a_spawn_recruitment_hockeystick

A pair of formula objects: $$ S = N_{al} p W_{al} $$ $$ R = R_0 \min{( S / B_{lim}, 1)} $$

\(N_{al}\)

Number of parent stock

\(W_{al}\)

Weight of parent stock

\(p\)

Proportion of parent stock spawning, from proportion_f

\(R_0\)

Argument r0 provided to function

\(B_{lim}\)

Argument blim provided to function

NB: This formula is differentiable, despite using min() in the definition above.

g3a_spawn

An action (i.e. list of formula objects) that will, for the given stock...

  1. Use proportion_f to calculate the total parent stock that will spawn

  2. Use recruitment_f to derive the total newly spawned stock

  3. Apply weightloss and mortality_f to the parent stock

... then, at recruitment stage ...

  1. Recruit evenly into output_stocks, using mean_f, stddev_f, alpha_f, beta_f as-per g3a_renewal_normalparam

Examples

ling_imm <- g3_stock(c('ling', maturity = 'imm'), seq(20, 156, 4)) |> g3s_age(0, 10)
ling_mat <- g3_stock(c('ling', maturity = 'mat'), seq(20, 156, 4)) |> g3s_age(3, 10)

actions <- list(
    g3a_time(1990, 1994, c(6, 6)),
    g3a_initialconditions_normalcv(ling_imm),
    g3a_initialconditions_normalcv(ling_mat),
    g3a_age(ling_imm),
    g3a_age(ling_mat),

    g3a_spawn(
        # Spawn from ling_mat
        ling_mat,
        # Use Ricker Recruitment Function to calculate # of recruits from total biomass
        recruitment_f = g3a_spawn_recruitment_ricker(),
        # Proportion of ling_mat spawning exponential relationship based on length
        proportion_f = g3_suitability_exponentiall50(
            alpha = g3_parameterized("spawn.prop.alpha", value = 4, scale = -1),
            l50 = g3_parameterized("spawn.prop.l50", value = 60)),
        # Proportion of ling_mat dying during spawning linear relationship to length
        mortality_f = g3_suitability_straightline(
            alpha = g3_parameterized("spawn.mort.alpha"),
            beta = g3_parameterized("spawn.mort.beta")),
        # Weight of spawning ling_imm should reduce by a fixed absolute amount (see g3a_weightloss)
        weightloss_args = list( abs_loss = g3_parameterized("spawn.weightloss", value = 0.1) ),
        # Spawn into ling_imm
        output_stocks = list(ling_imm),
        # Spawning should happen on the first step of every year
        run_f = ~cur_step==1 ),
    NULL )
model_fn <- g3_to_r(c(actions,
    g3a_report_detail(actions),
    g3a_report_history(actions, "__spawningnum$|__offspringnum$") ))

attr(model_fn, "parameter_template") |>
    # g3a_initialconditions_normalcv()
    g3_init_val("*.Linf", 50) |>
    g3_init_val("*.t0", -1.4) |>
    g3_init_val("*.walpha", 0.1) |>
    g3_init_val("*.wbeta", 1) |>
    # g3a_spawn_recruitment_ricker()
    g3_init_val("*.spawn.mu", 1e6) |>
    g3_init_val("*.spawn.lambda", 30) |>
    identity() -> params

r <- attributes(model_fn(params))
colSums(r$dstart_ling_imm__num)
#>        time
#> age     1990-01 1990-02 1991-01 1991-02 1992-01 1992-02 1993-01 1993-02 1994-01
#>   age0    10000   10000       0       0       0       0       0       0       0
#>   age1    10000   10000   10000   10000       0       0       0       0       0
#>   age2    10000   10000   10000   10000   10000   10000       0       0       0
#>   age3    10000   10000   10000   10000   10000   10000   10000   10000       0
#>   age4    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age5    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age6    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age7    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age8    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age9    10000   10000   10000   10000   10000   10000   10000   10000   10000
#>   age10   10000   10000   20000   20000   30000   30000   40000   40000   50000
#>        time
#> age     1994-02
#>   age0        0
#>   age1        0
#>   age2        0
#>   age3        0
#>   age4    10000
#>   age5    10000
#>   age6    10000
#>   age7    10000
#>   age8    10000
#>   age9    10000
#>   age10   50000
colSums(r$dstart_ling_mat__wgt)
#>        time
#> age     1990-01  1990-02  1991-01  1991-02  1992-01  1992-02  1993-01  1993-02
#>   age3      315 66.41122 66.41122  0.00000  0.00000  0.00000  0.00000  0.00000
#>   age4      315 67.75073 66.41122 64.24502  0.00000  0.00000  0.00000  0.00000
#>   age5      315 68.61062 67.75073 64.59135 64.24502 62.90600  0.00000  0.00000
#>   age6      315 68.97281 68.61062 65.18702 64.59135 62.84557 62.90600 61.74642
#>   age7      315 69.11103 68.97281 65.50674 65.18702 63.15756 62.84557 61.59535
#>   age8      315 69.16247 69.11103 65.63896 65.50674 63.37422 63.15756 61.74004
#>   age9      315 69.18147 69.16247 65.68958 65.63896 63.47144 63.37422 61.87163
#>   age10     315 69.18846 67.91383 66.01110 65.57397 64.29601 63.98007 62.88019
#>        time
#> age      1994-01  1994-02
#>   age3   0.00000  0.00000
#>   age4   0.00000  0.00000
#>   age5   0.00000  0.00000
#>   age6   0.00000  0.00000
#>   age7  61.74642 60.66980
#>   age8  61.59535 60.53587
#>   age9  61.74004 60.61440
#>   age10 62.64468 61.60524