Package 'Rceattle'

Description

Enforces parameter sharing by mapping parameters for fleets with a common Selectivity_index or Catchability_index to the same value as the initial index.

Usage

adjust_map_shared_params(map_list, data_list)

Arguments

Value

Updated map_list.

Description

Returns the model's derived population quantities in long form so that custom plots and post-processing don't have to walk the deeply nested quantities list or inherit the dimnames decisions in rename_output(). Two shapes are flattened into one tidy frame: species-by-year quantities (e.g. biomass, ssb, R, F_spp) and species-by-sex-by-age-by-year quantities (e.g. N_at_age, biomass_at_age, M_at_age). For the species-level shapes, sex and age are returned as NA. Cells of the 4D arrays that are padded out to max(nsex) / max(nages) for species with fewer sexes or ages are dropped.

Usage

## S3 method for class 'Rceattle'
as.data.frame(
  x,
  row.names = NULL,
  optional = FALSE,
  which = c("biomass", "ssb", "R", "biomass_depletion", "ssb_depletion", "F_spp"),
  ci_level = 0.95,
  ...
)

Arguments

Details

Standard errors (se) and confidence intervals (lwr, upr) are populated from the TMB sdreport for any quantity that was ADREPORT'd (biomass, ssb, R, exploitable_biomass, biomass_depletion, ssb_depletion); other quantities and fits produced with getsd = FALSE get NA for se / lwr / upr. Set ci_level to widen or narrow the band.

A strictly positive series takes its interval on the log scale, so it is right-skewed and cannot reach zero. This is the same construction plot_timeseries() draws, so the table and the figure agree.

Value

A data.frame with columns year, species, sex, age, quantity, value, se, lwr, upr. species is the character species name from data_list$spnames. Rows are sorted in the order which was given.

Description

A data list containing inputs for a single-species CEATTLE model fit to Aleutian Islands Atka mackerel data through 2022. See BS2017SS for format details.

Usage

Atka2022

Format

A list with the same structure as BS2017SS.

Description

A data list containing inputs for the three-species (walleye pollock, Pacific cod, arrowtooth flounder) multispecies CEATTLE model for the Eastern Bering Sea. See BS2017SS for format details.

Usage

BS2017MS

Format

A list with the same structure as BS2017SS.

Description

A data list of inputs for the single-species Eastern Bering Sea CEATTLE model (walleye pollock, Pacific cod, and arrowtooth flounder), covering 1979-2017.

Usage

BS2017SS

Format

Control

nspp

Number of species included in CEATTLE

styr

Start year of the hindcast

endyr

End year of the hindcast

projyr

End year of the forecast

nsex

Number of sexes to model in the population (1 = sexes combined, 2 = females and males modeled separately)

spawn_month

Spawning month of the population to adjust the numbers spawning

R_sexr

Percent of recruitment that is female (ignored if nsex = 1)

nages

Number of ages of each species included in the hindcast

minage

Minimum age for each population (i.e. does recruitment correspond to age 0, 1, 2?)

nlengths

Number of lengths of each species included in the hindcast

pop_wt_index

Weight-at-age (weight) index to use for calculation of each species population derived quantities (SSB, Consumption/Ration, Suitability, etc)

ssb_wt_index

Weight-at-age index used to compute each species' female spawning-stock biomass (SSB)

pop_age_transition_index

Age transition matrix (e.g. growth trajectory) index to use for derived quantities of the population to convert age to length (also used in length-based predation estimation)

sigma_rec

Standard deviation of the log-scale recruitment deviations (fixed, or an initial value if estimated)

other_food

Other food in the ecosystem for each species, in metric tons: it is summed into the available-food denominator beside numbers-at-age (thousands) times weight-at-age (kg)

estDynamics

Estimate or fix numbers-at-age. Accepts integer codes or the equivalent readable strings: 0 / "Estimated" = estimate dynamics; 1 / "Fixed" = use input numbers-at-age in NByageFixed; 2 / "FixedScaled" = multiply input numbers-at-age (NByageFixed) by a single scaling coefficient, estimated under predation (msmMode > 0) and fixed at 1 otherwise.

M1_model

Estimate residual (multi-species mode) or total natural mortality (single-species mode). 0 = use fixed natural mortality from M1_base, 1 = estimate sex- and age-invariant M1, 2 = sex-specific (two-sex model), age-invariant M1, 3 = estimate sex- and age-specific M1.

fleet_control

Survey and fishery data specifications

index_data

Survey index in weight (mt) or numbers (thousands of fish) data. The Observation column is on the NATURAL scale in the units fleet_control$Observation_units names, never a log: the likelihood takes the log itself, and a fitted observation must be positive. A negative Year predicts the row without fitting it.

index_cov

Optional named list of survey-index variance-covariance matrices, keyed by Fleet_name (or Fleet_code), used only by fleets whose fleet_control Index_distribution == "MVN". Each matrix must be square and symmetric with dimension equal to the number of fitted survey observations for that fleet (Year in styr:endyr, Observation > 0), ordered as in index_data. Inverted once internally to the precision matrix used in 0.5 * r' Sigma^-1 r. Leave unset for the default lognormal likelihood.

catch_data

Total catch in weight (mt) or numbers (thousands of fish) data

comp_data

Survey/fishery age or length composition data. Note if sex is 3, put female composition data then male composition data (similar to SS).

emp_sel

Empirical/fixed selectivity for surveys and fisheries (leave empty if not used)

age_trans_matrix

Age transition matrix (e.g. growth trajectory) used to convert age to length for length comp data. Can have multiple matrices for a species specified by Age_transition_index.

age_error

Ageing error matrices, one row per true age. One per species by default; give several and select them per fleet with the 'Ageing_error_index' column of both this table and fleet_control.

weight

Weight-at-age (weight) to use for calculation of derived quantities (SSB, Consumption/Ration, Suitability, Total Catch, Survey Biomass, etc). Can have multiple weight-at-age data-sets for each species, but must be full for all years of the hindcast.

maturity

Maturity-at-age for each species

sex_ratio

Percent female at age for each species

M1_base

Residual natural mortality for each species

aLW

Parameters for weight-at-length power function for each species. Used when estimating time-variant length-based gamma suitability (suitMode = 2) or time-variant length-based lognormal suitability (suitMode = 5)

Ceq

Which bioenergetics equation to use for each species for ft to scale max consumption: 1 = Exponential (Stewart et al 1983), 2 = Temperature-dependence for warm-water species (Kitchell et al 1977; sensu Holsman et al 2015), 3 = temperature dependence for cool and cold-water species (Thornton and Lessem 1979)

Cindex

Which environmental index in env_data to use to drive bioenergetics

Pvalue

This scales the maximum consumption used for ration for each species; Pvalue is in CmaxfTPvalue*ration_data

fday

Number of foraging days per year for each species

CA

Intercept of allometric mass function for calculating maximum consumption: CA * Weight ^ CB

CB

Slope of allometric mass function for calculating maximum consumption: CA * Weight ^ CB

Qc

Parameter for temperature scaling function of maximum consumption specified by Ceq

Tco

Parameter for temperature scaling function of maximum consumption specified by Ceq

Tcm

Parameter for temperature scaling function of maximum consumption specified by Ceq

Tcl

Parameter for temperature scaling function of maximum consumption specified by Ceq

CK1

Parameter for temperature scaling function of maximum consumption specified by Ceq

CK4

Parameter for temperature scaling function of maximum consumption specified by Ceq

env_data

Environmental indices such as bottom temperature data to incorporate into ration equation specified by Ceq and Cindex. Also used to drive catchability if Catchability = 5. Will use the mean for missing years. Temperature should be in Celsius.

ration_data

Annual relative foraging rate by age or input consumption at age. Multiplied by pvalue and fday to scale maximum consumption to the number of days in a year that foraging occurs.

UobsWtAge

Stomach proportion by weight for each predator, prey, predator age, prey age combination. Can also be year specific by including the column, "Year"

fleet_control: controls for survey and fisheries data

Fleet_name

Name of survey or fishery

Fleet_code

Index of survey/fishery ACROSS species

Fleet_type

0 or 'Off' = Do not estimate; 1 = 'Fishery'; 2 = 'Survey'

Species

Species number

Selectivity_index

Index used to give fleets the same selectivity (otherwise, same as Fleet_code). Fleets sharing a value share one selectivity parameter block, with its penalties and priors accumulated once on the group's first non-Off fleet. Sharing the parameters is not enough on its own: the columns that shape the curve are read per fleet and must agree across the group, or the fleets end up with different selectivities. data_check() reports any that differ. To mirror a fleet, copy its fleet_control row and change only the identity and catchability columns. See vignette("model-options-and-functionality").

Selectivity

Selectivity form: 0 = "Fixed"; 1 = "Logistic"; 2 = "NonParametric" (Ianelli et al. 2018); 3 = "DoubleLogistic"; 4 = "DescendingLogistic"; 5 = "Hake" non-parametric (Taylor et al. 2014); 6 = "2DAR1"; 7 = "3DAR1" (Cheng et al. 2024); 8 = "DoubleNormal"; 9 = "NonParametricPM" (non-parametric random walk with penalized increments, AMAK "pm"; not integrable under random_sel); 11 = "LogisticPM" (logistic with a free age-1 selectivity, AMAK "pm"); 13 = "NonParametricIntegrable" (the Ianelli base curve whose deviations have a proper density, so random_sel integrates them; Time_varying_sel picks the structure: "Off", "IID" or "RandomWalk"); 15 = "DoubleNormalSS3" (Stock Synthesis size pattern 24, six parameters on SS3's scales, varied through selectivity linkages). Whether a form is age- or length-based is set by Selectivity_dimension.

Selectivity_dimension

"Age" or "Length".

N_sel_bins

Number of bins (ages or length bins) to estimate for non-parametric selectivity.

Time_varying_sel

Whether a time-varying selectivity should be estimated. 0 = "Off"; 1 = "IID" penalized deviates given Time_varying_sel_sd, or random effects with that sd estimated under random_sel; 2 = "AR1" REMOVED in 5.16.0, it was scored as "IID", so data_check() errors on it; use "IID" for the same fit, or linkage_spec(~ ar1(1 | Year)) for a real one; 3 = "Block" time blocks with no penalty; 4 = "RandomWalk" following Dorn; 5 = "RandomWalkAscending" on the ascending limb of double logistic only.

Time_varying_sel_sd

The fixed or initial sd to use for time varying selectivity; always required by the 2DAR1 and 3DAR1 forms.

Sel_start_year

First year the fleet's time-varying selectivity deviations are estimated and penalized. Deviations before this year are fixed at 0: they have no data and no penalty (the selectivity penalties are anchored at this year), so estimating them would leave unidentified flat directions. Defaults to the earliest year of data across all fleets sharing the fleet's Selectivity_index (fleets sharing an index share one selectivity curve), bounded below by styr. Only used for time-varying selectivity; set explicitly to override.

Bin_first_selected

Age/length bin at which selectivity is non-zero

Sel_shape_sd, Sel_curvature_sd, Sel_devmag_sd

Selectivity penalties expressed as standard deviations (an intuitive alternative to the raw penalty weights in Sel_curve_pen1/2/3). Each penalty is a Gaussian SSQ, so the weight is 1 / (2 * sd^2): Sel_shape_sd is the shape/smoothness penalty (Sel_curve_pen1), Sel_curvature_sd the 2nd-difference curvature penalty (Sel_curve_pen2), Sel_devmag_sd the deviation-magnitude penalty (Sel_curve_pen3). Supplying an SD sets the corresponding Sel_curve_pen column; a legacy Sel_curve_pen value is never overwritten. Sel_shape_sd applies to NonParametric (2/9/13) and LogisticPM (11), Sel_devmag_sd to "NonParametricPM" (9) and "LogisticPM" (11) only, the slot being unread by "NonParametric" (2) and "NonParametricIntegrable" (13); Sel_curvature_sd is NonParametric-only (2/9/13; LogisticPM does not use Sel_curve_pen2). The 2DAR1/3DAR1 forms reuse Sel_curve_pen for logit-scale correlations and keep those columns.

Sel_shape_dir

Direction of the shape penalty when Sel_shape_sd is used: "Decreasing" (default, penalize decreasing selectivity) or "Increasing" (sets the sign of Sel_curve_pen1). Only "NonParametricPM" (9) with Sel_shape_mode = "Directional" branches on the sign; every other form, that one under "Smooth" included, applies Sel_curve_pen1 as written, so "Increasing" is refused there rather than negating it.

Sel_curve_pen1

Raw shape/smoothness penalty weight for non-parametric (type 2/9/13) and LogisticPM (11) selectivity. The intuitive alternative is Sel_shape_sd, which sets this to 1 / (2 * sd^2). The 2DAR1/3DAR1 forms reuse this column as a logit-scale correlation.

Sel_curve_pen2

Raw 2nd-difference curvature penalty weight for non-parametric selectivity (the intuitive alternative is Sel_curvature_sd). Reused as a logit-scale correlation by the 2DAR1/3DAR1 forms.

Sel_curve_pen3

Raw deviation-magnitude penalty weight, read only by "NonParametricPM" (9) and "LogisticPM" (11) (the intuitive alternative is Sel_devmag_sd).

Sel_pen_first_bin

First bin (age or length, on the fleet's Selectivity_dimension) for the non-parametric shape penalty. NA defaults to Bin_first_selected.

Sel_pen_last_bin

Last (left) bin of the shape-penalty pairs. NA defaults to nbins - 2.

Sel_shape_mode

Shape-penalty mode: "Directional" (default) or "Smooth" (two-sided second-difference penalty, RTMB).

Sel_avgsel_pen

Weight on the AMAK average-selectivity base-level penalty (type 9 only): 0 = off (default), 10 matches AMAK.

Sel_cap_bin

NonParametricPM selectivity bin cap. NA (default) applies no cap.

Sel_norm_bin

Where selectivity is normalized to 1. "Max" normalizes by the largest value, "Off" (or blank, the default) not at all, and "All" on a LogisticPM fleet covers its whole penalty range. A number is an absolute age for an age-based fleet (6 means age 6, not the 6th bin) or a 1-based length-bin ordinal for a length-based one, and must lie in the range the fleet is selected over; below it the curve is zeroed, so a value there is read as "Max". In a two-sex model this says only where the reference is taken, whether it is pooled across the sexes is Sel_norm_scope. See the sex-structure section of vignette("model-options-and-functionality").

Sel_norm_bin_upper

Optional upper age/length bin. "Off" (or blank, the default) normalizes at a single bin; a number normalizes by mean selectivity between Sel_norm_bin and this bin, on the same scale and within the same range.

Sel_norm_scope

Whether selectivity normalization pools its reference across sexes, orthogonal to Sel_norm_bin (which says where the reference is taken). "WithinSex" divides each sex by its own reference, so both reach 1 and only the shape differs by sex; "AcrossSexes" (default) pools one reference over both sexes, so the less-selected sex stays below 1 and relative sex-specific selectivity is retained. No effect on a one-sex species or where Sel_norm_bin is NA.

Observation_units

Units of the observation: 1 = weight (mt), 2 = numbers caught (thousands of fish). Drives both catch and index prediction.

Weight_index

Weight-at-age (weight) index to use for calculation of derived quantities

Age_transition_index

Age transition matrix (e.g. growth trajectory) index to use for derived quantities to convert age to length

Ageing_error_index

Ageing error matrix index, matching Ageing_error_index in age_error. Fleets sharing a value read the same matrix. Omit the column, or leave it NA, to use the fleet's own species, which is what a model with one matrix per species has always done.

Catchability_index

Index used to give fleets the same catchability (otherwise, same as Fleet_code). Fleets sharing a value share one q parameter, so the group holds only one answer to whether q is estimated: the group's first non-Off fleet decides for all of them, regardless of fleet type. Give a fleet its own value when it should be estimated independently. "Analytical" and "AnalyticalArith" are the exception, they solve q from each fleet's own index observations, so a group containing one does not share a catchability. data_check() reports that case, a Fixed lead that leaves fleets on different inits, and a Catchability or Time_varying_q that differs within a group. See vignette("model-options-and-functionality").

Catchability

Estimate catchability? 0 or "Fixed" = fixed at prior; 1 or "Estimated" = estimate a single parameter; 2 or "Estimated-with-prior" = estimate a single parameter with a prior; 3 or "Analytical" = geometric-mean analytical q (Walters and Ludwig 1994); 4 = power equation (not yet implemented); 5 or "Environmental" = linear equation log(q_y) = q_mu + beta * index_y; 6 or "AR1" = removed, and data_check() errors on it, the deviates were never estimated, so q came back constant; express the Rogers et al. (2024) form as a q linkage with ar1(1 | Year) and observe, passed to fit_mod(qFun = ); 7 or "AnalyticalArith" = arithmetic-mean analytical q = mean(obs)/mean(pred), the AMAK/ebswp form used with the MVN survey likelihood.

Catchability_init

Starting value (or fixed value) for catchability. When Catchability = "Estimated-with-prior" this also centers the q prior.

Catchability_prior_sd

Log-scale standard deviation of the q prior: dnorm(log_q, log(Catchability_init) - bias_adjust_proc * sd^2/2, sd). Used only when Catchability = "Estimated-with-prior".

Time_varying_q

Whether a time-varying q should be estimated. 0 = "Off"; 1 = "IID" penalized deviate or random effect; 2 = "AR1" REMOVED in 5.16.0, it was scored as "IID", so data_check() errors on it; use "IID" for the same fit, or linkage_spec(~ ar1(1 | Year)) for a real one; 3 = "Block" time blocks with no penalty; 4 = "RandomWalk" from the mean following Dorn 2018. Under Catchability = "Environmental" this column instead holds the env_data column index, not a mode.

Time_varying_q_sd

The sd to use for the random walk of time varying q if set to 1

Estimate_index_sd

Estimate survey/index observation SD. Accepts integer codes or readable strings: 0 / "Fixed" = use the fixed SD from the data CV; 1 / "Estimated" = estimate; 2 / "Analytical" = analytically estimate following Walters and Ludwig 1994.

Index_sd

Starting/fixed value for the index observation SD, used when Estimate_index_sd = "Estimated".

Estimate_catch_sd

Estimate fishery/catch observation SD. Accepts integer codes or readable strings: 0 / "Fixed" = use the fixed SD from the data CV; 1 / "Estimated" = estimate; 2 / "Analytical" = analytically estimate following Walters and Ludwig 1994.

Catch_sd

Starting/fixed value for the catch observation SD, used when Estimate_catch_sd = "Estimated".

Index_distribution

Survey/index biomass observation likelihood family. 0 or "Lognormal" (default) = independent lognormal on the log observation. 1 or "MVN" and 2 or "MVNORM" = multivariate normal on the natural-scale residual (obs - q*pred), with a user-supplied covariance in index_cov (e.g. a VAST-derived Sigma); "MVN" reports the bare quadratic form matching AMAK/ebswp and ADMB, "MVNORM" the properly normalized density, and the two give an identical fit. Pair either with Catchability = "AnalyticalArith". 3 or "Normal" = independent normal on the same residual, with an ABSOLUTE sd from Log_sd, not a log-scale CV; reproduces AMAK avo_like/cpue_like term for term. 4 or "TruncatedNormal" = that normal left-truncated at zero. Prefer it to "Normal" unless an exact ADMB comparison is needed: an index cannot be negative, so it is the only natural-scale family whose likelihood and simulator are the same distribution. See vignette("model-options-and-functionality").

Comp_distribution

Age/length composition data distribution. Accepts integer codes or readable strings: -1 or "MultinomialAFSC" (default) = AFSC multinomial; 0 or "Multinomial" = full multinomial; 1 or "DirichletMultinomial".

Comp_weights

Composition weight. Its scale follows Comp_distribution: under a multinomial it is the natural-scale multiplier on the input sample size, while under a Dirichlet-multinomial the model uses exp(Comp_weights), so the column holds the log of the starting weight, a value of 1 is a starting weight of $e$; use 0 for a weight of 1. It is read when a model is built from scratch, so a refit keeps whatever weight it is given. A fit reports the weight its residuals imply, a harmonic mean of the ratio of effective to input sample size (McAllister & Ianelli 1997), in Comp_weights_mcallister; reweight_comps tunes towards it.

CAAL_distribution

Conditional age-at-length composition distribution. Accepts integer codes or readable strings: -1 or "MultinomialAFSC" = AFSC multinomial; 0 or "Multinomial" (default) = full multinomial; 1 or "DirichletMultinomial".

CAAL_weights

Conditional age-at-length (CAAL) composition weight. Its scale follows CAAL_distribution exactly as Comp_weights follows Comp_distribution: a natural-scale sample-size multiplier under a multinomial, the log of the starting weight under a Dirichlet-multinomial.

Comp_accum_young

Young-tail composition accumulation bin (AFSC ac_yng): a 1-based ordinal on the fleet's composition dimension; age/length bins below it are folded into it before the likelihood. NA (default) or 1 = no young accumulation.

Comp_accum_old

Old-tail composition accumulation bin (AFSC ac_old): a 1-based ordinal on the fleet's composition dimension; age/length bins above it are folded into it before the likelihood. NA (default), 0, or a value at/above the number of bins = no old accumulation.

Month

Observation month for the fleet (0 = not specified).

Proj_F_proportion

The proportion of future fishing mortality assigned to this fleet.

Description

Function to build parameter bounds based on Holsman et al 2015 and Kinzey and Punt 2010

Usage

build_bounds(param_list = NULL, data_list)

Arguments

Value

List of upper and lower bounds

Description

Holds environmental linkages on survey/index catchability q. The effect of an env_data covariate is written as a formula and can hold priors, bounds, and an estimation phase like any other linkage.

Usage

build_catchability(linkages = NULL)

Arguments

Details

Catchability is the one process that accepts link = "exponential", Stock Synthesis's environmental link type 1: q^exp(beta * x), which MULTIPLIES ⁠log q⁠ where "log" shifts it. It needs an estimated base q and a lognormal or t index likelihood, and cannot share a fleet with a random-effect q linkage; each of those is refused with the reason. See vignette("environmental-linkages-and-priors").

Value

A list of catchability settings for fit_mod().

Examples

# One temperature effect on q per fleet. `by = ~ fleet` is the default, and a
# fleet without an estimated survey catchability is not linkable, so at fit time
# restrict the linkage to the estimated-q fleets (here fleet 7 of BS2017SS).
build_catchability(linkages = list(
  q = linkage_spec(~ temp, fleet = 7)))

# Restrict it to fleets 1 and 3, with a prior on the slope
build_catchability(linkages = list(
  q = linkage_spec(~ temp, fleet = c(1, 3),
                   priors = list(temp = prior_normal(0, 1)))))

Description

Holds priors on the Dirichlet-multinomial composition-weighting overdispersion. The DM weight (the "theta" that scales the effective sample size) is otherwise an unpenalized free parameter; a linkage lets you put a prior on it through the same grammar as every other parameter. The three parameters target the three DM likelihoods:

theta_comp

age / length composition weight (comp_weights), per fleet.

theta_caal

conditional-age-at-length weight (caal_weights), per fleet.

theta_diet

diet composition weight (diet_comp_weights), per predator.

This is a prior-only process: the DM weight is a scalar, not a year-varying quantity, so only intercept formulas (~ 1) with a priors entry are meaningful. The prior is placed on the natural-scale DM weight theta = exp(weight), so a gamma() prior reads naturally. A linkage on a fleet whose Comp_distribution (or CAAL_distribution) is not "DirichletMultinomial" errors, since the weight has no effect there.

Usage

build_composition(linkages = NULL)

Arguments

Value

A list of composition-weighting settings for fit_mod().

Examples

# A weak gamma prior on the DM overdispersion of every fleet's age comps
build_composition(linkages = list(
  theta_comp = linkage_spec(~ 1, by = ~ fleet,
                            priors = list(`(Intercept)` = gamma(2, 0.5)))))

Description

Assemble (or edit) an Rceattle data list in R code rather than from an xlsx workbook. Supply only the data blocks a model uses, dimensions, biology (weight, maturity, sex_ratio, M1_base), a fleet_control, and the observation tables (catch_data, index_data, comp_data, ...), and the optional blocks a single-species model does not need (caal_data, emp_sel, diet_data, ...) are default-filled by clean_data(). The result is the same bare list read_data() returns and round-trips through write_data() unchanged.

Usage

build_data(base = NULL, file = NULL, ..., .check = TRUE)

Arguments

Details

Three entry points, which may be combined:

• from blocks: pass the elements as named arguments: build_data(nspp = 1, styr = 1977, ..., fleet_control = fc, catch_data = catch).

• from a file: build_data(file = "model.xlsx", projyr = 2060) reads the workbook via read_data(), then applies the overrides.

• from an existing object: build_data(base = BS2017SS, projyr = 2060) starts from a data list (e.g. a bundled dataset or a combine_data() result) and overrides the named blocks. This is the common copy-and-edit / combine-and-restamp workflow.

Element names are checked against the recognized schema: an override whose name is a near-miss of a known element (e.g. maturty) errors with a suggestion, so a typo is caught here rather than surfacing much later in a fit. Legacy top-level names (fsh_biom, srv_biom, wt, pmature, Pyrs) are mapped to their canonical equivalents.

Full validation runs at fit time, inside fit_mod(). build_data() runs only a light presence pre-check (.check = TRUE) so a missing required block is reported at construction with a clear message. The pre-check reads an attached model_config(), so a configuration held on the object is accounted for. Requirements that depend on fit-time settings passed directly to fit_mod() and stored nowhere on the data list are not knowable here and are left to that fit-time validation; see data_requirements() to preview them.

Value

An Rceattle data list (a bare list, as from read_data()).

See Also

read_data(), clean_data(), data_requirements(), combine_data(), fit_mod().

Examples

# Copy-and-edit a bundled dataset.
dat <- build_data(base = BS2017SS, projyr = 2060)

# Preview what a configuration needs before building it.
data_requirements(msmMode = 0)

Description

Specify the growth model for Rceattle

Usage

build_growth(
  fun = "empirical",
  growth_age_L1 = NA,
  sd_plus_group = NA,
  linkages = NULL,
  sd_form = NA,
  plus_group_length = NA,
  plus_group_decay = NA,
  pop_lengths = NULL
)

Arguments

Details

Every switch above takes a scalar, recycled across species, or a length-nspp vector for per-species values. Each defaults to NA, which inherits the matching field from data_list when it is set – so a refit, or a model built by the SS3 converter, keeps the original choice – and otherwise uses the fallback named with it.

Plus-group mean length. Fish older than the oldest age are pooled, so the plus group's mean length lies between the growth curve at the oldest age, L_A, and L-infinity:

• "M1": mean of L_A + (a/n)(Linf - L_A) over a = 0..n (n = nages), weighted by survival at the oldest age's base M1.

• "none": L_A (SS3 Linf_decay = -998).

• "SS3.24": as "M1" with weights exp(-0.2 a) and a = 0..A, A the oldest age (SS3 Linf_decay = -999).

• "decay": L_A and 2A further ages, each one more year along the von Bertalanffy curve, weighted by exp(-d a) with d = plus_group_decay (SS3 Linf_decay = d).

Under "none", "SS3.24" and "decay" the plus group also grows within the year like every other age, as SS3 does; under "M1" it keeps its Jan-1 length through the year.

Maturity-at-length is set in the data, not here: the per-species control columns L50_mat_len and slope_mat_len.

Value

A list of switches defining the growth model.

Examples

# Sex-specific von Bertalanffy with temperature on K, by species + sex
build_growth(
  fun = "vonBertalanffy",   # or fun = 1
  linkages = list(
    K = linkage_spec(
      formula = ~ temp,
      by      = ~ species + sex,
      priors  = list(temp = normal(0, 1))
    )
  )
)

Description

Defines the harvest control rule and its associated reference points. The rules are tabulated, with which ones work under predation, in vignette("projections-and-reference-points"); the formulations of the three ramped rules are below.

Usage

build_hcr(
  HCR = 0,
  DynamicHCR = FALSE,
  Ftarget = 0.4,
  Flimit = 0.35,
  Ptarget = 0.4,
  Plimit = 0,
  Alpha = 0.05,
  Pstar = 0.45,
  Sigma = 0.5,
  Fmult = 1,
  HCRorder = 1
)

Arguments

Details

hcr = 5 or "NPFMC", the North Pacific Fishery Management Council Tier 3 spawner-per-recruit rule: Stock status: $SB > SB at Ftarget$

$Fofl = Flimit$

$Fuse = Ftarget$

Stock status: $Alpha < SB / SB at Ftarget \le 1$

$Fofl = Flimit * (SB / SB_{Ftarget} - Alpha) / (1 - Alpha)$

$Fuse = Ftarget * (SB / SB_{Ftarget} - Alpha) / (1 - Alpha)$

Stock status: $SB / SB at Ftarget \le Alpha$ or $SB < Plimit * SB0$

$Fofl = 0$

$Fuse = 0$

hcr = 6 or "PFMC", the Pacific Fishery Management Council category 1 40-10 annual catch limit rule (PFMC 2020), using Fspr in single species or the F reaching X% of SSB0 under predation. The uncertainty buffer is qnorm(Pstar, Flimit, Sigma), which equals Flimit + qnorm(Pstar, 0, Sigma) and so sits below Flimit when $Pstar < 0.5$. The taper runs between $SB0 * Plimit$ and $SB0 * Ptarget$, so the 40-10 shape needs Ptarget = 0.40 and Plimit = 0.10; the default Plimit = 0 gives a 40-0 rule. Stock status: $SB > SB0 * Ptarget$

$Fofl = Flimit$

$Fuse = qnorm(Pstar, Flimit, Sigma)$

Stock status: $SB0 * Plimit < SB \le SB0 * Ptarget$

$Fofl = Flimit$

$Fuse = qnorm(Pstar, Flimit, Sigma) * \frac{SB0 * Ptarget * (SB - SB0 * Plimit)}{SB * SB0 * (Ptarget - Plimit)}$

Stock status: $SB < SB0 * Plimit$

$Fofl = 0$

$Fuse = 0$

hcr = 7 or "SESSF", the Southern and Eastern Scalefish and Shark Fishery Tier 1 spawner-per-recruit rule, with F_Limit = F_(20%) and B_Limit = SB_20 (AFMA 2017): Stock status: $SB > SB0 * Ptarget$

$Fofl = Flimit$

$Fuse = Ftarget$

Stock status: $SB0 * Ptarget > SB > SB0 * Plimit$

$Fofl = Flimit * (SB / (SB0 * Plimit) - 1)$

$Fuse = Ftarget * (SB / (SB0 * Plimit) - 1)$

Stock status: $SB < SB0 * Plimit$

$Fofl = 0$

$Fuse = 0$

Value

A list containing the harvest control rule and associated biological reference points.

Examples

# Tier 3 NPFMC control rule: F40% as FABC, F35% as FOFL.
build_hcr(HCR = "NPFMC", DynamicHCR = FALSE, Ftarget = 0.4, Flimit = 0.35)

# No fishing: the projection a B0 reference point is read from.
build_hcr(HCR = "NoFishing")

Description

Reads a data list and map to update the map argument based on the HCR specified in build_hcr

Usage

build_hcr_map(
  data_list,
  map,
  debug = FALSE,
  all_params_on = FALSE,
  HCRiter = 1
)

Arguments

Value

a list of map arguments for each parameter

Description

Specify the residual natural mortality (M1) model for Rceattle

Usage

build_M1(
  M1_model = 0,
  M1_re = 0,
  updateM1 = FALSE,
  M1_use_prior = FALSE,
  M2_use_prior = FALSE,
  M_prior = 0.4,
  M_prior_sd = 0.35,
  M1_indices = NA,
  linkages = NULL
)

Arguments

Details

Which M1_re families work with which M1_model. Anything but 0 needs an M1_model of 0, 1, 2, or 3 with a by-year family (2 / 5); the by-age families are refused under 3 because an age-varying deviation and an age-specific level are not separately identifiable, and 4 and 5 are refused entirely (use a catchability-style linkage instead). M1_model = 2 needs two sexes under every M1_re family, because a single-sex species would take its deviations in the padding sex slot, free and scored by no density. fit_mod() downgrades M1_model to 1 there, so the distinction shows only in a direct build_map() call.

Value

A list of switches for defining the M1 model.

M1_base completeness

M1_base needs a finite, positive residual M (per year) for every age of every sex the species has; it is taken on the log scale, so a blank or a zero is refused, naming the offending ⁠(species, sex, age)⁠ cells. Columns past a species' own nages are ignored, so a row may stop at its last age.

Examples

# Sex/age-invariant M with a temperature linkage on M1
build_M1(
  M1_model = "sex_age_invariant",
  linkages = list(
    M1 = linkage_spec(
      formula = ~ temp,
      by      = ~ species,
      priors  = list(temp = normal(0, 0.5))
    )
  )
)

Description

Builds the TMB map, which tells TMB which parameters to estimate and which to hold fixed (a fixed parameter is mapped to NA), and which parameters share a single estimated value. One helper handles each process block (recruitment, M1, predation, selectivity, catchability, ...).

Usage

build_map(
  data_list,
  params,
  debug = FALSE,
  random_rec = FALSE,
  random_sel = FALSE,
  random_q = NULL
)

Arguments

Details

Fleets sharing a Selectivity_index or a Catchability_index have the lead fleet's map slice copied over the rest of the group, so a per-fleet setting that differs within a group is resolved here rather than reported. Those disagreements are checked in data_check(), which fit_mod() runs first; build_map() does not call it, so a caller invoking build_map() directly should run data_check() too.

Value

A list containing the factorized TMB map (mapFactor) and the original map matrix/array list (mapList).

Description

Maps catchability base parameters (index_log_q), time-varying deviations (index_q_dev), and environmental linkages (index_q_beta) for every fleet that holds fitted index_data, a fishery with a CPUE series as much as a survey. A fleet with no index rows gets none of them, whatever its Catchability says, since a q with no index to inform it is a flat direction. Sharing overrides this: adjust_map_shared_params() then copies the lead fleet's slice across each Catchability_index group.

Usage

build_map_catchability(map_list, data_list, nyrs_hind, random_q = FALSE)

Arguments

Value

Updated map_list.

Description

Sets all parameters in the map list to NA, except the dummy parameter, for use in debug or testing modes.

Usage

build_map_debug(map_list, debug)

Arguments

Value

Updated map_list.

Description

Maps fishing mortality parameters (log_F) and related targets. Also maps data weighting parameters (catch_log_sd, comp_weights, caal_weights).

Usage

build_map_f_and_data_weights(map_list, data_list, nyrs_hind)

Arguments

Value

Updated map_list.

Description

Turns off (sets to NA) most population and fleet parameters for species where the dynamics are fixed (estDynamics > 0).

Usage

build_map_fixed_natage(map_list, data_list)

Arguments

Value

Updated map_list.

Description

Maps the growth parameters (log_growth_pars, growth_log_sd, weight_length_pars) according to growth_model. Everything starts mapped off and each growth function turns on the parameters it uses.

Time-varying growth comes from the linkage grammar (build_growth(linkages = )), whose random effects hold their own density and map.

Usage

build_map_growth(map_list, data_list, nyrs_hind)

Arguments

Value

Updated map_list.

Description

Maps the fixed parameters (log_M1) and the random effects parameters (log_M1_dev, M1_dev_log_sd, M1_rho) based on M1_model and M1_re settings.

Usage

build_map_m1(map_list, data_list, nyrs_hind)

Arguments

Value

Updated map_list.

Description

Maps predation suitability parameters (log_gam_a, log_gam_b, log_phi) and diet weight parameters based on msmMode and suitMode.

Usage

build_map_predation(map_list, data_list)

Arguments

Value

Updated map_list.

Description

Maps the recruitment deviations (rec_dev, init_dev), stock-recruitment parameters (rec_pars) and Rec dev variances (R_log_sd)

see build_srr() for options,

Usage

build_map_recruitment(map_list, data_list, nyrs_hind, nyrs_proj, random_rec)

Arguments

Value

Updated map_list.

Description

Maps base selectivity parameters (log_sel_slp, sel_inf, sel_coff) and time-varying deviations, based on Selectivity and Time_varying_sel settings in fleet_control.

Usage

build_map_selectivity(map_list, data_list, nyrs_hind, random_sel)

Arguments

Details

Selectivity sets the curve's shape and Time_varying_sel whether it varies over time; for the allowed values of both, see "3. Selectivity (per fleet)" in vignette("model-options-and-functionality"), or rceattle_switch_maps for the code-to-name maps themselves. N_sel_bins is how many age or length bins the non-parametric and AR1 forms estimate (Selectivity 2, 5, 6, 7, 9 or 13); it is unused otherwise. random_sel in fit_mod() treats the deviations as random effects and estimates their variance.

The double-normal forms blend Gaussian ascending and descending limbs at the peak. DoubleNormal reads its four parameters from the shared arrays: sel_inf[1] = peak, sel_inf[2] = logit(right_floor) (the right-tail floor, SS3 P6/end_logit), log_sel_slp[1] = log(sigma_asc) and log_sel_slp[2] = log(sigma_desc). right_floor -> 0 is dome-shaped; right_floor -> 1 is logistic ascending only. DoubleNormalSS3 instead carries its own six-parameter array (sel_dn6).

Value

Updated map_list.

Description

Construct the TMB parameter list, with every parameter set to its starting value, for an Rceattle model.

Usage

build_params(data_list)

Arguments

Value

a named list of model parameters at their starting values

Description

Holds environmental linkages on selectivity parameters. An effect is written as a formula and composes additively with any Time_varying_sel process error on the same fleet, the two being separate mechanisms: a covariate effect against a deviation.

The linkages keys are the shape parameters of the parametric forms, so which apply depends on the fleet's Selectivity:

slp_asc, slp_desc

ascending / descending logistic slope (log scale); the double-normal widths, aliased sigma_asc / sigma_desc.

inf_asc, inf_desc

ascending / descending inflection age or length (natural scale); the double-normal peak and logit right-floor, aliased peak / right_floor.

coff

non-parametric selectivity-at-bin coefficients.

apical

a log-scale multiplier on one sex's whole curve, applied after the form and before normalization, so every estimated form takes it. See Details.

DoubleNormalSS3 (code 15) carries its own six under Stock Synthesis's names – dn_peak, top_logit, ascend_se, descend_se, start_logit, end_logit – and a linkage is the only way to vary that form over time, since Time_varying_sel is refused on it. A time block is ~ cut(Year, ...); annual deviations are a random-effect term.

Every parameter but apical takes link = "log" (multiplicative on the natural parameter) or link = "identity" (additive); apical is already a multiplier and takes "log" only.

Usage

build_selectivity(linkages = NULL)

Arguments

Details

The apical offset. Fishing mortality is one log_F per fleet and year shared by the sexes, and no form has a height parameter, so a sex difference in F can only come from selectivity. apical multiplies one sex's curve by exp(log_sel_apical) bin by bin, as Stock Synthesis's male-offset option does. Name the fleet and the sex that holds it (fleet = 3, sex = "male"; by defaults to ~ fleet + sex here) and the other sex is the reference. It equals the ratio of the sexes' peak heights only where their shapes peak equally, so for a sex-specific dome read it as the multiplier and take the peak ratio from fit$quantities$sel_at_age.

Only the contrast between sexes is identified, the common level being log_F, so one sex holds it and the fit is refused if both do, if no fleet or sex is named, if the species has one sex, on a Fixed, AR1 or Fleet_type = "Off" fleet, on a fleet sharing another's Selectivity_index, under link = "identity" (which could drive the multiplier negative), or where Sel_norm_scope = "WithinSex" would divide it straight back out. The contrast is informed only by joint composition (comp_data$Sex = 3); with single-sex compositions it rests on its prior, and fit_mod() warns when a fleet has neither. An intercept prior is on the multiplier's natural scale, so a lognormal() centred on 1 means no offset. See vignette("environmental-linkages-and-priors").

Priors on a selectivity parameter. An intercept-only formula (~ 1) with a priors entry puts a prior on the parameter itself, adding no year-to-year offset. Read it on the parameter's own scale: the slopes are log scale (lognormal()), the inflections natural scale (normal()). A prior targets one parameter, so in a two-sex model an unstratified ~ 1 prior constrains sex 1 only – use by = ~ sex. An init on a selectivity intercept has no effect, since the starting value comes from the data, and a prior on the double-normal right_floor is unsupported. Fleets sharing a Selectivity_index estimate one block, so place the prior on the group's lead fleet (its first that is not Off): a prior on a follower would penalize the shared block once per sharing fleet. On an Off fleet a prior is dropped with a warning and the fleet's linkage coefficients are held at their starting values, since the fleet has no data; on an Off member of a shared block it moves to the block's lead instead.

Like every selectivity linkage, a covariate acts in the hindcast years; projection years hold the last hindcast year's curve.

Value

A list of selectivity settings for fit_mod().

Examples

# A cold-pool effect on the ascending inflection of a logistic fleet
build_selectivity(linkages = list(
  inf_asc = linkage_spec(~ cold_pool, by = ~ fleet)))

# A normal prior on the ascending inflection (intercept-only formula)
build_selectivity(linkages = list(
  inf_asc = linkage_spec(~ 1, priors = list(`(Intercept)` = normal(0, 3)))))

Description

Sets the stock-recruit curve and how recruitment is estimated.

Relationships currently implemented:

• srr_fun = 0 or "mean": no stock-recruit relationship. Recruitment is $R0$ (log scale) plus annual deviates, and the reported steepness holds its placeholder 0.99 rather than being estimated.

  $R_y = exp(R0 + R_{dev,y})$

• srr_fun = 2 or "BevertonHolt":

  $R_y = \frac{\alpha_{srr} * SB_{y-minage}}{1+\beta_{srr} * SB_{y-minage}}$

• srr_fun = 4 or "Ricker":

  $R_y = \alpha_{srr} * SB_{y-minage} * exp(-\beta_{srr} * SB_{y-minage})$

Both curves are the deterministic mean; realized recruitment applies the annual log deviation, $R_y \cdot exp(R_{dev,y})$. For numerical stability the Ricker $\beta_{srr}$ is estimated divided by 1,000,000, so the fitted beta is 1e6 times the coefficient above. Bmsy_lim is on the natural scale: the template divides it by 1e6 before comparing, so supply it in metric tons.

With srr_pred_fun > 0 and srr_fun = 0, the hindcast estimates recruitment as in the mean form while the curve is fitted between srr_hat_styr and srr_hat_endyr as a penalty, then used in the projection.

Multispecies models. Spawning biomass per recruit is undefined once mortality includes predation, so under msmMode > 0 a curve fitted in the hindcast estimates its initial recruitment level (R0) rather than deriving it, and under the fished initModes (3, 4) that level trades off against the initial F. A steepness prior is refused and steepness is reported as 0; priors on alpha or beta go through linkages.

Two refusals worth knowing before you write a linkage: under a curve fitted in the hindcast a single-species R0 is derived from alpha and beta, so an R0 linkage is refused (a multispecies model accepts an intercept-only one), and a Ricker lognormal linkage prior on alpha duplicates the deprecated srr_est_mode = "LognormalPrior", so using both is refused.

Usage

build_srr(
  srr_fun = 0,
  srr_pred_fun = srr_fun,
  proj_mean_rec = TRUE,
  srr_mse_switchyr = NULL,
  srr_hat_styr = NULL,
  srr_hat_endyr = NULL,
  srr_est_mode = 1,
  srr_prior = 4,
  srr_prior_sd = 1,
  srr_alpha_init = NULL,
  srr_beta_init = NULL,
  srr_indices = NA,
  Bmsy_lim = NA,
  linkages = NULL
)

Arguments

Value

A list containing the stock recruitment relationship settings

Starting values

None of the starting values knows the stock's scale, and a bad one fails the fit rather than fitting badly. Mean recruitment (R0) starts at $e^9 = 8103$ thousand fish; alpha starts at srr_prior (default 4) wherever that is an alpha and at $e^3$ otherwise; beta starts at 3. Set them with srr_alpha_init / srr_beta_init or a linkage init.

$\beta$ sets the density dependence in $R = \alpha S / (1 + \beta S)$, so it must be on the order of $(\alpha - 1/\phi_0) / R_0$, typically $10^{-3}$ or smaller for a stock measured in tonnes. Starting three orders of magnitude away drives predicted recruitment to near zero and the optimizer returns NA/NaN gradient evaluation. From a steepness $h$ and unfished spawning biomass per recruit $\phi_0$:

$\alpha = \frac{4h}{\phi_0 (1 - h)}, \qquad \beta = \frac{\alpha - 1/\phi_0}{R_0}.$

Under predation, where $\phi_0$ is undefined, seed them from a curve fitted to an earlier fit's SSB-recruitment pairs.

Examples

# Mean recruitment: no stock-recruit relationship fitted.
build_srr(srr_fun = "mean")

# Beverton-Holt with a lognormal prior on alpha (mean 5, log-scale SD 0.5),
# species 1 only.
build_srr(srr_fun = "BevertonHolt",
          linkages = list(alpha = linkage_spec(~ 1, species = 1,
            priors = list(`(Intercept)` = prior_lognormal(log(5), 0.5)))))

# Mean recruitment with the curve as a penalty (Ianelli form).
build_srr(srr_fun = "mean", srr_pred_fun = "BevertonHolt")

# A prior on Beverton-Holt steepness, the one prior linkages cannot express
# (single-species only).
build_srr(srr_fun = "BevertonHolt", srr_est_mode = "LognormalPrior",
          srr_prior = 0.8, srr_prior_sd = 0.2)

Description

Check which saved MSE simulation files can be loaded

Usage

check_mse(dir = NULL, file = NULL)

Arguments

Value

list of MSE simulations/run

Description

Fills optional blocks with correctly-shaped empties, filters observations to the model's year range, extends catch to the projection years, and re-keys any survey-index covariance matrices. Whether a missing block is actually a problem is left to the fit-time validation, which knows the model configuration.

Usage

clean_data(data_list)

Arguments

Description

Extract estimated parameters from an Rceattle fit

Usage

## S3 method for class 'Rceattle'
coef(object, ...)

Arguments

Value

A named numeric vector of estimated fixed-effect parameters (the optimizer's par). Returns NULL if the model was not estimated.

Description

Concatenate two data lists into one. env_data is taken from data_list1, and diet data must be updated afterward.

Usage

combine_data(data_list1 = NULL, data_list2 = NULL)

Arguments

Description

Function to evaluate the simulation performance with regard to bias using the median relative error (MRE) and precision using the coefficient of variation.

Usage

compare_sim(operating_mod, simulation_mods, object = "quantities")

Arguments

Value

A data frame summarizing simulation performance metrics

Note

Compares against operating_mod only, so it is valid when the replicates redrew the observations alone. If they were produced with sim_mod(process = ) / self_test(process = ), the operating model's deviations are no longer what generated the data and the bias this reports is an artefact. Compare against attr(sims, "process_sim") in that case. Passing a self_test result that holds process draws warns, since the attribute makes that case detectable rather than only documented.

Description

Runs the post-fit convergence battery and returns a single structured object. Each check yields a record with a common schema (id, tier, severity, message, data); severity is one of "OK", "NOTE", "WARN", "FAIL". The object's status is the worst severity present.

Usage

convergence_diagnostics(object, ...)

Arguments

Details

fit_mod() runs this automatically and attaches the result as fit$convergence; call convergence_diagnostics() directly to re-run it on any fit. Checks cover the optimizer gradient, a requested sdreport that did not return, an sdreport that was never requested, Hessian positive-definiteness and conditioning, parameters on bounds, parameters past them, a deviation variance estimated to zero, phasing, parameter estimability, a numbers-at-age, Ricker-intercept or recruitment floor that was reached, and the stock-recruit curve.

Value

An object of class "Rceattle_convergence": a list with status (overall worst severity) and checks (named list of records).

Description

Introspects the Rceattle data-element catalogue and reports, for a given model configuration, which top-level data_list elements are Required, Optional (used if supplied, otherwise default-filled by clean_data()), or Ignored (not consulted because the feature that would use them is switched off). It answers "what do I actually need to supply for this model?" without having to read the validation code or the switch tables, the same conditions enforced at fit time (they share one declarative table).

Usage

data_requirements(
  data_list = NULL,
  msmMode = 0,
  growth_model = 0,
  estDynamics = 0,
  Selectivity = NULL,
  Index_distribution = NULL,
  Ceq = NULL
)

Arguments

Details

The configuration can be given either as an existing (possibly partial) data_list, its switches are normalized through clean_data() / switch_check() so the conditions evaluate against filled defaults, or, when no data_list is supplied, built from the convenience arguments.

Requirements are conditional: e.g. diet_data, ration_data and the bioenergetics scalars are Ignored under single-species (msmMode = 0) but Required under multispecies (msmMode > 0); caal_data is Optional unless growth_model > 0; NByageFixed is Ignored unless estDynamics > 0; emp_sel is Required only for fleets with Selectivity = "Fixed"; index_cov only for Index_distribution = "MVN".

Value

A data.frame with one row per data element and columns:

element

the data_list element name.

category

grouping (dimensions / biology / fishery / composition / predation / environment).

status

"Required", "Optional", or "Ignored".

condition

the condition under which the element is required ("always" for the core backbone).

default

for Optional elements, the clean_data() default used when the element is absent.

Rows are ordered Required, then Optional, then Ignored.

See Also

clean_data(), build_data().

Examples

# Single-species: predation/diet inputs are Ignored, comp_data Optional.
data_requirements(msmMode = 0)

# Multispecies: diet_data and bioenergetics become Required.
data_requirements(msmMode = 1)

# From a bundled dataset.
data_requirements(BS2017SS)

Description

A fitted Rceattle model object for the multispecies run of CEATTLE for the Eastern Bering Sea (walleye pollock, Pacific cod, arrowtooth flounder).

Usage

EBS_ms_run

Format

An object of class Rceattle.

Description

A fitted Rceattle model object for the single-species run of CEATTLE for the Eastern Bering Sea with natural mortality estimated.

Usage

EBS_ss_M_run

Format

An object of class Rceattle.

Description

A fitted Rceattle model object for the single-species run of CEATTLE for the Eastern Bering Sea (walleye pollock, Pacific cod, arrowtooth flounder) without estimating natural mortality.

Usage

EBS_ss_run

Format

An object of class Rceattle.

Description

Collects the optimizer and reporting settings into one object, so a call to fit_mod() stays about the model rather than how it is fit. Pass the result through fit_mod()'s fit_control argument, where it overrides the equivalent individual arguments, which are kept for back-compatibility.

Usage

fit_control(
  bias.correct = FALSE,
  getsd = TRUE,
  getJointPrecision = TRUE,
  getReportCovariance = FALSE,
  projection_uncertainty = FALSE,
  selectivity_se = FALSE,
  comp_offset = NULL,
  bias_adjust_obs = TRUE,
  bias_adjust_proc = TRUE,
  use_gradient = TRUE,
  rel_tol = 1,
  loopnum = 5,
  newtonsteps = 0,
  phase = FALSE,
  TMBfilename = NULL,
  verbose = 1,
  nlminb_control = list(eval.max = 1e+09, iter.max = 1e+09, trace = 0)
)

Arguments

Value

A list of class "Rceattle_fit_control".

Selectivity standard errors

selectivity_se needs getsd = TRUE, and reports the error belonging to whichever sdreport the fit ends on. Under estimateMode = "Estimate" with an HCR that re-optimizes, that is the projection fit, where every selectivity parameter is mapped off, so every error comes back exactly 0; "Projection" estimates no selectivity at all. Use "Hindcast", or "Estimate" with projection_uncertainty = TRUE, to get an error from a fit that estimated the curve. fit_mod() warns in each of these cases.

It is off by default because the delta method forms a Jacobian of every reported value against every parameter, so the cost is their product: on Atka2022, 1,012 values against 584 parameters. The error is on the log scale rather than the logit, because the non-parametric forms normalize to mean selectivity 1 instead of a maximum of 1 – 58% of Atka2022's sel_at_age exceeds 1, up to 3.06 – so a logit is undefined over most of the array.

Rows cover estimated, age-based lead fleets only, starting at each fleet's first selected bin. Four kinds of cell are a structural zero (a Fixed fleet's empirical curve, a length-based fleet's growth-matrix projection, a bin below Bin_first_selected, and array padding), and one log(0) = -Inf on the tape turns every quantity in the sdreport to NaN, biomass and SSB included. All four are identified from the data, so the reported set never depends on a parameter value, and a value that underflows to zero is floored so it cannot reintroduce the -Inf. See plot_selectivity(), which draws the interval.

Examples

# Quick-and-dirty fit: skip sdreport, single optimizer pass
ctl <- fit_control(getsd = FALSE, loopnum = 1)

# Production fit with bias correction and joint precision
ctl <- fit_control(bias.correct = TRUE, getJointPrecision = TRUE)

Description

Estimate CEATTLE population parameters by maximum likelihood, and optionally project the stock and apply a harvest control rule.

Usage

fit_mod(
  data_list = NULL,
  inits = NULL,
  map = NULL,
  bounds = NULL,
  file = NULL,
  estimateMode = 0,
  random_rec = FALSE,
  random_q = FALSE,
  random_sel = FALSE,
  HCR = build_hcr(),
  niter = 3,
  recFun = build_srr(),
  M1Fun = build_M1(),
  growthFun = build_growth(),
  qFun = build_catchability(),
  selFun = build_selectivity(),
  compFun = build_composition(),
  msmMode = 0,
  avgnMode = 0,
  initMode = "NonEquilibrium",
  suitMode = 0,
  suit_styr = NULL,
  suit_endyr = NULL,
  fit_control = NULL,
  config = NULL,
  quiet_data_check = FALSE,
  ...
)

Arguments

Details

CEATTLE is an age-structured population model, fit with or without predation mortality; msmMode selects which.

What config overlays, at two levels. The estimation controls (estimateMode, random_rec / random_q / random_sel, suit_styr / suit_endyr, fit_control) overlay only arguments you did not pass, so an explicit argument always wins. The stored model_config merges into the data object's field by field: only fields the config actually set are imposed, so a config built with model_config() keeps the linkages held on the data object. Where both set a field and disagree, the config wins and the difference is warned about by name. A config written before that field record existed counts as having set its non-default fields.

Value

A list of class "Rceattle" including:

• data_list: List of data inputs

• initial_params: List of starting parameters

• bounds: Parameter bounds used for estimation

• map: List of map used in TMB

• obj: TMB model object

• opt: Optimized model object from nlminb

• sdrep: Object of class sdreport from TMB, with standard errors

• estimated_params: List of estimated parameters

• quantities: Derived quantities from CEATTLE

• run_time: Model run time

Examples

data(BS2017SS)
ss_run <- fit_mod(
  data_list    = BS2017SS,
  estimateMode = 0,
  msmMode      = 0,
  fit_control  = fit_control(phase = FALSE, verbose = 0)
)

Description

A data list containing inputs for a multispecies CEATTLE model fit to three groundfish species on Georges Bank. See BS2017SS for format details.

Usage

GeorgesBank3spp

Format

A list with the same structure as BS2017SS.

Description

This function calculates a probability transition matrix that defines the probability of a fish of a given age belonging to specific length bins. It supports Von Bertalanffy and Richards growth models and includes a Stock Synthesis (SS) style plus-group correction.

Usage

get_growth_matrix_r(
  fracyr,
  nsex_sp,
  nages_sp,
  nlengths_sp,
  nyrs,
  lengths_sp,
  minage_sp,
  maxage_sp,
  growth_params_sp,
  growth_log_sd_sp,
  growth_model_sp,
  M_plus_sp
)

Arguments

Value

A 4D array of probabilities with dimensions (sex, age, length, year).

Description

Converts a growth matrix (length-at-age probabilities) into mean weight-at-age using a length-weight relationship (W = a * L^b).

Usage

get_weight_at_age_r(
  nsex_sp,
  nages_sp,
  nlengths_sp,
  nyrs,
  lengths_sp,
  length_at_age,
  growth_matrix,
  lw_params
)

Arguments

Details

The function calculates midpoints for length bins to avoid bias. For the first bin, it assumes the width is equal to the second bin's width. The final weight-at-age is the expected value across all length bins for that age.

Value

A 3D array of mean weights with dimensions (sex, age, year).

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska groundfish data through 2018. See BS2017SS for format details.

Usage

GOA2018SS

Format

A list with the same structure as BS2017SS.

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska arrowtooth flounder data. See BS2017SS for format details.

Usage

GOAatf

Format

A list with the same structure as BS2017SS.

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska arrowtooth flounder data through 2023. See BS2017SS for format details.

Usage

GOAatf2023

Format

A list with the same structure as BS2017SS.

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska Pacific cod data. See BS2017SS for format details.

Usage

GOAcod

Format

A list with the same structure as BS2017SS.

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska walleye pollock data. See BS2017SS for format details.

Usage

GOApollock

Format

A list with the same structure as BS2017SS.

Description

A list containing biomass, spawning stock biomass, and recruitment reference values from the 2018 Gulf of Alaska Stock Assessment and Fishery Evaluation (SAFE) report.

Usage

GOAsafe2018

Format

A list with components:

biomass

Total biomass time series

ssb

Spawning stock biomass time series

recruitment

Recruitment time series

Description

Refits the Rceattle model from starting values perturbed by N(0, sd) around the model's initial (pre-fit) parameters, to check convergence robustness.

Usage

jitter(
  object = NULL,
  njitter = 50,
  sd = 0.2,
  phase = FALSE,
  seed = 123,
  cores = NULL,
  getsd = NULL,
  timeout = Inf,
  fit_control = NULL,
  ...,
  Rceattle = NULL
)

Arguments

Value

a list of 1. Rceattle_list, the converged jitters, and 2. nll, their objective values. Non-converged (or timed-out) starts are dropped and reported in a message, so both can be shorter than njitter, and that count is itself the result, since the whole point is what fraction of random starts reach the same optimum.

Note

Attaching Rceattle masks jitter, because that function is not generic and this name is exported, so the mask cannot be removed by making this an S3 method. It is harmless for the spellings base callers use: a numeric first argument is handed straight to base::jitter(), so jitter(1:5) and jitter(1:5, factor = 2) behave as base does. jitter(x = 1:5) and a positional second argument do not. x is not a parameter here, so it is refused as an unused argument; a positional second argument binds njitter, so jitter(1:5, 2) jitters with base's defaults and silently drops the 2. Name base's arguments, or call base::jitter() directly.

Examples

data(BS2017SS)
ss_run <- fit_mod(data_list = BS2017SS,
    inits = NULL, file = NULL,
    estimateMode = 0, random_rec = FALSE,
    msmMode = 0, avgnMode = 0,
    phase = FALSE, verbose = 0)
jitters <- jitter(ss_run, njitter = 10)

Description

Capture a linkage specification

Usage

linkage_spec(
  formula,
  param = NULL,
  data = NULL,
  by = ~species,
  species = NULL,
  sex = NULL,
  fleet = NULL,
  link = "log",
  init = NULL,
  bounds = NULL,
  priors = NULL,
  est_phase = 1L,
  observe = NULL,
  obs_sd = NULL,
  obs_sd_est = FALSE,
  integrate = TRUE
)

Arguments

Details

The reserved keys sigma and rho in init / priors route the random-effect deviation SD and, for ar1, the correlation. sigma means different things by structure – the innovation (per-step) SD for rw(), the marginal (stationary) SD for ar1() – so the two are not comparable across structures; see vignette("environmental-linkages-and-priors").

integrate = FALSE leaves the deviations in the objective as a plain penalty and reports them with standard errors. That reproduces the ADMB/AMAK convention behind the legacy Time_varying_sel / Time_varying_q switches, which a Laplace-integrated rw() cannot match, because the marginal likelihood carries a log-determinant term the penalized form has no counterpart for. Permitted only with a fixed SD (init = list(sigma = ), no sigma prior, plus a fixed rho for ar1), since estimating deviations and their SD jointly as fixed effects is degenerate, and never with observe, which must stay integrated.

Value

An Rceattle_linkage_spec object.

Examples

# A covariate on a growth parameter, one coefficient per species. Attach to a
# parametric build_growth(); the default "empirical" form consumes no linkages.
linkage_spec(~ temp, param = "K")

# A covariate on catchability for one fleet, named rather than numbered.
# `by` defaults to ~ fleet for catchability, so it need not be given.
linkage_spec(~ temp, fleet = "Pollock_survey_1_shelikof_acoustic")

# IID annual deviations; rw() for a random walk, ar1() for AR1.
linkage_spec(~ (1 | Year))

# A random walk estimated as a penalized fixed effect, which requires a fixed
# SD: the ADMB/AMAK convention behind the legacy Time_varying_* switches.
linkage_spec(~ rw(1 | Year), init = list(sigma = 0.05), integrate = FALSE)

# Intercept-only: a prior on the base parameter itself.
linkage_spec(~ 1, priors = list(`(Intercept)` = normal(0, 1)))

Description

The inverse of save_config(): reads a run-config YAML file and rebuilds the model_config() structure, estimation controls, and fit_control() bundle (reconstructing linkage formulas and priors). Attach the result to a fit with fit_mod(data_list, config = load_config("run.yaml")), or read ⁠$model_config⁠ off it to attach to a data list.

Usage

load_config(file)

Arguments

Value

An Rceattle_run_config.

See Also

save_config(), fit_mod().

Examples

f <- file.path(tempdir(), "run.yaml")
save_config(model_config(msmMode = 1, initMode = "Equilibrium"), f)

cfg <- load_config(f)
cfg$model_config$msmMode
# Apply with fit_mod(data_list, config = cfg).
unlink(f)

Description

Load saved MSE simulation runs

Usage

load_mse(dir = NULL, file = NULL, exclude = NULL, include_em = TRUE)

Arguments

Value

list of MSE simulations/run

Description

Returns the joint negative-objective from nlminb, sign-flipped, as a "logLik" object. df is the number of estimated fixed-effect parameters. The nobs attribute is intentionally omitted: counting "observations" in a stock assessment likelihood (with composition cells, indices, catches, and priors) is not well-defined, so stats::AIC() works (uses df) while stats::BIC() does not.

Usage

## S3 method for class 'Rceattle'
logLik(object, ...)

Arguments

Value

An object of class "logLik", or NULL if the model was not estimated.

Description

Model average of derived quantities

Usage

model_average(
  object = NULL,
  weights = NULL,
  uncertainty = FALSE,
  nboot = 10000,
  Rceattle = NULL
)

Arguments

Value

an Rceattle object with derived quantities weighted by the specified weights. The length of the derived quantities spans the years which overlap across all models.

Examples

## Not run: 
# Equal weights across two candidate models.
model_average(list(fit1, fit2), weights = c(1, 1))
# Weighted by AIC, with bootstrapped uncertainty.
model_average(list(fit1, fit2), weights = c(0.7, 0.3), uncertainty = TRUE)

## End(Not run)

Description

Bundles the model-structure arguments of fit_mod(), predation mode, initialization, harvest control rule, and the ⁠build_*()⁠ process specifications, into a single validated object that can be stored on a data list (data_list$model_config, e.g. via build_data(..., model_config = model_config(...))). A configuration then travels with the data instead of living only in the fit_mod() call.

Usage

model_config(
  msmMode = 0,
  initMode = "NonEquilibrium",
  avgnMode = 0,
  suitMode = 0,
  niter = 3,
  HCR = build_hcr(),
  recFun = build_srr(),
  M1Fun = build_M1(),
  growthFun = build_growth(),
  qFun = build_catchability(),
  selFun = build_selectivity(),
  compFun = build_composition()
)

Arguments

Details

The defaults are exactly fit_mod()'s own argument defaults, so a data list holding model_config() fits identically to one with no slot at all. When a data list has a model_config, fit_mod() reads each field only for arguments the caller did not pass explicitly, an argument passed to fit_mod() always wins, even when passed at its default. Omit an argument to let the stored configuration take effect for that field.

Value

An object of class "Rceattle_model_config" (a named list) to be stored as data_list$model_config.

Persistence

The configuration is code-side model structure, not one of the workbook data sheets, so it is not written by write_data() and does not survive an xlsx round-trip, build_data(base = x, model_config = cfg) piped through write_data() then read_data() returns without the slot. Re-attach it in code, store it alongside the data, or persist it as a documented, git-diffable YAML with save_config() / load_config() and apply it with fit_mod(data_list, config = load_config("run.yaml")).

See Also

fit_mod(), build_data(), save_config(), load_config().

Examples

cfg <- model_config(msmMode = 1, initMode = "NonEquilibrium")
dat <- build_data(base = BS2017MS, model_config = cfg)
# fit_mod(dat) would then fit as multispecies without passing msmMode.

Description

Management strategy evaluation performance metric summary

Usage

mse_summary(mse, om_only = FALSE)

Arguments

Value

A named list, one element per entity dimension so each metric lives only where it applies (no NA padding):

• species, a data.frame with one row per species (keyed by Species) of the conservation / status metrics: per-species average catch, catch inter-annual variability (IAV), and P(Closed) = the probability the fishery is closed (catch ~ 0); the relative mean-squared error of terminal and average SSB; the estimation-model- and operating-model- perceived overfishing / overfished probabilities (via build_hcr) and the probability each status is misclassified (estimation model disagrees with the operating-model truth); and terminal biomass, SSB, dynamic SB0, SSB depletion (equilibrium and dynamic), average SSB depletion, and SSB-collapse counts.

• fleet, a data.frame with one row per fishery fleet (keyed by Fleet_code / Fleet_name) of average catch, catch IAV, and P(Closed).

• total, a named numeric of the across-fleet total average catch and catch IAV.

• meta, run provenance (nsim, nspp, nflts, HCR, projection-year range).

All metrics are averaged across projection years and simulations.

Metric columns hold syntactic names, so they can be typed without backticks: avg_catch, catch_iav, p_closed, ssb_rmse_avg, ssb_rmse_terminal, em_p_overfishing, em_p_overfished, om_p_overfishing, om_p_overfished, p_overfishing_false_pos, p_overfishing_false_neg, p_overfished_false_pos, p_overfished_false_neg, om_terminal_biomass, om_terminal_ssb, om_terminal_dynamic_sb0, om_terminal_depletion, om_terminal_depletion_dynamic, om_avg_depletion, om_sims_collapsed, om_no_f_sims_collapsed, om_sims_collapsed_from_f. An om_ prefix is the operating model's truth and em_ the estimation model's perception; the four ⁠*_false_pos⁠ / ⁠*_false_neg⁠ metrics are the probability the two disagree. The three ⁠*_sims_collapsed⁠ metrics are counts of simulations, not probabilities.

The three depletions divide by different references:

• om_terminal_depletion_dynamic: terminal SSB over the OM's own DynamicSB0, its history with no fishing on the stock-recruit curve with the realized recruitment deviations and, under predation, the suitability fitted in the hindcast. Defined for single- and multispecies OMs. A species with input numbers-at-age (estDynamics > 0) keeps those numbers in the dynamic run, so its dynamic depletion is near 1 and says nothing about fishing.

• om_terminal_depletion: terminal SSB over SB0, the equilibrium value for a single-species OM. Under msmMode > 0 it is MSSB0, SSB at the end of a no-fishing projection, which fit_mod() derives only when a harvest control rule is present; without one (HCR = "NoFishing") this is NA.

• om_avg_depletion: the OM's own ssb_depletion averaged over the projection. Under predation with no HCR that series divides by SSB in the last projection year, not by an unfished reference.

The collapse metrics compare the OM with OM_no_F, the OM refit with no fishing after the original operating model's endyr (remove_F()).

Each frame holds a "labels" attribute mapping those names to the long display strings (e.g. om_terminal_depletion_dynamic -> "OM: Terminal SSB Depletion (Dynamic)") for plots and tables: attr(summ$species, "labels").

Examples

## Not run: 
mse <- load_mse(dir = "mse_runs")
s <- mse_summary(mse)
s$species   # one row per species
s$fleet     # one row per fleet
s$total     # across the whole system
s$meta      # nsim, nspp, nflts, HCR

## End(Not run)

Description

A data list containing inputs for a single-species CEATTLE model fit to Gulf of Alaska northern rockfish data through 2022. See BS2017SS for format details.

Usage

NorthernRockfish2022

Format

A list with the same structure as BS2017SS.

Description

Computes the Stewart and Monnahan (2025) statistical diagnostics for a set of OSA residuals: the standard deviation of the normalized residuals (SDNR) and the lower/upper tail statistics, each with the 95% interval expected under the standard normal null hypothesis (so departures can be judged objectively rather than by eye). Computed per data source (type x fleet) and overall.

Under a correctly specified model OSA residuals are already iid standard normal, so the SDNR is simply their sample standard deviation. Its null interval follows the chi-square result for the sample standard deviation of n standard normals (Francis 2014). Each tail statistic is the r-th order statistic of the residuals, at index round(p * (n + 1)), and its null interval is exact: the r-th order statistic of n uniforms is Beta(r, n - r + 1). Nothing is simulated.

Usage

osa_diagnostics(osa, nsim = 10000, probs = c(0.025, 0.975), seed = 123)

Arguments

Value

A data frame (class "rceattle_osa_diagnostics", so it prints as a compact severity-tagged summary; every column is still there and $ works as before) with one row per data source plus an "all" row, with columns: group (the "<source> fleet <n>" label), source, fleet, n, sdnr, sdnr_lo, sdnr_hi, lower, lower_lo, lower_hi, upper, upper_lo, upper_hi, the order statistic each tail was read at (lower_r, upper_r) and its exact nominal probability r/(n + 1) (lower_p, upper_p, which differ from probs at small n), and the logical flags sdnr_ok, lower_ok, upper_ok (TRUE when the statistic is inside its null interval). On the "all" row source and fleet are NA.

References

Francis, R.I.C.C. 2014. Replacing the multinomial in stock assessment models: a first step. Fish. Res. 151:70-84.

Stewart, I.J., and Monnahan, C.C. 2025. Can. J. Fish. Aquat. Sci. 82:1-13.

See Also

osa_residuals()

Description

One-step-ahead residuals, also called forecast or quantile residuals (Thygesen et al. 2017), computed via TMB::oneStepPredict(). Unlike Pearson residuals they are iid standard normal under a correctly specified model even when observations are correlated through composition bins or the model holds random effects, so they support objective goodness-of-fit testing (Trijoulet et al. 2023; Stewart and Monnahan 2025). The residualization runs inside the assessment, so it also accounts for correlation induced by the model's own random effects.

Usage

osa_residuals(
  object = NULL,
  source = c("ecov", "index", "catch", "comp", "caal"),
  method = "oneStepGaussianOffMode",
  discrete = NULL,
  parallel = TRUE,
  seed = 123,
  trace = FALSE,
  ...,
  fit = NULL
)

Arguments

Details

Computed post hoc and expensive – TMB re-optimizes the random effects as each observation is added – so fit_mod() does not produce them, and the fit must have been optimized at estimateMode < 3. Composition data are decomposed into univariate conditional residuals (binomial / beta-binomial; Trijoulet et al. 2023), and the last bin of each composition is fixed by the sum-to-N constraint, so it has no residual (NA).

Choosing between the methods, what each costs, how they behave under random effects, the ceiling on ⁠|residual|⁠ under "cdf", and where compositions at scale defeat it are in vignette("model-diagnostics"), with the measurements behind each.

Value

A data frame of class rceattle_osa, one row per residualized observation, with columns source, fleet, fleet_name, species, sex, year, age_length_bin, accumulated (the bin folds neighbouring ages, so it stands for a range rather than the age named), length (the conditioning length bin for caal, NA otherwise), index_label, observed, predicted, sd and residual.

observed and predicted are on the residualization scale: log for lognormal catch and index, natural scale for "Normal" and "TruncatedNormal", whitened by the lower Cholesky of the survey covariance for "MVN"/"MVNORM", and bin counts for compositions. Three cases read differently, and none is an error:

• Under method = "cdf", predicted is NA on every row because the method forms no conditional mode. sd is NA under the default method too; only "oneStepGaussian" returns one.

• On a "TruncatedNormal" fleet under a Gaussian method, predicted is the truncated conditional mean E[x | x > 0], which sits above the fitted index (163.8 against a fitted 100.0 on the package fixture), and sd is NA. The residual is unaffected.

• A composition predicted can go slightly negative where a bin holds almost no fish, and the function warns, naming the count and the years. The bin's own count drives it, not the composition's sample size – median observed count 0.05 on those rows against 4.9 elsewhere – so the warning is not by itself evidence of thin data. Those residuals are biased positive.

Adding a "TruncatedNormal" fleet can move the other fleets' residuals on a random-effects model, because its rows interleave with the other index fleets by year instead of forming a contiguous block: a contiguous group reproduces a single call to 1.5e-14, while an interleaved split moved a residual by 5.8e-2. Both orderings are valid probability-integral-transform sequences, so neither is wrong. Fixed-effect models are unaffected.

Attributes carry seed and method – the string that was passed, or a named vector where a group of rows ran under its own, such as TruncatedNormal = "oneStepGeneric". Where composition types are present a "pearson" attribute holds the matching Pearson residuals, on proportions with the sample size in sample_size rather than the bin counts above, so the two are not comparable column by column.

References

Thygesen, U.H., et al. 2017. Validation of ecological state space models using the Laplace approximation. Environ. Ecol. Stat. 24:317-339.

Trijoulet, V., et al. 2023. Model validation for compositional data in stock assessment models. Fish. Res. 257:106487.

Stewart, I.J., and Monnahan, C.C. 2025. Diagnosing common sources of lack of fit to composition data using one-step-ahead residuals. Can. J. Fish. Aquat. Sci. 82:1-13.

See Also

osa_diagnostics(), plot.rceattle_osa(), process_residuals(), vignette("model-diagnostics")

Examples

## Not run: 
data(BS2017SS)
fit <- fit_mod(BS2017SS, estimateMode = "Hindcast")
osa <- osa_residuals(fit, source = c("index", "comp"))
plot(osa)

## End(Not run)

Description

The TMB parameter vector uses transformed, abbreviated names, log_M1, R_log_sd, index_log_q. This returns the table mapping each one to the quantity it represents on its natural scale, which process it belongs to, what it means, and its dimensions.

Usage

parameter_dictionary(internal = NULL, process = NULL)

Arguments

Details

A name beginning log_, or ending ⁠_log_sd⁠, is estimated on the log scale, so exp() turns it into the quantity named in the natural column. Most fitted values are already back-transformed in fitted$quantities, so this is mainly needed when reading estimated_params or an sdreport() directly.

Dimensions use the model's own notation: nspp species, nsex sexes, nages age bins, nyrs years including the projection, nyrs_hind hindcast years only, n_flt fleets and n_sel selectivity blocks. n_sel equals n_flt: fleets sharing a Selectivity_index share one block through the map, not through a narrower array.

Value

A data frame with one row per parameter and columns internal, natural, process, meaning and dims.

See Also

fit_mod() for the fitted object these names come from.

Examples

# The whole table
head(parameter_dictionary())

# What is R_log_sd?
parameter_dictionary("R_log_sd")

# Everything governing selectivity
parameter_dictionary(process = "selectivity")

# Search the meanings by keyword
dict <- parameter_dictionary()
dict[grep("catchability", dict$meaning, ignore.case = TRUE), ]

Description

One row per estimated (fixed-effect) parameter, giving the fleet, species, sex, age, bin, year or slot it occupies.

Usage

parameter_index(object)

Arguments

Details

TMB names each element of obj$par after its parameter block, so a convergence diagnostic can report that sel_coff_dev is non-identifiable but not which of its elements.

The mapping is recovered by pushing a tagged parameter vector back through TMB's own parList(), so it reflects build_map() exactly: parameters mapped off are absent, and fleets sharing a Selectivity_index or Catchability_index appear once, as the single parameter they are, with n_cells recording how many array cells they drive.

Ages are absolute. init_dev holds one deviation per age from minage + 1 upward, since age minage in the first year is recruitment.

Coordinate columns are NA where the axis does not apply to a block. The linkage and environmental-covariate blocks (beta_linkage, M1_beta, ...) are indexed by a linkage-table row rather than a model coordinate, so they hold no coordinates and are reported by element number.

A selectivity slot is named from the fleet's Selectivity: slot 2 of sel_inf is a descending inflection for the double-logistic family, the right-tail floor under DoubleNormal, and the age-1 selectivity under LogisticPM.

Under estimateMode = "Estimate" with any HCR but "NoFishing", object$obj is the projection object, whose only free parameters are log_Ftarget / log_Flimit, so the index describes those, not the hindcast. fit$convergence reads the hindcast index that fit_mod() captured before the projection remapped the model, and checks it against the vector it is labelling either way.

Value

A data frame with one row per estimated parameter and columns par_index (position in obj$par), block, species, fleet, sex, age, bin, year, slot, n_cells and label.

See Also

parameter_dictionary() for what each block means.

Description

Total biomass of each species consumed by its predators, summed over sex and age. For an MSE object the projection is summarized across simulations with a mean and a 95% interval.

Usage

plot_b_eaten(
  Rceattle,
  file = NULL,
  model_names = NULL,
  line_col = NULL,
  species = NULL,
  spnames = NULL,
  add_ci = FALSE,
  lwd = 3,
  save = FALSE,
  right_adj = 0,
  width = 7,
  height = 6.5,
  minyr = NULL,
  incl_proj = FALSE,
  mod_cex = 1,
  alpha = 0.4,
  mod_avg = rep(FALSE, length(Rceattle)),
  mse = FALSE,
  OM = TRUE,
  maxyr = NULL,
  lty = 1,
  incl_mean = FALSE,
  top_adj = 0.15
)

Arguments

Description

Function that plots the biomass consumed trends as estimated from Rceattle. Returns and saves a figure with the biomass eaten trajectory.

Usage

plot_b_eaten_prop(
  Rceattle,
  file = NULL,
  model_names = NULL,
  line_col = NULL,
  spnames = NULL,
  species = NULL,
  lwd = 3,
  right_adj = 0,
  top_adj = 0.15,
  minyr = NULL,
  mohns = NULL,
  width = 7,
  height = 6.5,
  incl_proj = FALSE,
  incl_mean = FALSE,
  add_ci = FALSE,
  mod_cex = 1,
  maxyr = NULL,
  lty = 1
)

Arguments

Details

Colour separates predators; line type separates models. Panels are prey species.

Description

Plots the mean minage+ biomass (million mt) and 95\

Usage

plot_biomass(
  Rceattle,
  file = NULL,
  model_names = NULL,
  line_col = NULL,
  species = NULL,
  spnames = NULL,
  add_ci = FALSE,
  lwd = 3,
  save = FALSE,
  right_adj = 0,
  legend.pos = "topright",
  width = 7,
  height = 6.5,
  minyr = NULL,
  maxyr = NULL,
  incl_proj = FALSE,
  mod_cex = 1,
  lty = rep(1, length(Rceattle)),
  alpha = 0.4,
  mod_avg = rep(FALSE, length(Rceattle)),
  mse = FALSE,
  OM = TRUE,
  reference = NULL,
  ylab = NULL
)

Arguments

Value

Returns and saves a figure with the biomass trajectory.

Description

Plots fitted fishery catch: observed points with lognormal 95% intervals and the model-predicted catch, faceted by fishery fleet. For MSE objects the projection period is summarized with 50% / 95% ribbons across simulations.

Usage

plot_catch(
  Rceattle,
  file = NULL,
  model_names = NULL,
  species = NULL,
  incl_proj = FALSE,
  width = 7,
  height = 6.5,
  error = TRUE,
  maxyr = NULL,
  mse = FALSE,
  fleets = NULL,
  line_col = NULL,
  right_adj = 0,
  top_adj = 1.2,
  single.plots = FALSE,
  alpha = 0.4,
  lwd = 2,
  ymax = NULL
)

Arguments

Details

The predicted line is the MEDIAN, as in plot_index(): catch_hat * exp(-bias_adjust_obs * sd^2 / 2). The catch likelihood is lognormal for every fishery, so this applies throughout.

Value

A ggplot object.

Description

Plots the fitted survey catchability q by year, faceted by fleet. Every fleet holding index_data is drawn, a fishery with a CPUE series as much as a survey.

Usage

plot_catchability(
  Rceattle,
  file = NULL,
  model_names = NULL,
  species = NULL,
  spnames = NULL,
  mse = FALSE,
  minyr = NULL,
  maxyr = NULL,
  incl_mean = FALSE,
  log = FALSE,
  line_col = NULL,
  lwd = 3,
  lty = 1,
  width = 7,
  height = 6.5,
  right_adj = 0,
  top_adj = 0.05,
  single.plots = FALSE
)

Arguments

Details

The value plotted is the realized catchability the model scaled the survey by, whichever route produced it: an estimated or fixed mean, annual deviations under Time_varying_q, an environmental regression (Catchability = "Environmental"), a q linkage such as rw(1 | Year) or ar1(1 | Year), or the closed-form value under "Analytical" / "AnalyticalArith". A fleet with a time-invariant catchability draws a flat line, which is the honest picture of a constant q rather than an empty panel.

Hindcast years only. The model applies catchability over the hindcast and does not project it, so the series stops at endyr regardless of the projection horizon. There is no incl_proj.

line_col and lty separate the models, so overlaying several fits compares their catchability series panel by panel.

Value

A ggplot object.

Examples

data(BS2017SS)
fit <- Rceattle::fit_mod(data_list = BS2017SS, msmMode = 0)
plot_catchability(fit)
plot_catchability(fit, log = TRUE)