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
)

Arguments

mods

(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.

do.table

T/F, produce table of AIC and/or Mohn's rho? Default = TRUE.

do.plot

T/F, produce plots? Default = TRUE.

fdir

character, path to directory to save table and/or plots. Default = getwd().

table.opts

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.

plot.opts

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".

Value

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)

Details

plot.opts$which specifies which plots to make:

1

3-panel of SSB (spawning stock biomass), F (fully-selected fishing mortality), and Recruitment

2

CV (coefficient of variation) for SSB, F, and Recruitment

3

Fleet selectivity (by block, averaged across years)

4

Index selectivity (by block, averaged across years)

5

Selectivity tile (fleets + indices, useful for time-varying random effects)

6

M time series (natural mortality, can specify which age with plot.opts$M.age)

7

M tile (useful for time-varying random effects)

8

3-panel of F X% SPR, SSB at F_X%SPR, and yield at F_X%SPR

9

2-panel of relative status (SSB / SSB at F_X%SPR and F / F_X%SPR)

10

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.

Examples

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)
}