R Under development (unstable) (2025-10-13 r88918 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2025 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> if (!interactive()) options(warn=2, error = function() { sink(stderr()) ; traceback(3) ; q(status = 1) })
> library(unittest)
> 
> library(gadget3)
> 
> st_imm <- g3_stock(c("stst", maturity = "imm"), c(10, 20, 30)) |> g3s_age(0, 5)
> st_mat <- g3_stock(c("stst", maturity = "mat"), c(10, 20, 30)) |> g3s_age(3, 15)
> stocks_st <- list(st_imm, st_mat)
> fl <- g3_fleet(c("fl", "surv"))
> 
> actions <- list(
+     g3a_time(1990, 1994, c(6,6)),
+     gadget3:::g3a_initialconditions_manual(st_imm,
+         quote( 100 + stock__minlen ),
+         quote( 1e4 + 0 * stock__minlen ) ),
+     gadget3:::g3a_initialconditions_manual(st_mat,
+         quote( 100 + stock__minlen ),
+         quote( 1e4 + 0 * stock__minlen ) ),
+     g3a_age(st_imm),
+     g3a_age(st_mat),
+ 
+     g3a_spawn(
+         st_mat,
+         g3a_spawn_recruitment_hockeystick(
+             r0 = g3_param_project(
+                 "rec",
+                 g3_param_project_logar1(),
+                 random = FALSE,
+                 scale = "rec.scalar",
+                 by_stock = stocks_st,
+                 by_step = FALSE )),
+             output_stocks = list(st_imm),
+             run_step = 1 ),
+ 
+     # NB: Only required for testing
+     gadget3:::g3l_test_dummy_likelihood() )
> full_actions <- c(actions, list(
+     g3a_report_detail(actions),
+     g3a_report_history(actions, 'proj_.*', out_prefix = NULL),
+     NULL))
> model_fn <- g3_to_r(full_actions)
> model_cpp <- g3_to_tmb(full_actions)
> 
> ok_group("No noise, fixed loglevel") ##########################################
# No noise, fixed loglevel
> 
> attr(model_cpp, 'parameter_template') |>
+     g3_init_val("stst.rec.#", rnorm(5, 1e5, 500)) |>
+     g3_init_val("stst.rec.proj.logar1.stddev", 1e-20) |>  # i.e. no noise
+     g3_init_val("stst.rec.proj.logar1.phi", 0.8) |>
+     g3_init_val("stst.rec.proj.logar1.level", round(runif(1, 10, 20))) |>
+     g3_init_val("stst_mat.spawn.blim", 1e2) |>  # blim too low to trigger
+     g3_init_val("proj_logar1_stst_weight.proj_logar1_rec_weight", 0) |>  # i.e. disable nll output (it'll be Inf)
+ 
+     g3_init_val("*.K", 0.3, lower = 0.04, upper = 1.2) |>
+     g3_init_val("*.Linf", max(g3_stock_def(st_imm, "midlen")), spread = 0.2) |>
+     g3_init_val("*.t0", g3_stock_def(st_imm, "minage") - 0.8, spread = 2) |>
+     g3_init_val("*.walpha", 0.01, optimise = FALSE) |>
+     g3_init_val("*.wbeta", 3, optimise = FALSE) |>
+ 
+     g3_init_val("project_years", 100) |>
+     identity() -> params.in
> nll <- model_fn(params.in) ; r <- attributes(nll) ; nll <- as.vector(nll)
> 
> ok(ut_cmp_equal(
+     as.vector(tail(r$proj_logar1_stst_rec__lvar, 10)),
+     rep(log(params.in["stst.rec.proj.logar1.level", "value"][[1]]), 10),
+     end = NULL ), "proj_logar1_stst_rec__lvar: Settles to loglevel in projection, regardless of initial value")
ok - proj_logar1_stst_rec__lvar: Settles to loglevel in projection, regardless of initial value
> ok(ut_cmp_equal(
+     as.vector( g3_array_agg(r$detail_stst_imm__spawnednum, c("year"), step = 1, age = 0) ),
+     as.vector( exp(r$proj_logar1_stst_rec__lvar) ),
+     end = NULL), "r$detail_stst_imm__spawnednum: projection variable used for recruitment")
ok - r$detail_stst_imm__spawnednum: projection variable used for recruitment
> 
> ok(ut_cmp_equal(nll, 0), "nll: 0, as we disabled it with weight")
ok - nll: 0, as we disabled it with weight
> ok(all(!is.finite(r$proj_logar1_stst_rec__nll)), "r$proj_logar1_stst_rec__nll: lstddev 0, so nll infinite")
ok - r$proj_logar1_stst_rec__nll: lstddev 0, so nll infinite
> 
> gadget3:::ut_tmb_r_compare2(model_fn, model_cpp, params.in)
# skip: not running TMB tests 0  <  1
NULL
> 
> ok_group("With noise, no projection") #########################################
# With noise, no projection
> 
> old_seed <- .Random.seed
> set.seed(1234)  # Fix seed so we always choose the same stst.rec.#
> attr(model_cpp, 'parameter_template') |>
+     g3_init_val("stst.rec.#", rnorm(5, 1e5, 500)) |>
+     g3_init_val("stst.rec.proj.logar1.stddev", 0.0003) |>
+     g3_init_val("stst.rec.proj.logar1.phi", 0.8) |>
+     g3_init_val("stst_mat.spawn.blim", 1e2) |>  # blim too low to trigger
+ 
+     g3_init_val("*.K", 0.3, lower = 0.04, upper = 1.2) |>
+     g3_init_val("*.Linf", max(g3_stock_def(st_imm, "midlen")), spread = 0.2) |>
+     g3_init_val("*.t0", g3_stock_def(st_imm, "minage") - 0.8, spread = 2) |>
+     g3_init_val("*.walpha", 0.01, optimise = FALSE) |>
+     g3_init_val("*.wbeta", 3, optimise = FALSE) |>
+ 
+     g3_init_val("project_years", 0) |>
+     identity() -> params.in
> .Random.seed <- old_seed
> nll <- model_fn(params.in) ; r <- attributes(nll) ; nll <- as.vector(nll)
> 
> ok(ut_cmp_equal(nll, sum(r$proj_logar1_stst_rec__nll)), "nll: Consistent with r$proj_logar1_stst_rec__nll (sums to same values)")
ok - nll: Consistent with r$proj_logar1_stst_rec__nll (sums to same values)
> ok(ut_cmp_equal(
+     as.vector(r$proj_logar1_stst_rec__nll),
+     as.vector(c(0, -dnorm(
+         tail(log(r$proj_logar1_stst_rec__lvar), -1) -
+         0.8 * head(log(r$proj_logar1_stst_rec__lvar), -1) -
+         0.2 * log(params.in["stst.rec.proj.logar1.level", "value"][[1]]),
+         0 - exp(2*-params.in["stst.rec.proj.logar1.stddev", "value"][[1]]) / 2,
+         params.in["stst.rec.proj.logar1.stddev", "value"][[1]],
+         1 ))),
+     tolerance = 1e7), "r$proj_logar1_stst_rec__nll: Consistent with proj_logar1_stst_rec__lvar, level not accounted for")
ok - r$proj_logar1_stst_rec__nll: Consistent with proj_logar1_stst_rec__lvar, level not accounted for
> 
> gadget3:::ut_tmb_r_compare2(model_fn, model_cpp, params.in)
# skip: not running TMB tests 0  <  1
NULL
> 
> ok_group("With noise") ########################################################
# With noise
> 
> attr(model_cpp, 'parameter_template') |>
+     g3_init_val("stst.rec.#", rnorm(5, 1e5, 500)) |>
+     g3_init_val("stst.rec.proj.logar1.stddev", exp(-8)) |>
+     g3_init_val("stst.rec.proj.logar1.phi", 0.8) |>
+     g3_init_val("stst.rec.proj.logar1.level", 1e5) |>
+     g3_init_val("stst_mat.spawn.blim", 1e2) |>  # blim too low to trigger
+ 
+     g3_init_val("*.K", 0.3, lower = 0.04, upper = 1.2) |>
+     g3_init_val("*.Linf", max(g3_stock_def(st_imm, "midlen")), spread = 0.2) |>
+     g3_init_val("*.t0", g3_stock_def(st_imm, "minage") - 0.8, spread = 2) |>
+     g3_init_val("*.walpha", 0.01, optimise = FALSE) |>
+     g3_init_val("*.wbeta", 3, optimise = FALSE) |>
+ 
+     g3_init_val("project_years", 100) |>
+     identity() -> params.in
> nll <- model_fn(params.in) ; r <- attributes(nll) ; nll <- as.vector(nll)
> 
> ok(ut_cmp_equal(
+     as.vector( g3_array_agg(r$detail_stst_imm__spawnednum, c("year"), step = 1, age = 0) ),
+     as.vector( exp(r$proj_logar1_stst_rec__lvar) ),
+     end = NULL), "r$detail_stst_imm__spawnednum: projection variable used for recruitment")
ok - r$detail_stst_imm__spawnednum: projection variable used for recruitment
> 
> ok(all(exp(r$proj_logar1_stst_rec__lvar) > 1e4), "r$proj_logar1_stst_rec__lvar: Noise not high enough for value to drop below 1e4")
ok - r$proj_logar1_stst_rec__lvar: Noise not high enough for value to drop below 1e4
> 
> ok(ut_cmp_equal(nll, sum(r$proj_logar1_stst_rec__nll)), "nll: Consistent with r$proj_logar1_stst_rec__nll (sums to same values)")
ok - nll: Consistent with r$proj_logar1_stst_rec__nll (sums to same values)
> ok(ut_cmp_equal(
+     as.vector(r$proj_logar1_stst_rec__nll),
+     as.vector(c(0, -dnorm(
+         tail(r$proj_logar1_stst_rec__lvar, -1) -
+         0.8 * head(r$proj_logar1_stst_rec__lvar, -1) -
+         0.2 * log(params.in["stst.rec.proj.logar1.level", "value"][[1]]),
+         0 - exp(2*-params.in["stst.rec.proj.logar1.stddev", "value"][[1]]) / 2,
+         params.in["stst.rec.proj.logar1.stddev", "value"][[1]],
+         1 ))),
+     tolerance = 1e7), "r$proj_logar1_stst_rec__nll: Consistent with proj_logar1_stst_rec__lvar")
ok - r$proj_logar1_stst_rec__nll: Consistent with proj_logar1_stst_rec__lvar
> 
> # plot(r$proj_logar1_stst_rec__lvar)
> 
> ok_group("lastx mode", local({ ################################################
+ 
+     st_imm <- g3_stock(c("stst", maturity = "imm"), c(10, 20, 30)) |> g3s_age(0, 5)
+     st_mat <- g3_stock(c("stst", maturity = "mat"), c(10, 20, 30)) |> g3s_age(3, 15)
+     stocks_st <- list(st_imm, st_mat)
+     fl <- g3_fleet(c("fl", "surv"))
+ 
+     actions <- list(
+         g3a_time(1990, 1994, c(6,6)),
+         gadget3:::g3a_initialconditions_manual(st_imm,
+             quote( 100 + stock__minlen ),
+             quote( 1e4 + 0 * stock__minlen ) ),
+         gadget3:::g3a_initialconditions_manual(st_mat,
+             quote( 100 + stock__minlen ),
+             quote( 1e4 + 0 * stock__minlen ) ),
+         g3a_age(st_imm),
+         g3a_age(st_mat),
+ 
+         g3a_spawn(
+             st_mat,
+             g3a_spawn_recruitment_hockeystick(
+                 r0 = g3_param_project(
+                     "rec",
+                     g3_param_project_logar1(lastx_f = g3_parameterized("lastx", value = 2)),
+                     random = FALSE,
+                     scale = "rec.scalar",
+                     by_stock = stocks_st,
+                     by_step = FALSE )),
+                 output_stocks = list(st_imm),
+                 run_step = 1 ),
+ 
+         # NB: Only required for testing
+         gadget3:::g3l_test_dummy_likelihood() )
+     full_actions <- c(actions, list(
+         g3a_report_detail(actions),
+         g3a_report_history(actions, 'proj_.*', out_prefix = NULL),
+         NULL))
+     model_fn <- g3_to_r(full_actions)
+     model_cpp <- g3_to_tmb(full_actions)
+ 
+     attr(model_cpp, 'parameter_template') |>
+         g3_init_val("stst.rec.#", rnorm(5, 1e5, 500)) |>
+         g3_init_val("stst.rec.proj.logar1.stddev", 1e-20) |>  # i.e. no noise
+         g3_init_val("stst.rec.proj.logar1.phi", 0) |>
+         g3_init_val("stst_mat.spawn.blim", 1e2) |>  # blim too low to trigger
+         g3_init_val("lastx", floor(runif(1, 1, 4))) |>
+         g3_init_val("proj_logar1_stst_weight.proj_logar1_rec_weight", 0) |>  # i.e. disable nll output (it'll be Inf)
+ 
+         g3_init_val("*.K", 0.3, lower = 0.04, upper = 1.2) |>
+         g3_init_val("*.Linf", max(g3_stock_def(st_imm, "midlen")), spread = 0.2) |>
+         g3_init_val("*.t0", g3_stock_def(st_imm, "minage") - 0.8, spread = 2) |>
+         g3_init_val("*.walpha", 0.01, optimise = FALSE) |>
+         g3_init_val("*.wbeta", 3, optimise = FALSE) |>
+ 
+         g3_init_val("project_years", 100) |>
+         identity() -> params.in
+     nll <- model_fn(params.in) ; r <- attributes(nll) ; nll <- as.vector(nll)
+ 
+     ok(!("stst.rec.proj.logar1.level" %in% names(params.in)), "level parameter disabled")
+ 
+     ok(ut_cmp_equal(nll, 0), "nll: 0, as we disabled it with weight")
+     ok(ut_cmp_equal(
+         as.vector( r$proj_logar1_stst_rec__lvar[length(r$proj_logar1_stst_rec__lvar)] ),
+         as.vector( mean(r$proj_logar1_stst_rec__lvar[as.character(seq(1994 - params.in["lastx", "value"][[1]] + 1, 1994))]) ),
+         tolerance = 1e-7 ), paste0("proj_logar1_stst_rec__lvar: Settled to mean of lastx (", params.in$lastx, ")"))
+ 
+     gadget3:::ut_tmb_r_compare2(model_fn, model_cpp, params.in)
+ }))
# lastx mode
ok - level parameter disabled
ok - nll: 0, as we disabled it with weight
ok - proj_logar1_stst_rec__lvar: Settled to mean of lastx ()
# skip: not running TMB tests 0  <  1
> 
> proc.time()
   user  system elapsed 
   2.70    0.12    2.81 
1..13
# Looks like you passed all 13 tests.