The multiStockassessment
package extends the
stockassessment
to allow linked
stocks (See e.g. Albertsen, C. M., Nielsen, A., and Thygesen, U. H.
(2018) Connecting single-stock assessment models through correlated
survival. ICES Journal of Marine Science, 75(1), 235-244. doi:
10.1093/icesjms/fsx114).
The package can be installed with
devtools
:
devtools::install_github("calbertsen/multi_SAM", subdir = "multiStockassessment")
To use the multiStockassessment
package, the individual stocks must
first be fitted with the stockassessment
package and combined to a
samset
using the c(...)
function.
The code below loads the North Sea cod data from the stockassessment package, fits it, simulates two new data sets and fit those individually.
library(stockassessment)
data(nscodData)
data(nscodConf)
data(nscodParameters)
set.seed(9876)
fit <- sam.fit(nscodData, nscodConf, nscodParameters,sim.condRE = FALSE)
datSim <- simulate(fit,nsim=2)
fitSim <- do.call("c",lapply(datSim,function(x)sam.fit(x,nscodConf, nscodParameters)))
library(multiStockassessment)
After the simulated data have been fitted and combined by the c
function, the multiStockassessment
package can be used to fit a model
with correlated survival.
The suggestCorStructure
function can be used to set up the correlation
matrix. The nAgeClose
argument can be used to construct the banded
structures used in the paper. The result is a logical symmetric matrix
indicating whether a correlation should be fixed at zero. The dimension
of the matrix is the number of stocks times the number of ages in the
data.
For instance, suggestCorStructure(fitSim,nAgeClose=1)
creates a matrix
where ages less than 1 appart are correlated between the stocks:
cs <- suggestCorStructure(fitSim,nAgeClose=1)
cs
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
## [1,] TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE
## [4,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE
## [5,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE
## [6,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [7,] FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [8,] TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [9,] TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [10,] TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [11,] TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE
## [12,] TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE
## [,12]
## [1,] TRUE
## [2,] TRUE
## [3,] TRUE
## [4,] TRUE
## [5,] TRUE
## [6,] FALSE
## [7,] TRUE
## [8,] TRUE
## [9,] TRUE
## [10,] TRUE
## [11,] TRUE
## [12,] TRUE
The suggestCorStructure
function has several options to configure the
correlation structure, and the result can be modified by hand for
complete freedom.
The model is fitted by the multisam.fit
function. The function
requires a samset
and a correlaiton structure. By default, the
correlation structure given models the partial correlations between
cohorts. This can be changed to regular correlations by setting
usePartialCors = FALSE
.
obj <- multisam.fit(fitSim,cs)
obj
## Multi-SAM model with 2 stocks: log likelihood is -253.2214. Convergence OK.
Most plotting and table functions from the stockassessment
package
have a corresponding version in multiStockassessment
.
The fitted correlation structure can be plotted by
corplot(obj)
Other plotting functions implemented from the stockassessment
package
are
fbarplot(obj)
ssbplot(obj)
ssbplot(obj)
tsbplot(obj)
recplot(obj)
catchplot(obj)
srplot(obj)
fitplot(obj,1)
Note that the fitplot
function above requires the stock to be
specified.
fitplot(obj,2)
Models can be compared with the modeltable
function. For instance, the
model above can be compared to individual assessments by
cs2 <- suggestCorStructure(fitSim,nAgeClose=0)
obj2 <- multisam.fit(fitSim,cs2)
modeltable(obj,obj2)
## log(L) #par AIC Pval( M1 -> M2 )
## M1 -253.2214 74 654.4428 NA
## M2 -255.0715 68 646.1429 0.7171872