SCA: Statistical catch-at-age (SCA) model

In tcarruth/MSEtool: Management Strategy Evaluation Toolkit

Description

A generic statistical catch-at-age model (single fleet, single season) that uses catch, index, and catch-at-age composition data. SCA parameterizes R0 and steepness as leading productivity parameters in the assessment model. Recruitment is estimated as deviations from the resulting stock-recruit relationship. In SCA2, the mean recruitment in the time series is estimated and recruitment deviations around this mean are estimated as penalized parameters (similar to Cadigan 2016). The standard deviation is set high so that the recruitment is almost like free parameters. Unfished and MSY reference points are inferred afterwards from the assessment output (SSB and recruitment estimates). SCA_Pope is a variant of SCA that fixes the expected catch to the observed catch, and Pope's approximation is used to calculate the annual harvest rate (U).

Usage

SCA(
  x = 1,
  Data,
  SR = c("BH", "Ricker"),
  vulnerability = c("logistic", "dome"),
  CAA_dist = c("multinomial", "lognormal"),
  CAA_multiplier = 50,
  I_type = c("B", "VB", "SSB"),
  rescale = "mean1",
  max_age = Data@MaxAge,
  start = NULL,
  fix_h = TRUE,
  fix_F_equilibrium = TRUE,
  fix_omega = TRUE,
  fix_sigma = FALSE,
  fix_tau = TRUE,
  early_dev = c("comp_onegen", "comp", "all"),
  late_dev = "comp50",
  integrate = FALSE,
  silent = TRUE,
  opt_hess = FALSE,
  n_restart = ifelse(opt_hess, 0, 1),
  control = list(iter.max = 2e+05, eval.max = 4e+05),
  inner.control = list(),
  ...
)

SCA2(
  x = 1,
  Data,
  SR = c("BH", "Ricker"),
  vulnerability = c("logistic", "dome"),
  CAA_dist = c("multinomial", "lognormal"),
  CAA_multiplier = 50,
  I_type = c("B", "VB", "SSB"),
  rescale = "mean1",
  max_age = Data@MaxAge,
  start = NULL,
  fix_h = TRUE,
  fix_F_equilibrium = TRUE,
  fix_omega = TRUE,
  fix_sigma = FALSE,
  fix_tau = TRUE,
  common_dev = "comp50",
  integrate = FALSE,
  silent = TRUE,
  opt_hess = FALSE,
  n_restart = ifelse(opt_hess, 0, 1),
  control = list(iter.max = 2e+05, eval.max = 4e+05),
  inner.control = list(),
  ...
)

SCA_Pope(
  x = 1,
  Data,
  SR = c("BH", "Ricker"),
  vulnerability = c("logistic", "dome"),
  CAA_dist = c("multinomial", "lognormal"),
  CAA_multiplier = 50,
  I_type = c("B", "VB", "SSB"),
  rescale = "mean1",
  max_age = Data@MaxAge,
  start = NULL,
  fix_h = TRUE,
  fix_U_equilibrium = TRUE,
  fix_sigma = FALSE,
  fix_tau = TRUE,
  early_dev = c("comp_onegen", "comp", "all"),
  late_dev = "comp50",
  integrate = FALSE,
  silent = TRUE,
  opt_hess = FALSE,
  n_restart = ifelse(opt_hess, 0, 1),
  control = list(iter.max = 2e+05, eval.max = 4e+05),
  inner.control = list(),
  ...
)

Arguments

x	A position in the Data object (by default, equal to one for assessments).
Data	An object of class Data
SR	Stock-recruit function (either "BH" for Beverton-Holt or "Ricker").
vulnerability	Whether estimated vulnerability is "logistic" or "dome" (double-normal). See details for parameterization.
CAA_dist	Whether a multinomial or lognormal distribution is used for likelihood of the catch-at-age matrix. See details.
CAA_multiplier	Numeric for data weighting of catch-at-age matrix if CAA_hist = "multinomial". Otherwise ignored. See details.
I_type	Whether the index surveys population biomass (B; this is the default in the DLMtool operating model), vulnerable biomass (VB), or spawning stock biomass (SSB).
rescale	A multiplicative factor that rescales the catch in the assessment model, which can improve convergence. By default, "mean1" scales the catch so that time series mean is 1, otherwise a numeric. Output is re-converted back to original units.
max_age	Integer, the maximum age (plus-group) in the model.
start	Optional list of starting values. Entries can be expressions that are evaluated in the function. See details.
fix_h	Logical, whether to fix steepness to value in Data@steep in the model for SCA. This only affects calculation of reference points for SCA2.
fix_F_equilibrium	Logical, whether the equilibrium fishing mortality prior to the first year of the model is estimated. If TRUE, F_equilibrium is fixed to value provided in start (if provided), otherwise, equal to zero (assumes unfished conditions).
fix_omega	Logical, whether the standard deviation of the catch is fixed. If TRUE, sigma is fixed to value provided in start (if provided), otherwise, value based on Data@CV_Cat.
fix_sigma	Logical, whether the standard deviation of the index is fixed. If TRUE, sigma is fixed to value provided in start (if provided), otherwise, value based on Data@CV_Ind.
fix_tau	Logical, the standard deviation of the recruitment deviations is fixed. If TRUE, tau is fixed to value provided in start (if provided), otherwise, value based on Data@sigmaR.
early_dev	Numeric or character string describing the years for which recruitment deviations are estimated in SCA. By default, equal to "comp_onegen", where rec devs are estimated one full generation prior to the first year when catch-at-age (CAA) data are available. With "comp", rec devs are estimated starting in the first year with CAA. With "all", rec devs start at the beginning of the model. If numeric, the number of years after the first year of the model for which to start estimating rec devs. Use negative numbers for years prior to the first year.
late_dev	Typically, a numeric for the number of most recent years in which recruitment deviations will not be estimated in SCA (recruitment in these years will be based on the mean predicted by stock-recruit relationship). By default, "comp50" uses the number of ages (smaller than the mode) for which the catch-at-age matrix has less than half the abundance than that at the mode.
integrate	Logical, whether the likelihood of the model integrates over the likelihood of the recruitment deviations (thus, treating it as a random effects/state-space variable). Otherwise, recruitment deviations are penalized parameters.
silent	Logical, passed to MakeADFun, whether TMB will print trace information during optimization. Used for dignostics for model convergence.
opt_hess	Logical, whether the hessian function will be passed to nlminb during optimization (this generally reduces the number of iterations to convergence, but is memory and time intensive and does not guarantee an increase in convergence rate). Ignored if integrate = TRUE.
n_restart	The number of restarts (calls to nlminb) in the optimization procedure, so long as the model hasn't converged. The optimization continues from the parameters from the previous (re)start.
control	A named list of agruments for optimization to be passed to nlminb.
inner.control	A named list of arguments for optimization of the random effects, which is passed on to newton.
...	Other arguments to be passed.
common_dev	Typically, a numeric for the number of most recent years in which a common recruitment deviation will be estimated (in SCA2, uninformative years will have a recruitment closer to the mean, which can be very misleading, especially near the end of the time series). By default, "comp50" uses the number of ages (smaller than the mode) for which the catch-at-age matrix has less than half the abundance than that at the mode.
fix_U_equilibrium	Logical, same as 'fix_F_equilibrium' for 'SCA_Pope'.

Details

The basic data inputs are catch (by weight), index (by weight/biomass), and catch-at-age matrix (by numbers). Annual F's are estimated parameters assuming continuous fishing over the year. Note: prior to version 1.2, catches were assumed to be known perfectly with an annual harvest rate from pulse fishing in SCA. That feature has now moved to SCA_Pope.

By default, steepness is fixed in the model to the value in Data@steep.

The annual sample sizes of the catch-at-age matrix is provided to the model (used in the likelihood for catch-at-age assuming a multinomial distribution), and is manipulated via argument CAA_multiplier. This argument is interpreted in two different ways depending on the value provided. If CAA_multiplier > 1, then this value will cap the annual sample sizes to that number. If CAA_multiplier <= 1, then all the annual samples sizes will be re-scaled by that number. By default, sample sizes are capped at 50.

Alternatively, a lognormal distribution with inverse proportion variance can be used for the catch at age (Punt and Kennedy, 1994, as cited by Maunder 2011).

For start (optional), a named list of starting values of estimates can be provided for:

• R0 Virgin recruitment, only for SCA.

• h Steepness, only for SCA. If not provided, the value in Data@steep is used.

• meanR Mean recruitment, only for SCA2.

• F_equilibrium Fishing mortality prior to the first year of model, e.g. zero means unfished conditions. Defaults to zero.

• vul_par Vulnerability parameters (length 2 vector for logistic or length 4 for dome, see below). Users should provide estimates of the parameters in normal space, e.g. vul_max between 0-1, and the function will perform the appropriate transformations for the model.

• F A vector of F's of length nyears, length(Data@Year). If not provided, defaults to 0.1.

• omega Standard deviation of catch. If not provided, the value based on Data@CV_Cat is used.

• sigma Standard deviation of index. If not provided, the value based on Data@CV_Ind is used.

• tau Standard deviation of recruitment deviations. If not provided, the value in Data@sigmaR is used.

Vulnerability can be specified to be either logistic or dome. If logistic, then the parameter vector vul_par is of length 2:

• vul_par[1]: a_95, the age of 95% vulnerability, via logit transformation to constrain a_95 to less than 75% of the maximum age: a_95 = 0.75 * max_age * plogis(vul_par[1]).

• vul_par[2]: a_50, the age of 50% vulnerability as an offset, i.e., a_50 = a_95 - exp(vul_par[2]).

A vague prior for vul_par[2] ~ N(0, sd = 3) is used to aid convergence, for example, when vulnerability >> 0.5 for the youngest age class.

With dome vulnerability, a double Gaussian parameterization is used, where vul_par is an estimated vector of length 4:

• vul_par[1]: a_asc, the first age of full vulnerability for the ascending limb, via logit transformation to constrain a_95 to less than 75% of the maximum age: a_asc = 0.75 * maxage * plogis(vul_par[1]).

• vul_par[2]: a_50, the age of 50% vulnerability for the ascending limb as an offset, i.e., a_50 = a_asc - exp(vul_par[2]).

• vul_par[3]: a_des, the last age of full vulnerability (where the descending limb starts) via logit transformation to constrain between a_asc and max_age, i.e., a_des = (max_age - a_asc) * plogis(vul_par[3]) + a_asc. By default, fixed to a small value so that the dome is effectively a three-parameter function.

• vul_par[4]: vul_max, the vulnerability (in logit space) at the maximum age.

Vague priors of vul_par[2] ~ N(0, sd = 3) and vul_par[3] ~ N(0, 3) are used to aid convergence, for example, when vulnerability >> 0.5 for the youngest age class.

Value

An object of class Assessment.

Required Data

• SCA, SCA_Pope, and SCA_Pope: Cat, Ind, Mort, L50, L95, CAA, vbK, vbLinf, vbt0, wla, wlb, MaxAge

Optional Data

• SCA: Rec, steep, sigmaR, CV_Ind, CV_Cat

• SC2: Rec, steep, CV_Ind, CV_Cat

• SCA_Pope: Rec, steep, sigmaR, CV_Ind

Author(s)

Q. Huynh

References

Cadigan, N.G. 2016. A state-space stock assessment model for northern cod, including under-reported catches and variable natural mortality rates. Canadian Journal of Fisheries and Aquatic Science 72:296-308.

Maunder, M.N. 2011. Review and evaluation of likelihood functions for composition data in stock-assessment models: Estimating the effective sample size. Fisheries Research 209:311-319.

Punt, A.E. and Kennedy, R.B. 1997. Population modelling of Tasmanian rock lobster, Jasus edwardsii, resources. Marine and Freshwater Research 48:967-980.

See Also

plot.Assessment summary.Assessment retrospective profile make_MP

Examples

res <- SCA(Data = DLMtool::SimulatedData)
res2 <- SCA2(Data = DLMtool::SimulatedData)

compare_models(res, res2)

SCA_assess <- SCA2(Data = DLMtool::Simulation_1)

## Not run: 
plot(res)

## End(Not run)

tcarruth/MSEtool documentation built on Oct. 19, 2020, 6:09 a.m.