In tcarruth/MSEtool: Management Strategy Evaluation Toolkit

knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
knitr::opts_chunk$set(dpi=85)
options(width = 650)

Selectivity and mortality

For fleets, selectivity is defined by blocks (indexed by $b$) which can then be assigned to any fleet $f$ and year $y$ to allow for time-varying selectivity. By default, each fleet is assigned to its own block for all years (no time-varying selectivity).

For flat-topped selectivity in block $b$, two parameters are used and expressed in terms of length units: the length of 5% selectivity ($L^5_b$) and the length of full selectivity $L^{\textrm{FS}}b$. For dome selectivity, a third parameter, the selectivity at $L{\infty}$, $V^{L_{\infty}}_b$ is also used. Length-based selectivity is converted to age-based selectivity in the age-structured model as:

$$ v_{y,a,b} = \begin{cases} 2^{-[(L_{y,a} - L^{\textrm{FS}}b)/(\sigma^{\textrm{asc}}_b)]^2} & \textrm{if } L{y,a} < L^{\textrm{FS}}b\ 1 & \textrm{if logistic and } L{y,a} \ge L^{\textrm{FS}}b,\ 2^{-[(L{y,a} - L^{\textrm{FS}}b)/(\sigma^{\textrm{des}}_b)]^2} & \textrm{if dome and } L{y,a} \ge L^{\textrm{FS}}_b \end{cases} $$

where $L_{y,a}$ is the mean length-at-age and $\sigma^{\textrm{asc}}b = (L^5_b - L^{\textrm{FS}}_b)/\sqrt{-\log_2(0.05)}$ and $\sigma^{\textrm{des}}_b = (L{\infty} - L^{\textrm{FS}}b)/\sqrt{-\log_2(V^{L{\infty}}_b)}$ control the shape of the ascending and descending limbs, respectively, of the selectivity function. In this parameterization, length-based selectivity is constant over time. The corresponding age-based selectivity within each block is constant so long as growth is not time-varying.

Selectivity can also be parameterized where $v_{y,a,b}$ are free independent parameters. In this case, selectivity does not vary among years.

Total mortality $Z$ in year $y$ and for age $a$ is the sum of fishing mortality $F$ from all fleets and natural mortality $M$,

$$ Z_{y,a} = M_{y,a} + \Sigma_f v_{y,a,f} F_{y,f},$$

where $v_{y,a,f}$ is the fleet selectivity after assigning blocks to fleets.

Survey selectivity

Survey selectivity is constant over time and is denoted as $v_{a,s}$ with either logistic, dome, or free parameterizations.

Initial population distribution

The population age distribution in the first year of the model $y=1$ is in equilibrium where $$ N_{1,a} = \begin{cases} R^{\textrm{eq}} \exp(-\Sigma_{i=1}^{a-1}Z^{\textrm{eq}}i) & a = 1, \ldots, A-1\ \dfrac{R^{\textrm{eq}} \exp(-\Sigma{i=1}^{a-1}Z^{\textrm{eq}}i)}{1 - \exp(-Z^{\textrm{eq}}_A)} & a = A, \end{cases} $$ where the $R^{\textrm{eq}}$ is the equilibrium recruitment and $Z^{\textrm{eq}}_a = M{1,a} + \Sigma_f v_{1,a,f} F^{\textrm{eq}}_f$ is the equilibrium total mortality rate. Unfished conditions are modeled by setting $F^{\textrm{eq}}_f = 0$. To estimate $F^{\textrm{eq}}_f$, the corresponding equilibrium catch in weight $\tilde{C}^{\textrm{eq}}_f$ prior to the first year of the model should be provided. In the equilibrium yield curve, $F^{\textrm{eq}}_f$ would be the fishing mortality corresponding to fishing at $F^{\textrm{eq}}_f$. Once $Z^{\textrm{eq}}_a$ is obtained, then the equilibrium recruitment is calculated as:

$$ R^{\textrm{eq}} = \begin{cases} \dfrac{\alpha^{\textrm{B}}\phi^{\textrm{eq}} - 1}{\beta^{\textrm{B}}\phi^{\textrm{eq}}} & \textrm{if Beverton-Holt stock-recruit relationship}\ \dfrac{\log(\alpha^{\textrm{R}}\phi^{\textrm{eq}})}{\beta^{\textrm{R}}\phi^{\textrm{eq}}} & \textrm{if Ricker stock-recruit relationship} \end{cases}, $$ where $\phi^{\textrm{eq}}$ is the spawners-per-recruit when the mortality is $Z^{\textrm{eq}}_a$. From steepness $h$, $\alpha^{\textrm{B}} = \frac{4h}{(1-h)\phi_0}$, $\beta^{\textrm{B}} = \frac{5h-1}{(1-h)B^S_0}$, $\alpha^{\textrm{R}} = \frac{(5h)^{1.25}}{\phi_0}$, $\beta^{\textrm{R}} = \frac{\log(5h)}{B^S_0}$, where $\phi_0$ and $B^S_0$ are unfished spawners-per-recruit and unfished spawning biomass, respectively.

Dynamics equations

After setting the equilibrium population age distribution in the first year of the model, the population abundance $N_{y,a}$ in subsequent years is $$ N_{y,a} = \begin{cases} R_y & a = 1\ N_{y-1,a-1} \exp(-Z_{y-1,a-1}) & a = 2, \ldots, A - 1,\ N_{y-1,a-1} \exp(-Z_{y-1,a-1}) + N_{y-1,a} \exp(-Z_{y-1,a}) & a = A \end{cases} $$ where $R_y$ is the recruitment and $A$ is the maximum-age as the plus-group. Recruitment is modelled as $$ R_y = \begin{cases} \dfrac{\alpha^{\textrm{BH}} B^S_{y-1}}{1 + \beta^{\textrm{BH}}B^S_{y-1}} \exp(\delta_y - 0.5 \tau^2) & \textrm{if Beverton-Holt stock-recruit relationship}\ \alpha^{\textrm{Ricker}} B^S_{y-1} \exp(-\beta^{\textrm{Ricker}} B^S_{y-1})\exp(\delta_y - 0.5 \tau^2) & \textrm{if Ricker stock-recruit relationship} \end{cases}, $$ where $\delta_y$ are recruitment deviates and $\tau$ is the standard deviation of the deviates.

The spawning biomass is $B^S_y$ is $$B^S_y = \sum_a w_{y,a} m_{y,a} N_{y,a},$$ where $m_{y,a}$ and $w_{y,a}$ are the maturity at age and weight at age, respectively.

Catch at age

The catch (in numbers) $C^N$ at age for fleet $f$ is $$ C^N_{y,a,f} = \dfrac{v_{y,a,f} F_{y,f}}{Z_{y,a}} N_{y,a} (1 - \exp[-Z_{y,a}]).$$

If the model is conditioned on catch, then $F_{y,f}$ can be estimated as parameters or solved iteratively to match the observed catch. If the model is conditioned on effort, then $$ F_{y,f} = q_f E_{y,f},$$ where $E_{y,f}$ is the observed effort and $q_f$ is the scaling coefficient.

The proportion of the catch-at-age is $$ p_{y,a,f} = \dfrac{C^N_{y,a,f}}{\sum_a C^N_{y,a,f}}.$$

Catch-at-length

The catch at length is calculated assuming a normally distributed length-at-age $P(\ell,a)$, where $$ C^N_{y,\ell,f} = \sum_a C^N_{y,a,f} P(\ell|a) $$ and

$$ P(\ell|a) = \begin{cases} \phi(L'{\ell+1}) & \ell = 1\ \phi(L'{\ell+1}) - \phi(L'\ell) & \ell = 2, \ldots, L - 1,\ 1 -\phi(L'\ell) & \ell = L \end{cases} $$ with $L'{\ell}$ as the length at the lower boundary of length bin $\ell$ and $\phi(L'{\ell})$ as the cumulative distribution function of a normal variable with mean $\tilde{L}{y,a}$ (the expected mean length at age $a$) and standard deviation $\tilde{L}{y,a} \times CV^L$ ($CV^L$ is the coefficient of variation in mean length at age).

The proportion of the catch-at-length is $$ p_{y,\ell,f} = \dfrac{C^N_{y,\ell,f}}{\sum_{\ell}C^N_{y,\ell,f}}.$$

Catch in weight

The catch in weight $\tilde{C}$ is $$ \tilde{C}{y,f} = \sum_a C^N{y,a,f} w_{y,a}.$$

Mean size

The mean length of the catch $\bar{L}{y,f}$ is $$ \bar{L}{y,f} = \dfrac{\sum_{\ell} L_{\ell} C^N_{y,\ell,f}}{\sum_{\ell} C^N_{y,\ell,f}},$$ where $L_\ell$ is the midpoint of the length bin $\ell$.

The mean weight of the catch $\bar{w}{y,f}$ is $$ \bar{w}{y,f} = \dfrac{\sum_a C^N_{y,a,f}w_{y,a}}{\sum_a C^N_{y,a,f}},$$

Survey

If the $s^{\textrm{th}}$ survey is biomass-based, then the survey value $I_{y,s}$ is calculated as $$ I_{y,s} = q_s \sum_a v_{y,a,s} N_{y,a} w_{y,a}, $$ where $q$ is the scaling coefficient and $s$ indexes survey.

If the survey is abundance-based, then $$ I_{y,s} = q_s \sum_a v_{y,a,s} N_{y,a} . $$

The proportion-at-age vulnerable to the survey is $$ p_{y,a,s} = \dfrac{v_{a,s}N_{y,a}}{\sum_a v_{a,s}N_{y,a}}.$$

The proportion-at-length vulnerable to the survey is $$ p_{y,\ell,s} = \dfrac{\sum_a v_{a,s} N_{y,a} P(\ell|a)}{\sum_{\ell} \sum_a v_{a,s} N_{y,a} P(\ell|a)}.$$

Likelihoods

If the model is conditioned on catch and fishing mortality rates are estimated parameters, then the log-likelihood component $\Lambda_1$ of the catch is $$\Lambda_1 = \sum_f \left[\lambda^{\tilde{C}}f \sum_y \left(-\log(0.01) - \dfrac{[\log(\tilde{C}^{\textrm{obs}}{y,f}) - \log(\tilde{C}^{\textrm{pred}}_{y,f})]^2}{2 \times 0.01^2}\right)\right],$$

where $\textrm{obs}$ and $\textrm{pred}$ indicate observed and predicted quantities, respectively, and $\lambda$ are likelihood weights. With a small standard deviation for the catch likelihood relative to the variance in other likelihood components, the predicted catch should match the observed catch.

The log-likelihood component $\Lambda_2$ of survey data is $$\Lambda_2 = \sum_s \left[ \lambda^I_s \sum_y \left(-\log(\sigma_{y,s}) - \dfrac{[\log(I^{\textrm{obs}}{y,s}) - \log(I^{\textrm{pred}}{y,s})]^2}{2\sigma_{y,s}^2}\right) \right].$$

The log-likelihood component $\Lambda_3$ of catch-at-age data is $$\Lambda_3 = \sum_f \lambda^A_f \left[\sum_y O^A_{y,f} \sum_a p^{\textrm{obs}}{y,a,f} \log(p^{\textrm{pred}}{y,a,f})\right],$$ where $O^A$ is the annual sample sizes for the age compositions.

The log-likelihood component $\Lambda_4$ of catch-at-length data is $$\Lambda_4 = \sum_f \lambda^L_f \left[ \sum_y O^L_{y,f} \sum_{\ell} p^{\textrm{obs}}{y,\ell,f} \log(p^{\textrm{pred}}{y,\ell,f})\right]$$ where $O^L$ is the annual sample sizes for the length compositions.

The log-likelihood component $\Lambda_5$ for the observed mean size in this catch is:

$$\Lambda_5 = \sum_f \lambda^{\bar{L}}f\left[ \sum_y \left(-\log(\omega{y,f}) - \dfrac{[\bar{L}^{\textrm{obs}}{y,f} - \bar{L}^{\textrm{pred}}{y,f}]^2}{2 \omega^2_{y,f}}\right)\right],$$ for mean lengths, or

$$\Lambda_5 = \sum_f \lambda^{\bar{w}}f\left[ \sum_y \left(-\log(\omega{y,f}) - \dfrac{[\bar{w}^{\textrm{obs}}{y,f} - \bar{w}^{\textrm{pred}}{y,f}]^2}{2 \omega^2_{y,f}}\right)\right],$$ for mean weights, where $\omega_{y,f}$ is the standard deviation of mean size. With constant cofficient of variation (CV), $\omega_{y,f} = \bar{w}^{\textrm{obs}}_{y,f} CV^{\omega}_f$.

The log-likelihood component $\Lambda_6$ of annual estimated recruitment deviates $\delta_y$ in log space is $$\Lambda_6 = \Sigma_y\left(-\log(\tau) - \dfrac{\delta_y^2}{2 \tau^2}\right),$$ where $\tau$ is the standard deviation of recruitment deviates.

The log-likelihood component $\Lambda_7$ of the equilibrium catch is $$\Lambda_7 = \sum_f \lambda^{\tilde{C}}_f \left(-\log(0.01) - \dfrac{[\log(\tilde{C}^{\textrm{eq,obs}}_f) - \log(\tilde{C}^{\textrm{eq,pred}}_f)]^2}{2 \times 0.01^2}\right),$$

The log-likelihood component $\Lambda_8$ of the survey proportion-at-age data is $$\Lambda_8 = \sum_s \lambda^{IA}s \left[\sum_y O^{IA}{y,s} \sum_a p^{\textrm{obs}}{y,a,s} \log(p^{\textrm{pred}}{y,a,s})\right],$$ where $O^{IA}$ is the annual sample sizes for the survey age compositions.

The log-likelihood component $\Lambda_9$ of the survey proportion-at-length data is $$\Lambda_9 = \sum_s \lambda^{IL}s \left[ \sum_y O^{IL}{y,s} \sum_{\ell} p^{\textrm{obs}}{y,\ell,s} \log(p^{\textrm{pred}}{y,\ell,s})\right]$$ where $O^{IL}$ is the annual sample sizes for the survey length compositions.

The total log-likelihood $\textrm{LL}$ to be maximized is $$\textrm{LL} = \sum_{i=1}^9\Lambda_i.$$

Estimated parameters

The estimated parameters, denoted in this section as $x$, are unconstrained over all real numbers and then transformed in order to constrain the corresponding model parameters. For optimization, the transformation is also designed to reduce the scale of all estimated parameters to within an order of magnitude.

Selectivity

For a fleet block $b$ for which selectivity is estimated, then parameters $x^{LFS}_b$ and $x^{L5}_b$ are estimated over all real numbers, where

$$ \begin{align} L^{\textrm{FS}}b &= 0.99 \times L{\infty} \times \textrm{logit}^{-1}(x^{LFS}_b)\ L^5_b &= L^{\textrm{FS}}_b - \exp(x^{L5}_b) \end{align}$$

If a third parameter $x^{V}b$is estimated for dome selectivity, then $$ V^{L{\infty}}_b = \textrm{logit}^{-1}(x^V_b)$$

If selectivity is parameterized as free parameters, then $$ v_{y,a,f} = \textrm{logit}^{-1}(x^v_{y,a,f}).$$

For surveys, parameterizations are identical except with indexing for survey $s$.

Fishing mortality

If $F_{y,f}$ are estimated parameters (condition = "catch"), then one parameter $x^F_f$ is the estimated $F$ in log-space in the middle of the time series is estimated and all others are subsequent deviations, represented as $x^{F_{dev}}_{y,f}$:

$$ F_{y,f} = \begin{cases} \exp(x^F_f) & y \textrm{ is midpoint of the time series}\ \exp(x^F_f) \times \exp(x^{F_{dev}}_{y,f}) & \textrm{otherwise}\ \end{cases} $$

If condition = "effort", then $q_f$ is estimated in log space, where $$F_{y,f} = q_f E_{y,f} = \exp(x^q_f) \times E_{y,f}$$

Index catchability

To scale biomass to index values, the index catchability $q_s$ is solved analytically in the model:

$$ q_s = \exp\left(\dfrac{\sum_y \log(I^{\textrm{obs}}{y,s}) - \sum_y \log(\sum_a v{y,a,s}N_{y,a,s})}{n_s}\right),$$ or $$ q_s = \exp\left(\dfrac{\sum_y \log(I^{\textrm{obs}}{y,s}) - \sum_y \log(\sum_a v{y,a,s}N_{y,a,s}w_{y,a})}{n_s}\right),$$ for an abundance or biomass based index, respectively, where $n_s$ is the number of years with index values and the summation is over those $n_s$ years.

Other parameters

Unfished recruitment is estimated in log-space, $R_0 = \dfrac{1}{z}\exp(x^{R_0})$ where $z$ is an optional rescaler, e.g. mean historical catch, to reduce the magnitude of the $x^{R_0}$ estimate. Recruitment deviations $\delta_y$ are directly estimated.

Additional resources

The help file for SRA_scope will provide more information on the possible inputs for the model. The help file for the SRA class (obtained by typing class?SRA into the R console) describes the outputs from the function. An additional vignette is available to describe set up of fleet and survey selectivity in the function call.

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