This tutorial starts off assuming that users are familiar with the setup with fishery selectivity.
Unlike fleet selectivity, index selectivity is unique to each series and no dummy fleets are used.
First, we consider the argument vector s_selectivity
,
which defines where the index selectivity is defined (the terms index
and survey are used interchangeably). Index selectivity may be identical
to fleet selectivity, i.e., vulnerable biomass, or could be related to
spawning or total biomass. If we have 4 indices with:
s_selectivity <- c("SSB", "B", 1, 2)
The first index is an index of spawning biomass as denoted by “SSB”
(maturity is configured in the OM
), the second index tracks
total biomass as denoted by “B” (selectivity = 1 for all ages), and the
third and fourth indices have the selectivity of the first and second
fleets, respectively. Integers for fleets refer to the true fleet and
not to selectivity blocks/dummy fleets.
In these cases, index selectivity is already defined somewhere else
in the model, and no further consideration of index selectivity is
needed. The RCM
function call can look like this:
output <- RCM(OM, RCMdata, selectivity = selectivity,
s_selectivity = s_selectivity,
...)
On the other hand, if index selectivity needs to be explicitly
defined, then the s_selectivity
vector can indicate the
functional form, using one of logistic
, dome
,
or free
. Let’s look at another situation with 5 indices
with this selectivity definition:
s_selectivity <- c(2, "SSB", "B", "dome", "free")
For the fourth and fifth index, the selectivity functions are dome-shaped and free parameters, respectively. We will need to explicitly consider what the parameters defining this functions are, either as starting values to be estimated or fixed values in the model.
Just as with the fleet selectivity parameters, the index selectivity
parameters by default use OM@LFS
, OM@L5
, and
OM@Vmaxlen
for start values when
s_selectivity = "logistic"
or "dome"
. Custom
start values are needed when index selectivity uses free parameters.
Custom start values are passed to the RCM in the
start$ivul_par
matrix with the same dimension and layout as
for start$vul_par
:
OM@maxage <- 5
s_selectivity <- c(2, "SSB", "B", "dome", "free")
ivul_par
## [,1] [,2] [,3] [,4] [,5]
## [1,] 0 0 0 55.0 1
## [2,] 0 0 0 40.0 0
## [3,] 0 0 0 0.5 0
## [4,] 0 0 0 0.0 0
## [5,] 0 0 0 0.0 0
## [6,] 0 0 0 0.0 0
The parameter slots in the first three columns (indices 1-3) are
ignored. Recall that selectivity is defined elsewhere (through fishery
selectivity, maturity, and total biomass, respectively) and the model
will ignore these columns in ivul_par
. Any value can be
placed here, but users must ensure that the columns in the corresponding
map
matrix are set to NA so that TMB will turn off
estimation of these parameters.
In column four, the first three rows are the start values for the three parameters of the dome function. The last three rows are not used and can be any value. In column five, we set up the selectivity only for age-0 animals, i.e., an index of recruits.
Recall that the map
argument tells TMB what to do to the
parameter of the corresponding row and column in ivul_par
.
Shared parameters are assigned a unique integer amongst themselves while
fixed or unused parameters are assigned NA
.
For the map
matrix, we do the following things:
NA
because index selectivity is defined
outside the ivul_par
matrix. TMB turns off estimation of
these parameters.ivul_par
are not used and the map should be set to
NA
.NA
because we fix the selectivity at
age schedule.map_ivul_par
## [,1] [,2] [,3] [,4] [,5]
## [1,] NA NA NA 1 NA
## [2,] NA NA NA 2 NA
## [3,] NA NA NA 3 NA
## [4,] NA NA NA NA NA
## [5,] NA NA NA NA NA
## [6,] NA NA NA NA NA
A function call utilizing this custom set-up for index selectivity would be:
output <- RCM(OM, data, selectivity = selectivity,
s_selectivity = s_selectivity,
start = list(ivul_par = ivul_par),
map = list(ivul_par = map_ivul_par),
...)
By default, each index series is time-invariant with respect to the selectivity and catchability (\(q\)). By default, catchability for index \(s\) is solved analytically in the model,
\[ q_s = \exp\left(\dfrac{\sum_y \log(I^{\textrm{obs}}_{y,s}) - \sum_y \log(x_{y,s})}{n_s}\right) \]
where \(I\) is the observed index, \(x\) is the biomass or abundance vulnerable to the index, and \(n_s\) is the number of data points.
This section goes through the three use cases if one or both of selectivity and catchability changes.
Simply break the index into two series and model them separately.
Break the index into two series and map the selectivity parameters so that the two series share the same parameters (see selectivity section above).
Break the index into two series. Then use the map argument to identify the index series which share a common \(q\). As an example with 5 index series,
RCM(..., map = list(q = c(1, 1, 2, 2, NA)))
This map argument sets up the model so that index series 1 & 2 share one \(q\) parameter and series 3 and 4 share another \(q\) parameter. When map is an integer, \(q\) is an estimated parameter while NA for the fifth index indicates that \(q\) will be calculated analytically.
RCM provides an option for fixing catchability to one:
RCMdata@abs_I <- c(0, 0, 0, 0, 1)
where \(q = 1\) for the fifth index.
There are two ways to fix \(q\) to an alternative value. First, rescale the index such that the new index will have \(q=1\). Second, set the q prior to with a tight prior:
prior <- list()
prior$q <- matrix(NA, 5, 2)
prior$q[5, ] <- c(0.5, 0.001)
RCM(..., prior = prior)
where the \(q\) for the fifth index has a prior mean of 0.5 and lognormal standard deviation of 0.001.
If two index series have offset catchability, e.g., the \(q\) of one index is twice of another index, simply re-scale one of them and share the \(q\) parameter between the two.