Overview of Model Options

SPoRC models are assembled through a pipeline of Setup_Mod_* functions (for the estimation model) and Setup_Sim_* functions (for the operating model / simulation). Each function appends to an input_list (or sim_list) that accumulates data, parameter starting values, and RTMB factor maps. This vignette catalogs every user-facing option across that pipeline. Where both an integer code and a string alias exist, either is accepted. Strings are generally preferred for readability. When n_sexes > 1, the first sex index is always female. When n_pop > 1, populations follow the order in which they are defined. Throughout this document, the notation [p × r × y × τ × a × s] is shorthand for array dimensions: population, region, year, seasons, age, sex.

Parameter sharing conventions

Many hyperparameters throughout the pipeline (process-error standard deviations, catchability, natural mortality, tag reporting rates, …) are controlled by a _spec string argument that follows the same convention: a single named dimension (e.g. region, season, fleet, sex) is either estimated independently for every level, shared across some subset of levels, or fixed. Concretely, for a hyperparameter varying across dimensions abbreviated r (region), seas (season), and f (fleet), the accepted strings are:

String	Meaning
"est_all"	Independent parameter for every combination of dimensions
"est_shared_r"	Shared across regions; independent per remaining dimension
"est_shared_seas"	Shared across seasons
"est_shared_f"	Shared across fleets
"est_shared_r_seas", "est_shared_r_f", "est_shared_seas_f"	Shared across the named pair of dimensions
"est_shared_r_seas_f"	A single parameter shared across everything
"fix"	Kept at its starting value (not estimated)

Which dimension abbreviations are valid, and in what order they can be combined, depends on the specific parameter (see each section below for its abbreviations, e.g. sigmaC_spec additionally has a y (year) dimension). sigmaF_spec, sigmaC_spec, sigmaR_spec, and Fdev_rho_spec all follow this convention via a single shared internal helper (build_pe_map/build_shared_spec_map), so they behave identically. Other _spec arguments using the same string vocabulary (e.g. sigmaC_pop_spec, and selectivity’s fish_fixed_sel_pars_spec, fishsel_pe_pars_spec, fish_sel_devs_spec, and their retention/survey equivalents) are implemented separately with their own hand-written logic, and some additionally support a fleet-sharing escape hatch, "est_shared_f_<n>", which copies the sharing structure of fleet <n> wholesale rather than collapsing a dimension.

---

How arguments are named

The setup functions take around five hundred arguments between them, but far fewer ideas: most names are one of a handful of stems assembled from the same pieces. Reading the grammar is faster than reading the list, and it lets you derive a name rather than search for it.

A data source’s three arrays

Every observed data source is supplied as a triple, and the names differ only in the prefix.

Prefix	Has
Obs	The observations
Use	A 0/1 flag per cell, saying which observations are fit
ISS_	Input sample size, for composition data sources

So the fishery age compositions are ObsFishAgeComps, UseFishAgeComps, ISS_FishAgeComps, and the survey twins swap Fish for Srv. An _SE suffix has observation error where a data source reports its own (ObsSrvIdx_SE).

The Use array is the switch. A data source with data present and its Use flags at zero is passed through the model and contributes nothing, which is how a data source is turned off without removing it.

How a data source is configured

Four suffixes appear on nearly every data source.

Suffix	Chooses
_Type	How the observation is split across regions and sexes: "agg", "spltRaggS", "aggRspltS", "spltRspltS" for at-age data sources; "agg", "spltRspltS", "spltRjntS", "none" for compositions
_LikeType	The density: "lognormal"/"normal" for indices, or the composition families
_sigma_form	Where observation error comes from: "none", "data", "est_additive", "est_quadrature"
_spec	How a parameter is shared or fixed across dimensions, per the sharing conventions above

A _Type may change part way through a series, and says so through a year and fleet specification: "agg_Year_1-terminal_Fleet_1". The years are positions in the year vector rather than calendar years, and the last may be "terminal". A vector of these describes a whole model, one entry per block, and every year of every fleet needs an entry covering it. Later entries overwrite earlier ones where they overlap, so a general setting can be given first and a period carved out of it afterwards.

The two families differ only in what a bare setting means. Composition and CAAL _Type strings must hold the range, and the bare setting on its own is rejected with an error naming the full form. At-age _Type strings accept either, so a bare "spltRaggS" sets the whole series, one value per fleet sets each fleet’s whole series, and the grammar is needed only when the setting changes.

Two infixes that multiply the surface

AA marks an age-disaggregated data source, fit age by age rather than as a total: ObsCatch against ObsCatchAA, with the configuration suffixes following (CatchAA_Type, CatchAA_LikeType, CatchAA_sigma_form, sigmaCAA_spec). A fleet fits one or the other, not both.

_pop marks the population-specific twin of a data source, fit for one population rather than summed across them. Over a hundred arguments are _pop variants, and every one of them is inert in a single-population model. ObsSrvIdx has ObsSrvIdx_pop; Wt_SrvIdx has Wt_SrvIdx_pop; rho_srv_idx_spec has rho_srv_idx_pop_spec.

The two compose, so a single data source can reach ObsCatchAA_pop_SE.

Worked expansion

Everything the survey age compositions accept, derived from the grammar rather than looked up:

Name	Role
ObsSrvAgeComps, UseSrvAgeComps, ISS_SrvAgeComps	The data triple
SrvAgeComps_Type, SrvAgeComps_LikeType	Structure and density
Wt_SrvAgeComps	Weight
ObsSrvAgeComps_pop, UseSrvAgeComps_pop, ISS_SrvAgeComps_pop	Population-specific triple
SrvAgeComps_pop_Type, SrvAgeComps_pop_LikeType, Wt_SrvAgeComps_pop	Its configuration and weight

The remaining stems follow the same shape: _key supplies a mapping matrix, _blocks a time-blocking string, rho_*_spec and AgeObsCorr_* the correlation across ages, and Wt_* the likelihood weight.

Where a setting names its own legal values, it does so in the error it raises when given something else, so passing a deliberate nonsense string is a quick way to see the options.

---

Spatial and Demographic Structure

Everything starts with Setup_Mod_Dim() (estimation model) or Setup_Sim_Dim() (operating model; simulation). These functions define the skeleton that every subsequent setup call validates against.

Core dimensions

Argument	Type	Description
years	integer vector	Calendar years included in the assessment (length determines $n_y$)
ages	integer vector	Modeled age classes; the final element is the plus group
lens	numeric vector or NULL	Length-bin midpoints; NULL disables length-based features
n_regions	integer	Spatial regions
n_sexes	integer	1 (aggregated) or 2 (sex-structured)
n_fish_fleets	integer	Fishery fleets
n_srv_fleets	integer	Survey fleets

Multi-population structure

SPoRC decouples biological populations from spatial regions. Populations have their own stock-recruit dynamics and natal assignment; regions define the arena through which all populations move.

Argument	Type	Description
n_pop	integer	Number of biological populations (default 1)
natal_region	integer vector (length n_pop)	Maps each population to its home region. Inferred automatically when n_pop == n_regions (one-to-one) or n_regions == 1 (all share region 1). Must be supplied explicitly otherwise

When n_pop > n_regions, multiple populations share a natal region, a contingent-population or mixed-stock structure. When n_pop < n_regions, the model represents a single population spread across more regions than there are spawning origins.

Seasonal sub-stepping

Argument	Type	Description
n_seas	integer	Seasons per year (default 1)
seasdur	numeric vector (length n_seas)	Duration of each season as a fraction of the year (must sum to 1). Defaults to equal-length seasons

Within each year, processes execute in sequence: recruitment → movement → mortality (with age advancement at the end of the final season).

Other dimension controls (EM only)

Argument	Description
n_proj_yrs_devs	How many projection-year slots to pre-allocate for deviation parameters (ln_RecDevs, move_devs, selectivity deviations). Default 0
store_config	If TRUE, archives all Setup_Mod_* arguments inside input_list$config for reproducibility

Additional simulation dimensions (OM only)

Argument	Description
n_sims	Number of Monte Carlo replicates
n_yrs	Projection horizon (years)
n_obs_ages	Number of observed age bins (can differ from n_ages when compositions pool ages differently)
run_feedback	TRUE for closed-loop MSE; FALSE (default) for open-loop simulation
feedback_start_yr	First year of feedback when run_feedback = TRUE

---

Population Initialization

Controlled via Setup_Mod_Rec() (argument init_age_strc) and Setup_Sim_Rec().

SPoRC offers four methods for deriving equilibrium numbers-at-age, plus a fifth, non-equilibrium option. The equilibrium methods project a constant recruitment ($R_0$, or a separate ln_rinit scalar when use_rinit = 1) through seasonal mortality, fishing, and (optionally) movement until the age structure stabilizes. After the equilibrium is computed, multiplicative log-scale initial deviations (ln_InitDevs) are applied to ages $2, \ldots, A$.

When use_rinit = 1, the regional recruitment scalar feeding the equilibrium calculation is bias-corrected: $\zeta_{p,r}\,\text{rinit}_p\exp(-\sigma_{\text{Rec,early}}^2/2)$, treating ln_rinit as the median of the assumed lognormal recruitment process and converting to the corresponding mean before it is used as a deterministic (non-stochastic) equilibrium seed. This mirrors the same correction applied to Setup_Sim_Rec()’s equivalent rinit_input pathway during simulation, so fitted and simulated equilibria remain on a consistent scale.

Equilibrium solution method (init_age_strc)

Value	String	How the equilibrium is solved
0	"iterative"	Brute-force iteration: runs the full seasonal cycle $n_\text{ages} \times 10$ times
1	"scalar_no_move"	Closed-form geometric series assuming no movement at any age. Plus group: $N_{A^+} = N_{A-1} e^{-Z_{A-1}} / (1 - e^{-Z_{A^+}})$
2	"matrix"	Builds seasonal transition matrices $\mathbf{T}_a = \prod_\tau \mathbf{M}_{a,\tau} \mathbf{S}_{a,\tau}$ that combine movement and survival, then solves $\mathbf{N}_{A^+} = (\mathbf{I} - \mathbf{T}_{A^+})^{-1} \mathbf{T}_{A-1} \mathbf{N}_{A-1}$. Default
3	"scalar_plus_only"	Hybrid: uses the matrix approach for ages below the plus group but switches to the scalar geometric-series for the plus group itself
4	"free"	No equilibrium at all: ln_InitDevs are the initial log numbers-at-age (ages 2+, apportioned by sex ratio), with age 1 still taken from recruitment

Recommendations. Use the default "matrix" (2) for spatial models; "scalar_no_move" (1) is fine, and marginally cheaper, for single-region models. Use "free" (4) when the initial age structure has no information about $R_0$ and should not be pulled toward an equilibrium, matching assessments in which initial numbers-at-age are freely estimated parameters. Note two consequences of "free": the initial condition becomes independent of init_F_par and of ln_rinit, and the deviations are on the scale of log-numbers rather than log-ratios about an equilibrium, so any penalty applied via equil_init_age_strc acts as a prior on log initial abundance. Pair "free" with InitDevs_pen_center = "own_mean" (see Recruitment below) to penalize only the roughness of the initial age structure rather than its level.

For no prior at all, note that equil_init_age_strc bundles two decisions that do not move together: "equil" (0) neither penalizes nor estimates the deviations, and "stoch_all" (2) both estimates and penalizes them. Under "free" you usually want the combination it does not offer, estimated and unpenalized, because the deviations are the numbers themselves. Either route reaches it with one extra line, and they give an identical fit:

# route A: no penalty from the setting, then re-free the deviations it fixed
input_list <- Setup_Mod_Rec(input_list, ..., init_age_strc = "free",
                            equil_init_age_strc = "equil")
input_list$map$ln_InitDevs <- factor(seq_along(input_list$par$ln_InitDevs))

# route B: estimated by the setting, then switch its penalty off
input_list <- Setup_Mod_Rec(input_list, ..., init_age_strc = "free",
                            equil_init_age_strc = "stoch_all")
input_list$data$init_devs_pen_use[] <- 0

Either line must come after the setup calls that touch ln_InitDevs, since a later call would overwrite an earlier assignment. Leaving it out under "equil" is quiet rather than loud: the deviations stay fixed at zero, which under "free" is one fish at every initial age, and the model still fits.

Initial deviation structure (equil_init_age_strc)

Value	String	Which ages receive initial deviations
0	"equil"	None, strict equilibrium
1	"stoch_no_plus"	All ages except the plus group. Default
2	"stoch_all"	All ages including the plus group
3	"stoch_shared_ages"	User-defined age sharing via init_age_devs_shared

init_F_par (array [n_regions × n_seas × n_fish_fleets]) optionally introduces fishing mortality into the equilibrium calculation, producing a fished initial condition. init_F_form sets what it means: "prop" treats it as a proportion of the mean fishing mortality (inverse-logit scale, bounded to (0,1)), so the initial age structure moves with mean F; "abs" treats it as an absolute rate (log scale) that is independent of mean F. init_F_spec ("fix" or "est") sets whether it is estimated, independently of the form. Prefer "abs" when the historical fishing mortality that shaped the initial condition is distinct from the mean F of the modeled period; under "prop" one parameter both depletes the initial age structure and scales the F series, and because catch constrains only their product the optimizer can fit catch equally well with a smaller, harder-fished stock. The older init_F_prop argument is still accepted and is converted to the "prop" form.

---

Recruitment

Controlled via Setup_Mod_Rec() (EM) or Setup_Sim_Rec() (OM).

Stock-recruit function (rec_model)

String	Description
"mean_rec"	Estimate mean $\ln R_0$ with annual deviations; no SSB feedback
"bh_rec"	Beverton-Holt: $R = 4hR_0 \cdot \text{SSB} / \left[(1-h)S_0 + (5h-1)\text{SSB}\right]$. Unfished spawning biomass per recruit ($S_0$) is computed internally by projecting a single recruit through all ages and seasons with movement
"ricker_rec"	Ricker in depletion form: $R = R_0 (S/S_0)\exp(\alpha(1 - S/S_0))$ with $\alpha = \log(4h/(1-h))$ (the “Dorn form” used by the EBS pollock assessment). The curve passes through $(S_0, R_0)$ and shares the Beverton-Holt’s compensation ratio at a given $h$, not the textbook $R(0.2S_0) = hR_0$ definition, so steepness values are not interchangeable between "bh_rec" and "ricker_rec", and a steepness prior calibrated for one should not be reused for the other

SR_ref_yr sets the year index whose inputs feed the unfished spawning-biomass-per-recruit calculation. That is every input to it, not only weight-at-age and maturity: natural mortality, movement, stray rate, sex ratio, and the selectivity, discard mortality and continuous-movement rate that enter through init_F all come from that year, in the operating model as well as the estimation model, and hence $S_0$ and the scale of the stock-recruit curve. Default 1 (first model year, the long-standing behavior); set length(years) to condition the curve on terminal biologicals, a convention several assessments use. With time-varying weight-at-age the two choices give different $S_0$; pick whichever the assessment you are bridging or the reference-point convention you follow uses, and keep it consistent with how reference points are computed.

SR_ref_yr fixes everything in the per-recruit part; the $R_0$ that multiplies spawning biomass per recruit to give $S_0$ is always the year’s own value, and it is the same $R_0$ the curve’s numerator uses. The two are one argument in the code and cannot be split: steepness is recruitment at $0.2\,S_0$ as a fraction of $R_0$, so drawing $S_0$ from a different year’s $R_0$ would leave the unfished state off the curve. With R0_blocks of 5 then 12 and $\phi_0 = 2.5$, building $S_0$ from the first block while the numerator uses the second returns 12.8 at the second block’s unfished biomass instead of 12. So: reference-year biologicals for the per-recruit part, current-year $R_0$ for the scale.

The consequence of that, when R0_blocks is used with a stock-recruit form, is that a block rescales the entire curve, so $S_0$, depletion, and any reference point derived from the curve step at the block boundary. That is the intended reading of “the curve’s scale is time-varying”, but it is a strong statement and worth being deliberate about. Under rec_model = "mean_rec" it does not arise. And a curve fitted as a penalty on mean-recruitment residuals (sr_penalty with sr_R0_spec = "shared") takes its scale from R0_ref_block, one fixed value, so the diagnostic curve stays put while the blocked $R_0$ drives recruitment.

Density dependence scope (rec_dd)

String	When to use
"local"	SSB and $S_0$ computed per population and/or per region. Required when n_pop > 1 with BH recruitment
"global"	SSB summed across all regions before entering the SRR. Single-population only

Spawning timing

Argument	Description
spawn_seas	Season index in which spawning occurs
t_spawn	Fraction of the spawning season elapsed before spawning occurs (0 = start, 0.5 = midpoint)
rec_lag	Delay (in seasons) between spawning and recruitment entry. 1 (default): recruitment driven by SSB from rec_lag seasons prior, may enter in any season. 0: age-0 recruitment, recruitment driven by that same year’s own SSB. Since that SSB isn’t known until spawn_seas is reached, recruits may only enter in spawn_seas itself or a later season in the same year (rec_seas_prop must be zero for every season before spawn_seas), and the recruit age class must have zero maturity everywhere

Recruitment allocation

Argument	Description
rec_region_prop_spec	How regional recruitment proportions are estimated/fixed
rec_seas_prop_spec	How seasonal recruitment proportions are estimated. Default "fix" (all recruitment enters in season 1)
sexratio_spec	Sex ratio at recruitment estimation. Default "fix" (equal or user-supplied)
sexratio_blocks	Time blocks for sex-ratio parameters, specified as "none_Pop_<p>_Region_<r>" or "Block_<b>_Year_<s>-<e>_Pop_<p>_Region_<r>"

Recruitment and Initial age deviations

Argument	Description
RecDevs_spec	Estimation structure for annual $\ln(\text{RecDevs})$
InitDevs_spec	Estimation structure for initial age deviations
InitDevs_sex_spec	"est_shared_s" (default): one initial age curve read by every sex. "est_all": one curve per sex, each penalized, with an "own_mean" InitDevs_pen_center pooled across sexes so the sexes share a single estimated level. This is how assessments with a common mean-log-initial and sex-specific initial deviations are written. Requires n_sexes > 1
sigmaR_spec	Recruitment variability: "est_all", "est_shared_r", "fix_early_est_late", "est_shared_all", "fix"
sigmaR_switch	Year index separating the early and late $\sigma_R$ periods
do_rec_bias_ramp	Methot & Taylor bias adjustment (0 = off, 1 = on)
bias_year, max_bias_ramp_fct	Bias ramp breakpoints and maximum correction factor
dont_est_recdev_last	Number of terminal-year recruitment deviations to fix at zero
equil_init_age_strc	Which initial age deviations are estimated and penalized: "equil" (neither), "stoch_no_plus", "stoch_all", "stoch_shared_ages", or "stoch_all_no_pen" (every age estimated, none penalized). The last is for init_age_strc = "free", where the deviations are the initial log numbers at age rather than departures from an equilibrium, so a penalty on them is a prior on initial abundance
dont_pen_recdev_first	Number of leading recruitment deviations that stay estimated but take no penalty, because the first years belong to the initial condition rather than to the recruitment process. 0 (default) penalizes every year
RecDevs_model	Process error on the recruitment deviations: "iid" (default), "rw", or "ar1". See Recruitment process error (RecDevs_model)
RecDevs_rho_spec	Sharing structure for the AR1 correlation RecDevs_rho; only active under RecDevs_model = "ar1"
RecDevs_rw_init_sigma	Standard deviation on the first estimated deviation of a random walk, which is what sets the level of the series. Default 5; NA starts the walk at zero under the estimated $\sigma_R$
ln_global_R0_spec	"est" (default) or "fix". "fix" maps ln_global_R0 off so the recruitment deviations hold log recruitment outright rather than as departures from a level, which is how SAM writes it. The recruitment counterpart of ln_F_mean_spec
init_age_devs_shared	Integer vector of length n_ages - 1 specifying age-sharing for ln_InitDevs. Positions with the same value share a single estimated parameter (e.g. c(1:42, rep(42, 9))). Required when equil_init_age_strc = 3; NULL (default) uses standard behavior
use_rinit	0 = population initialized using ln_global_R0 (default); 1 = separate ln_rinit used for initialization, with ln_global_R0 governing only the recruitment relationship.
R0_blocks	Time blocks for R0, one entry per population: "none_Pop_<p>" (default) or "Block_<b>_Year_<a>-<e>_Pop_<p>" in 1-based year indices with "terminal" allowed. Under rec_model = "mean_rec" R0 IS mean recruitment, so a block is a productivity regime; under a stock-recruit form it is the curve’s scale, so blocking makes the curve time-varying
R0_ref_block	Which block supplies R0 wherever a single value is needed rather than a year’s value: the initial age structure, the regional apportionment, the R0 prior, the ln_rinit penalty, and the stock-recruit scale under sr_R0_spec = "shared". Default 1

A simulation self-test is worth reading before relying on R0 blocks. With a doubling of productivity partway through a 30-year series, hand-placed blocks recovered R0 as 4.49 and 12.12 against a truth of 5 and 12 and cut the median SSB error from -10.3% to -6.3%. R0 and the recruitment deviations are partly confounded under mean recruitment, and the deviations absorb much of a regime a block is trying to capture, so place blocks on years you have an external reason to suspect. | Use_rinit_pen, rinit_pen_sd | 0/1 and a standard deviation. With use_rinit = 1 the two recruitment levels are linked only through the data, and the initial age structure is often the thinnest part of it. 1 places a normal on $\log(\text{rinit}/R_0)$, reported as rinit_nLL, saying how far the stock’s equilibrium recruitment may have sat from the level the modeled years show. An equilibrium recruitment averages over several year classes rather than one, so its spread is narrower than $\sigma_R$. A reasonable choice is $\sigma_R/(1/M - 0.5)$, the recruitment variability shrunk by the number of year classes the equilibrium averages over |

A deviation shared across cells of the parameter array through the map (est_shared_pop_r, est_shared_r, est_shared_s, init_age_devs_shared) is one parameter and has one penalty: each cell the penalty loop visits takes $1/n$ of it, $n$ the number of such cells holding the level. Under the bias ramp the initial age deviations are centered on the ramp read at each age’s own birth year, deviation index $1 - a$ for the deviation on age $a$, which is how a ramp defined on calendar years treats the years before the first model year; bias_year can therefore be negative.

Recruitment process error (RecDevs_model)

ln_RecDevs can follow the same three process-error structures as ln_F_devs can, specified via RecDevs_model:

String	Description
"iid"	Independent annual deviations about the mean (default)
"rw"	Random walk: each deviation is centered on the previous year’s
"ar1"	First-order autoregressive, reverting toward zero at rate RecDevs_rho

The independent form is the standard mean recruitment with lognormal deviations, and the walk is a state-space version of recruitment. The difference between iid and random walk or ar1 matters most in the terminal years and in data-poor stretches: independent deviations pull recruitment back to $R_0$ wherever the data are sparse, while a walk or ar1 holds it near the last year the data did inform.

Under "rw" and "ar1" the penalty is a statement about the change in recruitment rather than its level, so neither can be combined with do_rec_bias_ramp = 1 or RecDevs_pen_center = "own_mean", both of which assert a mean about zero; setup rejects both combinations. The level of the series is set instead by year one, which takes a normal of standard deviation RecDevs_rw_init_sigma (default 5, wide enough that $R_0$ and the data decide it), or by nothing at all under dont_pen_recdev_first = 1, which is the flat prior SAM places on its first year. NA starts the walk at zero under its own $\sigma_R$ instead, which pulls the first year toward mean recruitment. sigmaR_switch still applies, so a walk can take one standard deviation early and another late, and a step reads the standard deviation of the year it lands on.

Under rec_model = "mean_rec", log recruitment is ln_global_R0 plus a deviation, and the two are only both estimable when something reads the deviations’ level. An "iid" or "ar1" penalty does; a random walk does not, since it penalizes only the change from one deviation to the next. A walk combined with dont_pen_recdev_first >= 1, which removes the one remaining term that read the level, leaves the two exactly unidentified: the fit converges, the Hessian is singular, and every standard error comes back NA. Setup rejects that combination and points at ln_global_R0_spec = "fix". A walk with the first year still penalized is accepted with a warning, since the level is then held only by that one term and $R_0$’s standard error comes back near RecDevs_rw_init_sigma.

RecDevs_rho_spec controls the AR1 correlation parameter (RecDevs_rho, dimensioned as population by region) using the same "est_all" / "est_shared_<dims>" / "fix" convention as everywhere else, and is only used under "ar1"; the other two recruitment options map it to NA regardless of what is supplied.

The first year’s recruitment is the first year’s age one abundance. Under an equilibrium initialization that abundance belongs to the initial condition, and giving it a recruitment penalty as well prices it twice. WHAM keeps it as a separate initial numbers at age parameter with no process error, and dont_pen_recdev_first = 1 is the same statement in SPoRC: the deviation is still estimated, so the data set the first year’s recruitment freely, but it takes no prior from ln_sigmaR.

Mapping the deviation off through RecDevs_spec is a different thing. That fixes it at its starting value, which is not the same as leaving it free with no prior.

Every year has a recruitment, so a year left out of the penalty is a year the walk passes through rather than a gap in the series. The step into year $y$ always reads year $y-1$, whether or not that year is penalized, and a year left out contributes no term of its own. Two consequences follow. dont_pen_recdev_first frees the level of the series without breaking it, since year two still anchors on year one. And dont_est_recdev_last removes the terminal terms without pulling the last penalized deviation toward zero.

Deviation penalty centering (RecDevs_pen_center, InitDevs_pen_center, Fdev_pen_center)

The recruitment, initial age, and (iid) fishing mortality deviation penalties can each be centered in one of two ways:

String	Penalty mean	What it constrains
"fixed"	The asserted prior mean: zero, or the bias-corrected $-\sigma_R^2/2$ for recruitment under the bias ramp. Default	Both the level and the spread of the deviations
"own_mean"	The mean of the estimated deviations themselves	Only their spread; the level is left free (a sum of squares about the mean, matching assessments whose deviation vectors sum to zero)

Recommendations. "fixed" is the statistically coherent choice when $\sigma_R$ is estimated or the deviations are treated as random effects: the penalty is then a genuine distributional assumption, and the bias ramp routines depends on the mean being asserted. "own_mean" is primarily a bridging device: it reproduces assessments where the mean parameter (ln_global_R0, ln_F_mean) holds the level and the deviations have only shape, so the level is not penalized twice. Two cautions: under "own_mean" the deviations’ level must be fixed elsewhere (an $R_0$ prior, a fixed deviation, or informative data) or the likelihood is flat along it, and RecDevs_pen_center = "own_mean" cannot be combined with do_rec_bias_ramp = 1 (the $-\sigma^2/2$ offset is meaningless once the mean is estimated rather than asserted; setup errors out).

Recruitment level penalty

Separate from the deviation penalty, an optional penalty on the log recruitment series itself:

Argument	Description
Use_rec_level_pen	0 (default) / 1 toggle
rec_level_pen_sigma	Standard deviation of the penalty. A sum of squares with weight $w$ corresponds to $\sigma = 1/\sqrt{2w}$. Default 1
rec_level_pen_center	"own_mean" (default; penalizes only the series’ variability) or "fixed" (centers on zero)
rec_level_pen_yrs	Calendar years the penalty applies over; NULL (default) = all years

Rationale. Under a stock-recruit relationship the deviations are residuals about the predicted curve, so a model that also wants the realized recruitment series to stay regular has nowhere else to say so. This penalty is that second, independent statement, and reproduces the recruitment regularity penalties several existing assessments have (e.g., a penalty on log recruitment variability). Leave it off unless recruitment in data-poor years is wandering unreasonably; it is a tuning penalty, not a probability model, and it will shrink genuine recruitment variability if over-weighted (hence the rec_level_pen_sigma is left as data rather than an estimated parameter).

Initial age deviations, tie between sexes

When the sexes have their own initial age curves (InitDevs_sex_spec = "est_all"), a second statement can be made about them: not how variable each curve is (that is the initial-age penalty) but how far apart the two curves may sit,

$$-\ln L^{\text{InitSex}} = \sum_{s > 1}\sum_{i} -\ln \phi\left(\epsilon_{p,r,i,s}^{\text{Init}} - \epsilon_{p,r,i,1}^{\text{Init}};\ 0,\ \sigma_{\text{InitSex}}\right)$$

Argument	Description
Use_init_sex_pen	0 (default) / 1 toggle. Requires n_sexes > 1 and InitDevs_sex_spec = "est_all"; under "est_shared_s" the difference is identically zero and setup refuses it
init_sex_pen_sigma	$\sigma_{\text{InitSex}}$. A sum of squares with weight $w$ corresponds to $\sigma = 1/\sqrt{2w}$. Default 1

It runs over the same ages the initial-age penalty covers and enters the objective unweighted (reported as Init_Sex_nLL). Its natural use is the one flatfish assessments put it to: male and female initial numbers are estimated separately, the first few years of compositions are all that informs them, and the tie keeps the two together where those data are thin. Like the recruitment level penalty it is a tuning prior rather than a probability model, so its sigma is kept as data.

Stock-recruit curve as a penalty

Under rec_model = "mean_rec" a stock-recruit curve can be evaluated alongside the dynamics without ever advancing them, its residual penalized instead. Recruitment stays a mean with deviations, $R_y = R_0^{\text{mean}}\exp(\varepsilon_y)$, and the curve enters the objective only through the log residual $\xi^{\text{SR}}_y = \ln R_y - \ln \hat{R}_y$:

$$-\ln L^{\text{SR}} = \sum_{y \in \mathcal{Y}_{\text{SR}}} -\ln \phi\left(\xi^{\text{SR}}_y;\ 0,\ \sigma_{\text{SR}}\right)$$

Argument	Description
sr_penalty	"none" (default), "bh", or "ricker", naming the curve $\hat{R}_y$ follows. Valid only with rec_model = "mean_rec"; under "bh_rec" or "ricker_rec" the curve already generates recruitment, which makes the deviation the residual by construction, so penalizing the residual as well would count it twice and setup rejects the combination
sr_pen_sigma	$\sigma_{\text{SR}}$, the standard deviation of the residual. Kept as data rather than estimated. Default 1
sr_pen_yrs	Calendar years $\mathcal{Y}_{\text{SR}}$ the penalty is summed over; NULL (default) is every year that has a lagged spawning biomass, i.e. all but the first rec_lag. Years outside it keep their recruitment deviation estimated and contribute nothing, which is how a restricted stock-recruit window is expressed. Naming a year with no lagged spawning biomass is an error rather than a silent fallback: the curve’s prediction there is the equilibrium, so the residual would not be a stock-recruit residual
sr_R0_spec	Where the curve’s scale $R_0^{\text{SR}}$ comes from (see below)

String	The curve’s scale $R_0^{\text{SR}}$
"shared"	ln_global_R0, i.e. $R_0^{\text{mean}}$ itself, which under mean recruitment is the mean recruitment level. One parameter and well posed, but the curve is anchored on mean recruitment rather than on an unfished level. Default
"est"	The curve’s own estimated ln_sr_R0, identified by the curve fit alone. Reproduces templates with separate mean-recruitment and unfished-recruitment parameters, though nothing ties the two together, so the pair can slide
"rinit"	ln_rinit, the initial-equilibrium recruitment, so one parameter sets both the unfished age structure and the curve. Requires use_rinit = 1 and errors otherwise. This is the usual ADMB arrangement

Steepness is reachable in this mode: h_spec, Use_h_prior, and h_prior act on the penalty curve exactly as they do on a curve that generates recruitment.

What this is, and what it is not. It is a penalty on a derived quantity, a statement about how far realized recruitment sits from a curve, and not a data-generating process. The recruitment deviations remain free and keep their own penalty, so a weakly determined relationship informs the recruitment series rather than dictating it. The simulation side does not see it at all: Setup_Sim_Rec() and the operating model generate recruitment from the mean and its deviations under "mean_rec", and sr_penalty enters the estimation model’s objective alone. The curve’s prediction is reported as SR_pred and its penalty as SR_pen_nLL; see vignette("d_model_report"). MSY is the maximum of equilibrium yield over a curve that drives the dynamics, and this curve does not, so Get_Reference_Points() refuses an MSY request on a fit with the penalty and directs it to the SPR variants.

Recruitment penalty weighting (Setup_Mod_Weighting())

Argument	Description
Wt_Rec	Weight on the recruitment deviation penalty. Scalar (default 1), or an array [n_pop × n_regions × n_est_rec_devs] for per-deviation weighting; note the third dimension follows ln_RecDevs (moved by dont_est_recdev_last and n_proj_yrs_devs), not the number of years
Wt_Init_Rec	Weight on the initial age deviation penalty, [n_pop × n_regions × (n_ages - 1)] or scalar. NULL (default) inherits a scalar Wt_Rec; must be supplied explicitly when Wt_Rec is an array, since the two penalties are dimensioned differently

A per-deviation weight of zero removes that deviation from the penalty while it remains estimated, which is how a stock-recruit relationship is fit over a chosen window of years while recruitment stays effectively free elsewhere. That is distinct from dont_est_recdev_last, which removes the deviations themselves (recruitment reverts to the deterministic prediction), and from mapping a deviation off (fixed and unpenalized; see Which deviations are penalized). Beware the $\sigma_R$ interaction: deviations excluded from the penalty no longer inform an estimated ln_sigmaR, so weight-based windows are safest with $\sigma_R$ fixed.

Composition and index zero-guards (Setup_Mod_Weighting())

Argument	Description
addtocomp	Small constant guarding $\log(0)$ in composition likelihoods (default 1e-3)
comp_const_obs	0/1. Whether addtocomp is also added to the observed proportions used as multinomial weights. 1 (default) is the long-standing SPoRC behavior; 0 adds it only inside the logarithms, a convention several existing assessments use. The difference slightly reweights every bin, so set 0 when bridging such assessments and otherwise leave the default. The Dirichlet-multinomial also honors it, and there the constant is not neutral: every bin with no observed and no expected mass contributes $\log(\theta/(1+\theta))$, so compositions with many structurally empty bins (conditional age-at-length above all) bias $\theta$ upward under 1; use 0 for those. Setup_Mod_Weighting warns if any conditional age-at-length fleet uses the Dirichlet-Multinomial while this is 1
addtofishidx, addtosrvidx	Small constants guarding $\log(0)$ in the fishery/survey index likelihoods (default 1e-4 each)
addtotag	Small constant guarding zero tag releases/recoveries (default 1e-10)

addtocomp/comp_const_obs/addtofishidx/addtosrvidx/addtotag used to be arguments of Setup_Mod_Biologicals(). They still work there, forwarded with a message, but new code should pass them here instead.

Steepness priors

When a stock-recruit curve is present, either because it generates recruitment (rec_model = "bh_rec" or "ricker_rec") or because it enters as a penalty (rec_model = "mean_rec" with sr_penalty), steepness ($h$) can be penalized with a Beta distribution scaled to [0.2, 1] by default:

Argument	Description
Use_h_prior	0 = no prior, 1 = apply prior
h_prior	Data frame with columns pop, region, mu, sd, and optional lb, ub giving the beta’s support (default 0.2 and 1)
h_spec	Estimation structure for steepness parameters

The support matters, not just the mean and SD: a beta on $(0,1)$ is a different function of $h$ than one on $(0.2,1)$ (it has $\log(h)$ where the rescaled form has $\log(h - 0.2)$), and no choice of shape parameters reconciles the two. When bridging an assessment, match its prior’s support via lb/ub rather than approximating with the default.

R0 priors

A lognormal prior on $R_0$ can be applied per population:

Argument	Description
use_r0_prior	0 = no prior (default), 1 = apply lognormal prior on $\ln R_0$
r0_prior	Data frame with columns pop (population index), mu (prior mean on natural scale), and sd (prior SD on log scale). Required when use_r0_prior = 1

Population straying

When n_pop > 1, a fraction of recruits produced by population $p$ can “stray” and recruit into regions associated with other populations.

Argument	Description
stray_rate_spec	"fix" (default), "est_all", etc.
stray_rate_blocks	Time blocks: "none_Pop_<p>" or "Block_<b>_Year_<s>-<e>_Pop_<p>"
use_stray_rate_prior	0/1 toggle for Beta priors on stray rates
stray_rate_prior	Data frame with columns pop, block, mu (in (0,1)), sd

Spawning movement (single-season, multi-population)

When n_pop > 1 and n_seas == 1, individuals cannot physically move to their natal region within the seasonal cycle, so a separate spawning-movement matrix (sgl_seas_spawning_movement) routes SSB back to natal regions for the SRR calculation. If not supplied, SPoRC defaults to 100% natal homing.

---

Selectivity

Selectivity configuration is shared across fishery fleets (Setup_Mod_Fishsel_and_Q()), survey fleets (Setup_Mod_Srvsel_and_Q()), and retention curves (ret_sel_model within Setup_Mod_Fishsel_and_Q()). All three accept the same functional forms and time-varying structures.

Functional forms

Specified as character strings following the pattern "<form>_Fleet_<f>" (constant across all time blocks) or "<form>_Fleet_<f>_Block_<b>" (block-specific). The available forms are:

String	Functional form	Free parameters
"logist1"	Ascending logistic: $1 / (1 + \exp(-k(\text{bin} - a_{50})))$	2 ($\ln a_{50}$, $\ln k$)
"logist2"	Ascending logistic using $a_{50}$ and $a_{95}$: $1 / (1 + 19^{(a_{50} - \text{bin})/a_{95}})$	2 ($\ln a_{50}$, $\ln a_{95}$)
"gamma"	Dome-shaped gamma: $(b/a_\text{max})^{a_\text{max}/p} \exp((a_\text{max} - b)/p)$	2 ($\ln a_\text{max}$, $\ln \delta$)
"exponential"	Descending power: $1/\text{bin}^\beta$	1 ($\ln \beta$)
"dbnrml"	Double-normal with ascending and descending widths, plateau, and endpoint control	6
"nonpar"	Non-parametric: one logit-scale parameter per bin, transformed via $\text{logit}^{-1}$	$n_\text{bins}$
"nonparlog"	Non-parametric on the log scale, standardized so each year’s selectivity averages to 1 over a window of bins (every bin by default; see Normalization and length-based selectivity)	$n_\text{bins}$
"nonparfree"	Non-parametric on the log scale with no standardization at all: $\exp(\theta_\text{bin})$, so the values hold the height of the curve as well as its shape	$n_\text{bins}$
"asymplogist1"	Logistic with asymptote $\alpha \in (0,1)$: $\alpha / (1 + \exp(-k(\text{bin} - a_{50})))$	3 ($\text{logit}(\alpha)$, $\ln a_{50}$, $\ln k$)
"asymplogist2"	Logistic with asymptote, $a_{50}/a_{95}$ parameterization	3 ($\text{logit}(\alpha)$, $\ln a_{50}$, $\ln a_{95}$)
"bicubic"	Bicubic natural-cubic-spline surface over a bin-node $\times$ year-node grid	$n_\text{bin\_nodes} \times n_\text{yr\_nodes}$ (see below)

The three non-parametric forms differ only in what happens after the transform, and that difference is what decides whether selectivity can have a scale. Write $\theta_{a,b}$ for the estimated parameter at bin $a$ in block $b$, $b(y)$ for the block year $y$ falls in, and $s_{y,a}$ for the selectivity the model actually uses. Then

$$\texttt{"nonpar"}: \quad s_{y,a} = \frac{\operatorname{logit}^{-1}\!\left(\theta_{a,b(y)}\right)} {\dfrac{1}{n_\text{yrs} n_\text{bins}}\sum_{y'}\sum_{a'} \operatorname{logit}^{-1}\!\left(\theta_{a',b(y')}\right)}$$

$$\texttt{"nonparlog"}: \quad s_{y,a} = \frac{\exp\!\left(\theta_{a,b(y)}\right)} {\dfrac{1}{|B|}\sum_{a' \in B} \exp\!\left(\theta_{a',b(y)}\right)}$$

$$\texttt{"nonparfree"}: \quad s_{y,a} = \exp\!\left(\theta_{a,b(y)}\right)$$

where $B$ is the set of bins named by *_sel_norm_bins (every bin by default). The inverse logit bounds each raw value below one before anything else happens; the exponential does not.

Read the denominators. "nonpar" divides by a single number computed over years and bins together, so it fixes the grand mean of the whole surface at one: one degree of freedom is removed from the surface as a whole. "nonparlog" divides by a number computed within each year, so it fixes every year’s mean at one: one degree of freedom is removed per year, and only within-year contrasts survive. "nonparfree" divides by nothing, so every $\theta_{a,b}$ is a selectivity in its own right and the level is estimated along with the shape.

A standardized form cannot hold a level, so the level then sits in catchability or in the fishing mortality mean. A free form can hold it, which is why an index fit age by age uses one: the age multiplier is the whole quantity being estimated. It also makes it possible to estimate a level twice by accident: if selectivity has a level and the mean it multiplies is free, the two trade a common factor, the Hessian stops being positive definite, and the selectivity standard errors inflate while the fitted trajectory looks fine. Pin one bin, or one block’s reference group, by leaving it out of *_sel_nonpar_est_bins.

Choosing among the three non-parametric forms. "nonpar" bounds every raw value below one (logistic transform) and mean-standardizes jointly over years and bins; "nonparlog" leaves the scale free and centers within each year. Under "nonparlog" only the within-year differences among parameters are identified, since the level is absorbed by catchability or fishing mortality, so pair it with the selectivity parameter centering penalty (below) or fix a bin group via *_sel_nonpar_est_bins to pin the scale. Prefer "nonpar" when you want selectivity interpretable as a proportion without further constraints.

"nonparfree" does not standardize at all, so its values hold the level as well as the shape. An index fit age by age needs exactly that: a free catchability per age and a selectivity estimated at age are one quantity written two ways, so an at-age index puts the whole age multiplier in selectivity and leaves catchability alone. Group ages the way a catchability key would through *_sel_nonpar_est_bins, and leave an age out of the list entirely when the fleet never observes it, which maps it off rather than leaving it free and uninformed. Because nothing is canceled by standardization, "nonparfree" also accepts a "scale" sex offset, which "nonpar" and "nonparlog" refuse. Do not pair it with the centering penalty.

Bicubic spline selectivity ("bicubic")

Unlike the other functional forms above, "bicubic" selectivity is specified with its own extended syntax:

"bicubic_Bin_<n_bin_nodes>_Yr_<n_yr_nodes>_Fleet_<f>[_Block_<b>][_SelStyr_<year>][_NSelBins_<n>]"

A smooth 2-dimensional selectivity-at-bin-and-year surface is built from a small grid of $n_\text{bin\_nodes} \times n_\text{yr\_nodes}$ freely estimated log-scale node values (fish_fixed_sel_pars / srv_fixed_sel_pars, flattened column-major into a [yr_node × bin_node] matrix). Two natural-cubic-spline weight matrices are precomputed once at setup:

• $\mathbf{W}^{\text{bin}}$ ($n_\text{bins} \times n_\text{bin\_nodes}$), mapping bin-node values onto every bin,
• $\mathbf{W}^{\text{yr}}$ ($n_\text{yrs} \times n_\text{yr\_nodes}$), mapping year-node values onto every year,

and combined via a two-pass tensor product (bin-node values spline-interpolated across bins for every year-node, then those curves spline-interpolated across years) to give the full log-selectivity surface, which is then exponentiated. Setting n_yr_nodes = 1 collapses the surface to a time-invariant bin-only spline; combining n_yr_nodes = 1 with fish_sel_blocks/srv_sel_blocks re-fits an independent bin-only spline within each block.

Two optional suffixes restrict the fitted region of the surface, edge-holding (flat-lining) outside it:

• _SelStyr_<year>: a calendar year within the block. Only years from SelStyr through the block’s end are actually spline-fit across the year dimension; years within the block before SelStyr are kept constant at the SelStyr year’s fitted curve (“previous years are filled”).
• _NSelBins_<n>: restricts the spline fit to the first n bins (ages or lengths, per fish_selex_type/srv_selex_type). Bins beyond n are kept constant at the last fitted bin’s value (a plateau), rather than continuing the spline extrapolation.

Both suffixes only change which years/bins are considered part of the fitted surface; they do not change the total number of estimated node parameters.

Temporal variation

SPoRC provides two mutually exclusive mechanisms for time-varying selectivity within a fleet. You can use discrete blocks or continuous deviations, but not both on the same fleet.

Discrete blocks (fish_sel_blocks / srv_sel_blocks / ret_sel_blocks): defined as "Block_<b>_Year_<s>-<e>_Fleet_<f>" or "Block_<b>_Year_<s>-terminal_Fleet_<f>". Blocks must be non-overlapping and collectively span all model years. Each block gets its own set of fixed-effect selectivity parameters (and can even use a different functional form).

Continuous deviations (cont_tv_fish_sel / cont_tv_srv_sel / cont_tv_ret_sel): specified as "<type>_Fleet_<f>". For parametric forms ("logist1", "gamma", etc.), IID and random-walk deviations act multiplicatively on the transformed base parameters. For semi-parametric forms ("3dmarg", "3dcond", "2dar1"), deviations act multiplicatively on the selectivity curve at the bin level.

String	Description
"none"	Time-invariant (default)
"iid"	Independent annual deviations on selectivity parameters
"rw"	Random walk on selectivity parameters
"3dmarg"	3D Gaussian Markov random field, marginal variance parameterization
"3dcond"	3D GMRF, conditional variance parameterization
"2dar1"	Separable 2D AR(1) over bin × year

Ancillary controls for continuous time-variation include fishsel_pe_pars_spec (hyperparameter estimation), fish_sel_devs_spec (deviation estimation structure), fishsel_devs_shared_bins (bin grouping for shared deviations), and corr_opt_semipar (which correlation components to suppress in semi-parametric forms). Bin-sharing specs ("est_shared_b" variants) require the deviations to be indexed by bin: any GMRF/AR1 form, or "iid"/"rw" on a non-parametric fleet ("nonpar"/"nonparlog", where each deviation already belongs to one bin); "iid"/"rw" on a parametric form indexes deviations by parameter and is rejected. The hyperparameters follow whatever the deviations are shared over, since only one series per shared group is fit. Under "iid"/"rw" the process error standard deviations are indexed by bin too, so a bin group leaves one standard deviation to estimate at the group’s lowest bin; sharing deviations across sexes leaves one set for the first sex under any form. The rest are fixed, so the number of hyperparameters reported as estimated is the number the fit can actually move.

fishsel_pe_pars_spec takes the same four _b spellings, where they mean something related but distinct: one process error standard deviation across every bin the fleet reads, rather than one deviation series per bin group. That is the key matrix ICES assessments write for the fishing mortality process, with a single group. It needs a non-parametric form, since on a parametric one the slots are curve parameters on unrelated scales, and it needs "iid" or "rw", since the GMRF and AR1 forms hold correlations and a single standard deviation rather than one per bin. Grouped standard deviations, where some ages share one and others share another, are not available.

Three further controls tune the process-error penalty itself:

Argument	Description
fishsel_pe_wt / retsel_pe_wt / srvsel_pe_wt	Per-fleet multiplier on the selectivity process-error likelihood (default 1). 0 removes the distributional penalty entirely while the deviations remain estimated; use it when deviations should float subject only to explicit smoothness penalties, as several existing assessments do. Anything other than 0 or 1 makes an estimated PE sigma reinterpretable, so prefer 0/1 unless deliberately down-weighting
fishsel_rw_init_sigma / retsel_rw_init_sigma / srvsel_rw_init_sigma	Standard deviation on the first year of an "rw" deviation series. Default 5 (first year effectively free). NA starts the walk at zero under the walk’s own estimated sigma, making the first year as smooth as every later step; appropriate when the base parametric curve already describes the first year well
*_sel_bin_dev_bins, cont_tv_*sel_bin_devs	Bin-override deviations; see below

Bin-override selectivity deviations

Individual bins can be cut loose from the functional form entirely: bins named in fish_sel_bin_dev_bins / ret_sel_bin_dev_bins / srv_sel_bin_dev_bins (a list with one element per fleet, NULL for fleets with no overrides) take a freely estimated annual value $\exp(\epsilon_{y,b})$ in place of whatever the form produced, applied after every other transformation including standardization. The rest of the curve keeps its parametric shape.

Argument	Description
fish_sel_bin_dev_bins etc.	Which bins each fleet overrides (e.g., list(1, NULL) frees bin 1 of fleet 1 only)
cont_tv_fishsel_bin_devs etc.	Process error on the override deviations, per fleet: "none", "iid", or "rw" (with its own estimated sigma per bin and *sel_bin_devs_rw_init_sigma for the first year)

When to use. The canonical case is a gear whose curve is well described by a parametric form except for one bin governed by availability rather than the gear (e.g., age-1 availability to a trawl survey varying with year-class strength). Overriding that bin gives it free annual variation without abandoning the parametric form (or paying for a full semi-parametric surface) elsewhere. Give the override "rw" process error unless the bin jumps independently between years; with "none", each year’s value is informed only by that year’s compositions and can be poorly determined in sparse years.

Selectivity parameter centering penalty

Argument	Description
Use_fish_selex_penalty / Use_ret_selex_penalty / Use_srv_selex_penalty	0 (default) / 1 toggle
fish_selex_penalty etc.	Data frame with columns region, fleet, block, sex, par (a single index or list-column of integer vectors naming a set), wt. Each row penalizes $w\,[\log(\overline{\exp(\theta)})]^2$ over the named set

This pushes the average selectivity of a parameter set toward one, pinning the scalar of a "nonparlog" curve that catchability or fishing mortality would otherwise absorb. It is a softer alternative to fixing a bin outright, and the standard companion to "nonparlog" selectivity (weights of 10 to 100 are typical starting points; the penalty only needs to break a ridge, not dominate the fit). Because the expression averages on the natural scale it is only meaningful for log-scale parameter sets, so do not apply it to "nonpar" (logit-scale) or to logit-scale asymptote parameters.

Parameter sharing and fixing

The fish_fixed_sel_pars_spec / srv_fixed_sel_pars_spec argument controls how base selectivity parameters are estimated:

String	Meaning
"est_all"	Fully region-, sex-, and fleet-specific
"est_shared_r"	Shared across regions
"est_shared_s"	Shared across sexes
"est_shared_r_s"	Shared across both regions and sexes
"fix"	Fixed at starting values

Selectivity plateau (_NSelBins_<n>)

Any functional form’s model string accepts an optional _NSelBins_<n> suffix (e.g. "logist1_Fleet_1_NSelBins_17", composable with _Block_<b> in either order): bins beyond n are kept at bin n’s computed value rather than evaluated through the form. This is the plateau convention many existing assessments have, where selectivity is the functional form up to a chosen bin and constant after it. The bicubic form already had this suffix and builds its plateau into its spline weights; for every other form the hold is applied after the form and its parameter deviations, before bin-level semi-parametric deviations and bin overrides. The plateau changes the model whenever the curve has not saturated by bin n, so treat it as part of the selectivity specification rather than a display convention.

Sex offsets (fish_sel_sex_offset / ret_sel_sex_offset / srv_sel_sex_offset)

Between fully independent sexes ("est_all") and identical ones ("est_shared_s") sits a family of models where the sexes are linked: one sex’s curve is the other’s, shifted or scaled by a small number of estimated offsets. These per-fleet options express that directly for total fishery, retention, and survey selectivity alike (requires n_sexes > 1):

String	Meaning
"none"	Each sex’s stored parameters are its own. Default
"par"	The stored parameter slots of every sex beyond the first hold additive offsets on the first sex’s transformed-scale parameters, $\theta_s = \theta_1 + \delta_s$. For log-scale parameters the second sex’s natural value is the first sex’s times $e^{\delta}$, e.g. $k_M = k_F e^{\delta_k}$, $a^{50}_M = a^{50}_F e^{\delta_{50}}$
"scale"	Each sex keeps its own parameters, and the realized curve of every sex beyond the first is multiplied by a constant $e^{\gamma_s}$ (ln_fishsel_sex_scale / ln_retsel_sex_scale / ln_srvsel_sex_scale, estimated per region, block, and sex). The scaled curve may exceed one
"par_scale"	A par offset and a scale offset
"apical"	Double normal only. Each sex keeps its own parameters, and for every sex beyond the first the ascending and descending limbs are built up to $A_s = e^{\gamma_s}$ rather than to one, with the plateau sitting there too. The selectivity at the first and last bins stays where that sex’s own parameters put it
"par_apical"	A par offset and an apical offset

Under "par", an offset fixed at zero reproduces sex-shared parameters, so the option nests "est_shared_s" while letting the link be estimated; no new parameters are involved, the existing sex slots simply change meaning. Under "scale", the offset is a log-scale constant on the whole curve, which is how several existing assessments have a male selectivity offset; because a constant multiplier is canceled by mean standardization, setup refuses "scale" for the non-parametric forms ("nonpar", "nonparlog") and the semi-parametric time-varying structures.

Choosing between "scale" and "apical". A scale offset multiplies the finished curve, so selectivity at the first and last bins moves with it. An apical offset is the height the two limbs are built up to, so those bins stay where $p_5$ and $p_6$ put them and only the middle of the curve moves. Either can draw any single curve, so a plot will not separate them, but their derivatives differ and a fit responds to that. Use "apical" where the offset is meant as the sex’s maximum selectivity. It requires the double normal; setup refuses it elsewhere and points to "scale".

The anchor bin. By default the ascending limb is anchored at the first bin, so its selectivity-at-first-bin parameter means what it says. When the compositions begin above the population’s first length bin, that parameter is describing a bin the data never see. fish_sel_dbnrml_startbin and srv_sel_dbnrml_startbin (integer vectors, one bin index per fleet) anchor the limb at the first data bin instead, so the parameter is the selectivity at a bin the data inform, and every bin below it takes $(b/b_{\text{start}})^{2}$ times the selectivity there. Set it to the first bin the compositions have whenever those start above the population’s first bin.

Unanchored limbs. fish_sel_dbnrml_raw and srv_sel_dbnrml_raw (0/1 matrices, one row per fleet, columns for the ascending and descending limbs) leave a limb as the raw Gaussian $\exp(-(x - \text{peak})^2 / \text{width})$ built up to the apical value, with no rescaling to hit an endpoint; the default anchors both limbs. The two forms differ by the Gaussian’s value at the end bin, which is negligible when the peak sits many widths from that bin and not otherwise. Leave a limb raw when its end-bin parameter is not something the data can speak to, and anchored when the selectivity at the end bin is a quantity worth estimating.

Time-varying deviations are untouched by any of these options and still apply per sex to the effective parameters.

Normalization and length-based selectivity

Selectivity can be normalized relative to a specific bin, to the maximum, or to the mean across a bin range. When fit_lengths = 1, selectivity operates on length bins and is mapped to age-space via a user-supplied size-age transition matrix (SizeAgeTrans).

The "nonparlog" form has its own standardization: each year’s curve is divided by the mean of $\exp(\theta)$ over a window of bins $\mathcal{B}$,

$$\text{sel}_b = \exp(\theta_b) \Big/ \frac{1}{|\mathcal{B}|}\sum_{b' \in \mathcal{B}} \exp(\theta_{b'})$$

and that window is an input:

Argument	Description
fish_sel_norm_bins / srv_sel_norm_bins	The bins $\mathcal{B}$ each fleet standardizes over, supplied to Setup_Mod_Fishsel_and_Q() and Setup_Mod_Srvsel_and_Q() respectively. A list with one element per fleet naming that fleet’s bins, or NULL (the default) for every bin; a fleet’s element may also be NULL on its own. Stored as a 0/1 array [n_bins × n_fleets] and read only by "nonparlog" fleets

Why the window matters. A gear whose catchability is defined against only part of the bin range standardizes over that part. Two windows differ by a constant multiplier on the whole curve, which catchability would absorb completely if $q$ were free. Under an informative prior on $q$ it is not free, and the mismatch is instead levered onto the selectivity parameters by the factor $\partial(\text{index nLL}) / \partial \log q$, which in the BSAI Atka mackerel bridge is $-14.055$ (see vignette("ab_bsai_atka_mackerel_case_study")). A $q$ prior is also stated against a particular standardization, so the window and the prior mean have to follow the same convention or the same prior statement lands on a different number. With $q$ free and unpenalized the choice makes no difference to the fit.

Selectivity priors

Lognormal priors on selectivity parameters are toggled via Use_fish_selex_prior / Use_srv_selex_prior, with hyperparameters supplied in a data frame (fish_selex_prior / srv_selex_prior) containing region, fleet, block, sex, par, mu, and sd.

Selectivity smoothness penalty weights

All selectivity smoothness/regularization penalty weights are configured in a single place, Setup_Mod_Weighting(), via the fish_sel_pen_wts, ret_sel_pen_wts, and srv_sel_pen_wts arguments (one list per selectivity surface). Each is a named list/vector; any name not supplied defaults to 0 (off). All six terms are evaluated directly on a fleet’s realized selectivity-at-bin-at-year surface, so they apply to any selectivity functional form and any fleet, regardless of block/deviation structure.

Name	Description
smooth_dome	Hinge penalty discouraging decreases across adjacent bins (dome-shape control)
smooth_bin_curve	Second-difference penalty across bins, normalized by the number of fitted bins
smooth_bin_diff	Unconditional first-difference penalty across bins (both increases and decreases contribute, unlike smooth_dome), normalized by the number of fitted bins
smooth_yr_diff	First-difference penalty across years, normalized by the number of fitted years
smooth_yr_curve	Second-difference penalty across years, normalized by the number of fitted years
smooth_mean_center	Penalizes the per-year mean of log-selectivity away from zero; resolves the scale indeterminacy of the bicubic surface (a uniform per-year shift in log-selectivity otherwise trades off exactly against that year’s fishing mortality)

Per-fleet specifications. A single named specification is shared by every fleet; alternatively, pass an unnamed list with one named specification per fleet (use list() or NULL for fleets with no penalties) so, e.g., two surveys can have different smoothing.

Per-year weights, bin ranges, and normalization. Each weight may also be a vector with one value per model year (0 skips that year), so a penalty can act only in years where selectivity changes, or with year-specific strength. A specification may additionally have:

Name	Description
bin_range	Length-two vector c(first, last) restricting the bins the penalties act over; or a named list giving each term its own range. Confine a shape penalty (dome, curvature) to the older ages where the curve should flatten without constraining the ascending limb
normalize	TRUE (default) divides bin-wise weights by the number of penalized bins and year-wise weights by the number of years; settable per term. Turn off when weights are calibrated as explicit variances
yr_diff_ref	Reference log-selectivity vector for smooth_yr_diff’s first penalized year, which otherwise has no predecessor and goes unpenalized. Anchors an otherwise free series to a known selectivity before the data begin

A useful identity: with normalize = FALSE, a per-year smooth_yr_diff weight of $1/(2\sigma_y^2)$ is exactly the negative log-kernel of a random walk on the realized log-selectivity with year-specific standard deviation $\sigma_y$. This is how selectivity random walks with tabulated per-year sigmas are reproduced as penalties on the curve rather than as deviation parameters.

Recommendations. Prefer the distributional process-error forms (cont_tv_*) when you want an estimated, interpretable sigma; prefer these penalties when bridging assessments whose selectivity smoothing is a tuned penalty, or when regularizing a "bicubic"/"nonparlog" surface. smooth_mean_center (or the centering penalty above) should accompany any form with a free scale. Weights are on the scale of squared log-selectivity differences; start small (1 to 10), inspect the realized surfaces, and remember that normalized and unnormalized weights differ by a factor of $n_\text{bins}$ or $n_\text{yrs}$.

---

Catchability

Configured alongside selectivity in Setup_Mod_Fishsel_and_Q() / Setup_Mod_Srvsel_and_Q().

Estimation structure (fish_q_spec / srv_q_spec)

String	Description
"est_all"	Free parameter for each region × fleet combination
"est_shared_r"	One parameter shared across regions (per fleet)
"fix"	Kept at user-supplied value

Time blocks (fish_q_blocks / srv_q_blocks)

Same syntax as selectivity blocks: "none_Fleet_<f>" for a single block, or "Block_<b>_Year_<s>-<e>_Fleet_<f>" for structured change points.

Analytic catchability (srv_q_type / fish_q_type)

Per fleet, catchability can be concentrated out of the likelihood rather than estimated. srv_q_type is set in Setup_Mod_Srvsel_and_Q() and fish_q_type in Setup_Mod_Fishsel_and_Q(); the two take the same strings and behave identically, applied to whichever index that fleet fits:

String	Description
"est"	Estimate ln_srv_q / ln_fish_q (default)
"arith"	Solve analytically as the ratio of mean observed to mean predicted index
"geo"	Solve analytically on the log scale, $\hat{q} = \exp(\overline{\log \text{obs} - \log \text{pred}})$

Both analytic forms use only the years with observations and automatically fix that fleet’s catchability parameter regardless of srv_q_spec / fish_q_spec; they are incompatible with catchability covariates and priors on that fleet.

The solve is done within each catchability time block. A fleet with a single block, the default, gets one $q$ for the whole series; a fleet with srv_q_blocks = c("Block_1_Year_1-30_Fleet_1", "Block_2_Year_31-terminal_Fleet_1") gets two, each from the observations its own block owns. Writing $w_b(i)$ for the weight block $b$ has on observation $i$, which is one for the observations in that block and zero elsewhere:

$$\hat q_b^{\,\text{arith}}=\frac{\sum_i w_b(i)\,o_i}{\sum_i w_b(i)\,p_i}, \qquad \log \hat q_b^{\,\text{geo}}=\frac{\sum_i w_b(i)\left(\log o_i-\log p_i\right)}{\sum_i w_b(i)}$$

A block that owns no observations has nothing to solve from and takes the pooled value over the whole series.

Recommendations. Use "geo" with a lognormal index likelihood: it is the exact maximum-likelihood solution under a shared SE, removes one ridge-prone parameter per fleet, and typically speeds and stabilizes optimization. "arith" exists to match assessments that use the arithmetic ratio. Blocks are the way to let a solved $q$ shift: they cost nothing and the solve stays exact within each block. Stay with "est" whenever you need a $q$ prior, covariates, or when the index’s absolute scale is informative (e.g., a swept-area survey with a strong prior near 1, where concentrating $q$ out would discard that information).

Priors

Argument	Description
Use_fish_q_prior / Use_srv_q_prior	0/1 toggle
fish_q_prior / srv_q_prior	Data frame with region, fleet, block, mu (natural scale), sd (log scale). Penalty: $\text{Normal}(\ln(\mu), \sigma)$

---

Natural Mortality

Controlled via Setup_Mod_Biologicals().

Estimation vs. fixing (M_spec)

String	Description
"est_ln_M"	Estimate $\ln M$ across the defined block structure
"fix"	Fix at values supplied via Fixed_natmort (array [p × r × y × a × s], or [p × r × y × seas × a × s] when the rate differs between seasons)

What the value means

$M$ is a rate per year, and mortality within a season is that rate times the season’s duration. Four equal seasons at $M = 0.2$ removes $0.05$ per season and $0.2$ over the year, same as a single season model. Adding the season axis changed none of that.

So a rate is not a share. Two half-year seasons at $0.2$ and $0.4$ accumulate $0.1$ and $0.2$, an annual $0.3$; they don’t split $0.3$ between them. To get a stock whose season-1 share of annual $M$ is $0.38$, with two half-year seasons and annual $M$ of $0.3$, set season 1 to $0.38 \times 0.3 / 0.5 = 0.228$ and season 2 to $0.62 \times 0.3 / 0.5 = 0.372$.

Block structure

Natural mortality parameters can be shared or made specific along any combination of six axes. Each argument takes either 'constant' (all levels pooled) or a list of integer vectors defining blocks:

Argument	Axis	Example
M_popblk_spec	Population	list(c(1,2), 3) → pops 1-2 share $M$, pop 3 is separate
M_regionblk_spec	Region	list(1:3, 4:5)
M_yearblk_spec	Year	list(1:20, 21:40)
M_seasblk_spec	Season	list(1, 2) → a separate rate in each of two seasons; list(1:2, 3:4) splits a four season year in half
M_ageblk_spec	Age	list(1:5, 6:30) → young vs. old
M_sexblk_spec	Sex	list(1, 2) → sex-specific $M$

All block specifications are crossed to produce unique $M$ parameters. For instance, two age blocks × two sex blocks = four estimated $\ln M$ values (assuming everything else is 'constant').

Priors

Lognormal priors on $M$ are activated via Use_M_prior = 1, supplying a data frame with columns popblk, regionblk, yearblk, ageblk, sexblk, mu, and sd.

---

State-Space Numbers at Age

Set in Setup_Mod_Biologicals(). By default (NAA_re = "none") the numbers at age past the first are whatever mortality and ageing produce, and only recruitment and the initial age structure are stochastic. Turning the state on makes the log numbers themselves parameters for ages two and older, with the deterministic step as the prediction they are compared against. Further details can be found in the model equations vignette.

Turning it on

Argument	Description
NAA_re	"none" (default), or the structure over the age-year grid: "iid" (independent across both ages and years), "1dar1_a" (autoregression over ages, years independent), "1dar1_y" (over years, ages independent), "2dar1" (separable over both), "3dcond" and "3dmarg" (three-dimensional field over age, year and cohort, on the conditional or marginal variance)
NAA_re_ages, NAA_re_years	Ages and calendar years the state covers. NULL (default) uses everything from the second onward. Each must be a contiguous run, because the penalty compares the likelihood for one rectangular matrix
NAA_re_where	Population by region matrix, 1 where the state runs and 0 where a population never occupies that region. NULL (default) gives every cell a state
NAA_re_seasons	Seasons the state covers. "annual" (default) puts a state at season one only, so the numbers within a year stay deterministic and the innovation is purely annual; "all" puts one at the start of every season; an integer vector selects specific seasons and need not be contiguous
NAA_sigma_spec	"est" (default) or "fix", whether the process error standard deviations are estimated. The states themselves are always estimated

A population that never reaches a region

Under natal homing a population can be absent from a region entirely: the southern black sea bass component never enters the north, so its numbers there are zero in every year and age. A lognormal state on such a cell is undefined, and the penalty would take the logarithm of zero and return an infinite likelihood. NAA_re_where says which population and region cells hold a state:

# two populations, two regions, the second homing to region two and never leaving it
NAA_re_where = matrix(c(1, 0,
                        1, 1), nrow = 2, byrow = TRUE)

Cells set to 0 are dropped from the map, so they are not estimated, and from the penalty, so they contribute nothing. They are refused when the state correlates across regions or populations, since a dropped cell would otherwise sit inside a joint density it is being left out of. Leave the argument alone for a single population, or wherever every population can reach every region.

A cell of the state is the log numbers at the start of its season, so season one is the year boundary, after ageing and the plus group, and later seasons sit on the within-year survival and movement step where no ageing happens. The annual state is therefore the special case of the seasonal one, not a different model, and "annual" reproduces what the package did before seasons were an option.

Whether to go past "annual" is a question about data, not about routines. Within a year nothing structurally separates one season from the next except the observations assigned to them, so a season with no data returns its prior as its posterior and the seasons that do have data trade against it. NAA_re_seasons takes an arbitrary set for this reason: put states only in the seasons that have observations. It is the same weak-information story as a seasonal split of $M$, and the same advice follows.

The state covers the assessment years only. Do_Population_Projection() does not have it, so projected numbers at age advance deterministically from the terminal year with recruitment the only stochastic element, and a forecast from a state-space fit omits this process error. The closed loop operating model is a separate path and does have it forward, extending the active years over the projection period.

The split of what varies by season is deliberate. NAA_pe_pars has no season dim, so every active season shares the age, year and cohort correlations within a population, region and sex; the standard deviation is what varies, through NAA_sigma_seasblk_spec. A long season and a short one differ in how much process error accumulates, not in the shape of the correlation across ages.

Correlation across populations, regions, seasons and sexes

Argument	Description
NAA_re_pop, NAA_re_region, NAA_re_sex, NAA_re_season	"iid" (default) or "us", an unstructured correlation over that dim, $n(n-1)/2$ parameters. Composed with whatever NAA_re gives over the age-year grid, so a correlation over one dim says nothing about the others. The season correlation runs over the active seasons, so it needs more than one of them
NAA_re_region_spec, NAA_re_season_spec	How the region and season correlations are shared: "est_all" (default), "est_shared_p", "est_shared_s", "est_shared_p_s", "fix"
NAA_pe_spec	How the age, year and cohort correlations from NAA_re are shared across these same dims: "est_all" (default, a free set per population, region and sex), "est_shared_p", "est_shared_r", "est_shared_s", "est_shared_p_r", "est_shared_p_s", "est_shared_r_s", "est_shared_p_r_s" (one set for the whole model), "fix"

Unstructured is the only non-independent option on these dims because none of them has an ordering: regions are not in a line, so a correlation decaying with distance has nothing to decay along.

Sharing a correlation and correlating the innovations are different operations, and NAA_pe_spec does the first. "est_shared_r" gives every region the same $\rho$ while their shocks stay independent, which cuts the count: a three-region "2dar1" model goes from six correlations to two. NAA_re_region = "us" instead makes the shocks themselves covary, and adds parameters. Neither touches which cells of the state are estimated, only how many hyperparameters they draw on. Worth reaching for when a dim is thin, since each region otherwise estimates its own age correlation from however many ages NAA_re_ages leaves it.

Standard deviations

Argument	Description
NAA_sigma_popblk_spec, NAA_sigma_regionblk_spec, NAA_sigma_yearblk_spec, NAA_sigma_seasblk_spec, NAA_sigma_ageblk_spec, NAA_sigma_sexblk_spec	"constant" (default) or a list of integer vectors assigning indices to blocks, exactly as the M_*blk_spec arguments do

Blocking shares a standard deviation; it never removes a cell from the state. Only NAA_re = "iid" admits a standard deviation that varies over years or ages: every other form is separable or Markov in one of those dims, and a per-cell variance is neither. The setup errors rather than silently ignoring a year or age block under a correlated form.

The season dim is the exception, because it is whitened before the age and year density is reached rather than inside it. A season-varying standard deviation therefore works under every NAA_re, which is usually what a model with unequal season durations wants, and the one thing that rules it out is NAA_re_season = "us": a correlation has one scale across the dim it spans.

Movement

Configured via Setup_Mod_Movement() (EM) or Setup_Sim_Movement() (OM).

SPoRC implements two movement parameterizations, plus a fixed-matrix escape hatch.

Movement type (move_type)

Value	Description
0	Unstructured Markov. Transition probabilities estimated via multinomial logit (softmax), with region 1 as the implicit reference category. Supports blocking along population, age, year, season, and sex dimensions
1	CTMC. A continuous-time Markov chain builds an instantaneous rate matrix $Q = D + Z$ from diffusion ($D$, isotropic dispersal scaled by region area) and taxis ($Z$, directional preference). Transition probabilities are obtained by matrix exponentiation: $\mathbf{P} = \exp(Q)$

Fixed movement (use_fixed_movement = 1)

Bypasses estimation entirely. Supply Fixed_Movement as an array [p × r_\text{from} × r_\text{to} × y × \tau × a × s].

Unstructured Markov block structure (move_type = 0)

Each dimension can be pooled ('constant') or blocked:

Argument	Description
Movement_popblk_spec	Population blocks
Movement_ageblk_spec	Age blocks (e.g., list(1:5, 6:30) for age-dependent movement)
Movement_yearblk_spec	Year blocks
Movement_seasblk_spec	Season blocks
Movement_sexblk_spec	Sex blocks

CTMC configuration (move_type = 1)

Argument	Description
ctmc_move_dat	Data frame of covariates with required columns pop, regions, years, seas, ages, sexes, plus any variables referenced in the formulas. For projection years, covariate lookups are capped at the last historical year unless extended rows are supplied
diffusion_formula	R formula for diffusion (e.g., ~ 1, ~ bs(years, df = 4))
preference_formula	R formula for taxis/preference
adjacency_mat	Square [n_regions × n_regions] binary connectivity matrix, with a zero diagonal: residency falls out of the generator and is not specified here. A fully connected matrix is 1 - diag(n_regions); note that diag(1, n_regions) is the identity, which connects nothing and is rejected
area_r	Numeric vector of region areas (scales diffusion rates)
ctmc_diffusion_bounds	How the off-diagonal generator entries are kept non-negative, evaluated on the adjacency edges only. With d the preference gradient along an edge and theta the diffusion rate out of its origin: "none" (or 0) leaves theta + d unbounded; "softplus" (or 1) smooths its positive part with width ctmc_diffusion_eps; "upwind" (or 2) is the discontinuous Galerkin (finite volume) flux theta + max(d, 0)
	"upwind" has diffusion whole and hinges only the taxis, so positivity in the rate estiamtes don’t depend on the two cancelling
ctmc_diffusion_eps	Positive numeric (default 0.1): softplus width for ctmc_diffusion_bounds = "softplus", read by no other form. An edge where taxis exactly cancels diffusion has eps * log(2), so eps is a minimum exchange rate as well as a smoothing width and the fit moves with it; 0.001 behaves like a hard clamp and leaves shut edges gradient-dead
ctmc_scale_by_seasdur	0/1 (default 1): if 1, treats the generator as an annual rate and exponentiates Q * seasdur[s] each season, so movement and mortality share time units. 0 results in Q * 1 each season.
move_expm_nsub	Integer (default 0): how matrix exponentials of the generator are evaluated. 0 is exact (Matrix::expm); a power of two n >= 1 substitutes n implicit backward Euler substeps. See below

Movement and mortality sequencing (move_timing)

Controls how movement and mortality are ordered within a season. The three options coincide exactly when total mortality is constant across regions, when mortality is zero, and when movement is absent, so single-region models are unaffected by the setting.

Value	Sequencing
0 (default)	Movement then mortality. Historical SPoRC behavior
1	Mortality then movement
2	Continuous: movement and mortality act simultaneously, via expm(Q * seasdur - diag(Z))

move_timing = 2 requires an estimated CTMC generator (move_type = 1 with use_fixed_movement = 0); a discrete movement matrix has no guaranteed real generator, so it is rejected rather than approximated.

The setting is not confined to the projection step. Anything observed partway through a season has to be evaluated consistently with it:

Quantity	0	1	2
Spawning biomass	Post-movement location	Pre-movement location	Partially redistributed, expm(A * t_spawn)
Catch, discards	Baranov, post-movement	Baranov, pre-movement	Season-integrated abundance (spatial Baranov)
Fishery index	Post-movement N	Pre-movement N	Season-integrated abundance
Survey index	N * exp(-t_srv * Z)	N * exp(-t_srv * Z)	Snapshot, expm(A * t_srv) N
Tag recaptures	Baranov on tag cohort	Baranov on tag cohort	Season-integrated tag abundance

where A = t(Q) * seasdur - diag(Z). Equilibrium initialization (including the matrix plus-group solution), per-recruit reference points, and projections all use the same seasonal operator, so they inherit the setting as well. Note the survey/fishery contrast: a survey is a snapshot at an instant within the season and so uses a partial propagation, whereas catch and the fishery index accumulate over the season and so use the integral of that propagation. Full equations are in the Description of Model Equations vignette.

Evaluating the matrix exponential (move_expm_nsub)

Both places SPoRC turns a generator into fractions go through one routine: Get_Movement, converting Q into the Movement array, and every move_timing = 2 seasonal operator. move_expm_nsub selects how that exponential is evaluated.

Value	Scheme
0 (default)	Exact, Matrix::expm(A)
n >= 1, a power of two	Implicit (backward Euler), (I - A/n)^-n

n = 1 is plain solve(I - A). Larger n applies the same step over n subintervals, reached by squaring the substep operator: each squaring doubles the exponent, so n costs one linear solve plus log2(n) matrix products. That is why n is restricted to powers of two, which squaring reaches exactly. The scheme is first order either way, so the finer ladder a general integer would allow buys nothing.

The reason to want this is the cost of the gradient. Under reverse-mode AD the derivative of a linear solve is another solve, whereas the derivative of a matrix exponential is much more expensive than its forward evaluation. Benchmarking the exponential on its own over 600 six-region strata, the gradient runs about 18x faster at n = 1, 4x at n = 8 and 2.5x at n = 64, against a tape that is slightly more expensive to build.

Inside a real fit the gain is far smaller, because the exponential is only a slice of the gradient and binary powering’s log2(n) matrix products claw back the rest as n grows. On a three-region, one-season, 25-year estimation model the measured whole-model speedups were 2.7x at n = 1, 1.3x at n = 64 and 1.1x at n = 512 under move_timing = 2, and about 1.2x at any n under move_timing = 0. Since the accurate settings are the slow ones, the option is worth adopting only where the exponential is a large share of the work: cost grows with the cube of n_regions and with the number of (pop, year, season, age, sex) strata, so measure on the model you intend to fit rather than assuming the isolated figure.

What it costs is accuracy, and not by a small dim. The scheme is first order in 1/n, so it is a different discretization of the movement process rather than a cheap route to the same numbers. On a six-region chain generator with Z between 0.25 and 0.42:

n	max abs error in the operator	error in survival
1	1.3e-1	4.6e-2
8	2.1e-2	7.0e-3
64	2.8e-3	9.0e-4
1024	1.8e-4	5.6e-5

At n = 1 survival is 1 / (1 + Z) where the exponential gives exp(-Z). Since substeps are nearly free, n in the hundreds is the sensible setting if the intent is to approximate the exponential; small n is only defensible if backward Euler is being adopted as the movement process itself.

Two structural properties survive at any n, which is what makes the option safe to use, not just fast. the inverse of I - A/n is non-negative, so movement fractions cannot go negative (it is a non-singular M-matrix); with no mortality the result is exactly column-stochastic, so fractions still sum to one without renormalization. And because the seasonal operator and the catch integral come from the same Van Loan block, 1'T = 1' - z'Integral holds algebraically, so catch plus survivors equals the starting abundance to machine precision even where the operator itself is a coarse approximation.

Continuous movement deviations (cont_vary_movement)

Origin-destination deviations applied multiplicatively to off-diagonal rates (CTMC) or additively on the logit scale (unstructured).

String	Deviation structure
"none"	No deviations (default)
"iid_y"	Year only
"iid_a"	Age only
"iid_y_a"	Year × age
"iid_y_a_s"	Year × age × sex
"iid_y_seas_a_s"	Year × season × age × sex
"iid_p_y"	Population × year
"iid_p_a"	Population × age
"iid_p_y_a"	Population × year × age
"iid_p_y_a_s"	Population × year × age × sex
"iid_p_y_seas_a_s"	Population × year × season × age × sex

Additional movement controls

Argument	Description
do_recruits_move	0 = age-1 fish do not move; 1 = recruits follow the movement matrix
Use_Movement_Prior	0/1 toggle for Dirichlet movement priors (unstructured)
Movement_prior	Data frame with pop, region_from, year, seas, age, sex, alpha (Dirichlet concentration vector of length n_regions)
Movement_cont_pe_pars_spec	Estimation structure for process-error hyperparameters: "none", "fix", "est_all", "est_shared"

---

Catch, Fishing Mortality, and Discards

Configured via Setup_Mod_Catch_and_F().

Catch conditioning

Argument	Description
ObsCatch	Observed aggregate catch [r × y × τ × f]
ObsCatch_pop	Population-specific catch [p × r × y × τ × f]
catch_units	Per-fleet: 0 = abundance, 1 = biomass
Use_F_pen	0/1 toggle for the fishing mortality deviation penalty (Fmort_nLL); see Which deviations are penalized
sigmaC_spec / sigmaC_pop_spec	Catch observation-error estimation

Fishing mortality process error

ln_F_devs (annual log-scale deviations about ln_F_mean) can follow one of three process-error structures, set via Fdev_model:

String	Description
"iid"	Independent annual deviations (default)
"rw"	Random walk
"ar1"	First-order autoregressive

Fdev_pen_center ("fixed", the default, or "own_mean") sets where the iid deviation penalty is centered; see the deviation penalty centering discussion in the Recruitment section for the shared rationale. Under the mean-plus-deviations parameterization the level of $F$ is already set by ln_F_mean, so "own_mean" penalizes only the deviations’ spread and avoids constraining the level twice, matching assessments whose deviation vectors sum to zero. It leaves ln_F_mean and the deviations’ level mutually unidentified unless one is fixed, which ln_F_mean_spec = "fix" does; setup warns when it is combined with an estimated mean in configurations where nothing else reads the mean. Note also that the default iid penalty with "fixed" centering acts as a hidden “biomass tracks catch” prior when catch trends strongly; Fdev_model = "rw" is usually the better remedy for that than re-centering.

fishsel_dont_est_dev_first, and its retsel_ and srvsel_ siblings, start a selectivity walk in year two so the fixed selectivity parameters hold year one. It matters for the non-parametric forms ("nonpar", "nonparlog", "nonparfree"), which have one free base parameter per bin: year one’s deviation is then that same value written twice, and the walk’s first-year normal, *_rw_init_sigma, is the only thing separating them. That normal is a prior on a level usually meant to be free, and it costs $n_{\text{bins}} \log(\sigma_0\sqrt{2\pi})$ at the optimum, where the deviation sits at zero and the base parameter has absorbed the level. Default 0 for every fleet, since turning it on shifts a model’s reported likelihood by that constant. Refused for the GMRF and 2D AR1 forms, whose deviations are a field over years and bins rather than a walk anchored at year one.

ln_F_mean_spec ("est", the default, or "fix") chooses the fishing mortality parameterization. "est" is the mean-plus-deviations form. "fix" maps ln_F_mean off at its starting value (zero unless supplied), so the deviations are free annual log-F outright: $F = \exp(\epsilon)$. It must be paired with Fdev_pen_center = "own_mean", Fdev_model = "rw", or Use_F_pen = 0 — an "iid" or "ar1" penalty centered on the fixed zero mean would shrink the deviations toward $F = 1$, and setup rejects that combination.

sigmaF_spec controls sharing/fixing of the process-error standard deviation (ln_sigmaF) across region × season × fleet, using the same "est_all" / "est_shared_<dims>" / "fix" convention described in Parameter sharing conventions. Fdev_rho_spec analogously controls the AR1 correlation parameter (Fdev_rho) and is only active when Fdev_model = "ar1"; for "iid" or "rw", Fdev_rho is unused and mapped to NA regardless of what is supplied.

Catch-active years (UseCatch == 1 or any UseCatch_pop == 1) do not need to be contiguous under "rw" or "ar1", a fishery may close for several years and reopen later. The transition between two active years separated by a gap of $d$ closed years is taken directly over the elapsed gap (the same marginal transition obtained by estimating deviations for the closed years and integrating them out, without actually estimating them). See [Fishing Mortality Deviations] in vignette("c_model_equations") for the exact equations, and Which deviations are penalized for how the active sequence is determined.

Discard and retention framework

SPoRC decomposes total fishing mortality at age into retained and dead-discard components:

$$F_{p,r,\tau,a,s} = \sum_f F_{r,\tau,f} \left[\underbrace{\text{sel}_{a,f} \cdot \text{ret}_{a,f}}_{\text{retained}} + \underbrace{\text{sel}_{a,f} \cdot (1 - \text{ret}_{a,f}) \cdot \text{dmr}_f}_{\text{dead discards}}\right]$$

Argument	Description
ret_sel_model	Retention selectivity functional form (same options as fish_sel_model)
ret_sel_blocks	Retention selectivity time blocks
cont_tv_ret_sel	Continuous time-varying retention selectivity
dmr_mean_spec	Estimation structure for dead discard mortality rate means
dmr_dev_spec	Estimation structure for DMR deviations: "fix" (default) or "est_all"
Use_dmr_pen	0/1 toggle for the DMR deviation penalty (dmr_nLL); must be 1 when dmr_dev_spec = "est_all"
discard_units	Per-fleet: 0 = abundance, 1 = biomass, 2 = abundance fraction, 3 = biomass fraction
ObsDiscard / ObsDiscard_pop	Observed aggregate and population-specific discards
sigmaD_spec / sigmaD_pop_spec	Discard observation-error estimation

Estimating discard mortality rates

dmr_dev_spec = "est_all" estimates an annual deviation in every fished region × year × season × fleet cell, not only in cells with a discard observation. This is deliberate. From the decomposition above, dmr scales the dead-discard component of total mortality, so it enters $Z$ and therefore propagates into predicted retained catch, indices, and compositions:

$$C_{a} = \frac{F^{\text{ret}}_{a}}{Z_{a}} N_{a} \left(1 - e^{-Z_{a}}\right)$$

A cell is therefore informative about dmr whenever it is fished and retention is less than one, with or without discard data. Discard observations are not the dividing line, dmr in fact cancels out of the predicted discard itself, which divides the dead discards back through by it, leaving only the same $Z$ dependence.

The cells excluded are true closures, matching the condition under which the objective pins dmr to zero: no aggregate or population-specific catch is fit, and the aggregate catch observation is a recorded zero rather than missing (NA). A missing observation is treated as a fished year, on the assumption that fishing continued and we simply lack a value to fit.

Note that a fleet with full retention (ret_sel = 1, the package default) contributes no dead discards at all, so dmr drops out of the objective entirely and its deviations have a flat gradient. dmr_dev_spec and dmr_mean_spec are only meaningful once retention is actually modeled, see vignette("s_discard_retention").

Which deviations are penalized

The ln_F_devs, logit_dmr_devs and ln_RecDevs penalties are evaluated on exactly the deviations that are estimated, read off the map mirrors map_ln_F_devs, map_logit_dmr_devs and map_ln_RecDevs in the data list rather than recomputed from the catch indicators or the deviation index ranges. fit_model() refreshes those mirrors from the map immediately before building the model, so a deviation mapped off by hand after setup is neither estimated nor penalized:

# drop dmr deviations for a stretch of years, penalty included
mp <- array(input_list$map$logit_dmr_devs, dim = dim(input_list$par$logit_dmr_devs))
mp[, 10:20, , ] <- NA
input_list$map$logit_dmr_devs <- factor(mp)

Two consequences. Mapping a deviation off fixes it at whatever sits in $par, which is 0 by default, so dmr falls back on plogis(logit_dmr_mean), but a starting value supplied through ... is kept, not reset to zero. And under Fdev_model = "rw" or "ar1", a fixed deviation is dropped from the active sequence, widening the gap $d$ between the deviations either side of it rather than being treated as an actual year.

For ln_RecDevs this matters most when ln_sigmaR is estimated. A deviation fixed at zero still sits at the penalty’s mean, so were it penalized it would add to the sum of squares’ denominator without adding any spread, and $\sigma_R$ would be pulled low by roughly $\sqrt{n_{est} / n_{total}}$. Holding the first 80 of 255 quarterly deviations, for instance, biases $\sigma_R$ down by about 17%.

With ln_sigmaR fixed the excluded terms are constants: the objective shifts by a fixed amount and its gradient is unchanged everywhere, so the likelihood surface and the location of any optimum are untouched, and Rec_nLL becomes comparable across runs that hold different numbers of deviations. Note that this is not the same as the fitted values being unchanged. nlminb tests convergence on relative changes in the objective’s value, so on a model that stops short of convergence (a flat ridge, or a "singular convergence" or "false convergence" message), shifting the objective by a constant can move where the optimizer halts. If excluding these terms visibly changes your estimates, the identifiability of the model is the thing to look at, not the penalty.

---

Biological Inputs

Supplied via Setup_Mod_Biologicals() (EM) or Setup_Sim_Biologicals() (OM).

Input	Dimensions	Description
WAA	[p × r × y × τ × a × s]	Spawning weight-at-age (used for SSB)
WAA_fish	[p × r × y × τ × a × s × f]	Fishery-specific weight-at-age. Defaults to WAA if NULL
WAA_srv	[p × r × y × τ × a × s × f]	Survey-specific weight-at-age. Defaults to WAA if NULL
MatAA	[p × r × y × τ × a × s]	Maturity-at-age (proportions in [0, 1])
SizeAgeTrans	[p × r × y × τ × l × a × s]	Size-age transition matrix, shared by every fleet unless SizeAgeTrans_fish/SizeAgeTrans_srv overrides it. Required when fit_lengths = 1
SizeAgeTrans_fish, SizeAgeTrans_srv	NULL or [p × r × y × τ × l × a × s × f]	Optional per-fleet size-age keys, SizeAgeTrans with a trailing fleet dimension appended. NULL (default) reads every fleet from the shared SizeAgeTrans. Only valid under growth_model = "none"; a growth model already derives one key per fleet at that fleet’s own timing and rejects these to avoid two sources for the same quantity. The fixed-data counterpart of WAA_fish/WAA_srv overriding the shared WAA, and of Setup_Sim_Biologicals’s SizeAgeTrans_fish_input/SizeAgeTrans_srv_input on the operating-model side
fit_lengths	0/1	Toggle for fitting length compositions
do_caal	0/1	Toggle for reporting the joint arrays at length and age (Fish_caal, Fish_caal_discard, Srv_caal), dimensioned [p × r × y × τ × l × a × s × f] like CAA and CAL. Requires fit_lengths = 1, and is switched on by itself when conditional age-at-length data are supplied, since their expected rows are read from these arrays. Each array has a length and an age dimension at once, so they are the largest objects the model builds and every cell lands on the AD tape; leave off unless they are wanted
growth_model	"none" / "vb_schnute" / "richards"	"vb_schnute" builds each fleet’s size-age transition inside the model from von Bertalanffy parameters in Schnute’s form (L1 at growth_A1, L2 at growth_A2, K, CV1, CV2): linear growth from growth_L0 at age zero up to L1, CV interpolated on length between the references, the plus-group mean size adjusted for fish older than the accumulator age, and a binned normal key with the tails in the end bins, read at each fleet’s own timing (t_fish, t_srv). "richards" is the same curve applied to lengths raised to a sixth estimated parameter rho, which nests the von Bertalanffy form at rho = 1. Requires fit_lengths = 1
ln_growth_pars (via ...)	[p × r × s × n]	Starting values for L1, L2, K, CV1, CV2 (then rho under "richards"), log scale, passed through starting_values as every other parameter’s are. Defaults to the ends of the length bins with K = 0.15 and CVs of 0.1, so supply your own for any real model
growth_spec	"est_all" / "est_shared_r" / "est_shared_s" / "est_shared_r_s" / "fix"	How the growth parameters are shared across regions and sexes. growth_fix holds individual parameters at their starting values
growth_A1, growth_A2	scalars	Reference ages for L1 and L2; growth_A2 = "Linf" makes L2 the asymptote itself, with no second reference age to solve it from
growth_len_lower	[l]	Lower edges of the length bins; lens in Setup_Mod_Dim are midpoints
growth_cv_type, growth_sd_type, growth_dist	"len" / "age"; "cv" / "sd"; "normal" / "lognormal"	Whether the CV interpolates on mean length or on age between the reference ages; whether CV1, CV2 scale the mean or are standard deviations; the distribution of length at age behind the key
growth_plus_group	"mixture" / "curve"	"mixture" (default) takes the plus group’s mean length as the survivorship-weighted mixture of the ages it holds, numbers declining at an assumed 0.2 per year and length rising from the curve at the accumulator age to the asymptote; "curve" reads the curve at the accumulator age
growth_tv_model	NULL or per parameter	Time variation of the growth parameters. Either a vector of "none"/"iid"/"rw" in parameter order, or named by parameter (L1, L2, K, CV1, CV2, rho) with the rest constant. A varying parameter gains a deviation series ln_growth_devs[p, r, y, k, s] and a log sigma in the first data source of growth_pe_pars; its realized value by year is reported as growth_pars_y
growth_tv_years	NULL, a vector, or a list	Calendar years each varying parameter’s deviations are active in. NULL uses every year, a vector applies to all of them, a list names them per parameter. Deviations outside the range are kept at zero
growth_tv_link	"log" / "logit"	The scale a deviation enters on. "log" multiplies the parameter by exp(dev); "logit" keeps it strictly inside growth_par_bounds however large the deviation, which every growth parameter is
growth_par_bounds	[n × 2]	Lower and upper bounds per growth parameter, natural scale; required under the logit link
growth_tv_sigma_spec	"fix" / "est"	Whether the process error standard deviations of the deviations are estimated. Their starting values are the first data source of growth_pe_pars, one slot per growth parameter, supplied through starting_values and defaulting to log(0.1). growth_rw_init_sigma gives the first year of a random walk its own standard deviation
growth_tv_spec	"est_all" / "est_shared_r" / "est_shared_s" / "est_shared_r_s"	How the deviation series are shared across regions and sexes
growth_tv_type	"curve" / "cohort"	How the deviated parameters reach size at age. "curve" (default) reads each year’s sizes off that year’s own curve, so a fish has no memory of the years it lived through. "cohort" has size at age forward: each cohort grows by the increment the current year’s parameters imply from the size it reached, ages still in the linear phase keep the length at growth_A1 their birth year gave them, the first age past growth_A1 sits on the current year’s curve, and the plus group blends the cohort entering it with the fish already there by their numbers at age. The CV at age is then kept at the first year’s. Because the plus group reads abundance, growth under this option is evaluated inside the population dynamics year loop
growth_semipar	"none" / "iid" / "rw" / "2dar1" / "3dmarg" / "3dcond"	Semi-parametric growth: a year-by-age surface ln_growth_semipar_devs[p, r, y, a, s] multiplying the parametric mean length at age, so the curve stays the parametric part and the deviations hold departures from it. The structures are the same process errors the semi-parametric selectivity forms use; the correlated ones are what make the surface estimable, since the deviations are not identified year by year and age by age from length data alone. The spread follows the deviated mean, so under growth_cv_type = "len" a deviation also moves the fish along the CV ramp and under "age" the spread at age is untouched
growth_semipar_spec	"fix" / "est"	Whether the surface’s process error parameters are estimated. Their starting values are the second data source of growth_pe_pars, in the same slots the selectivity forms use, supplied through starting_values and defaulting to a scale of 0.05 with correlations of 0.3
growth_semipar_ages, growth_semipar_years	NULL or vectors	Ages and calendar years the surface is estimated over; outside them the deviations are kept at zero, which is how a surface is restricted to the ages the length data actually inform
LenBinMap	NULL or [l × l_\text{obs}]	Maps the model’s length bins onto the bins the compositions are recorded on, each row summing to one, for compositions on coarser bins than the population has. Expected compositions are mapped through it inside the likelihood exactly as the ageing error matrix maps ages, and the observed arrays are then dimensioned by the observed bins
AgeingError	NULL, [a × a_\text{obs}] or [y × a × a_\text{obs}]	Maps the model’s ages onto the ages the compositions are recorded on: genuine misclassification (a true age 7 read as 6, 7 or 8), a collapse onto coarser observed bins, or both. Each row sums to one, or to zero to drop a model age from the observations. The age-axis twin of LenBinMap, applied in the same position inside the likelihood and validated by the same check_bin_map()
AgeingError_fish, AgeingError_srv	NULL, [a × a_\text{obs} × n_f] or [y × a × a_\text{obs} × n_f]	Per-fleet ageing error, for when a fishery reading otoliths and a survey reading scales do not misclassify the same way. NULL (default) gives every fleet the shared AgeingError, which is what a model written before these existed does. Every fleet must land on the same observed age bins, since the observed composition arrays have one age dimension shared across fleets. Wired through the operating model too, as AgeingError_fish_input / AgeingError_srv_input in Setup_Sim_Biologicals
waa_model, wt_len_pars	"data" / "wt_len"; 2 values or [p × r × s × 2]	"data" (default) reads weight at age from the WAA* arguments. "wt_len" derives WAA (at spawning time) and each fleet’s WAA_fish and WAA_srv (at that fleet’s t_fish or t_srv) from its key and $W = aL^b$ at the bin midpoints instead of reading them as data. Get_Reference_Points reads the derived rep$WAA when it is present, falling back to data$WAA otherwise; that fallback is a zero placeholder under "wt_len", so pass rep from a fitted model rather than calling it on data alone

Growth is not seasonal. The curve depends only on real elapsed time since A1, so splitting a year into more seasons changes nothing about size at age: a von Bertalanffy or Richards curve with a single K gives the same length whether the time between two points is one whole year or several sub-annual steps that sum to it. There is no option to give K (or any other growth parameter) a within-year seasonal pattern, e.g. a multiplier per season the way Stock Synthesis’s MGparm_seas_effects does. Every fleet’s timing (t_fish, t_srv, t_spawn) still reads its own point on that single annual curve, so within-year growth differences across fleets are entirely a function of when in the year they sample, not of the growth rate itself varying by season.

---

Observation Model: Indices and Compositions

Indices and compositions are configured together for fishery fleets (Setup_Mod_FishIdx_and_Comps()) and survey fleets (Setup_Mod_SrvIdx_and_Comps()).

Indices of abundance

Argument	Description
ObsFishIdx / ObsSrvIdx	Observed index values
ObsFishIdx_SE / ObsSrvIdx_SE	Log-scale standard errors
fish_idx_type / srv_idx_type	"abd" abundance, "biom" biomass, "none". Survey fleets also take "recdev", an index of the recruitment deviations themselves
ObsFishIdx_pop / ObsSrvIdx_pop	Population-specific indices (separate likelihood contribution)
ObsFishIdx_pop_SE / ObsSrvIdx_pop_SE	Population-specific index SEs

Reporting once a year in a seasonal model

A season is a step in the population dynamics, not necessarily a reporting period. Many seasonal assessments run the dynamics in months or quarters but land one catch figure and one age composition per fleet per year, and an index is often a snapshot taken at one point in a season. Every data source therefore has a _seas_Type saying which of the two it is.

Value	What it does
"spltSeas"	Fit the observation against the prediction for the season it sits in. The default, and what every data source did before this setting existed
"aggSeas"	Sum the prediction over every season of the year and fit it against one observation

The observation still lives in whichever season it was placed in, and exactly one season per region and year may be turned on in the matching Use array. More than one is an error, because each would be fit against the same year total. The likelihood, the one-step-ahead residual and the reported negative log likelihood all land in that season, and the other seasons stay at zero.

The settings are per fleet, one for each data source:

Function	Arguments
Setup_Mod_Catch_and_F()	Catch_seas_Type, Discard_seas_Type, CatchAA_seas_Type, DiscardAA_seas_Type, and their _pop counterparts
Setup_Mod_FishIdx_and_Comps()	FishIdx_seas_Type, FishAgeComps_seas_Type, FishLenComps_seas_Type, the discard versions, and their _pop counterparts
Setup_Mod_SrvIdx_and_Comps()	SrvIdx_seas_Type, SrvIdxAA_seas_Type, SrvAgeComps_seas_Type, SrvLenComps_seas_Type, and their _pop counterparts

A composition set to "aggSeas" is built from the numbers at age summed over the year before they are turned into proportions, which is not the same as averaging the seasonal compositions: a season with more catch in it counts for more.

# a fleet that lands catch all year but reports one annual total and one annual
# age composition, both placed in season 1
input_list <- Setup_Mod_Catch_and_F(
  input_list,
  ObsCatch = ObsCatch,          # the year's total sitting in season 1, zero elsewhere
  UseCatch = UseCatch,          # 1 in season 1 only
  Catch_seas_Type = "aggSeas",
  ...
)

Two things follow from the fleet still fishing all year. Fishing mortality stays estimated season by season, so one annual observation cannot tell the seasons apart: share the deviations across seasons, or fix the seasonal pattern, or the split rides on whatever else in the model happens to constrain it. And the fleet is treated as open in every season of a year it reports in, so the initial fishing mortality and the closure logic both read the year rather than the single season the observation sits in.

An index measured at a point in time belongs in its own season with its own t_srv, not aggregated. Use "aggSeas" for a data source that accumulates across the year, such as a fishery index built from a whole year of effort.

An operating model reports the same way through Setup_Sim_Fishing() and Setup_Sim_Survey(), which take the same argument names. There the annual total is written into season one by convention, with the observation error applied once to that total rather than to each season and then added up. A multivariate normal index cannot be aggregated, since it is drawn once over a covariance the season layout defines.

Age-disaggregated observations

SPoRC fits catch as an aggregate with a composition, and an index the same way. ICES age-structured assessments generally do neither: they fit catch at age and index at age directly, every age its own observation with its own catchability and its own variance.

The two are not interchangeable. Splitting an at-age observation into a total and a composition is exact for Poisson and multinomial, where the total is a genuine count sum. However, both families are lognormal, and the sum of lognormals is not lognormal, so no composition likelihood recovers the at-age statement! A fleet uses one form or the other. Supplying both for the same fleet is an error, since it is the same information stated twice.

Argument	Data source
ObsCatchAA / UseCatchAA	Retained catch at age
ObsDiscardAA / UseDiscardAA	Discards at age
ObsSrvIdxAA / UseSrvIdxAA	Survey index at age

Each has a _pop counterpart with a leading population dimension. All six evaluate the same at-age likelihood through one function, differing only in which array supplies the prediction.

How a data source is reported

Every at-age array is dimensioned [n_regions, n_years, n_seas, n_ages, n_sexes, n_fleets], the layout of the prediction arrays the likelihood reads, so the two line up dim for dim. A fleet says which of regions and sexes it reports separately through a Type, using the same vocabulary the compositions use:

String	Regions	Sexes
"agg"	summed	summed
"spltRaggS"	separate	summed. Default
"aggRspltS"	summed	separate
"spltRspltS"	separate	separate

CatchAA_Type, DiscardAA_Type and SrvIdxAA_Type take one setting for every fleet, one per fleet, or the year and fleet grammar, and each has a _pop counterpart. A dim the fleet sums over has its observation in slot one of that dim, and a use flag anywhere else on it is refused rather than quietly ignored, since the prediction has already added those cells in. That check runs per year, so a fleet reporting sexes separately early and combined later is kept to each rule in the years where it applies.

A fleet whose reporting changed part way through the series states it directly:

CatchAA_Type = c("spltRspltS_Year_1-20_Fleet_1",
                 "spltRaggS_Year_21-terminal_Fleet_1")

One combination is refused. A fleet fitting "2dar1" cannot change its aggregation between years, because a separable correlation is defined over the whole block of years by ages, and a setting that changes inside that block has nothing to mean there. Either hold the aggregation constant for that fleet, or give each period its own fleet so that each block is its own observation.

The sex dim is required whatever the fleet reports. A data source summed over sexes still has its observation in the full array, in sex slot one, and an array one dimension short is refused rather than promoted: the alternative is a silent reinterpretation of what you passed.

Density and where the observation error comes from

An at-age data source has the same choices the aggregated index has, per fleet:

Argument	Options
*AA_LikeType	"lognormal" (default) or "normal"
*AA_sigma_form	"none" (default), "data", "est_additive", "est_quadrature"
Obs*AA_SE	Reported standard errors, shaped like the observations

"none" uses the estimated parameter alone, which is what a data source with no reported errors means. "data" uses the reported standard errors alone and holds the parameter fixed, since nothing reads it. The other two combine them, added or in quadrature. Disaggregating a survey index by age therefore no longer costs you its survey-design errors.

Key matrices

Parameters are coupled through integer key arrays [n_ages, n_sexes, n_fleets], in which equal entries share a parameter and NA excludes one. This is the key matrix convention ICES assessments use for coupling variances and catchability, with NA marking an excluded parameter, and the fleet last because that is where SPoRC puts it. One structure covers every sharing pattern that would otherwise need its own argument:

Intent	Key column for a fleet and sex
One parameter per age	1 2 3 4 5
By age group	1 1 2 2 2
One for the fleet	1 1 1 1 1
Age not observed	NA

The sex dim is required here too; a key coupling the sexes says so by repeating its entries across them. An age and sex a fleet never observes has no parameter whatever the key says, which is what holds the unused sexes of a sex-summed data source out.

sigmaCAA_key, sigmaDAA_key and sigmaSrvIdxAA_key hold the observation error. Catchability at age is not keyed here: an index fit age by age has its age shape in selectivity, through the "nonparfree" form (see Selectivity above), since a free catchability per age and a selectivity estimated at age are the same statement twice. A parameter informed by fewer than two observations is refused, and fewer than five warns: an observation error standard deviation with a single observation drives the likelihood to negative infinity rather than failing outright, and the optimizer reports convergence either way.

A season needing its own observation error is a separate fleet, the same way anything else needing its own selectivity or catchability is.

Correlation across ages

Each data source is set where its data are set, following how the package separates fishery from survey elsewhere, and the population-specific form has its own setting rather than borrowing the aggregated one:

Argument	Set in
AgeObsCorr_catch, AgeObsCorr_discard (and _pop)	Setup_Mod_Catch_and_F()
AgeObsCorr_srv_idx (and _pop)	Setup_Mod_SrvIdx_and_Comps()

String	Meaning	Parameters per fleet
"iid"	Ages independent within a cell. Default	none
"1dar1"	AR(1) in age distance	1
"us"	Unstructured across ages	$n_a(n_a-1)/2$
"2dar1"	Separable AR(1) over ages and years	2

Each takes one setting for every fleet or one per fleet, so a model can leave one fleet independent and correlate another. The three ICES obsCorStruct options are "iid", "1dar1" and "us"; "2dar1" adds the year dimension.

"1dar1" is a statement about age distance, not about position in the observed vector. A fleet observing ages 2, 3, 5 and 6 has a gap, and lag one across that gap is not lag one, so the covariance is built from the ages themselves. Where the observed ages are consecutive the autoregressive recursion is used instead, which is the same density more cheaply.

"us" places no shape on how ages covary, at the cost of a parameter per pair. The parameters are unconstrained: they fill the strict lower triangle of a matrix whose rows are then normalized to unit length, which makes it a Cholesky factor, so the result is a valid correlation matrix for any values and any principal submatrix is one too. A fleet asking for more correlations than it has observed cells is refused, and fewer than three cells per parameter warns.

"2dar1" correlates over ages and years jointly, with the covariance the Kronecker product of an AR(1) over each. It is defined on a complete grid, so a fleet’s observed ages and years must form one; a gap is refused at setup rather than filled in. Use "1dar1" or "us" for a ragged data source. The year correlation is trans_rho_<data source>_year.

How the correlations are shared is a spec string, the same one the rest of the package uses for sigmaF_spec and Fdev_rho_spec, rather than a structure of its own: rho_catch_spec, rho_discard_spec and rho_srv_idx_spec, each with a _pop counterpart.

The correlations sit over region, sex and fleet, with a leading population dim for the population-specific data sources, so the abbreviations are r, s, f and p:

rho_*_spec	Meaning
NULL	"est_shared_r_s", or "est_shared_p_r_s" for a _pop data source. Default
"est_shared_r_s"	One per fleet
"est_shared_s"	One per region and fleet
"est_shared_r"	One per sex and fleet
"est_shared_r_s_f"	A single correlation
"est_all"	One per region, sex and fleet
"fix"	Kept at their starting values

One spec governs the data source’s across-age correlation, its across-year correlation and its unstructured matrix together, so two fleets that share under "us" share a whole matrix rather than one entry of it. A region, sex or population a fleet never observes has no correlation parameter whatever the spec says, which is what holds the unused slots of a summed dim out: a data source summed over regions gets one correlation, in region one, however free the spec.

A cell with a single observed age has no correlation to describe and falls back to independent. Be careful estimating an AR(1) on a fleet observing only two or three ages: the correlation is barely identified and will sit on its bound.

A correlated cell contributes its whole density to the first age present, leaving the remaining ages at zero, so the reported likelihood array still sums to the cell’s contribution but is no longer readable age by age.

Units

catch_units and discard_units apply to the at-age data sources as they do to the aggregated ones: "abd" reads the numbers at age, "biom" multiplies them by the fleet’s weight at age.

For a lognormal data source with a weight at age supplied as data, the two are the same statement. The residual is $\log O_{a} - \log(\hat{C}_{a} W_{a}) = \log(O_{a}/W_{a}) - \log \hat{C}_{a}$, and $W_{a}$ is a constant, so it cancels and the likelihood is identical. Catch at age arrives in numbers, so "abd" is the natural choice.

"biom" becomes a different statement in two places. Under *AA_LikeType = "normal" the residual is $O_{a} - \hat{C}_{a} W_{a}$ on the arithmetic scale, so the weight rescales it against a fixed $\sigma_{a}$. And where weight at age is derived from estimated growth rather than supplied, $W_{a}$ has parameters, so a biomass observation at age informs growth and a numbers one does not.

Discards at age

Discards at age are the total discarded, not the dead fraction the population dynamics track, which is the same statement the aggregated discard data source makes: the prediction is $D_{a}/\text{dmr}$. They are reported in numbers or weight, so a fleet fitting discards at age needs discard_units of "abd" or "biom"; the fraction units are a property of the catch as a whole rather than of one age and are refused.

Identifiability

An estimated observation error needs the thing predicting the mean to be more constrained than the data. Fishing mortality with free annual deviations can reproduce every catch at age exactly, at which point the residuals vanish and the estimated standard deviation collapses towards zero, which the likelihood rewards without limit. Effort-driven fishing mortality, or a blocked or penalized structure, is what makes an estimated catch standard deviation meaningful.

Fishing mortality must be constrained in the same way for the at-age catch likelihood itself to be identifiable. Four ages by forty years by two seasons of free fishing mortality is roughly the number of observations, and even under a separable structure a ridge opens in which each age’s fishing mortality absorbs its own abundance.

Estimated index observation error (sigmaFishIdx_spec / sigmaSrvIdx_spec)

An index has a standard error from its own survey design, and an assessment may additionally estimate a component covering everything that design does not. this is often called an extra standard deviation. By default SPoRC uses the reported standard errors as they stand, which is "fix".

One parameter is estimated per fleet, ln_sigmaFishIdx or ln_sigmaSrvIdx. The aggregated and population-specific data sources are separate observations with separate error structures, so they have separate parameters, in the same way ln_sigmaC and ln_sigmaC_pop do. The population data sources are configured through sigmaFishIdx_pop_spec and sigmaSrvIdx_pop_spec.

String	Total standard deviation
"fix"	The reported standard error, unchanged. Default
"est_additive"	Reported standard error plus the estimated component. The additive convention, and what the bridges here reproduce
"est_quadrature"	Reported standard error and estimated component added in quadrature, treating them as independent variances
"est_replace"	The estimated standard deviation alone, discarding the reported errors. What several ICES assessments do

sigmaFishIdx_map and sigmaSrvIdx_map, with _pop variants, give the estimation groups as an integer vector of length n_fish_fleets or n_srv_fleets. Fleets sharing a value share a parameter and NA holds a fleet at its starting value, which is how a fleet a reference assessment pinned at a bound is kept without pretending it was estimated.

An estimated component is confounded with a likelihood weight, because a weight on a normal likelihood is the same statement as dividing the variance by that weight, so the estimate will absorb the weight; Setup_Mod_Weighting() warns when both are used. A fleet with a multivariate normal index likelihood takes its scale from the supplied covariance and ignores a standard deviation entirely, so requesting one there is an error rather than a silently unidentified parameter.

The simulator draws index observations from whatever total the fitted model used, so simulation_self_test() exercises the estimated component without further configuration.

Index error structure (FishIdx_LikeType / SrvIdx_LikeType)

Per fleet. A fleet’s population-specific index data source follows the same choice for "lognormal" and "normal", but stays lognormal under "mvn", whose covariance describes the region-aggregated series only. The simulator (Setup_Sim_Fishing / Setup_Sim_Survey, and simulation_self_test automatically) draws index observations under the same structures, with an "mvn" fleet drawn through a one-factor decomposition of its covariance:

String	Description
"lognormal"	Default. SEs are log-scale standard deviations
"normal"	Normal on the arithmetic scale; SEs are arithmetic standard deviations
"mvn"	Multivariate normal on the arithmetic scale with a fixed covariance supplied via FishIdx_Cov / SrvIdx_Cov (one matrix per "mvn" fleet, square with one row per fitted observation, ordered as observations appear scanning the fleet’s use flags in array order). Validated at setup for symmetry and positive definiteness

Recommendations. "lognormal" is recommended for almost all abundance indices (positive, multiplicative errors). Use "mvn" when the index comes with a covariance across years (e.g., a model-based index from VAST or sdmTMB whose inter-annual correlations are real information a diagonal likelihood would double-count). Use "normal" only for series that can legitimately go near zero or negative (e.g., recruitment indices), or when bridging an assessment that fits on the arithmetic scale.

Index timing and age restriction

Argument	Description
t_fish	[r × τ × f], fraction of the season elapsed when the fishery index is observed; numbers decay by $e^{-t Z}$ before the index is formed, mirroring t_srv. Default 0 (start of season, the historical behavior). Set 0.5 for a mid-season CPUE snapshot
fish_idx_ages / srv_idx_ages	Which ages contribute to each fleet’s index total: a list (one vector of ages per fleet, NULL = all) or an [n_ages × n_fleets] 0/1 array. The restriction applies to the index sum only; selectivity, catch, and the fleet’s compositions are untouched. Restricting a survey to one age turns it into an index of that age alone (e.g., an age-1 acoustic recruitment index), without needing a knife-edge selectivity that would corrupt the compositions

Selecting fish at length rather than at age

Fish of one age are spread over a range of lengths, and a length-based gear takes the long ones more often than the short ones. Selecting at age first averages that away: it gives the age one selectivity value and then treats every fish of the age as equally catchable. The options below decide whether that averaging happens, once for the compositions and once for the catch weight.

All of them require length-based selectivity, and the two weight options also require waa_model = "wt_len" so the model knows what a fish of each length weighs. srv_waa_selected only applies to an index in weight (srv_idx_type = "biom"), since an index in numbers has no weight at age in it.

If selectivity is flat across the lengths an age covers, every option below is the same as the default. The question to ask is how steep the selectivity curve is over one age’s length range.

Argument	Description
FishLenComps_sel / SrvLenComps_sel	"age" (default) or "length", per fleet. "age" selects the catch or index at age and spreads it over lengths afterwards, so the length composition within an age is the size-age key’s own. "length" spreads the fish at each age over the key first and selects them length by length, $C_{l} = s(l)\sum_{a}P(l \mid a)\,N_{a}\!\left( 1 - e^{-Z_{a}} \right)F/Z_{a}$, so the long fish of an age are taken more often. The key is the fleet’s own, at t_fish or t_srv. Use "length" when selectivity is length based and the length compositions are what inform it. The two are different expected compositions, not two roundings of the same one: in an EBS Pacific cod model the fishery length likelihood is 119.36 under one and 198.55 under the other
fish_waa_selected / srv_waa_selected	0/1 per fishery or survey fleet. 1 makes the fleet’s catch biomass, or a survey’s index in weight, use the mean weight of the fish it takes at each age, $\sum_{l}P(l \mid a)s(l)w(l)/\sum_{l}P(l \mid a)s(l)$, instead of the population mean weight at that age. Use 1 when the gear selects strongly within an age, so that a caught fish is larger or smaller than an average fish of its age. In an EBS Pacific cod model the caught two-year-olds weigh 40 percent more than the average two-year-old, which moves the catch by 0.44 percent. With flat or age-based selectivity the two are the same

Composition bin restriction (*_bins)

Which observed bins a fleet’s compositions are fitted over. Observed and expected compositions are both subset to the named bins and renormalized within them, so bins outside are left out of the likelihood rather than forced to be explained. Use for gears that never resolve part of the range (a fishery that never catches the youngest ages, whose structural zeros would otherwise have information); leave alone when the zeros are informative sampling zeros.

Every composition data source takes one, and they all share the same format: a list with one element per fleet, each a vector of bin indices or NULL for all bins, or an [n_obs_bins x n_fleets] array of 0/1 weights. Indices are into observed bins, that is after any AgeingError or LenBinMap. Default NULL fits every bin.

Data source	Fishery	Survey
Marginal age comps	FishAgeComps_bins	SrvAgeComps_bins
Marginal length comps	FishLenComps_bins	SrvLenComps_bins
Conditional age-at-length	Fish_caal_bins	Srv_caal_bins
Population-specific age comps	FishAgeComps_pop_bins	SrvAgeComps_pop_bins
Population-specific length comps	FishLenComps_pop_bins	SrvLenComps_pop_bins
Discard age comps	FishAgeComps_discard_bins	–
Discard length comps	FishLenComps_discard_bins	–
Population-specific discard age comps	FishAgeComps_discard_pop_bins	–
Population-specific discard length comps	FishLenComps_discard_pop_bins	–

# the acoustic survey resolves ages 2-15 only; every other fleet fits all ages
SrvAgeComps_bins = list(NULL, 2:15, NULL, NULL)

Every fleet must keep at least two bins. A composition fitted over one bin says nothing, since the proportion in a lone bin is one whatever the model predicts, and the one-step-ahead routines degenerates around it: the logistic-normal families spend a bin as the additive log-ratio reference and have none left, and the discrete families mark their only bin as the multinomial’s determined cell, leaving get_osa() nothing to keep. Setup refuses it and names the fleet.

The restriction applies whatever the composition type. For sex-joint comps (spltRjntS) the named bins are dropped from each sex’s block of the [bin x sex] stack, so the sex ratio the joint comps have becomes the ratio within the fitted bins, which is usually what you want when the excluded bins are unsampled rather than absent. Logistic-normal covariances are built across all observed bins and then cut down to the fitted ones, so a gap in the named bins still counts towards the AR1 lag between the bins either side of it. One-step-ahead residuals are packed and evaluated over the same restricted bins, so get_osa() and the likelihood always agree on what was fitted.

Restricting bins versus remapping them

*_bins and the bin maps (AgeingError, LenBinMap) are easy to confuse and do different things.

	What it does	When you need it
AgeingError / LenBinMap	Remaps. An [n_model_bins x n_obs_bins] matrix the expected composition is multiplied through, so it lands on the bins the data were recorded on. Mass is redistributed and, wherever a row sums to one, preserved	The data are on different bins from the model: 1 cm model bins fit to 5 cm length compositions, or true ages smeared across observed ages by reading error
*_bins	Drops. Names a subset of the observed bins; the rest are excluded from the likelihood and the remainder renormalized. Bin definitions do not change, mass in the dropped bins is discarded	The data are on the model’s bins but part of the range is not informative: a gear that never ages its youngest fish, or length bins nothing is ever sampled from

They compose, in that order: remap model bins onto observed bins, then restrict to the observed bins being fitted. Neither substitutes for the other. No subset of fifty 1 cm bins is ever ten 5 cm bins, so a restriction cannot do a map’s job; and a map that zeroed a bin out would leave an observed bin with a structural-zero expectation the composition likelihood cannot fit, so a map cannot do a restriction’s job. check_bin_map() enforces both halves of that: map rows must sum to one (or to zero, to drop a model bin from the observations entirely), and no observed bin may be left with nothing mapped into it.

Worked on six model bins with expected proportions 0.05, 0.10, 0.25, 0.30, 0.20, 0.10:

	Result	Bins	What happened to the first two bins’ mass
LenBinMap collapsing pairs	0.15, 0.55, 0.30	6 → 3	Summed into the first observed bin. Still there
*_bins = 3:6	0.29, 0.35, 0.24, 0.12	6 → 4	Discarded. The remaining four renormalize to one

Both shrink the vector, which is why they get confused, but only the map holds the first two bins’ fish somewhere rather than throwing them away. If bins 1 and 2 are 10 cm fish your survey does catch but records in a coarser bin, you want the first row. If they are 10 cm fish your survey never sees, so their observed proportion is a sampling artifact rather than a measurement, you want the second. Choosing the map when you meant the restriction quietly asks the model to explain a bin that has no information; choosing the restriction when you meant the map throws away fish that were really there.

A restriction leaves a visible trace at setup, so you can confirm it took:

SrvAgeComps_bins for fleet 2 is fitted over observed bins: 2, 3, 4, ..., 15

Composition likelihood families

Age and length compositions each accept one of eight likelihood families:

Value	Name	Likelihood
0	"Multinomial"	Multinomial
1	"Dirichlet-Multinomial"	Dirichlet-multinomial (overdispersion parameter $\theta$ estimated per fleet)
2	"iid-Logistic-Normal"	Logistic-normal, independent bins
3	"1d-Logistic-Normal"	Logistic-normal with AR(1) correlation across bins
4	"2d-Logistic-Normal"	Logistic-normal with AR(1) bin correlation and constant cross-sex correlation
5	"iid-Logistic-Normal-miss0"	Logistic-normal on the bins that were seen, independent
6	"1d-Logistic-Normal-miss0"	The same, correlated across bins
7	"2d-Logistic-Normal-miss0"	The same, correlated across bins and sexes

These are set per data source via comp_fishage_like, comp_fishlen_like, comp_srvage_like, comp_srvlen_like, and their _pop and _discard variants.

Two ways to handle an empty bin

A logistic normal takes a log ratio, so a bin with nothing in it has to be dealt with before the likelihood can be written. Forms 2 to 4 add addtocomp to every bin and keep the composition at its full length. Forms 5 to 7 drop the empty bins instead and renormalize the expected proportions over the bins that remain, so the density is the one the part of the composition that was actually seen has. This is the treatment WHAM calls logistic-normal-miss0, and it is the one to use when bridging an assessment written that way.

Three things follow from it:

• The standard deviation is divided by the square root of the input sample size, $\sigma^2_{y} = \exp(2\theta) / N_{y}$, so $\theta$ is a per fish quantity and a year sampled harder is fit more tightly. The other logistic normal forms read $\theta$ as the standard deviation outright and ignore the sample size.
• The change of variables from the log ratio is subtracted, so the result is a density on the composition rather than on its transform. That does not move the maximum, since it does not involve any parameter, but it does change the value, so do not compare an AIC across the two treatments.
• Forms 6 and 7 put their correlations through the logistic function, so they are positive. Forms 3 and 4 allow either sign, and also scale the marginal variance by $1/(1-\rho^2)$, which forms 6 and 7 do not. The lag is measured in bins, so a composition fit over a restricted or gapped bin range spaces the autoregression by the bins themselves rather than by position in the vector.
• Form 7 correlates bins with sexes through the same separable structure form 4 uses, so it needs a composition joint across sexes. Setup refuses it with "agg" or "spltRspltS", which have no bin by sex stack to correlate over, and refuses form 4 there too rather than silently returning nothing.

A cell with fewer than two bins seen has no log ratio to take and contributes nothing. One-step-ahead residuals are not available for these two forms, because the number of observations in a cell changes with the number of empty bins while the residual packing needs a fixed length; setup refuses the combination rather than computing something else.

Composition structure types

Each composition data source says how its observations are aggregated through a _Type, given as a year and fleet specification:

Value	Structure
"agg"	Aggregated across sexes and regions
"spltRspltS"	Split by sex and region, with no implicit sex-ratio information
"spltRjntS"	Joint across sexes, split by region, preserving sex-ratio information
"none"	No data for this fleet and year

A fishery that began sexing its samples in the eighth year of the series reads:

FishAgeComps_Type = c("agg_Year_1-7_Fleet_1",
                      "spltRjntS_Year_8-terminal_Fleet_1")

"agg" is refused where the likelihood is the 2d logistic normal, which needs a composition split by region and sex to correlate over.

Data sources

SPoRC supports a full matrix of composition data sources, each independently configurable:

Category	Aggregate	Population-specific
Fishery age	comp_fishage_like	comp_fishage_pop_like
Fishery length	comp_fishlen_like	comp_fishlen_pop_like
Fishery discard age	comp_fishage_discard_like	comp_fishage_discard_pop_like
Fishery discard length	comp_fishlen_discard_like	comp_fishlen_discard_pop_like
Survey age	comp_srvage_like	comp_srvage_pop_like
Survey length	comp_srvlen_like	comp_srvlen_pop_like

Each data source has its own ISS arrays, $\theta$ parameters (log-scale overdispersion), and correlation parameters.

Conditional age-at-length

Supplied on Setup_Mod_FishIdx_and_Comps() (ObsFish_caal, UseFish_caal, ISS_Fish_caal, Fish_caal_LikeType, Fish_caal_Type) and Setup_Mod_SrvIdx_and_Comps() (the Srv_ twins), with the weights Wt_Fish_caal and Wt_Srv_caal on Setup_Mod_Weighting(). A CAAL observation is the age composition of the fish aged from one length bin, so the arrays have a length dimension ahead of the age one ([r × y × τ × l × a × s × f]), the use flag and input sample size are per length bin ([r × y × τ × l × f], [r × y × τ × l × s × f]), and the sample size is the number aged from the bin, not the number measured. The expected row is the joint catch at length and age normalized across ages, which do_caal = 1 builds and the data switch on by themselves.

Argument	Options
Fish_caal_LikeType, Srv_caal_LikeType	"none", "Multinomial", "Dirichlet-Multinomial". The logistic-normal families are not offered: a bin’s age sample is small and mostly zeros
Fish_caal_Type, Srv_caal_Type	The same "CompType_Year_x-y_Fleet_z" strings as the marginal compositions; spltRspltS when each row holds one sex’s otoliths
ln_Fish_caal_theta, ln_Srv_caal_theta	Dirichlet-multinomial overdispersion, one per region and sex shared across length bins

Fit CAAL with the marginal length compositions and without the marginal age compositions of the same fleet, which would count the same fish twice. Use comp_const_obs = 0 with the Dirichlet-multinomial, since every structurally empty bin otherwise contributes to $ heta$. One-step-ahead residuals come from get_osa(comp_source = "Fish_caal") or "Srv_caal", with a len column for the bin each row was conditioned on.

---

Likelihood Weights

Set on Setup_Mod_Weighting(). Every weight multiplies its likelihood component in the joint negative log likelihood, so 0 drops a data source from the fit without removing it from the model or from the reported fits.

Two shapes appear. The scalar weights apply one number to the whole component. The composition weights are arrays [n_regions × n_years × n_seas × n_sexes × n_fleets], so a year, a region or a sex can be down-weighted on its own; a scalar is recycled.

Scalar weights

Argument	Component
Wt_Catch	Aggregated catch
Wt_Discard	Aggregated discards
Wt_FishIdx	Fishery indices
Wt_SrvIdx	Survey indices
Wt_F	Fishing mortality deviation penalty
Wt_D	Discard mortality rate deviation penalty
Wt_Rec	Recruitment deviation penalty (see Recruitment; may also be an array)
Wt_Init_Rec	Initial age deviation penalty
Wt_Tagging	Tag recapture likelihood

Composition weights

Argument	Component
Wt_FishAgeComps, Wt_FishLenComps	Fishery age and length compositions
Wt_SrvAgeComps, Wt_SrvLenComps	Survey age and length compositions
Wt_FishAgeComps_discard, Wt_FishLenComps_discard	Discard compositions
Wt_Fish_caal, Wt_Srv_caal	Conditional age-at-length

Population-specific weights

Every data source fit separately by population has a _pop twin, weighted independently of its aggregated counterpart: Wt_Catch_pop, Wt_Discard_pop, Wt_FishIdx_pop, Wt_SrvIdx_pop, and the four composition arrays Wt_FishAgeComps_pop, Wt_FishLenComps_pop, Wt_SrvAgeComps_pop, Wt_SrvLenComps_pop, plus their discard twins Wt_FishAgeComps_discard_pop and Wt_FishLenComps_discard_pop. These have no effect in a single-population model.

Recommendations. Use weights to switch a data source off (0) or to reproduce another assessment’s weighting, not to tune fit. Where compositions are over-weighted relative to their real information content, the Dirichlet-multinomial or a Francis iteration estimates that down-weighting rather than asserting it; run_francis() returns the multipliers to apply here. The at-age data sources take the same weight as their aggregated counterpart, so a fleet moved from aggregated to at-age keeps whatever weight it already had.

---

Tagging

Configured via Setup_Mod_Tagging() (EM) or Setup_Sim_Tagging() (OM).

SPoRC implements a conventional mark-recapture framework (Brownie-type likelihood) where tagged individuals follow the full seasonal population dynamics, movement, mortality, and fishing after release. Tags are tracked as release cohorts and recaptured by fishery fleets.

Switching tagging on

use_conv_fish_tagging is an integer vector of length n_fish_fleets, 1 for fleets that recapture and report tags. All zeros (the default) disables the whole module: the parameters collapse to length-1 placeholders and the likelihood contributes nothing.

Defining release cohorts

Argument	Description
conv_tag_release_indicator	Data frame of release cohorts, one row each. Required columns regions, tag_years, tag_seas. Every other tagging array is dimensioned by its row count
conv_tag_release_platform	Character matrix [n_cohorts × 2], columns platform and fleet, aligned row for row with the release indicator. "population" releases into the population directly (fleet is NA); "fishery" or "survey" releases through that fleet’s selectivity
conv_tag_max_liberty	Years at liberty tracked per cohort. Sets the recapture container size, so it is a direct cost in tape size. Default n_ages / 2

A cohort is a release event, not a release group. Splitting the same tags across more rows multiplies the recapture containers and the tape with them.

Observations

Argument	Description
conv_tagged_fish	Numbers released per cohort
obs_conv_tag_fish_recap	Observed recaptures, by cohort, year at liberty, and recapture cell

What a tag has (conv_fish_tag_attr)

Which biological dimensions are resolved at release and therefore retained at recapture. Region and fleet are always retained. Built from "p" (population), "a" (age) and "s" (sex) joined by underscores: "p_a_s", "a_s", "p_a", "p_s", "a", "s", "p", "none". Default "p_a_s".

A scalar applies to every event; a vector of length n_cohorts lets events differ, for a program where early releases recorded less than later ones.

The pooling arguments say how cells are aggregated in the likelihood, and must agree with what the attribution resolves:

Argument	Default	Use when
conv_tag_pop_pool	as.list(1:n_pop)	list(1:n_pop) when "p" is not attended
conv_tag_age_pool	as.list(1:n_ages)	list(1:n_ages) when "a" is not attended; custom groups such as list(1:5, 6:10) are allowed when it is
conv_tag_sex_pool	as.list(1:n_sexes)	list(1:n_sexes) when "s" is not attended

A pooling structure inconsistent with the attribution warns and is overridden to a single group, so the model runs but not with the structure that was asked for. Set both together.

Timing

Argument	Description
conv_tag_mixing_period	Years (or seasons in a seasonal model) after release before recaptures enter the likelihood. Lets tags mix before they inform movement. Default 1
conv_tag_t_tagging	Fraction of the season remaining at release, in [0, 1]: 1 at the start, 0.5 mid-season, 0 at the end. Scalar or one value per cohort. Default 1

Recapture likelihood (conv_fish_tag_like)

String	Description
"Poisson"	Poisson on recapture counts (default)
"NegBin"	Negative binomial, overdispersion from ln_conv_fish_tag_theta
"Multinomial_Release"	Multinomial conditioned on the release cohort
"Multinomial_Recapture"	Multinomial conditioned on the recapture cell
"Dirichlet-Multinomial_Release"	Dirichlet-multinomial on the release conditioning
"Dirichlet-Multinomial_Recapture"	Dirichlet-multinomial on the recapture conditioning

ln_conv_fish_tag_theta is read by the negative binomial and both Dirichlet-multinomial forms and ignored otherwise.

Tag-loss parameters

Both are log scale, scalar or one value per release event, and a scalar is recycled to every event. Both default to -1000, which is effectively zero.

Parameter	Spec argument	Sharing options
ln_init_conv_tag_mort	init_conv_tag_mort_spec	"fix", "est_shared"
ln_conv_tag_shed	conv_tag_shed_spec	"fix", "est_shared"

Leaving a spec at NULL keeps the parameter fixed at its starting value.

Reporting rates

Argument	Description
conv_tagrep_spec	Sharing across region, fleet and block: "est_all", "est_shared_r", "est_shared_f", "est_shared_r_f", "fix". Default "fix"
conv_tag_fish_reporting_blocks	Time blocking for reporting, same syntax as the selectivity and catchability blocks
use_conv_tag_fishrep_prior	0/1, whether reporting priors apply. Default 0
conv_tag_fishrep_prior	Data frame with columns region, block, fleet, mu, sd, type. Ignored unless the switch above is 1

"est_shared_r" requires identical block structure across regions within a fleet, and "est_shared_f" the same across fleets within a region.

Recommendations. Reporting rate and tag loss trade off against each other and against fishing mortality, so estimating all three at once rarely identifies. Fix reporting from a tag-seeding study where one exists, or put a prior on it, and estimate shedding only with a double-tagging component. Keep conv_tag_max_liberty no longer than the recoveries actually inform: it is the argument that most directly sets how long the model takes to build.

Weight

Wt_Tagging on Setup_Mod_Weighting().

---

Reference Points

Computed post-estimation via Get_Reference_Points().

The function accepts a type argument for spatial structure and a what argument for the reference-point method. All methods project a single recruit through the full age, season, and spatial structure to compute SBPR and YPR.

Available methods

type	what	Description	Multi-pop?
"single_region"	"SPR"	$F_{\text{SPR}_x}$, no movement	,
"single_region"	"BH_MSY"	Beverton-Holt $F_\text{MSY}$, no movement	,
"multi_region"	"independent_SPR"	Per-region $F_{\text{SPR}_x}$ ignoring movement	✓
"multi_region"	"independent_BH_MSY"	Per-region $F_\text{MSY}$ ignoring movement	✓
"multi_region"	"global_SPR"	Single $F_{\text{SPR}_x}$ applied uniformly across regions, with movement	✓
"multi_region"	"global_BH_MSY"	Single $F_\text{MSY}$ with movement (single-pop only)	,
"multi_region"	"local_BH_MSY"	Region-specific $F_\text{MSY}$ values jointly maximizing total yield under movement. Uses Newton-Raphson to solve equilibrium recruitment by origin	✓

Controls

Argument	Description
SPR_x	Target SPR fraction (e.g., 0.4 for $B_{40\%}$)
n_avg_yrs	Terminal years to average demographic rates (selectivity, $M$, WAA, maturity, movement)
calc_rec_st_yr	First year for computing mean historical recruitment
rec_age	Recruitment lag used to exclude most-recent years from the mean
is_discard_fleet	Integer vector (per fleet) flagging discard-only fleets to exclude from landed yield in MSY calculations
local_bh_msy_newton_steps	Newton iterations for the local $F_\text{MSY}$ solve (default 6)

---

Simulation and Closed-Loop MSE

The operating-model side uses a parallel set of Setup_Sim_* functions that mirror the estimation-model pipeline but populate a sim_list rather than input_list. Key simulation-specific functions:

Function	Purpose
Simulate_Pop_Static()	Forward-project the OM without feedback (open-loop)
condition_closed_loop_simulations()	Closed-loop MSE: periodically re-fits the EM, applies an HCR, and updates $F$
simulation_data_to_SPoRC()	Converts OM output (with observation error) to input_list format for the EM
simulation_self_test()	Simulation-estimation test: fits the EM back to OM-generated data
get_closed_loop_reference_points()	Computes reference points inside the MSE feedback loop

Key closed-loop arguments in condition_closed_loop_simulations():

Argument	Description
closed_loop_yrs	Year range over which feedback is active
assessment_period	Frequency of EM re-fitting (e.g., every 2 years)
use_true_values	If TRUE, the HCR uses OM-truth instead of EM estimates (perfect-information benchmark)

---

Estimation and Optimization

fit_model() constructs the RTMB automatic-differentiation function, optimizes via nlminb, and refines with Newton steps.

Argument	Default	Description
data	,	input_list$data
parameters	,	input_list$par
mapping	,	input_list$map
random	NULL	Parameter names to marginalize as random effects (Laplace approximation)
newton_loops	3	Post-convergence Newton steps ($\Delta\theta = -H^{-1}g$) to reduce residual gradients
do_optim	TRUE	FALSE returns the un-optimized MakeADFun object for debugging
nlminb_control	list(iter.max = 1e5, eval.max = 1e5, rel.tol = 1e-15)	Passed to stats::nlminb

---

Diagnostics

Function	What it does
do_retrospective()	Sequentially peels terminal years and re-fits; reports Mohn’s $\rho$ per quantity
do_jitter()	Refits from perturbed starting values to test convergence stability
do_likelihood_profile()	Profiles the likelihood surface over user-specified parameters
do_francis_reweighting()	Computes Francis TA1.8 weights for composition data
run_francis()	Iterative Francis reweighting loop (re-fits after each adjustment)
get_osa()	One-step-ahead residuals, compositions, conventional tagging, and Catch/Discard/FishIdx/SrvIdx indices (external post-hoc, or internal model-based via do_internal_comp_osa/do_internal_conv_tag_osa). See vignette("u_osa_residuals")
do_runs_test()	Runs test for serial correlation in residuals
get_model_rep_from_mcmc()	Extracts model report quantities across MCMC posterior draws (compatible with adnuts / tmbstan)
marg_AIC()	Marginal AIC for models with random effects

---

Plotting

Function	Output
plot_all_basic()	Multi-panel diagnostic overview
get_ts_plot()	Time series: SSB, recruitment, $F$, catch, depletion
get_idx_fits_plot()	Index fits (observed vs. predicted)
get_catch_fits_plot()	Catch fits
get_comp_prop()	Composition fits (bubble plots / proportion plots)
get_selex_plot()	Selectivity-at-age or selectivity-at-length
get_biological_plot()	WAA, maturity, natural mortality
get_nLL_plot()	Likelihood component breakdown
get_retrospective_plot()	Retrospective trajectories
get_data_fitted_plot()	Data-availability timeline
plot_resids()	Residual diagnostics
get_retrospective_relative_difference()	Relative difference plots for retrospective analyses

All plotting functions accept either a single fitted model or a list of models for side-by-side comparison.