| Title: | State-Space Assessment Model |
|---|---|
| Description: | Fitting SAM... |
| Authors: | Anders Nielsen [aut, cre], Casper Berg [aut], Christoffer Moesgaard Albertsen [aut], Kasper Kristensen [aut], Mollie Brooks [aut], Vanessa Trijoulet [aut], Olav Nikolai Breivik [aut] |
| Maintainer: | Anders Nielsen <[email protected]> |
| License: | GPL-2 |
| Version: | 0.12.0 |
| Built: | 2024-11-22 10:24:58 UTC |
| Source: | https://github.com/fishfollower/SAM |
SAM add forecasts
addforecast( fit, what, dotcol = "black", dotpch = 19, dotcex = 1.5, intervalcol = gray(0.5, alpha = 0.5), ... ) ## S3 method for class 'samforecast' addforecast( fit, what, dotcol = "black", dotpch = 19, dotcex = 1.5, intervalcol = gray(0.5, alpha = 0.5), ... )addforecast( fit, what, dotcol = "black", dotpch = 19, dotcex = 1.5, intervalcol = gray(0.5, alpha = 0.5), ... ) ## S3 method for class 'samforecast' addforecast( fit, what, dotcol = "black", dotpch = 19, dotcex = 1.5, intervalcol = gray(0.5, alpha = 0.5), ... )
fit |
the object returned from sam.fit |
what |
what to plot |
dotcol |
color for dot |
dotpch |
pch for dot |
dotcex |
cex for dot |
intervalcol |
color for interval |
... |
extra arguments not currently used |
internal plotting fun
Add stock-recruitment curve to srplot
addRecruitmentCurve( fit, CI = TRUE, col = rgb(0.6, 0, 0), cicol = rgb(0.6, 0, 0, 0.3), plot = TRUE, PI = FALSE, picol = rgb(0.6, 0, 0), pilty = 2, ... ) ## S3 method for class 'sam' addRecruitmentCurve( fit, CI = TRUE, col = rgb(0.6, 0, 0), cicol = rgb(0.6, 0, 0, 0.3), plot = TRUE, PI = FALSE, picol = rgb(0.6, 0, 0), pilty = 2, year = NA_real_, lastR = NA_real_, ... )addRecruitmentCurve( fit, CI = TRUE, col = rgb(0.6, 0, 0), cicol = rgb(0.6, 0, 0, 0.3), plot = TRUE, PI = FALSE, picol = rgb(0.6, 0, 0), pilty = 2, ... ) ## S3 method for class 'sam' addRecruitmentCurve( fit, CI = TRUE, col = rgb(0.6, 0, 0), cicol = rgb(0.6, 0, 0, 0.3), plot = TRUE, PI = FALSE, picol = rgb(0.6, 0, 0), pilty = 2, year = NA_real_, lastR = NA_real_, ... )
fit |
Object to show SR-curve for |
CI |
Add confidence intervals? |
col |
Color of fitted line |
cicol |
Color of confidence intervals |
plot |
Add the curve to a plot? |
PI |
Add prediction intervals? |
picol |
Color of prediction interval line |
pilty |
Line type of prediction interval line |
... |
not used |
year |
Show recruitment calculated conditional on this year (for recruitment functions that depend on year) |
lastR |
Show recruitment calculated conditional on this previous recruitment (for recruitment functions that depend on recruitment the previous year) |
srplot
SAM equilibrium biomass in the absence of fishing plot
b0plot(fit, ...) ## Default S3 method: b0plot(fit, ...) ## S3 method for class 'samforecast' b0plot(fit, ...) ## S3 method for class 'hcr' b0plot(fit, ...)b0plot(fit, ...) ## Default S3 method: b0plot(fit, ...) ## S3 method for class 'samforecast' b0plot(fit, ...) ## S3 method for class 'hcr' b0plot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of deterministic equilibrium biomass in the absence of fishing assuming biological parameters and selectivity for that year remains unchanged in the future.
B0 biomass table
b0table(fit, ...) ## Default S3 method: b0table(fit, ...)b0table(fit, ...) ## Default S3 method: b0table(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Spline basis for use with formula interface
bc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)bc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)
x |
Points to evaluate the basis in |
df |
Degrees of freedom |
knots |
Internal knots. If NULL, they are selected from quantiles of x. |
Boundary.knots |
Boundary knots. Defaults to range of x |
intercept |
Include an intercept in basis? |
A spline basis
Collect sam objects
## S3 method for class 'sam' c(...)## S3 method for class 'sam' c(...)
... |
one or more sam fits (as returned from the |
...
SAM catchbyfleet plot
catchbyfleetplot(fit, obs.show = FALSE, ...)catchbyfleetplot(fit, obs.show = FALSE, ...)
fit |
the object returned from sam.fit |
obs.show |
if observations are to be shown also |
... |
extra arguments transferred to plot |
Plot of estimated (and optionally observed) total catch in weight
CatchByFleet table
catchbyfleettable(fit, obs.show = FALSE)catchbyfleettable(fit, obs.show = FALSE)
fit |
object returned from sam.fit |
obs.show |
logical add a column with catch sum of product rowsums(C*W) |
...
SAM catch plot
catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'sam' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'samset' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'samforecast' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'hcr' catchplot(fit, obs.show = TRUE, drop = NULL, ...)catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'sam' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'samset' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'samforecast' catchplot(fit, obs.show = TRUE, drop = NULL, ...) ## S3 method for class 'hcr' catchplot(fit, obs.show = TRUE, drop = NULL, ...)
fit |
the object returned from sam.fit |
obs.show |
if observations are to be shown also |
drop |
number of years to be left unplotted at the end. Default (NULL) is to not show years at the end with no catch information |
... |
extra arguments transferred to plot including the following: |
Plot of estimated (and optionally observed) total catch in weight
Catch table
catchtable(fit, obs.show = FALSE, ...) ## S3 method for class 'sam' catchtable(fit, obs.show = FALSE, ...)catchtable(fit, obs.show = FALSE, ...) ## S3 method for class 'sam' catchtable(fit, obs.show = FALSE, ...)
fit |
object returned from sam.fit |
obs.show |
logical add a column with catch sum of product rowsums(C*W) |
... |
extra arguments not currently used |
...
Catch-at-age in numbers table
caytable(fit, fleet = which(fit$data$fleetTypes == 0))caytable(fit, fleet = which(fit$data$fleetTypes == 0))
fit |
a fitted object of class 'sam' as returned from sam.fit |
fleet |
the fleet number(s) to return catch summed for (default is to return the sum of all residual fleets). |
...
remove void catches
clean.void.catches(dat, conf)clean.void.catches(dat, conf)
dat |
data for the sam model as returned from the setup.sam.data function |
conf |
model configuration which can be set up using the |
an updated dataset without the catches where F is fixed to zero
Extract fixed coefficients of sam object
## S3 method for class 'sam' coef(object, ...)## S3 method for class 'sam' coef(object, ...)
object |
sam fitted object as returned from the |
... |
extra arguments |
fixed coefficients of sam object
Area plot of spawning components
componentplot(fit, ...) ## S3 method for class 'sam' componentplot( fit, onlyComponentYears = FALSE, ylab = "Composition", colSet = c("#332288", "#88CCEE", "#44AA99", "#117733", "#999933", "#DDCC77", "#661100", "#CC6677", "#882255", "#AA4499"), legend.pos = "bottom", bg = "white", ncol = length(cf), ... )componentplot(fit, ...) ## S3 method for class 'sam' componentplot( fit, onlyComponentYears = FALSE, ylab = "Composition", colSet = c("#332288", "#88CCEE", "#44AA99", "#117733", "#999933", "#DDCC77", "#661100", "#CC6677", "#882255", "#AA4499"), legend.pos = "bottom", bg = "white", ncol = length(cf), ... )
fit |
sam fit |
... |
passed to legend |
onlyComponentYears |
If true, x axis is limited to the range with spawning component data. Otherwise, the model years are used. |
ylab |
Label for y axis |
colSet |
Colors |
legend.pos |
Legend position. See ?legend |
bg |
Background of legend. See ?legend |
ncol |
Number of columns in legend. See ?legend |
Nothing
Christoffer Moesgaard Albertsen
Plots between-age correlations by fleet, either estimated or empirical using residuals.
corplot(x, ...) ## S3 method for class 'sam' corplot(x, ...) ## S3 method for class 'samres' corplot(x, ...)corplot(x, ...) ## S3 method for class 'sam' corplot(x, ...) ## S3 method for class 'samres' corplot(x, ...)
x |
Either a sam fit as returned by sam.fit OR the object returned from residuals.sam |
... |
extra arguments to plot |
Common function for plotting correlation matrices.
corplotcommon(x, fn, ...)corplotcommon(x, fn, ...)
x |
a list of correlation matrices |
fn |
a vector of fleet names |
... |
extra arguments to plotcorr |
SAM Data plot
dataplot(fit, col = NULL, fleet_type = NULL, fleet_names = NULL) ## S3 method for class 'sam' dataplot(fit, col = NULL, fleet_type = NULL, fleet_names = NULL)dataplot(fit, col = NULL, fleet_type = NULL, fleet_names = NULL) ## S3 method for class 'sam' dataplot(fit, col = NULL, fleet_type = NULL, fleet_names = NULL)
fit |
the object returned from sam.fit |
col |
color to use for each fleet, default is two sequential colors |
fleet_type |
character vector giving the type of data per fleet. The default uses fit$data$fleetTypes as follows: |
fleet_names |
character vector giving fleet names. The default is given by attr(fit$data,"fleetNames") |
Plot data available for the stock
Setup basic minimal configuration for sam assessment
defcon(dat, level = 1)defcon(dat, level = 1)
dat |
sam data object |
level |
1 or 2 (1 most basic configuration, 2 configuration with AR correlation structure on surveys) |
The configuration returned by defcon is intended as a help to set up a syntactically correct configuration for the sam model. The dimensions are set from the data (years, age-classes, and fleet types available). The configuration is intended to be fairly simplistic in the hope that the model configured will at least converge (not guaranteed). Most importantly: No model validation has been performed, so it should not be assumed that the returned model configuration will result in a sensible assessment of the stock. The actual model configuration is the responsibility of the user.
a list containing the elements needed to configure a sam model (e.g. minAge, maxAge, maxAgePlusGroup, keyLogFsta, ...).
Setup initial values for all model parameters and random effects.
defpar(dat, conf, spinoutyear = 10)defpar(dat, conf, spinoutyear = 10)
dat |
sam data object as returned from the function |
conf |
sam configuration list, which could be read from a configuration file via the |
spinoutyear |
Technical setting only used for biological parameter process models to insure equilibrium distribution in final edge year |
The model parameters and random effects are not initialized in any clever way - most are simply set to zero. If convergence problems occour different initial values can be tested, but it is more likely a problem with the model configuration.
a list containing initial values for all model parameters and random effects in the model.
The function estimates reference points based on deterministic per-recruit calculations with no process variance. The following reference points are implemented:
F fixed to x, e.g., "F=0.3"
F in the last year of the assessment
F in the y years before the last in the assessment, e.g., "StatusQuo-1"
F that maximizes yield
Fs that gives 0.x*100% of MSY, e.g., "0.95MSY"
F that maximizes yield per recruit
F such that the derivative of yield per recruit is 0.x times the derivative at F=0, e.g., "0.1dYPR"
F such that spawners per recruit is 0.x times spawners per recruit at F=0, e.g., "0.35SPR"
F such that biomass is 0.x times the biomass at F=0, e.g., "0.2B0"
deterministicReferencepoints(fit, referencepoints, ...)deterministicReferencepoints(fit, referencepoints, ...)
fit |
A fitted SAM model |
referencepoints |
list of reference points to calculate (See details) |
... |
other arguments not used |
List of estimated reference points
## Not run: deterministicReferencepoints(fit, c("MSY","0.95MSY","Max","0.35SPR","0.1dYPR","StatusQuo-3")) ## End(Not run)## Not run: deterministicReferencepoints(fit, c("MSY","0.95MSY","Max","0.35SPR","0.1dYPR","StatusQuo-3")) ## End(Not run)
Function to calculate reference points for the embedded deterministic model of a SAM fit
## S3 method for class 'sam' deterministicReferencepoints( fit, referencepoints, catchType = "catch", nYears = 300, Fsequence = seq(0, 2, len = 50), aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), biasCorrect = FALSE, newton.control = list(), run = TRUE, ... )## S3 method for class 'sam' deterministicReferencepoints( fit, referencepoints, catchType = "catch", nYears = 300, Fsequence = seq(0, 2, len = 50), aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), biasCorrect = FALSE, newton.control = list(), run = TRUE, ... )
fit |
A fitted SAM model |
referencepoints |
list of reference points to calculate (See details) |
catchType |
landing, catch, or discard |
nYears |
Number of years in per-recruit calculations |
Fsequence |
Sequence of F values for plotting and starting values |
aveYears |
Years to average over for biological input |
selYears |
Years to average over for selectivity |
biasCorrect |
Should bias correction be used in sdreport? |
newton.control |
Control arguments passed to the newton optimizer (See newton) |
run |
Run estimation? If false, a list of arguments to MakeADFun is returned. |
... |
other arguments not used |
List of estimated reference points
Plots the residual between-age correlation matrices by fleet.
empirobscorrplot(res, ...) ## S3 method for class 'samres' empirobscorrplot(res, ...)empirobscorrplot(res, ...) ## S3 method for class 'samres' empirobscorrplot(res, ...)
res |
the object returned from residuals.sam |
... |
extra arguments to plot |
SAM equilibrium biomass plot
equilibriumbiomassplot(fit, ...) ## Default S3 method: equilibriumbiomassplot(fit, ...) ## S3 method for class 'samforecast' equilibriumbiomassplot(fit, ...) ## S3 method for class 'hcr' equilibriumbiomassplot(fit, ...)equilibriumbiomassplot(fit, ...) ## Default S3 method: equilibriumbiomassplot(fit, ...) ## S3 method for class 'samforecast' equilibriumbiomassplot(fit, ...) ## S3 method for class 'hcr' equilibriumbiomassplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of deterministic equilibrium spawners per recruit assuming biological parameters and selectivity for that year remains unchanged in the future.
equilibrium biomass table
equilibriumbiomasstable(fit, ...) ## Default S3 method: equilibriumbiomasstable(fit, ...)equilibriumbiomasstable(fit, ...) ## Default S3 method: equilibriumbiomasstable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
F-at-age table
faytable(fit, ...) ## S3 method for class 'sam' faytable(fit, fleet = which(fit$data$fleetTypes == 0), ...)faytable(fit, ...) ## S3 method for class 'sam' faytable(fit, fleet = which(fit$data$fleetTypes == 0), ...)
fit |
a fitted object of class 'sam' as returned from sam.fit |
... |
extra arguments not currently used |
fleet |
the fleet number(s) to return F summed for (default is to return the sum of all residual fleets). |
...
SAM Fbar plot
fbarplot(fit, ...) ## S3 method for class 'sam' fbarplot( fit, partial = TRUE, drop = NULL, pcol = "lightblue", page = NULL, plot = TRUE, ... ) ## S3 method for class 'samset' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... ) ## S3 method for class 'samforecast' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... ) ## S3 method for class 'hcr' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... )fbarplot(fit, ...) ## S3 method for class 'sam' fbarplot( fit, partial = TRUE, drop = NULL, pcol = "lightblue", page = NULL, plot = TRUE, ... ) ## S3 method for class 'samset' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... ) ## S3 method for class 'samforecast' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... ) ## S3 method for class 'hcr' fbarplot( fit, partial = FALSE, drop = NULL, pcol = "lightblue", page = NULL, ... )
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
partial |
true if included partial F's are to be plotted |
drop |
number of years to be left unplotted at the end. Default (NULL) is to not show years at the end with no catch information |
pcol |
color of partial lines |
page |
partial ages to plot |
plot |
true if fbar should be plotted |
Plot the defined fbar.
Fbar table
fbartable(fit, ...) ## Default S3 method: fbartable(fit, ...)fbartable(fit, ...) ## Default S3 method: fbartable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Read a fitted model from stockassessment.org
fitfromweb(stockname, character.only = FALSE, return.all = FALSE)fitfromweb(stockname, character.only = FALSE, return.all = FALSE)
stockname |
The short-form name of a stock on stockassessment.org. This will (currently?) not work for stocks defined via the AD Model builder version of SAM. |
character.only |
a logical indicating whether 'stockname' can be assumed to be a character string |
return.all |
a logical indicating whether everything from model.RData should be returned in an environment |
...
Plots fit to data
fitplot(fit, log = TRUE, ...) ## S3 method for class 'sam' fitplot(fit, log = TRUE, fleets = unique(fit$data$aux[, "fleet"]), ...)fitplot(fit, log = TRUE, ...) ## S3 method for class 'sam' fitplot(fit, log = TRUE, fleets = unique(fit$data$aux[, "fleet"]), ...)
fit |
the object returned from sam.fit |
log |
should the plot be against log-obs |
... |
extra arguments to plot |
fleets |
an integer vector of fleets to plot. Default is all of them |
forecast function to do shortterm
forecast( fit, fscale = NULL, catchval = NULL, catchval.exact = NULL, fval = NULL, nextssb = NULL, landval = NULL, cwF = NULL, nosim = 1000, year.base = max(fit$data$years), ave.years = max(fit$data$years) + (-4:0), rec.years = max(fit$data$years) + (-9:0), label = NULL, overwriteSelYears = NULL, deterministic = FALSE, processNoiseF = TRUE, customWeights = NULL, customSel = NULL, lagR = FALSE, splitLD = FALSE, addTSB = FALSE, useSWmodel = (fit$conf$stockWeightModel >= 1), useCWmodel = (fit$conf$catchWeightModel >= 1), useMOmodel = (fit$conf$matureModel >= 1), useNMmodel = (fit$conf$mortalityModel >= 1), savesim = FALSE, cf.cv.keep.cv = matrix(NA, ncol = 2 * sum(fit$data$fleetTypes == 0), nrow = length(catchval)), cf.cv.keep.fv = matrix(NA, ncol = 2 * sum(fit$data$fleetTypes == 0), nrow = length(catchval)), cf.keep.fv.offset = matrix(0, ncol = sum(fit$data$fleetTypes == 0), nrow = length(catchval)), estimate = median )forecast( fit, fscale = NULL, catchval = NULL, catchval.exact = NULL, fval = NULL, nextssb = NULL, landval = NULL, cwF = NULL, nosim = 1000, year.base = max(fit$data$years), ave.years = max(fit$data$years) + (-4:0), rec.years = max(fit$data$years) + (-9:0), label = NULL, overwriteSelYears = NULL, deterministic = FALSE, processNoiseF = TRUE, customWeights = NULL, customSel = NULL, lagR = FALSE, splitLD = FALSE, addTSB = FALSE, useSWmodel = (fit$conf$stockWeightModel >= 1), useCWmodel = (fit$conf$catchWeightModel >= 1), useMOmodel = (fit$conf$matureModel >= 1), useNMmodel = (fit$conf$mortalityModel >= 1), savesim = FALSE, cf.cv.keep.cv = matrix(NA, ncol = 2 * sum(fit$data$fleetTypes == 0), nrow = length(catchval)), cf.cv.keep.fv = matrix(NA, ncol = 2 * sum(fit$data$fleetTypes == 0), nrow = length(catchval)), cf.keep.fv.offset = matrix(0, ncol = sum(fit$data$fleetTypes == 0), nrow = length(catchval)), estimate = median )
fit |
an assessment object of type sam, as returned from the function sam.fit |
fscale |
a vector of f-scales. See details. |
catchval |
a vector of target catches. See details. |
catchval.exact |
a vector of target catches which will be met without noise. See details. |
fval |
a vector of target f values. See details. |
nextssb |
a vector target SSB values the following year. See details |
landval |
a vector of target catches. See details. |
cwF |
a vector target custom weighted F values. customWeights must also be specified |
nosim |
number of simulations default is 1000 |
year.base |
starting year default last year in assessment. Currently it is only supported to use last assessment year or the year before |
ave.years |
vector of years to average for weights, maturity, M and such |
rec.years |
vector of years to use to resample recruitment from |
label |
optional label to appear in short table |
overwriteSelYears |
if a vector of years is specified, then the average selectivity of those years is used (not recommended) |
deterministic |
option to turn all process noise off (not recommended, as it will likely cause bias) |
processNoiseF |
option to turn off process noise in F |
customWeights |
a vector of same length as number of age groups giving custom weights (currently only used for weighted average of F calculation) |
customSel |
supply a custom selection vector that will then be used as fixed selection in all years after the final assessment year (not recommended) |
lagR |
if the second youngest age should be reported as recruits |
splitLD |
if TRUE the result is split in landing and discards |
addTSB |
if TRUE the total stock biomass (TSB) is added |
useSWmodel |
if TRUE the catch mean weight predicted from the assessment model is used (can only be used for configurations supporting this) |
useCWmodel |
if TRUE the catch mean weight predicted from the assessment model is used (can only be used for configurations supporting this) |
useMOmodel |
if TRUE the proportion mature predicted from the assessment model is used (can only be used for configurations supporting this) |
useNMmodel |
if TRUE the natural mortality predicted from the assessment model is used (can only be used for configurations supporting this) |
savesim |
save the individual simulations |
cf.cv.keep.cv |
exotic option |
cf.cv.keep.fv |
exotic option |
cf.keep.fv.offset |
exotic option |
estimate |
the summary function used (typically mean or median) |
There are three ways to specify a scenario. If e.g. four F values are specified (e.g. fval=c(.1,.2,.3,4)), then the first value is used in the last assessment year (base.year), and the three following in the three following years. Alternatively F's can be specified by a scale, or a target catch. Only one option can be used per year. So for instance to set a catch in the first year and an F-scale in the following one would write catchval=c(10000,NA,NA,NA), fscale=c(NA,1,1,1). The length of the vector specifies how many years forward the scenarios run.
an object of type samforecast
Estimating Fmsy
forecastMSY( fit, nYears = 100, nlminb.control = list(eval.max = 2000, iter.max = 2000), rec.years = c(), ave.years = max(fit$data$years) + (-9:0), processNoiseF = FALSE, ... ) ## S3 method for class 'sam' forecastMSY( fit, nYears = 100, nlminb.control = list(eval.max = 2000, iter.max = 2000, trace = 1), rec.years = c(), ave.years = max(fit$data$years) + (-9:0), processNoiseF = FALSE, jacobianHScale = 0.5, nCatchAverageYears = 20, ... )forecastMSY( fit, nYears = 100, nlminb.control = list(eval.max = 2000, iter.max = 2000), rec.years = c(), ave.years = max(fit$data$years) + (-9:0), processNoiseF = FALSE, ... ) ## S3 method for class 'sam' forecastMSY( fit, nYears = 100, nlminb.control = list(eval.max = 2000, iter.max = 2000, trace = 1), rec.years = c(), ave.years = max(fit$data$years) + (-9:0), processNoiseF = FALSE, jacobianHScale = 0.5, nCatchAverageYears = 20, ... )
fit |
a SAM fit |
nYears |
Number of years to forecast |
nlminb.control |
list of control variables for nlminb |
rec.years |
Numeric vector of years to use (to calculate mean and standard deviation) for recruitment. An empty vector will use the recruitment model. |
ave.years |
vector of years to average for weights, maturity, M and such. Following ICES guidelines, the default is the last 10 years. |
processNoiseF |
Should random walk process noise be used for F? |
... |
other arguments passed to forecast |
jacobianHScale |
Scale step size in jacobian calculation |
nCatchAverageYears |
Number of years to average catch over for finding MSY |
Albertsen, C. M. and Trijoulet, V. (2020) Model-based estimates of reference points in an age-based state-space stock assessment model. Fisheries Research, 230, 105618. doi: 10.1016/j.fishres.2020.105618
SAM F-selectivity plot
fselectivityplot(fit, cexAge = 1, ...) ## S3 method for class 'sam' fselectivityplot(fit, cexAge = 1, ...)fselectivityplot(fit, cexAge = 1, ...) ## S3 method for class 'sam' fselectivityplot(fit, cexAge = 1, ...)
fit |
An object returned from sam.fit |
cexAge |
cex variable giving the size of the age numbers |
... |
extra arguments transferred to barplot and text |
Plots selectivity in F.
SAM generation length plot
generationlengthplot(fit, ...) ## Default S3 method: generationlengthplot(fit, ...) ## S3 method for class 'samforecast' generationlengthplot(fit, ...) ## S3 method for class 'hcr' generationlengthplot(fit, ...)generationlengthplot(fit, ...) ## Default S3 method: generationlengthplot(fit, ...) ## S3 method for class 'samforecast' generationlengthplot(fit, ...) ## S3 method for class 'hcr' generationlengthplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of life expectancy
Generation length table
generationlengthtable(fit, ...) ## Default S3 method: generationlengthtable(fit, ...)generationlengthtable(fit, ...) ## Default S3 method: generationlengthtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Update sam fit with additional derived values
getAllDerivedValues(fit)getAllDerivedValues(fit)
fit |
sam fit returned by sam.fit |
Updated sam fit
Extract a fleet observed or predicted value from a fitted object
getFleet(fit, fleet, pred = "FALSE")getFleet(fit, fleet, pred = "FALSE")
fit |
A fitted object as returned from sam.fit |
fleet |
The fleet number |
pred |
Should it be predicted value, default is observed |
Extract for example the observed or predicted catch at age of fleet "fleet"
A matrix of observed or predicted values for fleet "fleet"
Bounds
getLowerBounds(parameters, conf)getLowerBounds(parameters, conf)
parameters |
initial values for the model in a format similar to what is returned from the defpar function |
conf |
model configuration in a format similar to what is returned from the defcon function |
a named list
Extract a list of catch fleets
getResidualFleets(fit, pred = "FALSE")getResidualFleets(fit, pred = "FALSE")
fit |
A fitted object as returned from sam.fit |
pred |
Should it be predicted value, default is observed |
A list of matrices of observed or predicted values for catch fleets
Bounds
getUpperBounds(parameters, conf)getUpperBounds(parameters, conf)
parameters |
initial values for the model in a format similar to what is returned from the defpar function |
conf |
model configuration in a format similar to what is returned from the defcon function |
a named list
Calculate gradient of a function
grad( func, x, h = abs(1e-04 * x) + 1e-04 * (abs(x) < sqrt(.Machine$double.eps/7e-07)), ... )grad( func, x, h = abs(1e-04 * x) + 1e-04 * (abs(x) < sqrt(.Machine$double.eps/7e-07)), ... )
func |
function |
x |
parameter values |
h |
step size |
... |
passed to func |
gradient vector
The formula below is used to determine a new F based on the previous SSB.
If and , is always returned.
hcr(fit, ...) ## S3 method for class 'sam' hcr( fit, nYears = 20, Ftarget, Btrigger, Forigin = 0, Borigin = 0, Fcap = 0, Bcap = 0, nosim = 10000, ave.years = max(fit$data$years) + (-4:0), rec.years = numeric(0), preForecast = list(), currentSSB = FALSE, ... )hcr(fit, ...) ## S3 method for class 'sam' hcr( fit, nYears = 20, Ftarget, Btrigger, Forigin = 0, Borigin = 0, Fcap = 0, Bcap = 0, nosim = 10000, ave.years = max(fit$data$years) + (-4:0), rec.years = numeric(0), preForecast = list(), currentSSB = FALSE, ... )
fit |
A SAM fit |
... |
additional arguments passed to modelforecast |
nYears |
Number of years to forecast |
Ftarget |
Target F for high SSB |
Btrigger |
SSB that triggers the control rule |
Forigin |
F used for SSB = Borigin |
Borigin |
Between Blim and Btrigger, F values are selected based on linear interpolation from Forigin to Ftarget |
Fcap |
F for SSB < Bcap |
Bcap |
SSB for which Fcap is used below |
nosim |
Number of simulations to do. If NULL a model forecast based on the Laplace approximation is used |
ave.years |
vector of years to average for weights, maturity, M and such |
rec.years |
vector of years to use to resample recruitment from. If an empty vector is given, recruitment is based on the fitted model. |
preForecast |
list of forecast parameters (i.e., fval, fscale, catchval, landval, or nextssb) to use before the HCR |
currentSSB |
if TRUE, SSB at the begining of the control rule year is used. If FALSE, SSB at the begining of the previous year is used. |
model forecast using a harvest control rule
hcr model forecast object
modelforecast
Integrated spline basis for use with formula interface
ibc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)ibc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)
x |
Points to evaluate the basis in |
df |
Degrees of freedom |
knots |
Internal knots. If NULL, they are selected from quantiles of x. |
Boundary.knots |
Boundary knots. Defaults to range of x |
intercept |
Include an intercept in basis? |
A spline basis
Forecast with an ICES advice rule
icesAdviceRule( x, Fmsy, MSYBtrigger, Blim, nosim = 10000, ave.years = max(x$data$years) + (-4:0), rec.years = numeric(0), preForecast = list(), currentSSB = FALSE, ... )icesAdviceRule( x, Fmsy, MSYBtrigger, Blim, nosim = 10000, ave.years = max(x$data$years) + (-4:0), rec.years = numeric(0), preForecast = list(), currentSSB = FALSE, ... )
x |
Fitted assessment model |
Fmsy |
ICES Fmsy which is used as target F |
MSYBtrigger |
ICES MSYBtrigger below which F is reduced |
Blim |
ICES Blim below which F is set to zero. |
nosim |
Number of simulations to do. If NULL a model forecast based on the Laplace approximation is used |
ave.years |
vector of years to average for weights, maturity, M and such |
rec.years |
vector of years to use to resample recruitment from. If an empty vector is given, recruitment is based on the fitted model. |
preForecast |
list of forecast parameters (i.e., fval, fscale, catchval, landval, or nextssb) to use before the HCR |
currentSSB |
if TRUE, SSB at the begining of the control rule year is used. If FALSE, SSB at the begining of the previous year is used. |
... |
Other arguments passes to hcr |
hcr object
The function does not make a short term forecast to see if fishing can continue below Blim.
ICES (2021) Advice on fishing opportunities. DOI: 10.17895/ices.advice.7720
Double integrated spline basis for use with formula interface
iibc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)iibc(x, df = 3L, knots = NULL, Boundary.knots = range(x), intercept = FALSE)
x |
Points to evaluate the basis in |
df |
Degrees of freedom |
knots |
Internal knots. If NULL, they are selected from quantiles of x. |
Boundary.knots |
Boundary knots. Defaults to range of x |
intercept |
Include an intercept in basis? |
A spline basis
Function to test if x is ...
is.whole.positive.number(x, tol = .Machine$double.eps^0.5)is.whole.positive.number(x, tol = .Machine$double.eps^0.5)
x |
number |
tol |
precision |
...
Calculate jacobian of a function
jacobian( func, x, h = abs(1e-04 * x) + 1e-04 * (abs(x) < sqrt(.Machine$double.eps/7e-07)), ... )jacobian( func, x, h = abs(1e-04 * x) + 1e-04 * (abs(x) < sqrt(.Machine$double.eps/7e-07)), ... )
func |
function |
x |
parameter values |
h |
step size |
... |
passed to func |
jacobian matrix
Jitter runs
jit( fit, nojit = 10, par = defpar(fit$data, fit$conf), sd = 0.25, ncores = detectCores() )jit( fit, nojit = 10, par = defpar(fit$data, fit$conf), sd = 0.25, ncores = detectCores() )
fit |
a fitted model object as returned from sam.fit |
nojit |
a list of vectors. Each element in the list specifies a run where the fleets mentioned are omitted |
par |
initial values to jitter around. The defaule ones are returned from the defpar function |
sd |
the standard deviation used to jitter the initial values (most parameters are on a log scale, so similar to cv) |
ncores |
the number of cores to attemp to use |
...
A "samset" object, which is basically a list of sam fits
leaveout run
leaveout( fit, fleet = as.list(2:fit$data$noFleets), ncores = detectCores(), ... )leaveout( fit, fleet = as.list(2:fit$data$noFleets), ncores = detectCores(), ... )
fit |
a fitted model object as returned from sam.fit |
fleet |
a list of vectors. Each element in the list specifies a run where the fleets mentioned are omitted |
ncores |
the number of cores to attemp to use |
... |
extra arguments to |
...
SAM life expectancy plot
lifeexpectancyplot(fit, atRecruit = TRUE, ...) ## Default S3 method: lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...) ## S3 method for class 'samforecast' lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...) ## S3 method for class 'hcr' lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...)lifeexpectancyplot(fit, atRecruit = TRUE, ...) ## Default S3 method: lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...) ## S3 method for class 'samforecast' lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...) ## S3 method for class 'hcr' lifeexpectancyplot(fit, atRecruit = TRUE, ylimAdd = fit$conf$maxAge, ...)
fit |
the object returned from sam.fit |
atRecruit |
If true, show life expectancy given survival until minAge, otherwise show life expectancy at birth |
... |
extra arguments transferred to plot including the following: |
ylimAdd |
values to add when calculating ylim for the plot |
Plot of life expectancy
Life expectancy table
lifeexpectancytable(fit, atRecruit = TRUE, ...) ## Default S3 method: lifeexpectancytable(fit, atRecruit = TRUE, ...)lifeexpectancytable(fit, atRecruit = TRUE, ...) ## Default S3 method: lifeexpectancytable(fit, atRecruit = TRUE, ...)
fit |
... |
atRecruit |
If true, show life expectancy given survival until minAge, otherwise show life expectancy at birth |
... |
extra arguments not currently used |
...
Loads a model configuration from a file
loadConf(dat, file, patch = TRUE)loadConf(dat, file, patch = TRUE)
dat |
sam data list as returned from the function setup.sam.data |
file |
the file to read the configuration from |
patch |
logical if TRUE missing entries will be automatically filled with default values |
function useful loading a model configuration. Such a configuration can be saved via the saveConf function
Log likelihood of sam object
## S3 method for class 'sam' logLik(object, ...)## S3 method for class 'sam' logLik(object, ...)
object |
sam fitted object as returned from the |
... |
extra arguments |
log likelihood of a sam model run
Description of model
modelDescription(fit, ...)modelDescription(fit, ...)
fit |
returned object from sam.fit |
... |
Additional parameters to be passed to ... |
...
Model based forecast function
Model based forecast function
modelforecast(fit, ...) ## S3 method for class 'sam' modelforecast( fit, constraints = NULL, fscale = NULL, catchval = NULL, fval = NULL, nextssb = NULL, landval = NULL, findMSY = NULL, hcr = NULL, nosim = 0, year.base = max(fit$data$years), ave.years = c(), rec.years = c(), label = NULL, overwriteSelYears = NULL, deterministicF = FALSE, processNoiseF = FALSE, resampleFirst = !is.null(nosim) && nosim > 0, customSel = NULL, lagR = FALSE, splitLD = FALSE, addTSB = FALSE, biasCorrect = FALSE, returnAllYears = FALSE, nCatchAverageYears = 1, returnObj = FALSE, hcrConf = numeric(0), hcrCurrentSSB = 0, progress = TRUE, estimate = median, silent = TRUE, newton_config = NULL, ... )modelforecast(fit, ...) ## S3 method for class 'sam' modelforecast( fit, constraints = NULL, fscale = NULL, catchval = NULL, fval = NULL, nextssb = NULL, landval = NULL, findMSY = NULL, hcr = NULL, nosim = 0, year.base = max(fit$data$years), ave.years = c(), rec.years = c(), label = NULL, overwriteSelYears = NULL, deterministicF = FALSE, processNoiseF = FALSE, resampleFirst = !is.null(nosim) && nosim > 0, customSel = NULL, lagR = FALSE, splitLD = FALSE, addTSB = FALSE, biasCorrect = FALSE, returnAllYears = FALSE, nCatchAverageYears = 1, returnObj = FALSE, hcrConf = numeric(0), hcrCurrentSSB = 0, progress = TRUE, estimate = median, silent = TRUE, newton_config = NULL, ... )
fit |
an assessment object of type sam, as returned from the function sam.fit |
... |
other variables used by the methods |
constraints |
a character vector of forecast constraint specifications |
fscale |
a vector of f-scales. See details. |
catchval |
a vector of target catches. See details. |
fval |
a vector of target f values. See details. |
nextssb |
a vector target SSB values the following year. See details. |
landval |
a vector of target catches. See details. |
findMSY |
Should not be used. See forecastMSY. |
hcr |
Should not be used. See hcr. |
nosim |
number of simulations. If 0, the Laplace approximation is used for forecasting. |
year.base |
starting year default last year in assessment. Currently it is only supported to use last assessment year or the year before |
ave.years |
vector of years to average for weights, maturity, M and such |
rec.years |
vector of years to use to resample recruitment from. If the vector is empty, the stock recruitment model is used. |
label |
optional label to appear in short table |
overwriteSelYears |
if a vector of years is specified, then the average selectivity of those years is used (not recommended) |
deterministicF |
option to set F variance to (almost) zero (not recommended) |
processNoiseF |
option to turn off process noise in F |
resampleFirst |
Resample base year when nosim > 0? |
customSel |
supply a custom selection vector that will then be used as fixed selection in all years after the final assessment year (not recommended) |
lagR |
if the second youngest age should be reported as recruits |
splitLD |
if TRUE the result is split in landing and discards |
addTSB |
if TRUE the total stock biomass (TSB) is added |
biasCorrect |
Do bias correction of reported variables. Can be turned off to reduce running time (not recommended). |
returnAllYears |
If TRUE, all years are bias corrected. Otherwise, only forecast years are corrected. |
nCatchAverageYears |
Should not be used. See forecastMSY. |
returnObj |
Only return TMB object? |
hcrConf |
Should not be used. See hcr. |
hcrCurrentSSB |
Should not be used. See hcr. |
progress |
Show progress bar for simulations? |
estimate |
the summary function used (typically mean or median) for simulations |
silent |
Passed to MakeADFun. Should the TMB object be silent? |
newton_config |
Configuration for newton optimizer to find F values. See ?TMB::newton for details. Use NULL for TMB defaults. |
Function to forecast the model under specified catch constraints. In the forecast, catch constraints are used to set the mean of the process. Therefore, catch constraints are not matched exactly in individual simulations (unlike the forecast function simulations). Likewise, the summary of a specific set of simulations will not match exactly due to random variability.
By default, recruitment is forecasted using the estimated recruitment model. If a vector of recruitment years is given, recruitment is forecasted using a log-normal distribution with the same mean and variance as the recruitment in the years given. This is different from the forecast function, which samples from the recruitment estimates.
Catch scenarios are specified by a vector of target constraints. The first value determines F in the year after the base year.
an object of type samforecast
Forecasts for F values are specified by the format "F[f,a0-a1]=x" where f is the residual catch fleet and a0-a1 is an age range. For example, "F[2,2-4]=0.3" specifies that the average F for the second fleet over ages 2-4 should be 0.3. If an "*" is added to the target value, the target will be relative to the year before. For example, "F[2,2-4]=0.9*" specifies that the average F for the second fleet over ages 2-4 should be 90 If the fleet is omitted (e.g., F[2-4]), the target is for the total F. If the age range is omitted (e.g., F[2]), the fbar range of the model is used. Likewise, both fleet and age range can be omited (e.g., F=0.3) to specify a value for total F with the range used in the model.
Forecasts for catch and landing values are specified by the format "C[f,a0-a1]=x" for catch and "L[f,a0-a1]" for landings. If the age range is omitted, all modelled ages are used. Otherwise, the format is similar to F based scenarios. If an "*" is added to the target value, the target will be relative to the year before. Further, the catch target for a fleet can be relative to the total by adding "*C" or to another fleet by adding "*C[f]" where f is the fleet number. The same age range will always be used. Likewise, relative landing targets can be specified using "*", "*L", or "*L[f]" for targets relative to last year, the total, or fleet f, respectively.
Forecasts for spawning stock biomass (SSB) and total stock biomass (TSB) values are specified by the format "SSB[a0-a1]=x" for SSB and "TSB[a0-a1]" for TSB. The format is similar to catch/landing based scenarios. However, fleets have no effect. If an "*" is added to the target value, the target will be relative to the year before. Note that SSB amd TSB used for catch constraints are at the beginning of the next year to avoid dependence on future F values. Therefore, the values will differ from the output SSB and TSB estimates if propM or propF are not zero.
Constraints for different fleets can be combined by "&". For example, "F[2-4]=0.5 & C[2]=10000" specifies that total Fbar over ages 2-4 should be 0.5 while the catch for the second residual catch fleet should be 10,000t. The constraints cannot affect within-fleet selectivity. Therefore, a fleet can at most have one constraint per year, and the total number of constraints cannot exceed the number of catch fleets. That is, if a constraint is given for the sum of fleets, there must be at least one fleet without any constraints. For fleets where no constraints are given, a constraint is set to keep their relative Fs constant.
Catch constraints specified as specific values are inherently different from catch constraints specified as relative values, even if they lead to the same F. Catch constraints specified as relative values will propagate the uncertainty in, e.g, F from previous years whereas constraints specified as specific values will not. This is different from the forecast function where, for example, a forecast using fval is the same as a forecast using fscale, if they lead to the same F.
It is also possible to specify forecast constraints in a way similar to the forecast function. There are four ways to specify a scenario. If e.g. four F values are specified (e.g. fval=c(.1,.2,.3,4)), then the first value is used in the year after the last assessment year (base.year + 1), and the three following in the three following years. Alternatively F's can be specified by a scale, or a target catch. Only one option can be used per year. So for instance to set a catch in the first year and an F-scale in the following one would write catchval=c(10000,NA,NA,NA), fscale=c(NA,1,1,1). If only NA's are specified in a year, the F model is used for forecasting. The length of the vector specifies how many years forward the scenarios run. Unlike the forecast function, no value should be given for the base year. Internally, the old specification is translated such that "fval=x" becomes "F=x", "fscale=x" becomes "F=x*", "catchval=x" becomes "C=x", "nextssb=x" becomes "SSB=x", and "landval=x" becomes "L=x".
NA
Long term forecasts with random walk recruitment can lead to unstable behaviour and difficulties finding suitable F values for the constraints. If no suitable F value can be found, an error message will be shown, and F values will be NA or NaN. Likewise, forecasts leading to high F values in some years may cause problems for the optimization as they will be used as starting values for the next years. Since the model works on log space, all target values should be strictly positive. Values too close to zero may cause problems.
forecast
model table
modeltable(fits, ...) ## S3 method for class 'sam' modeltable(fits, ...) ## S3 method for class 'samset' modeltable(fits, ...)modeltable(fits, ...) ## S3 method for class 'sam' modeltable(fits, ...) ## S3 method for class 'samset' modeltable(fits, ...)
fits |
A sam fit as returned from the sam.fit function, or a collection c(fit1, fit2, ...) of such fits |
... |
extra arguments not currently used |
...
Description of model
modelVersionInfo(fit, ...)modelVersionInfo(fit, ...)
fit |
returned object from sam.fit |
... |
Additional parameters to be passed to ... |
Writes a string to install the version of the package which was used to run the model.
Mohn's rho calculation
mohn(fits, what = NULL, lag = 0, ...) ## S3 method for class 'samset' mohn(fits, what = NULL, lag = 0, ...)mohn(fits, what = NULL, lag = 0, ...) ## S3 method for class 'samset' mohn(fits, what = NULL, lag = 0, ...)
fits |
a samset object as returned from the retro function. |
what |
a function computing the quantity to calculate Mohn's rho for (default NULL computes Fbar, SSB, and R). |
lag |
lag applied to fits and reference fit. |
... |
extra arguments |
...
Extract number of observations from sam object
## S3 method for class 'sam' nobs(object, ...)## S3 method for class 'sam' nobs(object, ...)
object |
sam fitted object as returned from the |
... |
extra arguments |
...
nscodConf
data("nscodConf")data("nscodConf")
The format is: $ minAge $ maxAge $ maxAgePlusGroup $ keyLogFsta $ corFlag $ keyLogFpar $ keyQpow $ keyVarF $ keyVarLogN $ keyVarObs $ stockRecruitmentModelCode $ noScaledYears $ keyScaledYears $ keyParScaledYA $ fbarRange
...
...
...
data(nscodConf) ## maybe str(nscodConf) ; plot(nscodConf) ...data(nscodConf) ## maybe str(nscodConf) ; plot(nscodConf) ...
nscodData
data("nscodData")data("nscodData")
The format is: $ noFleets $ fleetTypes $ sampleTimes $ noYears $ years $ nobs $ idx1 $ idx2 $ aux $ logobs $ propMat $ stockMeanWeight $ catchMeanWeight $ natMor $ landFrac $ disMeanWeight $ landMeanWeight $ propF $ propM
...
...
...
data(nscodData) ## maybe str(nscodData) ; plot(nscodData) ...data(nscodData) ## maybe str(nscodData) ; plot(nscodData) ...
nscodParameters
data("nscodParameters")data("nscodParameters")
The format is: List of 14 $ logFpar $ logQpow $ logSdLogFsta $ logSdLogN $ logSdLogObs $ rec_loga $ rec_logb $ itrans_rho $ logScale $ logScaleSSB $ logPowSSB $ logSdSSB $ logF $ logN
...
...
...
data(nscodParameters) ## maybe str(nscodParameters) ; plot(nscodParameters) ...data(nscodParameters) ## maybe str(nscodParameters) ; plot(nscodParameters) ...
N table
ntable(fit, ...) ## S3 method for class 'sam' ntable(fit, ...)ntable(fit, ...) ## S3 method for class 'sam' ntable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Plots the estimated correlation matrices by fleet.
obscorrplot(fit, ...) ## S3 method for class 'sam' obscorrplot(fit, ...)obscorrplot(fit, ...) ## S3 method for class 'sam' obscorrplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments to plot |
Extract observation covariance matrices from a SAM fit
obscov(fit, corr = FALSE, ...) ## S3 method for class 'sam' obscov(fit, corr = FALSE, ...) ## S3 method for class 'samset' obscov(fit, corr = FALSE, ...)obscov(fit, corr = FALSE, ...) ## S3 method for class 'sam' obscov(fit, corr = FALSE, ...) ## S3 method for class 'samset' obscov(fit, corr = FALSE, ...)
fit |
the object returned from sam.fit |
corr |
if TRUE return correlation matrices rather than covariances |
... |
extra arguments not currently used |
a list of matrices
SAM parameter plot
parplot(fit, cor.report.limit = 0.95, ...) ## S3 method for class 'sam' parplot(fit, cor.report.limit = 0.95, ...) ## S3 method for class 'samset' parplot(fit, cor.report.limit = 0.95, ...)parplot(fit, cor.report.limit = 0.95, ...) ## S3 method for class 'sam' parplot(fit, cor.report.limit = 0.95, ...) ## S3 method for class 'samset' parplot(fit, cor.report.limit = 0.95, ...)
fit |
the object returned from sam.fit |
cor.report.limit |
correlations with absolute value > this number is reported in the plot |
... |
extra arguments transferred to plot |
Plot of all estimated model parameters (fixed effects). Shown with confidence interval.
parameter table
partable(fit, ...) ## S3 method for class 'sam' partable(fit, ...)partable(fit, ...) ## S3 method for class 'sam' partable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Plot hcr object
## S3 method for class 'hcr' plot(x, ...)## S3 method for class 'hcr' plot(x, ...)
x |
output from the |
... |
extra arguments |
...
Plot sam object
## S3 method for class 'sam' plot(x, ...)## S3 method for class 'sam' plot(x, ...)
x |
fitted object as returned from the |
... |
extra arguments (not possible to use add=TRUE — instead collect to a list of fits using e.g the c(...), and then plot that collected object). |
gives a 3 plot overview plot og ssb, fbar, and recruits. These plots are available individually via the functions ssbplot, fbarplot, and recplot.
Plot samforecast object
## S3 method for class 'samforecast' plot(x, ...)## S3 method for class 'samforecast' plot(x, ...)
x |
fitted object as returned from the |
... |
extra arguments |
...
Plot sam residuals
## S3 method for class 'samres' plot(x, type = "bubble", ...)## S3 method for class 'samres' plot(x, type = "bubble", ...)
x |
an object of type 'samres' as returned from functions |
type |
either "bubble" (default) or "summary" |
... |
extra arguments |
In the "bubble" type red indicate negative residuals and blue positive. The area of the circles scales with the absolute size of the residuals.
## Not run: data(nscodData) data(nscodConf) data(nscodParameters) fit <- sam.fit(nscodData, nscodConf, nscodParameters) par(ask=FALSE) plot(residuals(fit)) ## End(Not run)## Not run: data(nscodData) data(nscodConf) data(nscodParameters) fit <- sam.fit(nscodData, nscodConf, nscodParameters) par(ask=FALSE) plot(residuals(fit)) ## End(Not run)
Plot sam object
## S3 method for class 'samset' plot(x, ...)## S3 method for class 'samset' plot(x, ...)
x |
fitted object as returned from the |
... |
extra arguments |
...
Plot sam object
## S3 method for class 'samypr' plot(x, ...)## S3 method for class 'samypr' plot(x, ...)
x |
... |
... |
extra arguments |
...
Plot by one or two
plotby( x = NULL, y = NULL, z = NULL, x.line = NULL, y.line = NULL, z.line = NULL, by = NULL, bubblescale = 1, x.common = !is.null(x), y.common = !is.null(y), z.common = !is.null(z), xlab = NULL, ylab = NULL, xlim = NULL, ylim = NULL, zmax = NULL, axes = TRUE, ... )plotby( x = NULL, y = NULL, z = NULL, x.line = NULL, y.line = NULL, z.line = NULL, by = NULL, bubblescale = 1, x.common = !is.null(x), y.common = !is.null(y), z.common = !is.null(z), xlab = NULL, ylab = NULL, xlim = NULL, ylim = NULL, zmax = NULL, axes = TRUE, ... )
x |
numeric vector of points to be plotted |
y |
numeric vector of points to be plotted |
z |
numeric vector of points to be plotted |
x.line |
numeric vector of points of line to be added |
y.line |
numeric vector of points of line to be added |
z.line |
numeric vector of points of line to be added |
by |
vector or two column matrix to create sub sets from |
bubblescale |
scaling of bubble size |
x.common |
logical: use same x-axis for all plots |
y.common |
logical: use same y-axis for all plots |
z.common |
logical: use same z-axis for all plots |
xlab |
normal graphical parameter |
ylab |
normal graphical parameter |
xlim |
normal graphical parameter |
ylim |
normal graphical parameter |
zmax |
internally used to scale bubbles similarly |
axes |
normal graphical parameter |
... |
additional graphical parameters |
Function used for splitting plots e.g. used to plot residuals
exdat<-expand.grid(age=1:5, year=1950:2016, fleet=1:3) exdat$perfectres<-rnorm(nrow(exdat)) attach(exdat) par(ask=FALSE) plotby(year,age,perfectres, by=fleet) detach(exdat)exdat<-expand.grid(age=1:5, year=1950:2016, fleet=1:3) exdat$perfectres<-rnorm(nrow(exdat)) attach(exdat) par(ask=FALSE) plotby(year,age,perfectres, by=fleet) detach(exdat)
Plot helper
plotit(fit, what, ...) ## S3 method for class 'sam' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ylim = NULL, ylimAdd = NA, ... ) ## S3 method for class 'samset' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = rep(FALSE, length(fit)), drop = 0, unnamed.basename = "current", xlim = NULL, ... ) ## S3 method for class 'samforecast' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ylim = NULL, ... ) ## S3 method for class 'hcr' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ... )plotit(fit, what, ...) ## S3 method for class 'sam' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ylim = NULL, ylimAdd = NA, ... ) ## S3 method for class 'samset' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = rep(FALSE, length(fit)), drop = 0, unnamed.basename = "current", xlim = NULL, ... ) ## S3 method for class 'samforecast' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ylim = NULL, ... ) ## S3 method for class 'hcr' plotit( fit, what, x = fit$data$years, ylab = what, xlab = "Years", ex = numeric(0), trans = function(x) x, add = FALSE, ci = TRUE, cicol = gray(0.5, alpha = 0.5), addCI = NA, drop = 0, unnamed.basename = "current", xlim = NULL, ... )
fit |
the fitted object from sam.fit of a set of such fits c(fit1,fit2) |
what |
quoted name of object to extract |
... |
extra arguments transferred to plot |
x |
x-values |
ylab |
label on y-axis |
xlab |
label on x-axis |
ex |
extra y's to make room for |
trans |
function to transform values by |
add |
logical, plotting is to be added on existing plot |
ci |
logical, confidence intervals should be plotted |
cicol |
color to plot the confidence polygon |
addCI |
A logical vector indicating if confidence intervals should be plotted for the added fits. |
drop |
number of years to be left unplotted at the end. |
unnamed.basename |
the name to assign an unnamed basefit |
xlim |
xlim for the plot |
ylim |
ylim for the plot |
ylimAdd |
values to add when calculating ylim for the plot |
The basic plotting used bu many of the plotting functions (e.g. ssbplot, fbarplot ...)
Prediction-standard deviation plot
predstdplot( fit, fleet, age = NULL, type = "log", ylim = NULL, ylab = "Standard deviation", xlab = "Prediction", main = "Pred-std relation", ... )predstdplot( fit, fleet, age = NULL, type = "log", ylim = NULL, ylab = "Standard deviation", xlab = "Prediction", main = "Pred-std relation", ... )
fit |
A sam fit object returned from sam.fit. |
fleet |
Fleet number to plot relation between prediction and standard deviation. |
age |
Relation at age. Only used in cases with more than one relation within the same fleet. |
type |
Either 'log' or 'natural': relation for observations on a log or natural scale. |
ylim |
Optional, sent to plot |
ylab |
Optional, sent to plot |
xlab |
Optional, sent to plot |
main |
Optional, sent to plot |
... |
Sent to plot |
Plot the relation between observation prediction and standard deviation.
Print hcr object
## S3 method for class 'hcr' print(x, ...)## S3 method for class 'hcr' print(x, ...)
x |
a sam hcr object as returned by |
... |
extra arguments |
prints the HCR forecast
Print sam object
## S3 method for class 'sam' print(x, ...)## S3 method for class 'sam' print(x, ...)
x |
the fitted object as returned from the |
... |
extra arguments |
prints the log-likelihood and the main convergence criteria
Print referencepoint object
## S3 method for class 'sam_referencepoints' print(x, tables = c("F", "Biomass", "Yield"), digits = 4, format = "f", ...)## S3 method for class 'sam_referencepoints' print(x, tables = c("F", "Biomass", "Yield"), digits = 4, format = "f", ...)
x |
a sam referencepoint object as returned by |
tables |
tables to print |
digits |
number of digits to print |
format |
printing format for numbers |
... |
extra arguments |
prints the F reference point table
Print samcoef object
## S3 method for class 'samcoef' print(x, ...)## S3 method for class 'samcoef' print(x, ...)
x |
samcoef object as returned from the |
... |
extra arguments |
Print samforecast object
## S3 method for class 'samforecast' print(x, ...)## S3 method for class 'samforecast' print(x, ...)
x |
an object as returned from the forecast function |
... |
extra arguments |
...
Print samres object
## S3 method for class 'samres' print(x, ...)## S3 method for class 'samres' print(x, ...)
x |
a sam residual object as returned from either |
... |
extra arguments |
prints the residuals as a data.frame
Print samset object
## S3 method for class 'samset' print(x, ...)## S3 method for class 'samset' print(x, ...)
x |
a list of sam models |
... |
extra arguments |
...
Print samypr object
## S3 method for class 'samypr' print(x, ...)## S3 method for class 'samypr' print(x, ...)
x |
an object as returned from the ypr function |
... |
extra arguments |
...
Compute process residuals (single joint sample)
procres(fit, map = fit$obj$env$map, ...)procres(fit, map = fit$obj$env$map, ...)
fit |
the fitted object as returned from the |
map |
map from original fit |
... |
extra arguments (not currently used) |
Single joint sample residuals of log(F) and log(N)
an object of class samres
table of survey catchabilities
qtable(fit, ...)qtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
table of survey catchabilities
## S3 method for class 'sam' qtable(fit, ...)## S3 method for class 'sam' qtable(fit, ...)
fit |
A sam fit as returned from the sam.fit function |
... |
extra arguments not currently used |
plot survey catchabilities
plot survey catchabilities
qtableplot(qt, exp = FALSE) ## S3 method for class 'samqtable' qtableplot(qt, exp = FALSE)qtableplot(qt, exp = FALSE) ## S3 method for class 'samqtable' qtableplot(qt, exp = FALSE)
qt |
An object of class 'samqtable' as returned from qtable |
exp |
if true return on natural scale rather than log |
Read all standard data SAM files and return a list as created by 'setup.sam.data'
read.data.files(dir = ".")read.data.files(dir = ".")
dir |
Directory to read from |
Read all standard SAM data files
list (as created by 'setup.sam.data')
Function to read ICES/CEFAS data files and validate if input makes sense
read.ices(filen)read.ices(filen)
filen |
The filename |
First two lines are ignored and can be used for comments. Can read formats 1 full, 2 row, 3 scalar, and 5 column
Tests: Formatcode is valid, years and ages are pos. integers minimum <= maximum for years and ages number of rows and coulmns match year and age ranges data contains only numbers.
Returns: A validated data matrix.
Function to read ices survey format
read.surveys(filen)read.surveys(filen)
filen |
the file |
...
Function to supress incomplete final line warning
read.table.nowarn(...)read.table.nowarn(...)
... |
arguments |
...
SAM Recruits plot
recplot(fit, lagR = FALSE, ...) ## S3 method for class 'sam' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'samset' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'samforecast' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'hcr' recplot(fit, lagR = FALSE, ...)recplot(fit, lagR = FALSE, ...) ## S3 method for class 'sam' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'samset' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'samforecast' recplot(fit, lagR = FALSE, ...) ## S3 method for class 'hcr' recplot(fit, lagR = FALSE, ...)
fit |
the object returned from sam.fit |
lagR |
use the age after the youngest as R |
... |
extra arguments transferred to plot including the following: |
Plot of numbers of recruits (youngest age class)
Recruit table
rectable(fit, lagR = FALSE, ...) ## Default S3 method: rectable(fit, lagR = FALSE, ...)rectable(fit, lagR = FALSE, ...) ## Default S3 method: rectable(fit, lagR = FALSE, ...)
fit |
... |
lagR |
use the age after the youngest as R |
... |
extra arguments not currently used |
...
reduce helper function to reduce data
reduce(data, year = NULL, fleet = NULL, age = NULL, conf = NULL)reduce(data, year = NULL, fleet = NULL, age = NULL, conf = NULL)
data |
a data object as returned by the function setup.sam.data |
year |
a vector of years to be excluded. |
fleet |
a vector of fleets to be excluded. |
age |
a vector of ages fleets to be excluded. |
conf |
an optional corresponding configuration to be modified along with the data change. Modified is returned as attribute "conf" |
When more than one vector is supplied they need to be of same length, as only the pairs are excluded
Estimate reference points
referencepoints( fit, nYears, Fsequence, aveYears, selYears, SPRpercent, catchType, MSYreduction, newtonSteps = 3, optN = 100, jacobianHScale = 0.5, ... ) ## S3 method for class 'sam' referencepoints( fit, nYears = 100, Fsequence = seq(0, 4, len = 200), aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), SPRpercent = c(0.35), dYPRpercent = c(0.1), B0percent = c(0.2), catchType = "catch", MSYreduction = c(0.05), newtonSteps = 3, optN = 20, jacobianHScale = 0.5, fixRP = c(), RecCorrection = "median", biasCorrect = FALSE, nlminb.control = list(eval.max = 1000, iter.max = 1000), ... )referencepoints( fit, nYears, Fsequence, aveYears, selYears, SPRpercent, catchType, MSYreduction, newtonSteps = 3, optN = 100, jacobianHScale = 0.5, ... ) ## S3 method for class 'sam' referencepoints( fit, nYears = 100, Fsequence = seq(0, 4, len = 200), aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), SPRpercent = c(0.35), dYPRpercent = c(0.1), B0percent = c(0.2), catchType = "catch", MSYreduction = c(0.05), newtonSteps = 3, optN = 20, jacobianHScale = 0.5, fixRP = c(), RecCorrection = "median", biasCorrect = FALSE, nlminb.control = list(eval.max = 1000, iter.max = 1000), ... )
fit |
an object to calculate reference points for |
nYears |
Number of years to use in per-recruit calculations |
Fsequence |
Sequence of F values used to report per-recruit and equilibrium values |
aveYears |
Vector of year indices used to calculate average natural mortality, weights, etc. Following ICES guidelines, the default is the last 10 years (starting at 0) |
selYears |
Vector of year indices used to calculate selectivity (starting at 0) |
SPRpercent |
Vector of x values for F[x * 100%] reference points. Default is 0.35. |
catchType |
Catch type used: (total) catch, landings, discard. |
MSYreduction |
Vector of proportions for MSY ranges. Default is 0.05 giving an MSY range corresponding to no more than a 5% yield reduction. |
newtonSteps |
Number of additional Newton steps at the end of the reference point optimization. |
optN |
N used for numerical optimizers to find equilibrium biomass |
jacobianHScale |
Scale step size in jacobian calculation |
... |
not used |
dYPRpercent |
Defunct |
B0percent |
Defunct |
fixRP |
Defunct |
RecCorrection |
Defunct |
biasCorrect |
Defunct |
nlminb.control |
Defunct |
a sam_referencepoints fit
Albertsen, C. M. and Trijoulet, V. (2020) Model-based estimates of reference points in an age-based state-space stock assessment model. Fisheries Research, 230, 105618. doi: 10.1016/j.fishres.2020.105618
Re-fit a model from stockassessment.org
refit(fit, newConf, startingValues, ...)refit(fit, newConf, startingValues, ...)
fit |
a sam fit or the name of a fit from stockassessment.org |
newConf |
list changes to the configuration |
startingValues |
list of parameter values to use as starting values |
... |
Arguments passed to sam.fit |
A new sam fit
Extract residuals from sam object
## S3 method for class 'sam' residuals(object, discrete = FALSE, ...)## S3 method for class 'sam' residuals(object, discrete = FALSE, ...)
object |
sam fitted object as returned from the |
discrete |
logical if model contain discrete observations |
... |
extra arguments for TMB's oneStepPredict |
one-observation-ahead quantile residuals are calculated
...
retro run
retro(fit, year = NULL, ncores = detectCores(), ...) ## S3 method for class 'sam' retro(fit, year = NULL, ncores = detectCores(), ...)retro(fit, year = NULL, ncores = detectCores(), ...) ## S3 method for class 'sam' retro(fit, year = NULL, ncores = detectCores(), ...)
fit |
a fitted model object as returned from sam.fit |
year |
either 1) a single integer n in which case runs where all fleets are reduced by 1, 2, ..., n are returned, 2) a vector of years in which case runs where years from and later are excluded for all fleets, and 3 a matrix of years were each column is a fleet and each column corresponds to a run where the years and later are excluded. |
ncores |
the number of cores to attempt to use |
... |
extra arguments to |
...
SAM rmax plot
rmaxplot(fit, ...) ## Default S3 method: rmaxplot(fit, ...) ## S3 method for class 'samforecast' rmaxplot(fit, ...) ## S3 method for class 'hcr' rmaxplot(fit, ...)rmaxplot(fit, ...) ## Default S3 method: rmaxplot(fit, ...) ## S3 method for class 'samforecast' rmaxplot(fit, ...) ## S3 method for class 'hcr' rmaxplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of life expectancy
rmax table
rmaxtable(fit, ...) ## Default S3 method: rmaxtable(fit, ...)rmaxtable(fit, ...) ## Default S3 method: rmaxtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
rmvnorm helper function to draw multivariate normal samples
rmvnorm(n = 1, mu, Sigma)rmvnorm(n = 1, mu, Sigma)
n |
the number of samples. |
mu |
the mean vector. |
Sigma |
a positive-definite symmetric matrix specifying the covariance matrix. |
Generates samples via the Cholesky decomposition, which is less platform dependent than eigenvalue decomposition.
If n = 1 a vector of the same length as mu, otherwise an n by length(mu) matrix with one sample in each row.
runwithout helper function
runwithout(fit, year, fleet, ...) ## S3 method for class 'sam' runwithout(fit, year = NULL, fleet = NULL, map = fit$obj$env$map, ...)runwithout(fit, year, fleet, ...) ## S3 method for class 'sam' runwithout(fit, year = NULL, fleet = NULL, map = fit$obj$env$map, ...)
fit |
a fitted model object as returned from sam.fit |
year |
a vector of years to be excluded. When both fleet and year are supplied they need to be of same length, as only the pairs are excluded |
fleet |
a vector of fleets to be excluded. When both fleet and year are supplied they need to be of same length, as only the pairs are excluded |
... |
extra arguments to sam.fit |
map |
map to use |
...
Fit SAM model
sam.fit( data, conf, parameters, newtonsteps = 3, rm.unidentified = FALSE, run = TRUE, lower = getLowerBounds(parameters, conf), upper = getUpperBounds(parameters, conf), sim.condRE = TRUE, ignore.parm.uncertainty = FALSE, rel.tol = 1e-10, penalizeSpline = FALSE, fullDerived = FALSE, ... )sam.fit( data, conf, parameters, newtonsteps = 3, rm.unidentified = FALSE, run = TRUE, lower = getLowerBounds(parameters, conf), upper = getUpperBounds(parameters, conf), sim.condRE = TRUE, ignore.parm.uncertainty = FALSE, rel.tol = 1e-10, penalizeSpline = FALSE, fullDerived = FALSE, ... )
data |
data for the sam model as returned from the setup.sam.data function |
conf |
model configuration which can be set up using the |
parameters |
initial values which can be set up using the |
newtonsteps |
optional extra true newton steps |
rm.unidentified |
option to eliminate unidentified model parameters based on gradient in initial value (somewhat experimental) |
run |
if FALSE return AD object without running the optimization |
lower |
named list with lower bounds for optimization (only met before extra newton steps) |
upper |
named list with upper bounds for optimization (only met before extra newton steps) |
sim.condRE |
logical with default |
ignore.parm.uncertainty |
option passed to TMB:::sdreport reported uncertainties will not include fixed effect parameter uncertainties |
rel.tol |
option passed to stats:::nlminb sets the convergence criteria |
penalizeSpline |
Add penalization to spline recruitment? |
fullDerived |
Report all derived values? |
... |
extra arguments to MakeADFun |
The model configuration object conf is a list of different objects defining different parts of the model. The different elements of the list are:
A single integer defining the the lowest age class in the assessment.
A single integer defining the the highest age class in the assessment.
Is last age group considered a plus group (1 yes, or 0 no).
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes. The matrix describes the coupling of the fishing mortality states (the first rows are the catch fleet without effort). '-1' is used for entries where no fishing mortality applies (e.g. age groups in survey fleets, or unobserved age groups). For the valid entries consecutive integers starting at zero must be used, because they are used as indices in the corresponding state vector. If the same number is used for two fleet-age combinations, then the fishing mortality for those are assumed equal (linked to the same state).
An integer vector to specify the correlation structure of log-scale of fishing mortality increments (0 independent, 1 compound symmetry, or 2 AR(1)). The length of the vector is equal to the number of catch fleets without effort information.
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes. The matrix describes the coupling of survey catchability parameters (so only used for survey fleets). '-1' is used for entries where catchability should not be specified (e.g. fleet - age groups combinations where fishing mortality is specified above, or unobserved fleet - age group combinations). For the valid entries consecutive integers starting at zero must be used, because they are used as indices in the corresponding parameter vector. If the same number is used for two age classes, then the catchability for those age classes are assumed equal (linked to the same parameter).
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes. The matrix describes the coupling of density dependent catchability power parameters. This can only be applied to fleets - age combinations where a catchability is defined. '-1' is used for entries where this cannot be applied (e.g. fleet - age groups combinations where fishing mortality is specified above, or unobserved fleet - age group combinations). '-1' is also used to specify that density dependent catchability power parameters is turned off (the most common setup). For entries where density dependent catchability power parameter is to be estimates entries consecutive integers starting at zero must be used. If the same number is used for two age classes, then the density dependent catchability power parameter for those age classes are assumed equal (linked to the same parameter).
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes. The matrix describes the coupling of variance parameters for the different states in the log-scale fishing mortality random walk process. '-1' should be used for entries where no fishing mortality state is defined in keyLogFsta above. For the valid entries consecutive integers starting at zero must be used, because they are used as indices in the corresponding parameter vector. If the same number is used for two age classes, then the catchability for those age classes are assumed equal (linked to the same parameter). ((a curiosity of this setup is that it is possible to set different variance parameter indices for F-states that are coupled in keyLogFsta. This is ignored and the index corresponding to the lowest F-state number is used)).
A vector of integers. The length of the vector is equal to the number of age classes. The vector describes the coupling of variance parameters for the log(N)-process. Consecutive integers starting at zero must be used, because they are used as indices in the corresponding parameter vector. If the same number is used for two age classes, then the catchability for those age classes are assumed equal. A typical setup is to use a unique index for the first age group, because that corresponds to the variance in the (stock-)recruitment, which is often not similar to the variance in the survival process from year to year.
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes. The matrix describes the coupling of observation variance parameters. '-1' should be used for entries where no observations are available. For the valid entries consecutive integers starting at zero must be used, because they are used as indices in the corresponding parameter vector. If the same number is used for two age classes, then the observation variance for those age classes are assumed equal (linked to the same parameter).
A factor specifying the covariance structure used across ages for each fleet. The length of the factor is equal to the number of fleets. The possible options are: ("ID" independent, "AR" AR(1), or "US" for unstructured).
A matrix of integers. The number of rows is equal to the number of fleets and the number of columns is equal to the number of age classes _minus_ _one_. The matrix describes the coupling AR correlations between age classes, and hence is only meaningful for fleets where the "AR" observation correlation structure is chosen. '-1' should be used for entries where no observations are available. Notice that the matrix has one column less than the number of age classes, which is because the correlation between age classes is described. Consecutive integers starting at zero must be used. If the same number is used for a given fleet it means that a normal AR(1) structure is used. If different numbers are used for a fleet it means that the correlation parameter changes where the numbers differ. If the "AR" structure is specified above, then the corresponding row in this matrix must have valid non-negative entries.
A single integer to specify the stock recruitment connection to use:
| Code | Model |
| 0 | plain random walk on log recruitment |
| 1 | Ricker |
| 2 | Beverton-Holt |
| 3 | piece-wise constant |
| 61 | segmented regression (hockey stick) |
| 62 | AR(1) on log-recruitment |
| 63 | bent hyperbola (smooth hockey stick) |
| 64 | power function with degree < 1 |
| 65 | power function with degree > 1 |
| 66 | Shepherd |
| 67 | Deriso/Hassel |
| 68 | Saila-Lorda |
| 69 | sigmoidal Beverton-Holt |
| 90 | CMP spline (Non-increasing spline on log(R/S)) |
| 91 | Integrated spline on log(R/S) |
| 92 | Spline on log(R/S) |
See Albertsen & Trijoulet (2020) for details.
A vector of years to determine piece-wise constant recruitment periods for recruitment model 3. A vector of knot placements on log-SSB for spline recruitment models (90, 91, 92).
Determines the smoothness of recruitment model 63. The smoothness is estimated if set to NA.
A single integer specifying the number of years where catch scaling is to be estimated (most often 0, as this is a somewhat exotic option).
A vector of the years where catch scaling is applied (length should match noScaledYears) (most often empty, as this is a somewhat exotic option).
A matrix of integers specifying the couplings of scale parameters (nrow = noScaledYears, ncols = no ages) (most often empty, as this is a somewhat exotic option).
An integer vector of length 2 specifying lowest and highest age included in Fbar (average fishing mortality summary).
A vector of integers with length equal to the number of fleets. '-1' should be used for entries where the corresponding fleet is not a mass index. A the corresponding fleet is a mass index, then three options are available (0 SSB index, 1 catch index, and 2 FSB index).
A factor specifying the type of likelihood to use for each fleet. The length of the factor is equal to the number of fleets. The possible options are: ("LN" for log-normal and "ALN" Additive logistic normal).
A single integer. If weight attribute is supplied for observations this option defines how it is treated (0 as relative weight, 1 as a fixed variance = weight).
an object of class sam
Albertsen, C. M. and Trijoulet, V. (2020) Model-based estimates of reference points in an age-based state-space stock assessment model. Fisheries Research, 230, 105618. doi:10.1016/j.fishres.2020.105618
data(nscodData) data(nscodConf) data(nscodParameters) fit <- sam.fit(nscodData, nscodConf, nscodParameters, silent = TRUE)data(nscodData) data(nscodConf) data(nscodParameters) fit <- sam.fit(nscodData, nscodConf, nscodParameters, silent = TRUE)
Saves a model configuration list to a file
saveConf(x, file = "", overwrite = FALSE)saveConf(x, file = "", overwrite = FALSE)
x |
sam configuration list as returned from defcon or loadConf |
file |
the file to save the configuration to |
overwrite |
logical if an existing file should be overwritten (FALSE by default) |
function useful for saving a model configuration. A saved configuration can be read back in via the loadConf function
Plots the sd of the log observations as estimated in SAM in increasing order
sdplot(fit, barcol = NULL, marg = NULL, ylim = NULL, ...) ## S3 method for class 'sam' sdplot(fit, barcol = NULL, marg = NULL, ylim = NULL, show.rel.w = FALSE, ...)sdplot(fit, barcol = NULL, marg = NULL, ylim = NULL, ...) ## S3 method for class 'sam' sdplot(fit, barcol = NULL, marg = NULL, ylim = NULL, show.rel.w = FALSE, ...)
fit |
the object returned from sam.fit |
barcol |
color for each fleet and age |
marg |
margin for plot (mar in par()) |
ylim |
bounds for y-axis |
... |
extra arguments to plot |
show.rel.w |
plots the relative weight of each observation rather than the sd, estimated as (1/sd^2)/max(1/sd^2) |
small helper function
setS(x)setS(x)
x |
vector if indices |
internal
small helper function
setSeq(min, max)setSeq(min, max)
min |
from |
max |
to |
internal
Combine the data sources to SAM readable object
setup.sam.data( fleets = NULL, surveys = NULL, residual.fleets = NULL, prop.mature = NULL, stock.mean.weight = NULL, catch.mean.weight = NULL, dis.mean.weight = NULL, land.mean.weight = NULL, natural.mortality = NULL, prop.f = NULL, prop.m = NULL, land.frac = NULL, recapture = NULL, sum.residual.fleets = NULL, keep.all.ages = FALSE )setup.sam.data( fleets = NULL, surveys = NULL, residual.fleets = NULL, prop.mature = NULL, stock.mean.weight = NULL, catch.mean.weight = NULL, dis.mean.weight = NULL, land.mean.weight = NULL, natural.mortality = NULL, prop.f = NULL, prop.m = NULL, land.frac = NULL, recapture = NULL, sum.residual.fleets = NULL, keep.all.ages = FALSE )
fleets |
comm fleets vith effort (currently unimplemented) |
surveys |
surveys |
residual.fleets |
fleet, or list of fleets without effort information |
prop.mature |
pm |
stock.mean.weight |
sw |
catch.mean.weight |
cw |
dis.mean.weight |
dw |
land.mean.weight |
lw |
natural.mortality |
nm |
prop.f |
... |
prop.m |
... |
land.frac |
... |
recapture |
... |
sum.residual.fleets |
... |
keep.all.ages |
... |
...
Simulate data from fitted model and re-estimate from each run
simstudy(fit, nsim, ncores = detectCores())simstudy(fit, nsim, ncores = detectCores())
fit |
a fitted model returned from sam.fit |
nsim |
number of simulations |
ncores |
number of cores to be used |
Simulate from a sam object
## S3 method for class 'sam' simulate( object, nsim = 1, seed = NULL, full.data = TRUE, keep.process = FALSE, ... )## S3 method for class 'sam' simulate( object, nsim = 1, seed = NULL, full.data = TRUE, keep.process = FALSE, ... )
object |
sam fitted object as returned from the |
nsim |
number of response lists to simulate. Defaults to 1. |
seed |
random number seed |
full.data |
logical, should each inner list contain a full list of data. Defaults to TRUE |
keep.process |
Keep logN and logF processes when full.data = TRUE? |
... |
extra arguments |
simulates data sets from the model fitted and conditioned on the random effects estimated
returns a list of lists. The outer list has length nsim. Each inner list contains simulated values of logF, logN, and obs with dimensions equal to those parameters.
SAM SPR plot
sprplot(fit, ...) ## Default S3 method: sprplot(fit, ...) ## S3 method for class 'samforecast' sprplot(fit, ...) ## S3 method for class 'hcr' sprplot(fit, ...)sprplot(fit, ...) ## Default S3 method: sprplot(fit, ...) ## S3 method for class 'samforecast' sprplot(fit, ...) ## S3 method for class 'hcr' sprplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of deterministic equilibrium spawners per recruit assuming biological parameters and selectivity for that year remains unchanged in the future.
SPR table
sprtable(fit, ...) ## Default S3 method: sprtable(fit, ...)sprtable(fit, ...) ## Default S3 method: sprtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Plots the stock recruitment
srplot(fit, ...) ## S3 method for class 'sam' srplot( fit, textcol = "red", years = TRUE, linetype = "l", linecol = "black", polycol = do.call("rgb", c(as.list(col2rgb("black")[, 1]), list(alpha = 0.1))), polyborder = do.call("rgb", c(as.list(col2rgb("black")[, 1]), list(alpha = 0.3))), polylty = 3, polylwd = 1, xlim, ylim, add = FALSE, CIlevel = 0.95, addCurve = TRUE, ... )srplot(fit, ...) ## S3 method for class 'sam' srplot( fit, textcol = "red", years = TRUE, linetype = "l", linecol = "black", polycol = do.call("rgb", c(as.list(col2rgb("black")[, 1]), list(alpha = 0.1))), polyborder = do.call("rgb", c(as.list(col2rgb("black")[, 1]), list(alpha = 0.3))), polylty = 3, polylwd = 1, xlim, ylim, add = FALSE, CIlevel = 0.95, addCurve = TRUE, ... )
fit |
the object returned from sam.fit |
... |
extra arguments to plot |
textcol |
color of years on plot |
years |
the plotting symbols are the years |
linetype |
type for the plot (default line) |
linecol |
color of lines between points |
polycol |
Inner color of error ellipses |
polyborder |
Border color of error ellipses |
polylty |
Border line type of error ellipses |
polylwd |
Border line width of error ellipses |
xlim |
bounds for x-axis |
ylim |
bounds for y-axis |
add |
false if a new plot should be created |
CIlevel |
Confidence level for error ellipses on stock-recruitment pairs |
addCurve |
Call addRecruitmentCurve? |
SAM SSB plot
ssbplot(fit, ...) ## Default S3 method: ssbplot(fit, ...) ## S3 method for class 'samforecast' ssbplot(fit, ...) ## S3 method for class 'hcr' ssbplot(fit, ...)ssbplot(fit, ...) ## Default S3 method: ssbplot(fit, ...) ## S3 method for class 'samforecast' ssbplot(fit, ...) ## S3 method for class 'hcr' ssbplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of spawning stock biomass
SSB table
ssbtable(fit, ...) ## Default S3 method: ssbtable(fit, ...)ssbtable(fit, ...) ## Default S3 method: ssbtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
The function estimates reference points based on stochastic model forecasts. The following reference points are implemented:
F fixed to x, e.g., "F=0.3" (NOT IMPLEMENTED YET)
F in the last year of the assessment (NOT IMPLEMENTED YET)
F in the y years before the last in the assessment, e.g., "StatusQuo-1" (NOT IMPLEMENTED YET)
F that maximizes yield
Fs that gives 0.x*100% of MSY, e.g., "0.95MSY"
F that maximizes yield per recruit
F such that the derivative of yield per recruit is 0.x times the derivative at F=0, e.g., "0.1dYPR"
F such that spawners per recruit is 0.x times spawners per recruit at F=0, e.g., "0.35SPR"
F such that biomass is 0.x times the biomass at F=0, e.g., "0.2B0"
stochasticReferencepoints(fit, referencepoints, ...) ## S3 method for class 'sam' stochasticReferencepoints( fit, referencepoints, method = "Q0.5", catchType = "catch", nYears = 300, Frange = c(0, 2), nosim = 1000, aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), newton.control = list(), seed = .timeToSeed(), formula = ~ibc(F, 5), nosim_ci = 200, derivedSummarizer = median, ... )stochasticReferencepoints(fit, referencepoints, ...) ## S3 method for class 'sam' stochasticReferencepoints( fit, referencepoints, method = "Q0.5", catchType = "catch", nYears = 300, Frange = c(0, 2), nosim = 1000, aveYears = max(fit$data$years) + (-9:0), selYears = max(fit$data$years), newton.control = list(), seed = .timeToSeed(), formula = ~ibc(F, 5), nosim_ci = 200, derivedSummarizer = median, ... )
fit |
a sam fit |
referencepoints |
a character vector of reference points to estimate (see Details) |
... |
additional parameters that can be passed on |
method |
estimation method (See Details) |
catchType |
catch type: catch, landing, discard |
nYears |
Number of years to forecast |
Frange |
Range of F values to consider |
nosim |
Number of simulations for estimation |
aveYears |
Years to average over for biological input |
selYears |
Years to average over for selectivity |
newton.control |
List of control parameters for optimization |
seed |
Seed for simulations |
formula |
Formula to estimate optimization criteria as a function of F |
nosim_ci |
Number of simulations for bootstrap confidence intervals |
derivedSummarizer |
Function to summarize derived per-recruit values |
Reference points can be estimated using these methods:
Use least squares to estimate mean equilibrium values
Use quantile regression to estimate the 0.x quantile of equilibrium values
Identical to Q0.5
(NOT IMPLEMENTED YET)
To estimate median equilibrium yield, as required by ICES, the method "Q0.5" should be used.
sam reference point object
reference point object
## Not run: stochasticReferencepoints(fit, c("MSY","0.95MSY","Max","0.35SPR","0.1dYPR")) ## End(Not run)## Not run: stochasticReferencepoints(fit, c("MSY","0.95MSY","Max","0.35SPR","0.1dYPR")) ## End(Not run)
Deprecated and defunct functions
referencepointsFor referencepoints, use deterministicReferencepoints.
Summary of sam object
## S3 method for class 'sam' summary(object, ...)## S3 method for class 'sam' summary(object, ...)
object |
sam fitted object as returned from the |
... |
extra arguments |
summary table containing recruits, SSB, and Fbar
Table helper
tableit(fit, what, x = fit$data$years, trans = function(x) x, ...) ## S3 method for class 'sam' tableit(fit, what, x = fit$data$years, trans = function(x) x, ...) ## S3 method for class 'samforecast' tableit(fit, what, x = fit$data$years, trans = function(x) x, ...)tableit(fit, what, x = fit$data$years, trans = function(x) x, ...) ## S3 method for class 'sam' tableit(fit, what, x = fit$data$years, trans = function(x) x, ...) ## S3 method for class 'samforecast' tableit(fit, what, x = fit$data$years, trans = function(x) x, ...)
fit |
returned object from sam.fit |
what |
quoted name of what to extract |
x |
rownames of table |
trans |
function to be applied |
... |
extra arguments not currently used |
...
SAM TSB plot
tsbplot(fit, ...) ## Default S3 method: tsbplot(fit, ...)tsbplot(fit, ...) ## Default S3 method: tsbplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of total stock biomass
TSB table
tsbtable(fit, ...) ## Default S3 method: tsbtable(fit, ...)tsbtable(fit, ...) ## Default S3 method: tsbtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...
Write all data files from a list as created by 'setup.sam.data'
write.data.files(dat, dir = ".", writeToOne = TRUE, ...)write.data.files(dat, dir = ".", writeToOne = TRUE, ...)
dat |
A list as created by 'setup.sam.data' |
dir |
Directory where the files are written |
writeToOne |
Write multi fleet data to one file if data is equal for all fleets |
Write all data files from a list as created by 'setup.sam.data'
Write ICES/CEFAS data file from matrix
write.ices(x, fileout, writeToOne = TRUE, ...)write.ices(x, fileout, writeToOne = TRUE, ...)
x |
a matrix where rownames are taken as years and colnames are taken as ages |
fileout |
file name or connection |
writeToOne |
Write multi fleet data to one file if data is equal for all fleets |
... |
Arguments to be passed to write |
Takes the data and writes them in the ICES/CEFAS format. It is assumed that rows represent consecutive years and cols consecutive ages
Write surveys in ICES/CEFAS data file from a model object
write.surveys(fit, fileout, ...)write.surveys(fit, fileout, ...)
fit |
A fitted object as returned from sam.fit |
fileout |
file name or connection |
... |
Arguments to be passed to write |
Takes the survey data from the fitted object and writes them in the ICES/CEFAS format.
SAM years lost to fishing plot
yearslostplot(fit, cause, ...) ## Default S3 method: yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...) ## S3 method for class 'samforecast' yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...) ## S3 method for class 'hcr' yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...)yearslostplot(fit, cause, ...) ## Default S3 method: yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...) ## S3 method for class 'samforecast' yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...) ## S3 method for class 'hcr' yearslostplot(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...)
fit |
the object returned from sam.fit |
cause |
Fisning, Other, or LifeExpectancy |
... |
extra arguments transferred to plot including the following: |
Plot of years lost to fishing
Years Lost table
yearslosttable(fit, cause, ...) ## Default S3 method: yearslosttable(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...)yearslosttable(fit, cause, ...) ## Default S3 method: yearslosttable(fit, cause = c("Fishing", "Other", "LifeExpectancy"), ...)
fit |
... |
cause |
Fisning, Other, or LifeExpectancy |
... |
extra arguments not currently used |
...
Yield per recruit calculation
ypr( fit, Flimit = 2, Fdelta = 0.01, aveYears = min(15, length(fit$data$years)), sprProp = 0.35, ... ) ## S3 method for class 'sam' ypr( fit, Flimit = 2, Fdelta = 0.01, aveYears = min(15, length(fit$data$years)), sprProp = 0.35, ... )ypr( fit, Flimit = 2, Fdelta = 0.01, aveYears = min(15, length(fit$data$years)), sprProp = 0.35, ... ) ## S3 method for class 'sam' ypr( fit, Flimit = 2, Fdelta = 0.01, aveYears = min(15, length(fit$data$years)), sprProp = 0.35, ... )
fit |
the object returned from sam.fit |
Flimit |
Upper limit for Fbar |
Fdelta |
increments on the Fbar axis |
aveYears |
Number of years back to use when calculating averages (selection, weights, ...) |
sprProp |
Proportion of SPR at F=0, for example 0.35 if F0.35SPR |
... |
extra arguments not currently used |
SAM YPR plot
yprplot(fit, ...) ## Default S3 method: yprplot(fit, ...) ## S3 method for class 'samforecast' yprplot(fit, ...) ## S3 method for class 'hcr' yprplot(fit, ...)yprplot(fit, ...) ## Default S3 method: yprplot(fit, ...) ## S3 method for class 'samforecast' yprplot(fit, ...) ## S3 method for class 'hcr' yprplot(fit, ...)
fit |
the object returned from sam.fit |
... |
extra arguments transferred to plot including the following: |
Plot of deterministic equilibrium yield per recruit assuming biological parameters and selectivity for that year remains unchanged in the future.
YPR table
yprtable(fit, ...) ## Default S3 method: yprtable(fit, ...)yprtable(fit, ...) ## Default S3 method: yprtable(fit, ...)
fit |
... |
... |
extra arguments not currently used |
...