Build a self-contained data list for TAM

Prepares catch, survey, weight, and maturity data for the chosen model years and ages. Checks input coverage, combines older ages into a plus group when needed, and sets up mortality, catchability, and observation-error formulas. Optional projection settings append future years. fit_tam() calls this function automatically; call it directly to inspect a model specification.

Usage

prepare_tam(
  data,
  years = NULL,
  ages = NULL,
  N_settings = list(process = "iid", init = "exp", mu_form = NULL, sd_form = ~1),
  rec_settings = list(process = "rw", init = "free", mu_form = NULL, sd = NULL, phi = NULL),
  F_settings = list(process = "rw", mu_form = NULL, sd_form = ~1, age_breaks = NULL, cor_form = ~1,
        mean_ages = NULL),
  M_settings = list(process = "off", mu_form = NULL, mu_supplied = ~I(0.2), sd_form = ~1,
        age_breaks = NULL, first_dev_year = NULL, mean_ages = NULL),
  catch_settings = list(sd_form = ~1, sd_supplied = NULL),
  index_settings = list(sd_form = ~1, sd_supplied = NULL, q_form = ~q_block, q_link = "log"),
  proj_settings = NULL
)

# Optional projections:
# proj_settings = list(n_proj = 3, n_mean = 3, F_mult = 1, catch = NULL)

Arguments

data

A list of tidy observation data.frames: catch, index, weight, and maturity. See Details.

years

Consecutive historical model years (at least two). Inferred from observed data (non-projection) if NULL.

ages

Consecutive model ages (at least two). Inferred from observed data (non-projection) if NULL. If the ages in the data extend beyond max(ages), the "obs" column is summed for catch and index data within each year (and survey for indices). Older-age weight and maturity are retained for abundance-based aggregation inside nll_fun(), preserving biomass and mature biomass in the plus group.

N_settings

A list with elements:

• process: "off" for deterministic cohort survival, "iid" for independent cohort residuals, "rw" for residuals that accumulate through time, or "ar1" for residuals correlated between years and ages. See tinyAM-model.

• init: "exp" (default), "free", "latent", "mono" or "capped", independently of process. All use first-year recruitment log_r0 as the starting anchor. "exp" uses deterministic survivorship; "free" estimates fixed older-age log_n0 states; "latent" integrates them without an initial density. "mono" constrains abundance to decline from r0; "capped" bounds each older age by r0 without imposing age ordering. See Details.

• mu_form: ordinary fixed formula for proportional cohort gains/losses, default NULL. Covariates come from pop_covariates at the destination year and older age, after survivors are combined (including the plus group). Recruitment and initialization are unaffected. These effects do not identify immigration separately from mortality or observation changes. RW residuals require independent temporal mean contrasts; constant directions cancel and are rejected. See tinyAM-model.

• sd_form: ordinary formula for log SD of older-age survival residuals, default ~ 1. Use ~ factor(age) or a supplied age-sharing column in data$pop_covariates. Its design must be constant across years. Recruitment SD is controlled separately by rec_settings; initial abundance has no separate SD. Random-effect SD terms are not enabled.

rec_settings

A list controlling recruitment separately from cohort survival:

• process: "off" for recruitment exactly following its mean formula, or "rw" (default), "iid", or "ar1" for annual log-recruitment residuals. Boundary recruitment states are controlled by init.

• init: "free" (default) estimates first-year recruitment and additional pre-parent-SSB boundary states as fixed parameters; "latent" integrates them without an initial density or separate SD. Independent of process, including "off", and of older-age N_settings$init. Simulations condition on these states; joint uncertainty draws can vary them. See recruitment_formulas.

• mu_form: one-sided formula for the expected log recruitment. NULL is allowed for RW; off/IID/AR1 require an explicit formula, such as ~ 1 for a constant mean. Use bh() or ricker() for a median stock-recruit curve with no residuals or IID/AR1 residuals. Ordinary covariates come from data$pop_covariates at the youngest modeled age. See recruitment_formulas.

• sd: positive fixed recruitment SD, or NULL to estimate it. This is an innovation SD for RW/AR1 and a marginal SD for IID. Must be NULL when process is off.

• phi: fixed AR1 correlation in [0, 1), or NULL to estimate it. Must be NULL for other processes. Estimated correlations are positive.

F_settings

A list with elements:

• process: "off" makes historical F equal its explicit mean formula, without F latent states, residual SDs or correlations. "iid" gives independent departures from the mean log F, "rw" for departures that accumulate over time, or "ar1" for departures correlated between years and ages. "cor_rw" accumulates changes over time with AR1 correlation between annual increments at neighbouring states within each fleet. It estimates one signed age correlation per active correlation group; its SD is the marginal SD of annual log-F increments.

• sd_form: ordinary age-based formula for log process SD, default ~ 1, evaluated on data$catch. Sharing an SD does not share latent F states. IID SDs describe residuals; RW/correlated-RW SDs describe annual increments; stationary AR1 retains its innovation-scale convention.

• age_breaks: optional increasing integer interval boundaries for shared actual F. Intervals include their lower boundary and exclude their upper boundary; the final boundary must be max(ages) + 1. For ages 1:6, c(1, 2, 3, 4, 5, 7) shares F at ages 5 and 6. Ages below the first boundary have exactly zero F in all years, with no F state or process density. Positive catch observations at those ages are rejected. NULL retains separate F states at every modeled age. All mean components and process SDs must agree within a shared block; refine blocks for age-varying means. Process correlation uses ordered blocks, not replicated fish ages. With fleet labels, supply a common vector or a named list covering all fleets; NULL entries retain separate states by age within that fleet.

• cor_form: categorical correlation groups constant within each fleet. Default ~ 1 shares correlations; ~ fleet estimates separate values. For AR1, list(age = ~ fleet, year = ~ 1) allows different sharing by dimension. A single formula applies to both dimensions. Sharing parameters never correlates innovations between fleets. Unused groups are omitted.

• mu_form: an optional formula for mean-\(F\) (coefficients estimated as log-scale parameters log_mu_f), evaluated on data$catch. Required for "off", "iid" and "ar1"; use ~ 1 for one estimated mean level or ~ factor(age) for age-specific means. For "rw", NULL is allowed. Mean columns constant through time cancel from RW increments and are removed with a warning. Remaining columns must be identifiable from historical year-to-year differences. Ordinary F mean designs use historical support only; projection rows carry the terminal historical design by age. No future-year coefficients are fitted. With IID residuals, iid(), rw() and ar1() can add random mean effects, e.g. ~ factor(age) + rw(year). Temporal mean terms must have disjoint support, for example one per fleet; overlapping temporal mean terms and structured means with temporal residual processes are not enabled. mono() and logistic() can instead describe a rising mean age curve with no residual process or with IID/AR1 residuals, e.g. ~ logistic(age). Ordinary terms supply baseline F; realized F can depart from the rising mean. Static curves are not identifiable with RW/correlated-RW residuals. See catchability_curves. double_logistic() supplies WHAM's unnormalized four-parameter dome. Catch covariates must be available even where observations are missing. They must agree across biological ages pooled into a modeled plus group. The catch response obs cannot enter an F formula. With fleet labels, each fleet has its own F states; total F is their sum for survival and Z. for_fleet() gives one fleet its own formula, e.g. ~ 0 + for_fleet("Trawl", ~ factor(year) + logistic(age)) + for_fleet("Gillnet", ~ factor(year) + factor(age)). Ordinary terms outside scopes share their coefficients. Each inner formula includes an intercept unless removed explicitly; avoid redundant global and scoped baselines. smooth_formulas adds penalized s() and age-year t2() relationships. Smooth means support "off", "iid" and "ar1" F processes; the latter combinations are experimental and require mean/process separation checks.

• mean_ages: optional vector of ages to include in population weighted average F (F_bar) calculations. All ages used if absent.

M_settings

A list with elements:

• process: "off" uses only supplied/mean M. Otherwise "iid", "rw", and "ar1" have the same meanings as for F, on selected years/age blocks.

• sd_form: ordinary age-based formula for log process SD, default ~ 1. Evaluated on data$pop_covariates and required to be constant within each fitted M age block and across years. Use finer age_breaks for distinct state SDs.

• mu_form: optional formula for mean-\(M\) (applied on the log scale) built on data$pop_covariates, yielding coefficients mu_m. These enter the log-\(M\) surface directly and may therefore be positive or negative; they intentionally avoid a log_ prefix to prevent automatic back-transformation when tidied. If provided together with mu_supplied, the intercept in mu_form is dropped (warning) so supplied levels act as fixed offsets. iid(), rw() and ar1() support random mean effects with IID residuals or the residual process off. Use ~ 0 + rw(year) to vary supplied M through time. All mean components must be constant within fitted M age blocks. With AR1 residuals, one age-defined IID random slope is also supported, e.g. ~ 0 + capelin + iid(M_age_group, by = capelin). The grouping column must be constant across years at each age and the by column numeric. This experimental combination retains shared actual M within blocks; see formula_effects for the model and statistical cautions.

• mu_supplied: optional one-sided formula giving supplied (non-estimated) \(M\), e.g. ~ I(0.2) or a column reference such as ~ M_assumption stored in the data$pop_covariates data.frame.

• age_breaks: optional integer break points used by cut_ages() to define age_blocks sharing absolute latent M states across ages. Existing specifications retain their original interpretation. Alternatively, end at max(ages) + 1 to use explicit half-open intervals as for F: with ages 1:6, c(2, 5, 7) makes blocks 2–4 and 5–6. Ages outside M process blocks retain their existing mean M; they are not set to zero. When a narrower set of age_breaks than modeled ages is provided, process states are only estimated for ages within the specified range; ages outside this range are fixed to their mean or assumed levels. By default, blocks use cut_ages(ages[-1], unique(range(ages[-1]))), usually one block excluding the youngest age; two older ages form separate singleton blocks. This restricts how M can trade off against recruitment; it does not guarantee identifiability.

• first_dev_year: one historical modeled year at which M process states begin. Defaults to the second year if NULL. Earlier M stays at its mean or supplied value; earlier years are not coupled to the first latent state. Ignored when the M process is off.

• mean_ages: optional vector of ages to include in population weighted average M (M_bar) calculations. All ages used if absent.

catch_settings

A list with elements:

• sd_form: formula for log observation SD for catch-at-age data. Ordinary terms and Gaussian formula_effects are supported, e.g. ~ iid(age).

• sd_supplied: optional one-sided formula giving supplied SDs (on the natural scale of the log-observation residuals) for catch-at-age data. When provided, the intercept is removed from sd_form so supplied SDs act as offsets.

index_settings

A list with elements:

• sd_form: formula for log observation SD for index-at-age data. Ordinary terms and Gaussian formula_effects are supported, e.g. ~ survey + rw(age, by = survey).

• sd_supplied: optional one-sided formula giving supplied SDs (on the natural scale of the log-observation residuals) for index-at-age data. When provided, the intercept is removed from sd_form so supplied SDs act as offsets.

• q_form: formula for catchability, evaluated on the index table. Ordinary ~ q_block allows any pattern across blocks; ~ mono(q_block) is non-decreasing across ordered blocks. See mono() for independent non-decreasing curves by survey. Structured effects such as ~ survey + ar1(year, by = survey) allow shared annual changes across ages. ~ survey + logistic(age, by = survey) fits rising selectivity curves. ~ q_key + power(N, by = q_power_group) allows abundance-dependent catchability using internally calculated abundance at survey time. Missing power-group labels mean p = 1, not missing observations. Powers require the log q link. See power() and formula_effects. double_logistic() allows dome-shaped selectivity and multiplies baseline q after its inverse link. Penalized s() and t2(age, year) terms instead modify the selected q-link predictor. See smooth_formulas for supported bases, grouping and penalties. for_survey() applies a formula to one named survey, e.g. ~ 0 + for_survey("RV", ~ logistic(age)) + for_survey("Sentinel", ~ q_block). Effects are built on each survey's own reporting ages, then contribute zero to the predictor elsewhere.

• q_link: "log" (default) allows any positive catchability; "logit" restricts catchability to between zero and one. Use the logit link only when the survey index has an absolute abundance scale and this restriction is scientifically appropriate. Relative indices or surveys affected by herding may need the log link. The same formula and optional mono() terms work with either link. With "logit", coefficients and dq act on logit-q rather than log-q; covariate effects are no longer constant multipliers of q. This is a restriction, not a prior favouring a value near one. q can approach a boundary with large coefficients and SEs.

proj_settings

Optional list with elements:

• n_proj: number of years to project (default NULL disables projections).

• n_mean: number of recent historical years (default 3) used to average each weight type and maturity by age. Catch/index observations in projection years are missing. Other columns are copied from the last available year in each table.

• F_mult: multiplier of terminal historical F (default 1, status quo). Can be a value of length 1 or length = n_proj. When it is a vector of length one, that multiplier is recycled across all projection years.

• catch: total annual catch weight, in catch-number times mean-weight units (e.g. tonnes for thousand fish and kg). Supply one value recycled across years, or n_proj values. Specify either catch or F_mult, not both. catch = "observed" is restricted to fit_hindcast() evaluation; it is not a forecast input for an ordinary fit. See tinyAM-model for the catch solver and feasibility safeguards.

Value

A named list dat containing:

• years, ages — modeled ranges (years includes proj_years, if used)

• is_proj — logical vector identifying whether year is projected

• proj_years — integer vector of projection years, if used

• obs — per-type tables restricted to years x ages (including proj_years, if used)

• spawn_time, spawn_time_plus — timing matrices on modeled ages and, when retained, biological plus ages

• W, W_catch, P — stock weight, catch weight and proportion mature matrices (year x age). Optional weight$type distinguishes complete stock and catch weights; stock alone supplies both uses.

• plus_ages, W_plus_input, P_plus_input — retained biological ages and inputs, present only when biological data extend above max(ages). The last retained age represents a hidden plus group. Effective modeled W/P and hidden abundance N_plus can then be inspected in the likelihood report (fit$rep).

• obs_map — stack of obs$catch and obs$index mapping variables

• log_obs, is_missing, is_observed, observed - vector of log observations (with NA), logical vector indicating missing and observed values, and vector of non-missing values, respectively.

• design matrices: sd_catch_modmat, sd_index_modmat, q_modmat, and optionally F_modmat, M_modmat

• mean-level placeholders: log_mu_f and/or mu_m (or log_mu_supplied_m)

• process settings: N_settings, rec_settings, F_settings, M_settings, catch_settings, index_settings

• projection settings: proj_settings

• AR(1) parameter assumptions, logit_phi_*, if applicable

Details

Settings lists

Supply only the entries to change from the defaults shown in Usage. Formulas and vectors replace an entire entry. Optional NULL entries clear a formula, offset or age selection; sd_form = NULL restores ~ 1. Required entries cannot be NULL. update() retains unspecified entries from the original fit instead.

Inputs and preparation

See as_tam_data() for table keys, units, fleet-specific catch weights and covariate routing, and cod_obs for example data. Catch and survey zeros and missing values are excluded from the lognormal observation likelihood; complete predictions are still calculated.

Biological ages above the modeled plus age are retained in W_plus_input, W_catch_plus_input and P_plus_input. The population calculation derives effective plus-group weights and maturity from their hidden abundances; prepared terminal biological values are placeholders. Formula covariates use the supplied modeled-age rows, avoiding circular mortality calculations. The equations and spawning-time definition are in tinyAM-model.

Formula designs

Ordinary formula terms create design matrices. Structured formula_effects, mono() curves and smooth_formulas also store unique effect states and their grouping metadata. make_par() constructs the associated parameters. Supplied M and observation SDs act as offsets; their mean/SD formulas omit an intercept to avoid estimating another baseline. The settings below define available terms, restrictions and data sources.

Projections

proj_settings appends future years and marks them with is_proj in the prepared grid and tables. Stock/catch weights and maturity use recent-year averages by biological age; covariates and spawning timing retain terminal historical values. See tinyAM-model for forecasting and initial-state equations, and fit_hindcast() for held-out-data evaluation.

See also

fit_tam(), tinyAM-model, stats::model.matrix(), cut_ages()

Examples

dat <- prepare_tam(
  cod_obs,
  years = 1983:2024,
  ages = 2:14,
  N_settings = list(process = "iid", init = "exp"),
  F_settings = list(process = "rw", mu_form = NULL),
  M_settings = list(process = "off", mu_supplied = ~ I(0.3)),
  catch_settings = list(sd_form = ~ 1),
  index_settings = list(sd_form = ~ 1, q_form = ~ q_block)
)
names(dat)
#>  [1] "years"                 "ages"                  "N_settings"           
#>  [4] "rec_settings"          "F_settings"            "M_settings"           
#>  [7] "catch_settings"        "index_settings"        "proj_settings"        
#> [10] "obs"                   "is_proj"               "pop_data"             
#> [13] "has_fleet"             "fleets"                "F_data"               
#> [16] "F_columns"             "F_cor"                 "plus_ages"            
#> [19] "W_plus_input"          "W_catch_plus_input"    "P_plus_input"         
#> [22] "spawn_time_plus"       "spawn_time"            "P"                    
#> [25] "W_catch"               "W"                     "obs_map"              
#> [28] "log_obs"               "is_missing"            "is_observed"          
#> [31] "observed"              "sd_catch_modmat"       "log_sd_catch_supplied"
#> [34] "sd_index_modmat"       "log_sd_index_supplied" "q_modmat"             
#> [37] "N_modmat"              "M_terms"               "F_terms"              
#> [40] "log_mu_f"              "F_modmat"              "mu_m"                 
#> [43] "M_modmat"              "log_mu_supplied_m"     "logit_phi_n"          
#> [46] "process_sd"            "rec"                  

## With projection settings
dat <- prepare_tam(
  cod_obs,
  N_settings = list(process = "iid", init = "exp"),
  F_settings = list(process = "rw", mu_form = NULL),
  M_settings = list(process = "off", mu_supplied = ~ I(0.3)),
  catch_settings = list(sd_form = ~1),
  index_settings = list(sd_form = ~1, q_form = ~q_block),
  proj_settings = list(n_proj = 3, n_mean = 3, F_mult = 1)
)