After fitting multiple WHAM (or ASAP) models, compare_wham_models
produces plots
and a table of AIC and Mohn's rho to aid model comparison.
compare_wham_models(
mods,
do.table = TRUE,
do.plot = TRUE,
fdir = getwd(),
table.opts = NULL,
plot.opts = NULL,
fname = NULL,
sort = NULL,
calc.rho = NULL,
calc.aic = NULL,
do.print = NULL
)
(named) list of fit WHAM/ASAP models. To read in ASAP model output, use read_asap3_fit
. If no names are given, 'm1', 'm2', ...
will be used.
T/F, produce table of AIC and/or Mohn's rho? Default = TRUE.
T/F, produce plots? Default = TRUE.
character, path to directory to save table and/or plots. Default = getwd()
.
list of options for AIC/rho table:
$fname
character, filename to save CSV results table (.csv will be appended). Default = 'model_comparison'
.
$sort
T/F, sort by AIC? Default = TRUE.
$calc.rho
T/F, calculate Mohn's rho? Retrospective analysis must have been run for all modes. Default = TRUE.
$calc.aic
T/F, calculate AIC? Default = TRUE.
$print
T/F, print table to console? Default = TRUE.
$save.csv
T/F, save table as a CSV file? Default = TRUE.
list of options for plots:
$out.type
character, either 'pdf'
or 'png'
(default = 'png'
because I am not sure system('pdftk')
will work across platforms.)
$ci
vector of T/F, length = 1 (applied to all models) or number of models
$years
vector, which years to plot? Default = all (model and projection years).
$which
vector, which plots to make? Default = all. See details.
$relative.to
character, name of "base" model to plot differences relative to.
$alpha
scalar, (1-alpha)% confidence intervals will be plotted. Default = 0.05 for 95% CI.
$ages.lab
vector, overwrite model age labels.
$kobe.yr
integer, which year to use in Kobe plot (relative status). Default = terminal model year.
$M.age
integer, which age to use in M time-series plot. Default = max(data$which_F_age)
(max age of F to use for full total F).
$return.ggplot
T/F, return a list of ggplot2 objects for later modification? Default = TRUE.
$kobe.prob
T/F, print probabilities for each model in each quadrant of Kobe plot? Default = TRUE.
$refpt
"XSPR" or "MSY", which reference point to use. Default = "XSPR".
a list with the following components:
daic
Vector of delta-AIC by model (if do.table=T
and table.opts$calc.aic=T
)
aic
Vector of AIC by model (if do.table=T
and table.opts$calc.aic=T
)
rho
Matrix of Mohn's rho by model (if do.table=T
and table.opts$calc.rho=T
)
best
Name of best model (lowest AIC) (if do.table=T
and table.opts$calc.aic=T
)
tab
Results table of AIC and Mohn's rho (if do.table=T
)
g
List of ggplot2 objects for later modification (if do.plot=T
and plot.opts$return.ggplot=T
)
plot.opts$which
specifies which plots to make:
3-panel of SSB (spawning stock biomass), F (fully-selected fishing mortality), and Recruitment
CV (coefficient of variation) for SSB, F, and Recruitment
Fleet selectivity (by block, averaged across years)
Index selectivity (by block, averaged across years)
Selectivity tile (fleets + indices, useful for time-varying random effects)
M time series (natural mortality, can specify which age with plot.opts$M.age)
M tile (useful for time-varying random effects)
3-panel of F X% SPR, SSB at F_X%SPR, and yield at F_X%SPR
2-panel of relative status (SSB / SSB at F_X%SPR and F / F_X%SPR)
Kobe status (relative SSB vs. relative F)
If plot.opts$return.ggplot = TRUE
, a list g
is returned holding the above ggplot2 objects for later modification.
g[[i]]
holds the plot corresponding to i
above, e.g. g[[2]]
is the CV plot.
if (FALSE) {
base <- read_asap3_fit()
m1 <- fit_wham(input1)
m2 <- fit_wham(input2)
mods <- list(base=base, m1=m1, m2=m2)
res <- compare_wham_models(mods)
}